St. Petersburg Paradox

圣彼得堡悖论: 掷硬币,若第一次掷出正面,你就赚2元.若第一次掷出反面,那就要再掷一次,若第二次掷的是正面,你便赚4元.若第二次掷出反面,那就要掷第三次,若第三次掷的是正面,你便赚8元,如此类推,即可能掷一次游戏便结束,也可能反复掷没完没了.问题是:你最多肯付多少钱参加这个游戏?

很容易求出获利的期望值是正无穷大,那么结论是无论出多少钱都应该玩这个游戏.这个结论很显然是有问题的.伯努利提出的这个问题,他给的答案就是效用函数(utility function),现在已经是经济学里的一个最基础的概念了.

简单讲一个例子,如果有一个赌博游戏,要求你压上全部的家当,1%的可能性赚1000倍,而99%的可能性全输掉.

如果单纯从期望值来考虑,值得玩.但是我们考虑输掉全部家当和家当翻1000倍对我们的影响,也就是效用.

假设效用函数是对数函数ln的话,输掉全部家当相当于$ln(0)=-\infty$,而赚1000倍也只是个正值X,$E=0.99 \times -\infty + 0.01 \times X = -\infty$. 这样算就不划算了~

同理圣彼得堡悖论也可以用对数效用函数来解决.

不过这不是我要讲的重点啦~

我的思路呢,把这个游戏当作一个可重复游戏来看待.这样子的话,我们不考虑玩了多少局,而是考虑投了多少次硬币.

假设每局的费用是c,我们一共投了N次硬币,那么这N次硬币有N/2次是反面,那么我们一共开局了(N+1)/2次,也就是说投N次硬币,需要付费的期望值$E_{pay}(N) = (N + 1)/2 \times c$,这一部分比较简单啦~

然后我们考虑投了N次硬币后,我们已经赚到的期望钱数.这时候比较复杂一些~,需要分情况讨论一下,首先投完N次硬币的时候,根据最后一次投币情况,存在两种状态:1,如果最后一次投币是正面,游戏还未结束.2,如果最后一次投币是反面,游戏结束.那么根据这两种情况,我们可以建立从第1次投币到第N次投币的动态转移方程,也就是递推公式啦~

每个阶段有两种状态, $E_{now}(N)$代表第N次投币后还没有结束的情况下,未结束的游戏当前赌注的期望值. $E_{win}(N)$代表第N次投币后已经赚到的钱数的期望.

然后就人肉一下状态转移方程啦~

解得:

这样我们就得到了$E_{win}(N)$了~ 所以我们只需要解一下$E_{win}(N) > E_{pay}(N) \times c$,就可以知道需要投多少次就开始期望获利了!假设每局费用是100元的话,那么解出N>199的时候获益期望为正~说明我们只要投币次数大于199次,那么就开始获益了~

那么这个结论的重点在于,如果投币次数少于199,获益期望是负的,199次投币的局数期望是100局.可以粗略的认为,前100局内你的整体获益期望都是负的,也就是要做好先输100局的心理准备~

一小段暴力验证的代码,brute枚举2^n种情况,来强行计算$E_{now}(N), E_{win}(N), E_{pay}(N)$, solve就是直接使用的公式~

def brute(n) :
    m = 1 << n
    a, b, c = 0, 0, 0
    for i in range(m) :
        t = 2
        for j in range(n):
            if t == 2 :
                c += 1
            if (i >> j) & 1 :
                b += t
                t = 2
            else :
                t *= 2
        if t != 2 :
            a += t
    return a / m, b / m, c / m

def solve(n) :
    return n + 1.0, (n - 1) * (n + 4) / 4 + 1, (n + 1) / 2

最后~还没写程序来模拟验证一下我的思路,等有时间再写吧~很忙很忙很忙~

再感叹一下,这真是最近写的为数不多的靠谱blog啊~ _(:з」∠)_

Part 2

好基友质疑这篇blog,就帮我写了个模拟!然后妥妥被好基友啪啪啪打脸了!哈哈哈哈~

所以我也写了个模拟!

import random

def once():
    b = 1
    while(random.SystemRandom().randint(0,1)) :
        b += 1
    return b

def run(n, c):
    i = 0
    e = 0
    while(i < n):
        e -= c
        b = once()
        e += (1 << b)
        i += b
    return e

def sample(m, n, c):
    e = 0
    for i in range(m):
        e += run(n, c)
    return e / m

用法很简单啦,sample(m, n, c) m是取样次数,n是投币多少次,c是每局的代价. 测试c=100肯定会被打脸,感觉一是取样次数达不到要求,二是伪随机的原因,我肯定是迷信了!

用n=60,c=20做一个说明吧~当c=20的时候,n>40就可以期望获利,接下来我们看一下取样不同而导致的结果

>>> sample(10, 60, 20)
-285.8
>>> sample(100, 60, 20)
-232.3
>>> sample(1000, 60, 20)
-116.096
>>> sample(10000, 60, 20)
-67.4866
>>> sample(100000, 60, 20)
76.618

这是一个很有趣的现象,取样次数太少,期望极为可能是负值.随着取样次数的增多,期望值稳定上升.

我们从另一个角度来考虑这个问题,对于n=60取样10次,我们可以认为其实是相当于投了600次硬币.而取样100000次n=60,相当于投币6000000次.

这时候我们可以把问题转化一下,认为对于c=20,虽然投币600次和投币6000000次的获益期望都是正的,但投币600获益的可能性很小,而投币6000000次获益的可能性很大.

这说明我们不仅要考虑投币多少次后期望获益为正,还要考虑投币多少次后有很大的概率是获益的.所以,接下来的工作就是来求当多少次投币之后,我们获益的概率大于给定的p.

这个问题等价于投币N次之后,已获利为正的概率是多少,公式很简单,考虑第n次投币如果游戏未结束,获利和n-1一样.f(n)代表第n次结束获利的概率 $P(N) = (P(N-1) + f(n))/2$.

这里的问题在于f(n)很难计算,顺推f(n)的话,需要知道投完前n-1次所有局面的精确获益才能进行统计,复杂度是2^n.所以精确计算不太现实,正在研究怎么近似计算.

一个比较简单的思路:n次里面只要连投$O(\log n)$次正面就不会赔钱,所以我们统计log(n)次正面出现的概率来近似赚钱的概率.

未完待续~