摘要
mathe于2008年7月引用百度知道中一个抛硬币的概率问题
抛硬币100次,出现10次以上连续正面的概率是多少?
我们将它推广一下,抛硬币n次,其中出现过t次连续正面的概率是多少?
shshsh_0510发现这个问题
和广义t阶Fibonacci数列 有关系, 其中不出现t次连续正面的概率为 2nFn+2(t)。
或者说对于本题,就是要计算1−2100F102(10)。
最后我们得出广义t阶Fibonacci数列Fn(t)=⌊2t−(t+1)rtrt−1(1−rt)rtn⌉, 其中rt是方程 xt+1−2xt+1=0 的唯一一个大于1的实数根,⌊x⌉表示x四舍五入到最接近的整数后的结果。
具体内容
在mathe转载了内容 而shshsh_0510发现了wolfram网站上给出了问题和广义t阶Fibonacci数列的关系以后
mathe先通过矩阵方法验证 了这种关系。
广义Fibonacci数列
我们需要计算硬币抛n次过程中出现t次连续正面的概率。
我们对抛硬币过程中出现的状态进行分类:
其中s0表示还没有出现t次连续正面的情况,而且当前最后一次硬币出现状态为反面或初试状态。
其中s1表示还没有出现t次连续正面的情况,而且当前最后一次硬币出现状态为正面但是最后才连续1次正面
其中s2表示还没有出现t次连续正面的情况,而且当前最后一次硬币出现状态为正面而且最后正好连续2次正面
...
其中st−1表示还没有出现t次连续正面的情况,而且当前最后一次硬币出现状态为正面而且最后正好连续t-1次正面
其中st表示已经出现t次连续正面的情况。
所以开始抛硬币之前处在状态s0,然后抛硬币过程中,如果当前状态为sk,其中k<t,那么抛一次硬币有21的概率
转移到s0,还有21的概率转移到sk+1;而如果当前状态为st,那么状态永远不会再改变。
所以我们可以得到状态转移矩阵M为
⎣⎢⎢⎢⎢⎢⎢⎢⎡212100…0210210…0210021…0………………21000…210000…1⎦⎥⎥⎥⎥⎥⎥⎥⎤
其中我们知道对于初试状态,状态只可能为s0,对应概率分布为b=⎣⎢⎢⎢⎢⎢⎢⎢⎡100…00⎦⎥⎥⎥⎥⎥⎥⎥⎤
其中b是一个长度为t+1的列向量,b的第h个分量表示状态在sh−1的概率
而经过k步以后,状态分布的概率变化为Mkb
而我们需要的值为经过n步以后到达状态st的概率,也就是Mnb最后一个分量,或者我们记
u=(0,0,0,...,0,1)
那么所求的结果就是an=uMnb
为此我们可以先对矩阵M进行分析,由于矩阵2M更加简单,我们改为分析矩阵2M
而矩阵2M的特征多项式很简单,为
g(x)=(2−x)(1+x+...+xt−1−xt)
=x−1(xt+1−2xt+1)(2−x)
所以矩阵M的特征多项式为f(x)=g(2x)
所以我们知道如果计算结果为an,
我们知道数列一个特征根为1,而由于n趋向无穷时an→1,
可以假设an=1+u1(2x1)n+u2(2x2)n+...+ut(2xt)n
其中x1,x2,...,xt时方程xt+1−2xt+1除去1以后的n个根
我们可以假设bn=2n(1−an)
那么我们可以知道x−1xt+1−2xt+1是数列bn的特征多项式
也就是说,数列bn满足递推式bn+t+1=2bn+t−bn或者bn+t=bn+t−1+bn+t−2+...+bn,这个应该就是所谓的广义Fibonacci数列
数值计算结果
有了上面结论,我们就可以比较轻松数值计算题目中要求的结果了。
也就是为了对于给定的t,如果我们要计算对应通项公式,可以通过解方程xt+1−2xt+1=0的所有解,然后就可以有通项公式了。
谁来数值计算一下F102(10)和F1002(10)?
其中F0(10)=0,F1(10)=1,F2(10)=1,F3(10)=2,F4(10)=4,F5(10)=8,F6(10)=16,F7(10)=32,F8(10)=64,F9(10)=128,F10(10)=256,F11(10)=512,F12(10)=1023
而Fn+11(10)=2Fn+10(10)−Fn(10)
而抛100次,存在连续正面10次以上的概率为1−2100F102(10)。
无心人很快给出了数值结果
F102(10)=1211700015849788251502892752696
至于那个分式 = 0.9558627713595783
于是投100次硬币出现连续10次以上的概率为1-0.9558627713595783=0.0441372286404217。
后面mathe继续用后面分析的公式解进行数值计算给出了精度更高的结果
而对于t=10,使用计算器就可以算出rt=1.9990186327101011386634092391292
所以b(102)=1211700015849788251502892752685.8,b(1002)=6.5849587983847754109895686668823e+300
而投100次不出现连续10个正面的概率为2100b(102)=0.95586277135957831225932140641244,
而投1000次不出现连续10个正面的概率为21000b(1002)=0.61455024758751836408966755676318,
公式解分析
不过mathe显然不满足于对前面的数值解,此后花费了很多精力去分析其公式解。由于已经知道其特征多项式为x−1xt+1−2xt+1, 所以其公式解必然和这个多项式的根有关系。
他先做了如下的理论推导
我们现在看看能否给出t阶Fibonacci数列b(n)=Fn(t)的通项公式
我们采用初试值b(−t+2)=b(−t+3)=...=b(0)=0,b(1)=1,b(2)=1
为了方便起见,我们可以定义b(−t+1)=1,在这样定义以后b(1)的值也可以通过递推公式来递推得到了。
我们定义h(y)=yt+1−2y+1=k(y)(y−1),那么我们知道y−1h(y1)yt+1是数列的特征方程,所以h(y)的所有根是数列特征值的倒数。
对于复平面上的圆∣z∣=d其中32≤d<1,我们有
∣zt+1∣=dt+1≤2d−1≤∣2z−1∣
所以根据儒歇定理,我们知道h(y)在圆∣z∣≤d内部解的数目同方程−2z+1=0的数目相同。也就是在∣z∣<1内方程h(z)=0只有一个解。而如果∣z∣=1满足h(z)=0,那么容易得到z=1,而且根的重数为1。
由此我们得到h(z)有一个根在单位圆内部,t-1个根在单位圆外部,此外还有一个根1。
假设单位圆内部根为z1,单位圆外部的根为z2,z3,...,zt,于是我们可以假设数列b(n)的通项公式为
u1z1−n+u2z2−n+...+utzt−n
其中u1,u2,...,ut为待定系数。
采用类似随机游走中的概率问题的方法可以得到其中 us=k′(zs)1
而我们知道
k(z)=(z−z1)(z−z2)(z−z3)...(z−zt)
k′(zs)=∏h=s(zs−zh)
h′(zs)=(z−1)∏h=s(zs−zh)
所以k′(z−s)=zs−1h′(zs)
us=h′(zs)zs−1
而h′(zs)=(t+1)zst−2=(t+1)(2−zs−1)−2=2t−(t+1)zs−1
所以
us=2t−(t+1)zs−1zs−1
记rt=z1−1
而递推式中其他各式绝对值都远远小于,所以我们得到,对于n充分大
b(n)=⌊2t−(t+1)rtrt−1(1−rt)rtn⌉
而他怀疑这个表达式对所有的n都成立,不过这是还没有验证,但是如果仅仅数值计算,不需要精确值,由于∣r∣接近2,而其他项都递减,所以用上面公式精度已经非常高了。
继续对b(n)=2t−(t+1)rtrt−1(1−rt)rtn做误差分析可以得出
除了上面一项以外,实际上还有t-1项uszs−n,s=2,3,...,t.
其中对于n≥1,有∣uszs−n∣≤∣us∗zs∣≤∣2t−(t+1)∣1+∣zs−1∣<t−12
所以t-1项之和的绝对值小于t−12(t−1)=2
也就是说使用这个公式b(n)的误差必然小于2. 而则换成计算概率,第n项的计算误差小于2n−11,对于比较大的n仅使用第一项计算精度非常高,实际计算过程中可能r采用的精度都没有这么高。
后面他使用数值计算进行验算发现:
而对于公式b(n)=⌊u1rtn⌉,计算前面若干项我们可以发现猜想均成立:
n |
b(n) |
u1rtn |
1 |
1 |
0.50222005921650437539314900083299 |
2 |
1 |
1.0039472560945626042296717505261 |
3 |
2 |
2.0069092711912102894563100627129 |
4 |
4 |
4.0118490272698787619203081600742 |
5 |
8 |
8.0197609571323822998698197638132 |
6 |
16 |
16.031651583188626895454838284701 |
7 |
32 |
32.047570227910457178387829251568 |
8 |
64 |
64.063690018678506437470207581223 |
9 |
128 |
128.06451002750246161567818042825 |
10 |
256 |
256.00334173386700758850659443979 |
11 |
512 |
511.75545016205159804636690872342 |
12 |
1023 |
1023.0086802648866917173406684459 |
13 |
2045 |
2045.0134132736788208304516651412 |
14 |
4088 |
4088.0199172761664313714470202193 |
15 |
8172 |
8172.0279855250629839809737322779 |
然后他发现这个问题中特征方程xt+1−2xt+1是随机游走中的概率问题中进t退1的特殊情情况,最后可以得出我们使用上面的近似公式那么b(n)的误差为q(n)−q(n+1),其中q(n)是那个问题中5#中定义的当A在B后面n步时能够返回原点的概率.所以我们只要能够证明0<q(n)<21那么就证明了只使用一项然后四舍五入的方法是正确的。而且也可以解释上面例子中除了第一项以外,后面所有项的误差都是正的情况(因为q(n)应该是单调的)
而且显然对于固定的n,q(n)关于t是递减的,所以如果我们能够对于t=2的情况计算出q(1)≤21,那么对于所有这一类题目,都能够直接使用那个近似计算方法,然后四舍五入,而结果就会是精确的。
以此为基础,他又开始对上面猜测正确性的分析:
为了证明猜想公式的正确性,即b(n)=⌊(t+1)r−2tr−1(r−1)rn⌉
我们先证明
(t+1)r−2tr−1(r−1)<21
其中1<r<2而且rt+1−2rt+1=0,
首先我们容易得出(t+1)r−2t>0
这个是因为方程h(x)=xt+1−2xt+1中,显然h(t+12t)<0,h(2)>0,而且前面已经证明r是h(x)唯一一个绝对值大于1的解,所以t+12t<r<2
所以我们需要证明2(1−r−1)<(t+1)r−2t也就是(2−r)r<t+12
由于我们知道rt(2−r)=1,所以只要证明rt−1>2t+1也就是r>t−12t+1(其中t≥2)
而其中右边容易得出在t=2时有最大值1.5,而r>t+12t,所以在t≥4时都已经有不等式成立。
而对于t=3时,我们知道r3=1.83⋯>1.5,对于t=2,r2=1.618⋯>1.5
所以我们证明了(t+1)r−2tr−1(r−1)<21
而更加容易的,我们可以证明(t+1)r−2tr−1>21 (这个等价于r<2),
而又有前面的不等式得到(t+1)r−2tr−1<2r<1
现在我们定义序列e(n)=b(n)−2t−(t+1)rr−1(1−r)rn也就是这个误差,其中对于所有的n≥−t+2全部定义。
而其中b(n)在n≤0时定义为b(−t+2)=b(−t+3)=...=b(0)=0(这样所有这些项都可以用同一个递推公式推导了)
显然e(n)也满足递推公式e(n+t)=e(n+t−1)+e(n+t−2)+...+e(n),这个递推式对于所有的n≥−t+2都成立。
而且我们知道e(0)=b(0)−2t−(t+1)rr−1(1−r)=−2t−(t+1)rr−1(1−r),所以−21<e(0)<0
同样,我们知道−21<e(−1)<0,−21<e(−2)<0,…,−21<e(−t+2)<0
而e(1)=1−(t+1)r−2tr−1,所以0<e(1)<21
我们现在定义一个模型,有一个人在实数轴的整点上行走(开始坐标大于1),他每次分别以均等的概率随机后退t格或前进1格,知道到达坐标位置[−t+2,1]之间停止;如果他最后停留在坐标x,那么我们说他最终的得分为e(x),那么请问如果他的起始坐标为x,那么他最终得分的期望值是多少呢?(容易证明他最终会停止下来的概率为1)
非常有意思,对于任何一个位置n,他最终得分的期望值为e(n),由于每个最终得分都在−21和21之间,我们可以知道所有的e(n)也必然只能在−21和21之间。也就是这种方法可以证明上面的猜想。
这个是因为假设开始在位置n的得分为f(n),我们知道对于任意n≥2有f(n)=2f(n+1)+f(n−t)或者说f(n+1)=2f(n)−f(n−t),而且对于n≤1有边界条件f(n)=e(n)
根据f(n)的递推公式我们知道它对应的特征方程为xt+1−2xt+1,这个方程有一个根的绝对值大于1,还有一个根为1,其余t-1个根绝对值小于1,根据
随机游走中的概率问题中的结论,由于数列f(n)有界,那么绝对值大于1的根所对应的项的系数必然为0。也就是f(n)的通项公式可以由其余t个特征根决定,也就是前t项已经唯一决定了数列f(n),而由于e(n)也满足同样递推式,必然得出f(n)=e(n)
由此我们可以证明了前面猜想的正确性,即
b(n)=⌊(t+1)r−2tr−1(r−1)rn⌉
而抛硬币n次出现连续t次正面的概率为
p(n)=1−2n⌊(t+1)r−2tr−1rn+1⌉;r>1;rt+1−2rt+1=0
其它应用
广义Fibonacci数列还会以其它一些形式出现,比如倪举鹏提出过的寿命10年的兔子问题:
第一年有1对兔子,以后每一年活着的兔子都能繁殖出一倍兔子来(活着的兔子乘以2)。兔子寿命都是10年。第n年有多少对兔子.
实际上计算的就是9阶广义Fibonacci数列。