python3を使った、サンクトペテルブルグのパラドックスの一考察

サンクトペテルブルクのパラドックス

コインを投げて表が出るまで続けるゲームで、投げた回数を N とすると、2^(N-1) 円の賞金が手に入ります。期待値を計算すると、無限大に発散します。

W=Σk=1,∞(1/2^k・2^(k-1))=∞

これは直感に反するため、サンクトペテルブルクのパラドックスと呼ばれています。

まず、期待されるコイン投げ回数 S を求め、

次に期待される賞金金額 μ を計算してみます (方法 1)。…①

S=Σk=1,∞ (1/2^k)=2

S=2 なので、μ=2^(S-1)=2

したがって、答えは 2 です。

これは、意味論的には妥当な推論ですが、数学的には間違っています。1 回のゲームのコインを投げる回数が無限大になる可能性があるためです。

一般に、期待値については f(g(x))=g(f(x)) は成り立ちません。

ウィリアム フェラーもサンプリングでこれを計算しましたが、正解は期待値が無限大であるということです。

ゲームの数が有限であれば、期待値ははるかに小さい値に収束します。

ダニエル ベルヌーイはこのパラドックスを提示し、賞金の対数を取ることで効用という形で解決策を見つけましたが、賞金が増えるとお金の価値が下がると言うことは人間の主観を導入しているので、客観的な答えではありません。

方法 1 とシミュレーションを n 回繰り返したときの平均値を計算するプログラムを以下に示します。

ここでは、n=10000 です。

```ev0.py
#!/usr/bin/python3
import os
import binascii
import random
from getkey import getkey,keys

#Method 1
# Approximate value of S when ∞=10000
# S=Σk=1,10000 (1/2^k)
#
s=0
for i in range(10000):
s=s+1/2**i
print("Method 1:",2**(s-1))

print("Hit any key")
key=getkey()

# Calculate the average when the game is repeated 10000 times
N=10000
l=[]
r=[]
f=open("/dev/random",'rb') # Open /dev/random

# Create a list of results when simulating n times
for i in range(N):
k=0
while(True):
k=k+1
n=f.read(1)
r1=binascii.hexlify(n)
r2=int(r1,16)
if r2%2: break
r+=[k]
l+=[2**(k-1)]

# Take the average
a=sum(l)/N
b=sum(r)/N

print("Simulation")
print("Number of times:",r)
print("Payoff:",l)
print("Average number of coin tosses: ",b)
print("μ(%d)=%f"%(N,2**(b-1)))
print("Average payoff: ",a)
```

プログラム実行の結果

```
$ ev.py
方法 1:2.0
任意のキーを押す
シミュレーション
回数: [3, 1, 2, 1, 1, 2, 2, 1, 1, 1, 1, 1, 1, 2, 1, 3, 2, 1, 2, 1, 1, 2, 1, 2, 1, 2, 2, 4, 3, 1, 3, 4, 5, 3, 1, 3, 3, 1, 1, 1, 1, 3, 1, 1, 1, 1, 2, 1, 1, 3, 1, 1, 1, 2, 1, 1, 3, 1, 1, 3, 4, 1, 1, 5, 5, 4, 1, 1, 2, 5, 1, 1, 1, 1, 2、1、2、1、2、3、2、3、2、2、1、2、1、2、1、1、2、3、1、1、1、3、1、10、、長いので省略、8、1、7、1、1、1]
利益額: [4.0、1.0、2.0、1.0、1.0、2.0、2.0、1.0、1.0、1.0、2.0、1.0、1.0、1.0、1.0、1.0、1.0、1.0、2.0、1.0、1.0、1.0、2.0、1.0、4.0、2.0、1.0、2.0、1.0、1.0、1.0、2.0、2.0、8.0、4.0、 1.0、4.0、8.0、16.0、4.0、1.0、4.0、4.0、1.0、1.0、1.0、1.0、1.0、4.0、1.0、1.0、1.0、1.0、2.0、1.0、1.0、4.0、1.0、 1.0、4.0、8.0、1.0、1.0、16.0、16.0、8.0、1.0、1.0、2.0、16.0、1.0、1.0、1.0、1.0、2.0、1.0、2.0、1.0、2.0、4.0、2.0、 4.0、2.0、 2.0、1.0、2.0、1.0、1.0、2.0、4.0、1.0、1.0、1.0、4.0、1.0、512.0、長いので省略、128.0、1.0、64.0、1.0、1.0、1.0]
コイントスの平均回数: 1.9905
μ(10000)=1.974132271
平均ペイオフ: 6.5987
$
```

10000回のシミュレーションを実行すると、コイントスの平均回数は1.9905なので、
μ(10000回のシミュレーション)＝2^(1.9905-1)≒1.974132271

これは平均ペイオフ 6.5987 とは大きな差があります。

r の平均: (r[0]+r[1]+r[2]+r[3]+r[4]+...+r[9999])/10000=1.9905

r の平均は (r1+r2+r2+...+rn)/n と表されます。

r の回数がすべて 2 の場合、(∑(i=1,n){2^(2-1)})n=(∑(i=1,n){2})n=2 となります。

lの平均は、la=(l[0]+l[1]+l[2]+...+l[9999])/10000=6.5987 //シミュレーション
=(∑(i=1,n){l[i]})/n=(∑(i=1,n){2^(r[i]-1)})/n。

極限では、①でl=μ、r=Sとすると、l=2^(r-1)となる。

​​コンピュータシミュレーションによる解

ゲームを10億回シミュレーションして得られた解は以下のとおり。

プログラム

```ev.py
#!/usr/bin/python3
import os
import binascii
import random
from getkey import getkey,keys

N=1000000000
# ゲームを N 回繰り返したときの平均を計算します
l=0
r=0

f=open("/dev/random",'rb') # /dev/random を開きます

# n 回シミュレーションしたときの結果のリストを作成します

for i in range(N):
k=0
while(True):
k=k+1
n=f.read(1)
r1=binascii.hexlify(n) # 16 進文字列に変換します
r2=int(r1,16) # 整数に変換します
if r2%2:
break
print(i,chr(13),end='',flush=True)
r+=k
l+=2**(k-1)

f.close() # /dev/random を閉じる

# 平均を取る
a=l/N
b=r/N

print("シミュレーション回数:",N)
print("コイン投げの平均回数: ",b)
print("平均ペイオフ: ",a)
```

実行結果

```
シミュレーション回数: 1000000000
コイン投げの平均回数: 2.000005233
平均ペイオフ: 19.908359876
```

したがって、「現実的に」期待値はせいぜい 20 円です。このギャンブルはしない方がよいでしょう。

我々は、数学の新しい理論で期待値を再正規化することでこの問題を解決するエゼキエル・ベルヌーイの登場を待っている。

演算をしたら情報が消えることも、この問題を解けない原因の一つでしょう。

シミュレーションがN回で終了する場合、平均報酬H(N)は確率的未知数を含むNの関数になる。

​​結局、無限に行うと、lim N→∞ H(N)＝∞となる。