数字"1"在三角阵中的生存概率

KeyTo9_Fans

我们按照下面的方式构造一个由"0"、"1"组成的三角阵：

三角阵的前3行是：

    1
  1   1
0   1   0

接下来的每一行都比上一行多1个数，

每行的第1个数和最后1个数都是"0"，

其余的每个数都是按照下面的方式确定的：

如果上方的两个数都是"0"，那么这个数为"0"，

否则这个数有$p$的概率为"1"，$1-p$的概率为"0"。

令$N\rightarrow\infty$，

我们将【第$N$行至少有1个"1"的概率】称为【生存概率】。

问题$1$：当$p=75%$时，生存概率是多少？

问题$2$：要令生存概率为$0$，$p$的最大值是多少？

wayne

呵，概率型的元胞自动机啊
====
其实用0填充三角阵为一个N*N的网格的话， fans给的三角阵前两行可以去掉了，

wayne

先贴一个图，跟主题不完全相关，来热热场面：


看来30号规则，生生不息啊

wayne

有点意思，写了个程序，数值计算发现，第一问p=3/4的情况下， 貌似 N趋于无穷的时候， 三角阵中间区域的格点为1的概率都是2/3

1. p=.75;t=Table[f[n,k]=If[n<4,If[n==3&&k!=2,0,1],If[k==1||k==n,0,p(1-(1-f[n-1,k-1])(1-f[n-1,k]))]],{n,1000},{k,n}];t//Last

复制代码

==================
至于0到2/3之间的过渡到底 有多少项，仍然需要进一步挖掘。


不过，我可以先来分析一下，为啥是2/3
假设n,k都很大，则 第n行第k个坑 是1的概率 为 f(n,k)，于是存在关系：
$f(n,k) = p*(1-(1-f(n-1,k-1))(1-f(n-1,k)))$
当N很大，趋近于无穷时，第N行的各个f(n,k) 趋于一个稳定的值，设为x，则     x=p*(1-(1-x)^2)  ，解得 x= 2-1/p 。
正好解释了上面的程序跑的结果。


所以，反过来 自动解决了第二问。p最大为1/2

wayne

1. {4,-(-2+p) p}
   
2. {5,p^2 (4-3 p-p^2+p^3)}
   
3. {6,-p^3 (-8+8 p+p^2+2 p^3-7 p^4+3 p^5)}
   
4. {7,-p^4 (-16+20 p-p^2+4 p^3-p^4-26 p^5+21 p^6+3 p^7-6 p^8+p^9)}
   
5. {8,-p^5 (-32+48 p-10 p^2+9 p^3-5 p^4-2 p^5-100 p^6+143 p^7-55 p^8+42 p^9-80 p^10+49 p^11-5 p^12-4 p^13+p^14)}
   
6. {9,-p^6 (-64+112 p-40 p^2+23 p^3-17 p^4-5 p^5-14 p^6-328 p^7+648 p^8-322 p^9+22 p^10+126 p^11-559 p^12+644 p^13-175 p^14-91 p^15+4 p^16+62 p^17-32 p^18+5 p^19)}
   
7. {10,p^7 (128-256 p+128 p^2-64 p^3+51 p^4+8 p^5+29 p^6-47 p^7+1369 p^8-3356 p^9+2965 p^10-1925 p^11+1969 p^12-2125 p^13+4695 p^14-7889 p^15+6404 p^16-3593 p^17+3569 p^18-2861 p^19+143 p^20+1296 p^21-755 p^22+70 p^23+69 p^24-23 p^25+2 p^26)}
   
8. {11,-p^8 (-256+576 p-368 p^2+184 p^3-145 p^4-2 p^5-58 p^6+114 p^7-153 p^8-4573 p^9+13937 p^10-14671 p^11+7867 p^12-1753 p^13-7838 p^14+18205 p^15-38366 p^16+72345 p^17-77018 p^18+41147 p^19-4000 p^20-31503 p^21+57997 p^22-35700 p^23-12167 p^24+23861 p^25-4265 p^26-5709 p^27+1059 p^28+2617 p^29-1745 p^30+385 p^31+10 p^32-17 p^33+2 p^34)}
   
9. {12,p^9 (512-1280 p+992 p^2-528 p^3+402 p^4-49 p^5+117 p^6-273 p^7+347 p^8-750 p^9+20128 p^10-67707 p^11+92897 p^12-80154 p^13+69041 p^14-65115 p^15+99988 p^16-201794 p^17+407421 p^18-766966 p^19+1015304 p^20-883707 p^21+644523 p^22-594232 p^23+815567 p^24-1296237 p^25+1454424 p^26-860199 p^27+246183 p^28-172974 p^29+82579 p^30+285786 p^31-405006 p^32+137130 p^33+92974 p^34-93503 p^35+20715 p^36+6794 p^37-3160 p^38-744 p^39+708 p^40-167 p^41+14 p^42)}
   
10. {13,p^10 (1024-2816 p+2560 p^2-1488 p^3+1100 p^4-267 p^5+251 p^6-644 p^7+799 p^8-1678 p^9+3379 p^10+64725 p^11-271162 p^12+429833 p^13-380486 p^14+226164 p^15+18395 p^16-334500 p^17+765978 p^18-1620975 p^19+3401815 p^20-6704585 p^21+10052504 p^22-10303931 p^23+7530266 p^24-3546709 p^25-1942936 p^26+9849666 p^27-20911953 p^28+30315591 p^29-26709886 p^30+9978505 p^31+7125483 p^32-19117315 p^33+23240395 p^34-11591406 p^35-7821343 p^36+13149193 p^37-2892121 p^38-4613999 p^39+2477110 p^40+944577 p^41-701383 p^42-573968 p^43+661325 p^44-177479 p^45-54714 p^46+47159 p^47-9909 p^48-712 p^49+674 p^50-111 p^51+6 p^52)}
    
11. {14,p^11 (2048-6144 p+6400 p^2-4096 p^3+2984 p^4-1032 p^5+599 p^6-1502 p^7+1864 p^8-3751 p^9+7623 p^10-14309 p^11+300065 p^12-1278352 p^13+2379424 p^14-2715054 p^15+2551240 p^16-2319545 p^17+2620653 p^18-4476024 p^19+8702295 p^20-17450668 p^21+36202573 p^22-73052581 p^23+121944055 p^24-153237615 p^25+149696611 p^26-124715197 p^27+102325740 p^28-123116879 p^29+237851808 p^30-476010268 p^31+783674373 p^32-952920422 p^33+808632522 p^34-488449283 p^35+293955087 p^36-400423492 p^37+738078872 p^38-883483361 p^39+469108700 p^40+122035948 p^41-262875075 p^42+103736160 p^43-123219348 p^44+202840032 p^45-41966380 p^46-174199905 p^47+159823970 p^48-8395853 p^49-60866902 p^50+30613411 p^51+4177389 p^52-8015251 p^53+1245909 p^54+1217573 p^55-545015 p^56-35447 p^57+86139 p^58-26134 p^59+2630 p^60+270 p^61-87 p^62+6 p^63)}
    
12. {15,p^12 (4096-13312 p+15616 p^2-11008 p^3+8032 p^4-3468 p^5+1601 p^6-3488 p^7+4391 p^8-8394 p^9+17191 p^10-32251 p^11+62439 p^12+923874 p^13-4989144 p^14+10546139 p^15-12858106 p^16+10902609 p^17-5430054 p^18-3039023 p^19+13904360 p^20-32415871 p^21+70509213 p^22-149579106 p^23+318108553 p^24-657053538 p^25+1154067887 p^26-1549672787 p^27+1553666604 p^28-1156724625 p^29+544261096 p^30+100641903 p^31-927484333 p^32+2787375355 p^33-6857583300 p^34+13176165189 p^35-18662831116 p^36+18347890951 p^37-11086974108 p^38+1113352386 p^39+7972383321 p^40-17434936328 p^41+30322366508 p^42-42278630995 p^43+38889964361 p^44-14678697067 p^45-13328449150 p^46+27909111610 p^47-28628670288 p^48+15795330427 p^49+8833252853 p^50-24167995805 p^51+12726766381 p^52+7766262215 p^53-10795707399 p^54+108694893 p^55+4943736787 p^56-1235366086 p^57-1744743143 p^58+638896639 p^59+863605772 p^60-706634622 p^61+10530104 p^62+208312860 p^63-93169360 p^64-2651487 p^65+15012096 p^66-4381969 p^67-222905 p^68+359991 p^69-39093 p^70-20038 p^71+7434 p^72-1023 p^73+54 p^74)}
    
13. 
    

复制代码

针对第一问，代入p=3/4,得到

1. {{4,0.9375},{5,0.905273},{6,0.884537},{7,0.869776},{8,0.858622},{9,0.849854},{10,0.842762},{11,0.836902},{12,0.831978},{13,0.827784},{14,0.824172},{15,0.821033}}
   

复制代码

貌似 收敛到0.82

mathe

不同位置的数据不独立，结果不一定正确

mathe

比如p=1/2
那么第四行有四种可能，各自1/4概率并且独立没有问题
0 0 0 0
0 1 0 0
0 0 1 0
0 1 1 0
但是第五行有八种可能
0 0 0 0 0 :13/32
0 1 0 0 0 :3/32
0 0 1 0 0:5/32
0 0 0 1 0:3/32
0 1 1 0 0:3/32
0 1 0 1 0:1/32
0 0 1 1 0:3/32
0 1 1 1 0:1/32

于是我们看第二个位置，为1的概率是（3+3+1+1）/32=1/4，而如果我们指定第二格为0，那么第一格为1的概率是（3/32+1/32）/（13/32+3/32+3/32+1/32）=4/20=1/5
所以两者关系不是独立的

wayne

wayne 发表于 2013-8-1 00:57
次数如此之高的多项式是如何得到的

程序算出来的，由于每一种组合类型都不是等可能的，而每一种组合类型产生新的组合类型也没有通项公式，所以我只能写程序了。

程序的大致思路就是算每一种0,1组合情况的概率情况，然后根据规则产生新的组合情况，统计一下新的概率分布。

代码很乱，需要进一步优化的，不过，还是贴出来吧：

1. t={{{0,1,0},1}};
   
2. d=Table[{n+3,t=Table[{g[[1,1]],Factor[Total[g[[All,2]]]]},{g,GatherBy[Flatten[Table[Table[{{0,Apply[Sequence,k[[All,1]]],0},i[[2]](k[[All,2]]/.List->Times)},{k,Tuples[Table[If[i[[1,ii]]+i[[1,ii+1]]==0,{{0,1}},{{0,(1-p)},{1,p}}],{ii,Length[i[[1]]]-1}]]}],{i,t}],1],#[[1]]&]}];1-t[[1,2]]//Factor},{n,5}]
   

复制代码

gxqcn

1. @_@

复制代码

现在可以了。

顺带贴上 wayne 要贴的代码：

1. t={{{0,1,0},1}};
   
2. d=Table[{n+3,t=Table[{g[[1,1]],Factor@Total@g[[All,2]]},{g,GatherBy[Flatten[Table[Table[{{0,Sequence@@k[[All,1]],0},(k[[All,2]]/.List->Times)i[[2]]},{k,Tuples[Table[If[i[[1,ii]]+i[[1,ii+1]]==0,{{0,1}},{{0,(1-p)},{1,p}}],{ii,Length[i[[1]]]-1}]]}],{i,t}],1],#[[1]]&]}];
   
3. 1-t[[1,2]]//Factor},{n,5}]

复制代码

//==============================

关于这个问题（在【code】里面不能添加at字符，否则论坛的数据库报错），
我曾请教程序开发商，
答复居然是：请勿在后台点击用户表优化相关的设置。

mathe

wayne的这个数据是什么：
{{4,0.9375},{5,0.905273},{6,0.884537},{7,0.869776},{8,0.858622},{9,0.849854},{10,0.842762},{11,0.836902},{12,0.831978},{13,0.827784},{14,0.824172},{15,0.821033}}
是p=1/2时余下数目1的期望值？如果这样，说明p=1/2最终期望是趋向0。
实际上，考虑到对于每个孤零零的1，它对下一行的作用是将产生的期望是2p个1。但是两个相邻的1，由于1，1和1,0产生1的概率都相等，也就是相邻的1效用会降低。所以一行期望为h的行的效果不会比只有一个中有h的概率为1其余都为0的情况好。
也就是说，对于某个p,如果我们对初始状态(n=3)作用k次后，余下1的期望确定小于1，那么必然说明最终1必然会消亡。
由此，wayne在5#中给出的数据已经可以用来给出问题2不错的结论。我们只要计算出第n行余下数字1的期望（关于p的多项式），求这个多项式等于1时p的值就可以给出问题2的一个下界