プログラミングのための確率統計

Last-modified: 2011-07-09 (土) 15:45:23

Rubyで書かれたプログラムをPythonに変換
ついでに少し拡張したり

目次

today: ?
yesterday: ?
total: ?

更新履歴

2010/02/23
モンティホール問題を追加
（追加）ファイルを添付

2011/05/28
ビンゴ大会を追加

2011/0709
ランダムウォークを追加

モンティホール問題

有名な問題らしいです．

#><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><#
# -*-coding: utf-8 -*-
# monty.py
# Time-stamp: <2011-02-23 12:50:20 kato>
#
# プログラミングのための確率統計 [p. 25]
# モンティーホール問題のシミュレート
# ドアの数と試行回数は変更可能
#><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><#
import random

class Monty(object):
    def __init__(self):
        self.elem = {'ans':1, 'mine':1, 'yours':1}
        self.cont = {'keep':0.0, 'change':0.0}

    def define(self, doornum):
        dooriter = (random.randint(1, doornum)  for x in range(2))
        self.elem['ans']  = next(dooriter)
        self.elem['mine'] = next(dooriter)

    def select(self, doornum, command):
        self.define(doornum)
        if command is 'keep':
            self.cont['keep'] += 1 if self.elem['ans'] == self.elem['mine'] else 0

        elif command is 'change':
            while (self.elem['yours'] == self.elem['ans'] or
                   self.elem['yours'] == self.elem['mine']):
                self.elem['yours'] = random.randint(1, doornum)
            oldmine = self.elem['mine']
            while (self.elem['mine'] == self.elem['yours'] or
                   self.elem['mine'] == oldmine):
                self.elem['mine'] = random.randint(1, doornum)
            self.cont['change'] += 1 if self.elem['ans'] == self.elem['mine'] else 0

        else:
            pass

    def display(self, loop, num):
        correct   = num / loop * 100
        print(' O: %2.1f' %correct + ' %' + '   X: %2.1f' %(100 - correct) + ' %')
### class Monty end

if __name__ == '__main__':
    DOOR_NUM = 5
    LOOP_NUM = 5000
    mty = Monty()
    for i in range(LOOP_NUM):
        mty.select(DOOR_NUM, 'keep')
        mty.select(DOOR_NUM, 'change')

    print('Door num. = %s \nLoop num. = %s' %(DOOR_NUM, LOOP_NUM))
    print('I don\'t change my select.')
    mty.display(LOOP_NUM, mty.cont['keep'])
    print('I change my select.')
    mty.display(LOOP_NUM, mty.cont['change'])

サンプルコード

filemonty.py

ビンゴ大会

ビンゴのシミュレーションができます．
多分コメントアウト部分を読めばわかるかと．

## -*-coding: utf-8 -*-
## $Id: bingo.py 2011/05/23 17:11:36 kato Exp $
## $ python bingo.py
import random

class BingoGame(object):#==================================#
    def ask(self, bg):### ビンゴしたか確認する
        t = map(list, zip(*bg))### リストを転置する
        cross = [0, 0]
        for k, (i, j) in enumerate(zip(bg, t)):
            if(sum(i) == 0 or sum(j) == 0):### 横（縦）でビンゴしたか確認する
                return True
            cross[0] += i[k]
            cross[1] += i[4 - k]
        if(cross[0] == 0 or cross[1] == 0):### 斜めでビンゴしたか確認する
            return True
        return False

    def start(self, num):### ビンゴをスタートする
        card = []### ビンゴカードを用意する
        for i in range(1, 75, 15):
            card.append(random.sample(range(i, i+15), 5))
        card[2][2] = 0### 中央はFREE
        ball = random.sample(range(75), 75)### ビンゴの球を用意する
        for i in range(num):### 球を出した回数
            for j in range(5):### 球とカードを比べる
                if(ball[i] in card[j]):### 出た球とカードが一致した？
                    card[j][card[j].index(ball[i])] = 0### カードに穴（0）を開ける
                    if(self.ask(card)):### ビンゴしたか確認する
                        for k, e in enumerate(map(list, zip(*card))):
                            print('|%2d %2d %2d %2d %2d|'
                                  %(e[0], e[1], e[2], e[3], e[4]))
                        print('%2d times\n' %(i + 1))
                        return i + 1### ビンゴ！
        return False
#class

def main():
    bingoGame = BingoGame()
    count = []
    times = 10000### 試行回数
    num = 10### 球を出す回数を制限する（4 <= num <= 75）

    for i in range(times):
        check = bingoGame.start(num)
        if(check):
            count.append(float(check))

    print('%d回でビンゴする確率を求める' %(num))
    print('試行回数は%d回とする' %(times))

    if(len(count)):
        print('ビンゴ率 = %6.3f ％' %(float(len(count))*100/times))
        ave = sum(count)/float(len(count))
        print('期 待 値 = %6.3f 回' %(ave))
    else:
        print('誰もビンゴしませんでした')

if __name__ == '__main__':
    main()

サンプルコード

filebingo.py

ランダムウォーク

ランダムウォークのシミュレーションができます．
多分コメントアウト部分を読めばわかるかと．
一次元と二次元のランダムウォークです．

## -*-coding: utf-8 -*-
## $Id: randomWalk.py 2011/06/29 21:25:16 kato Exp $
## ランダムウォークの検証
import random

class RndWalk(object):# ====================================
    def make(self,num):### 一次元のランダムウォークを実行し，合計を計算する
        rnd = []
        while num > (len(rnd)-rnd.count(0)):### 1，-1の数を数える
            rnd.append(random.randint(-1,1))### ランダムに-1，0，1を出す
        return sum(rnd)### -1，1の合計を計算する

    def make2D(self,num):
        save = []
        for i in range(2):### 二次元なのでランダムウォークを2回実行
            rnd = []
            add = []
            while num > (len(rnd)-rnd.count(0)):### 1，-1の数を数える
                rnd.append(random.randint(-1,1))### ランダムに-1，0，1を出す
            [rnd.remove(0) for i in range(rnd.count(0))]### 0を消去
            ### 例) [1,-1,1,1,-1,1] -> [1,0,1,2,1,2]
            [add.append(sum(rnd[:i+1])) for i in range(len(rnd))]
            save.append(add)
        return map(list, zip(*save))### 行と列を入れ替えて返す
# class RndWalk(object) ====================================

def main():
    check = []
    num = 2**12### 試行回数の決定
    print('試行回数 %d 回\n計算開始'%(num))
    rw = RndWalk()
    for i in range(num):### RandomWalkを繰り返す
        if(i%100 == 0):
            print('%5.1f ％'%(i*100.0/num))
        check.append(rw.make(num*2))### num*2だけRandomWalkを実行する
    check.sort()### 降順に並び替え

    count = 0### RandomWalkの合計数のヒストグラムを計算する
    old = check[0]
    hist = []### ヒストグラムを格納する
    for i, elem in enumerate(check):
        if(count != 0 and elem != old):
            hist.append([old,count])
            count = 0### カウントリセット
        count += 1
        old = elem
    hist.append([old, count])### 最後のヒストグラムをここで格納する

    with open('walk.dat','w') as fp:### 一次元のRandomWalkの結果を保存
        [fp.writelines('%4d %2d\n' %(i[0],i[1])) for i in hist]

    with open('walk2D.dat','w') as fp2d:### 二次元のRandomWalkの結果を保存
        [fp2d.writelines('%4d %2d %2d\n' %(i, e[0],e[1]))
         for i,e in enumerate(rw.make2D(num*2))]

if __name__ == '__main__':
    main()

参考画像

plot.jpg

サンプルコード

filerandomWalk.py