请教大家一个复杂概率问题

673024707

.·.·. 发表于 2018-12-12 13:37
然而这是1%的唯一解释……
LZ并不是统计方向的
这意味着LZ可能会犯这种错误

您有兴趣就再研究一下吧，真正哪个对我也不敢保证

.·.·.

673024707 发表于 2018-12-12 15:21
您有兴趣就再研究一下吧，真正哪个对我也不敢保证

解方程x^5-0.6x^4-0.6*0.4x^3-0.6*0.4^2 x^2-0.6*0.4^3 x-0.6*0.4^4=0

结果0.5409

把中间步骤写出来然后那你就知道错在哪里了
因为如果你写完中间步骤，你可以依次验算
（1）前4次都未命中的概率
（2）第一次命中后四次未命中的概率
（3）第2,3,4次命中第5次未命中的概率
用容斥原理（1）+（2）-（3）能得到，五次中连续四次未命中的概率
感觉你的中间步骤里面，这两个概率至少会有一个对不上

kastin

10楼解答有些模糊，这里有个类似的递推方法：
设 `X_n` 为连续 `n` 次射击，其中最后出现连续0,1,2,3,4次未命中的概率为 `X_n`，易知 `X_5=(0.6,0.24,0.096,0.0384,0.01536)^\mathrm T`.
不难得到递推关系\[X_{n+1}=\begin{pmatrix}0.6 & 0.6 & 0.6 & 0.6 & 0.6\\0.4 & 0 & 0 & 0 & 0\\0 & 0.4 & 0 & 0 & 0\\0 & 0 & 0.4 & 0 & 0\\0 & 0 & 0 & 0.4 & 0\end{pmatrix}=AX_n\]那么 `X_{100}=A^{95}X_5=(0.326627465358228,0.1314849435936656,0.05292968970894266,0.02130701794528448,0.008577209052558868)^\mathrm T`，求和即得100次连续射击最后不出现连续4次以上未命中的情形的概率为0.5409263256586795.
注意，上面的递推关系并未约束连续 `n` 次射击中间是否出现“连续4次以上未命中”事件，它只是基于射击最后出现连续未命中情形的递推。也就是说，上述递推关系并不能保证10楼说的“且未出現連續4次以上不中”。因此怀疑10楼“且未出現連續4次以上不中”这个信息并未体现在递推关系式中（事实上A的特征方程与10楼特征方程是一样的）。
不过，如果中间未出现过“超过连续4次以上未命中”事件，那么根据递推关系可知，后面也不会出现这种情况。其中初始值恰好就是没有出现5次连续射击未命中情形，因此后续射击也不会出现连续5次及以上未命中事件。

如果上述解答有问题，到底问题出现在何处？

kastin

用蒙特卡洛模拟10万遍100次连续射击试验，得到概率为0.54
matlab代码

1. n = 10^5; % 模拟10^5遍
   
2. I = rand(100,n)<=0.6;  % 以0-1向量表示100次连续射击命中情况，命中概率为0.6
   
3. J = I==0;              % 标定未命中的序列
   
4. k = diff(I)~=0;
   
5. count = 0;
   
6. for iter = 1:n
   
7.     r = max(find(J(:,iter)&[k(:,iter);true])-find(J(:,iter)&[true;k(:,iter)])+1); % 最大连续不命中次数
   
8.     if r<5, count=count+1; end % 若不会出现连续4次以上没射中，频数加1
   
9. end
   
10. p = count/n;
    
11. disp(['100次连续射击中间不会出现连续4次以上没射中的概率为：' num2str(p)])

复制代码

.·.·.

kastin 发表于 2018-12-13 14:21
用蒙特卡洛模拟10万遍100次连续射击试验，得到概率为0.54
matlab代码

附一个表演性质的R代码

1. set.seed(1)
   
2. n=10^5
   
3. res=matrix(runif(n*100),n)>0.6
   
4. mismax=apply(res,1,function(x){
   
5.   max(diff(c(0,cumsum(c(x,0))[diff(c(0,x,0))<0])))
   
6. })
   
7. sum(mismax<=4)/n

复制代码

kastin

@.·.·.  
见仁见智了，没有完美的工具，只有最适合的工具。总的来说，失之东隅，收之桑榆。某个软件在一个地方便捷优秀，必然在另外一个地方变得不便。比如，有人抱怨Java语言很罗嗦，也有人抱怨C++太庞杂——但这不影响各种Java程序员和C++程序员用它们去做项目，开发软件、解决问题。如果真的java全面优秀于C++，那么为何还有人去学C++，还有那么多公司用C++？反之亦然。

开机速度，反正我这里不超过3秒，比起代码开发效率，这个时间没什么受不了的。与之相对的一个例子是，python没有这么多库预先加载，打开很快，也清爽，但每当写一个算法时都要import一大堆库，这时候多花的时间不止3秒吧？但是，有的人不在乎开机时间，在乎后面开发算法的高效率、不用为一些非算法层面细枝末节的东西所打断；有的人不在乎多写几行代码，在乎使用体验、在乎开源或者跨平台特性等等。所以，没有绝对的好坏，就看你的取舍和在乎什么而已。

说到论坛画风，我不知道你上的是哪个论坛，我经常上的是ilovematlab.cn这个论坛，这是matlab官方支持维护的，里面设计我没觉得挺好的。此外，我不懂论坛画风跟软件有啥关系，难道论坛画风影响软件的使用、或是影响编程效率？反倒是我觉得小木虫论坛界面够乱的，但还是有很多科研人员注册使用它，这没毛病吧。

关于收费问题……毕竟windows也收费，也没见所有人都去用linux上这个论坛，对吧？此外，如果在学校或科研单位，都有购买正版，也不用自己付费。所以，收费也不是主要问题。

R、SAS、SPSS都是统计学方面比较擅长的工具，所以数据统计、数据挖掘、统计可视化很方便。但如果涉及到贝叶斯统计方面，就需要更专业的OpenBUGS、WinBUGS了。可见，没有最好，只有最专业（而专精意味着丧失通用性）。而前面这些工具的功能，matlab都能很方便地实现。由于学新的工具（非指简单入门）需要不少时间，而在自己精通的通用工具中去解决问题，未必比专门学一个新工具再去解决问题要慢。

sheng_jianguo

这个概率问题很有意思。我以前分析过，现将我的计算方法总结如下：
注意我理解楼主问题是，“中间不会出现连续4次以上全射中的概率是多少？”也就是问在100次射击中，对整数j(1≤j≤96)，在j、j+1、j+2、j+3、j+4次射击中都不会全射中的概率（而不是第1-5次无一命中、第6-10次无一命中、……）。否则题目意义不大。

设连续射击k（1≤k≤100）次中，不会出现连续4次以上全射中的概率是A(k)，按条件不难得出：
A(1)=1,a(2)=1,A(3)=1,A(4)=1,A(5)=1-0.6^5
A(k)的递推关系为：当5＜k时
A(k)=0.4*A(k-1)+0.4*0.6*A(k-2)+0.4*0.6^2*A(k-3)+0.4*0.6^3*A(k-4)+0.4*0.6^4*A(k-4)
通过上式不断迭代求出需要的A(k)，但k很大时，若不用编程或不采用大数计算器，很难正确求得正确结果。
    下面想法找出一般计算公式求A(k)：
    上面方程的特征方程为x^5-0.4*x^4-0.4*0.6*x^3-0.4*0.6^2*x^2-0.4*0.6^3*x-0.4*0.6^4=0
    此特征方程的五个根（近似）为：
x(1)= 0.9623083691506949166212055929538172566539997867760030480467724403
x(2)= - 0.38711620508151738556331326335356-0.25991027362730615787739261049269i
x(3)= - 0.38711620508151738556331326335356+0.25991027362730615787739261049269i
x(4)= 0.105962020506169927252710466876648+ 0.48636637872922435807617136922825i
x(5)= 0.105962020506169927252710466876648- 0.48636637872922435807617136922825i
则 A(k)=C(1)*x(1)^k+C(2)*x(2)^k+C(3)*x(3)^k+C(4)*x(4)^k+C(5)*x(5)^k
将A(1)=1,a(2)=1,A(3)=1,A(4)=1,A(5)=1-0.6^5代入上式,得系数C(i)(i=1,2,3)的近似解
C(1)= 1.1263561382167118697347515675811
C(2)=- 0.12758025003881332929113666272279- 0.094771950328737674874730845239327i
C(3)=- 0.12758025003881332929113666272279+ 0.094771950328737674874730845239327i
C(2)= 0.064402180930457394423760878932225+ 0.15118506614264270936054641980588i
C(3)= 0.064402180930457394423760878932225- 0.15118506614264270936054641980588i
所以 ，所要求的概率
A(100)= C(1)*x(1)^100+C(2)*x(2)^100+C(3)*x(3)^100+C(4)*x(4)^100+C(5)*x(5)^100
≈0.024160282888513560545318451519189
对于此类问题，（当k大到一定时候）一般计算第一项就有足够精度了：
A(100) ≈C(1)*x(1)^100 ≈0.024160282888513560545318451519166

有人可能不相信所要求的概率会那么小，我们可以从以下两方面验证所要求的概率就是那么小：
1.如果将单次射击的命中率改为0.5，结果是多小呢？
Mathe在本论坛分析过，按他的计算方法，不难计算出
当p=0.5时，不会出现连续4次以上全射中的概率≈0.159621（p=0.6时应更小）。
参见：https://bbs.emath.ac.cn/thread-667-1-1.html
2.按14#kastin方法，用蒙特卡洛模拟10万遍100次连续射击试验，得到概率为0.024。
Kastin的代码中r计算和计算要求不同，我改了一下：
matlab代码

1. n = 10^5; % 模拟10^5遍
   
2. I = rand(100,n)<=0.6;  % 以0-1向量表示100次连续射击命中情况，命中概率为0.6
   
3. L(1,10^5)=false;
   
4. A=[L;I;L];
   
5. J = A==0;              % 标定未命中的序列
   
6. k = diff(A)~=0;
   
7. count = 0;
   
8. for iter = 1:n
   
9.     r = max(find(A(:,iter)&[k(:,iter);false])-find(J(:,iter)&[k(:,iter);false])); % 最大连续不命中次数
   
10.     if r<5, count=count+1; end % 若不会出现连续4次以上全射中，频数加1
    
11. end
    
12. p = count/n;
    
13. disp(['100次连续射击中间不会出现连续4次以上全射中的概率为：' num2str(p)])

复制代码

sheng_jianguo

仔细看了楼主的问题，我可能理解错了，不是“中间不会出现连续4次以上全射中的概率是多少？”而是“中间不会出现连续4次以上全射不中的概率是多少？”，如果是这样，只要将我解答中，0.6改为0.4，0.4改为0.6就可以了，以下是分析结果：
特征方程为x^5-0.6*x^4-0.6*0.4*x^3-0.6*0.4^2*x^2-0.6*0.4^3*x-0.6*0.4^4=0
此特征方程的五个根（近似）为：
x(1)= 0.9936573920360673876107579036146306291592568595162600008703777456
x(2)= -0.28323613467937871017316479044484+0.19312553907460499967725161950068i
x(3)= - 0.28323613467937871017316479044484-0.19312553907460499967725161950068i
x(4)= 0.086407438661345016367785838637527934+ 0.35223427711571730626911931828829i
x(5)= 0.086407438661345016367785838637527934- 0.35223427711571730626911931828829i
则 A(k)=C(1)*x(1)^k+C(2)*x(2)^k+C(3)*x(3)^k+C(4)*x(4)^k+C(5)*x(5)^k
将A(1)=1,a(2)=1,A(3)=1,A(4)=1,A(5)=1-0.4^5代入上式,得系数C(i)(i=1,2,3)的近似解
C(1)= 1.0220481302172501341780492330212
C(2)=- 0.045533280243984992535829105787221+ 0.038558922616239936849449802204675i
C(3)=- 0.045533280243984992535829105787221- 0.038558922616239936849449802204675i
C(2)= 0.034509215135359925446804489276606+ 0.0460337091708444654604231255451i
C(3)= 0.034509215135359925446804489276606- 0.0460337091708444654604231255451i
所以 ，所要求的概率
A(100)= C(1)*x(1)^100+C(2)*x(2)^100+C(3)*x(3)^100+C(4)*x(4)^100+C(5)*x(5)^100
≈C(1)*x(1)^100 ≈0.54092632565868016230376990688816

wayne

sheng_jianguo 发表于 2018-12-21 10:51
仔细看了楼主的问题，我可能理解错了，不是“中间不会出现连续4次以上全射中的概率是多少？”而是“中间不 ...

1. {n,p}={4, 3/5};
   
2. m=DiagonalMatrix[ConstantArray[1-p,n],-1,n+1];
   
3. m[[1]]=ConstantArray[p,n+1];m[[1]][[n+1]]=0;m[[n+1]][[n+1]]=1;
   
4. MatrixPower[m,100][[-1,1]]//N

复制代码

这道题跟你发出来的链接https://bbs.emath.ac.cn/thread-667-1-1.html其实是同一道题目。 得到的是 $0.800788$
修改一下代码，{n,p}={10,1/2}; 得到答案$0.0441372$,跟mathe一样。

mathe

连接中用到辅助函数
$1+x+...+x^{t-1}-x^t$
对应到本题应该相当于
$1+x+...+x^{t-1}-{1-p}/p x^t$
这个方程特征根应该没有那么好的性质了