如果炸弹每秒钟爆炸概率提高一点，从数学期望上来看最有可能在哪一秒爆炸？ 第1页

@小侯飞氘 用物理模型解决这个问题的思路很赞。然而可惜他算的不是题目问的数学期望，另外虽然有些其他回答里用了数值模拟计算了接近的结果，但是数值模拟本身可能也有因为方法导致的误差，并且这个题目凑巧是可以有解析解给出近似结果的。

这里我直接用数学方法，因为针对这类问题已经有专门的数学模型和很多讨论了 -- 见生存分析 | Wikiwand

问题描述

先泛化一下问题。一个对象（炸弹）存在两个状态（为了简单先假设二元），状态A和B，其中A可以转化到B，但是反之不可以。我们经常关心的是由A状态转换到B状态的时间T。同时又可以假设完成观察对象状态一次需要单位时间 , 而由A->B 的转换是一个随机事件，因此T显然是一个随机变量。关于T的一系列问题就是生存分析（survival analysis)要解决的问题。

这类分析既可以针对没生命的系统，比如一台机器发生故障，一台手机死机，也可以用来研究一个自然人的寿命等等。这里为了简便，可以将A状态定义为生存，B状态定义为死亡。

接下来用数学语言描述：

定义T的“生存”分布函数 ，表示从0时刻开始对象在 期间保持A状态的概率为 ，类似可以得到对象在 期间保持A状态的概率为 ，由概率定义可知

同样可以定义T的概率密度函数 。

一般来说，我们最关心从A状态转换到B状态的瞬间，假设现在为 ,在这一瞬间，对象切换状态的概率是 , 又因为该事件发生的前提条件是对象“幸存”到了时间x，因此显然该事件发生的概率是条件概率:

当 该概率为：

一般称该函数 为危险函数Hazard function,其意义为对象在“生存”到x时刻后在x+dx时刻突然“死亡”得概率。利用这个函数我们可以得到几个有趣得结论。

T得期望值 利用分部积分，可得

接下来有了理论武器，我们回到这个炸弹问题本身。

问题可以分成两问，既可以算概率密度最大的点，也可以看时间T的数学期望。

如果按照提问的字面意义来理解的话，其实大多数回答计算的概率密度最大点并不是炸弹爆炸时间T的数学期望，所以说很可惜大多数答案都是错误的，但是不妨我们依旧来研究一下概率密度函数

先来看密度函数最大值

按照题目的描述，显然 就是炸弹在x到x+1秒的时间段内爆炸的概率（假设一秒已经足够短），得到 ，根据前述公式，可以得到

同样可以带入公式得到 , 令

解方程得 和其他答主的解一致。

以上解答了炸弹在哪个时刻，爆炸的概率密度最大。但并不是提问者问的爆炸期望时间

现在来看爆炸时间T的数学期望

已知其公式为

带入已经算出的 和问题规定的时间上限100，得到

，显然这是一个超越函数，属于带参数的高斯函数，或者说是误差函数。当上下限为有限实数时这个函数本身没有解析解，但是不妨碍我们得到数值解，一下列出三个解法：

方法一，用数值积分

知乎竟然不支持R语言。。。吐槽一下

       myfn <- function(x, k) exp(-k * x^2) integrate(myfn, lower=0, upper=100,k=1/200) 12.53314 with absolute error < 0.00093      

近似结果为12.533，kid's stuff，不解释

方法二，近似积分

        x <- seq(1,150,0.5)  y <- myfn(x,0.005)  data <- data.frame(x,y)  library(plotly)  p <- plot_ly(data, x = ~x, y = ~y, type = 'scatter', mode = 'lines')  p     

可以看到这个函数衰减很快，在0到100区间内的积分数值上和0至正无穷的差可以忽略不计。因此可以用 代替0到100区间上的积分。而那个积分是有解析解的：

带入 约等于12.533，和上一个方法的结果一致

方法三、用离散期望间接求解

假设函数 内服从均一分布（其实也是题目隐含的假设），可以得到

由定义可知 ,而第二项可以简单求和求解

       sum(myfn(1:100,0.005)) [1] 12.03314     

加上0.5，正好和其他两种方法得到的结果一致，voila!齐活

ps.想不到最后码了那么多字，其实也就是很简单的小题目，如果有哪里不小心用错符号或者打错字了欢迎指正。

ps又ps，有兴趣的同学可以算一下这个 也就是当炸弹在x时刻前都还未爆炸，那它在接下哪个时刻估计会爆炸（数学期望）？

另外大家如果有兴趣，我可以再写一点用这个理论分析人口数据的问题。

反正正确答案都已经被给出了，我随便凑个热闹说点相关的吧。

先看直观一点的情况。

情况1

如果0s爆炸概率是0，50s爆炸概率是0.5，100s爆炸概率是1。

只在三个时间点观察炸弹有没有爆炸，一共进行1000次实验。

0s时炸弹存活1000次，

50s时炸弹存活500次，

100s时炸弹存活0次。

情况2

如果0s爆炸概率是0，25s爆炸概率0.25，50s爆炸概率是0.5，75s爆炸概率0.75，100s爆炸概率是1。

只在这五个时间点观察炸弹情况：

0s时1000个存活

25s时750个存活

50s时375个存活

75s时93个存活

100s时0个存活

可以看到在第50s时存活的炸弹数两种情况不一致。

若如题设，将100s分割为100份，每过1s，爆炸概率提高1%。

则炸弹在第N-1秒存活而在第N秒之前爆炸的概率为：

如果将100s分割为1000份，每过0.1秒，爆炸概率提高0.1%。

则炸弹在第N-1秒存活而在第N秒之前爆炸的概率为：

其中 表示存活到N-1秒的概率， 表示在第N-1秒到第N秒之间爆炸的概率。

这个概率显然和上者不同。

可以证明，对时间的取样约细，最概然的爆炸时刻越靠前。

对时间的取样数为m份时，炸弹在第N-1秒存活而在第N秒爆炸的概率为：

（好复杂的式子啊，我为什么要把它展开啊。。。好好写连乘不好吗。。。好在这式子好像没有错，真佩服我自己~）

这是一个随m取值变化的函数。

同样的1s，同一时间段，初始条件变化时，爆炸的概率也会随之变化。

为什么会这样呢？

因为假定的独立事件的数量发生了变化，10个发生概率为0.1的事件同时不发生的概率，和100个发生概率为0.01的事件同时不发生的概率是不相同的。

因此如题所述，每秒增加1%概率的情况，和每0.1秒增加0.1%的情况是完全不同的，计算机模拟的结果也显示，两种已知条件其实是完全不同的，不仅仅是已知条件的精确度的问题。

每秒增加1%爆炸概率，最概然点在第10秒处；而每0.1s增加0.1%爆炸概率，则最概然点会提前到第3.2秒；如果将100秒分成m份，则最概然点在 秒的那一秒。

实际操作的时候，大概还是讲分布的时候多一些，因为包含的信息量更大，也更直观；某一状态下发生某一事件的概率其实并不直观，有时候甚至会给人误导，尤其是这种一旦发生就没有然后的事情。

比如与HIV携带者无保护啪啪啪，如果啪100次感染的概率才会积累到100%，那么大多数被感染者真的可能就在前10次就中了，而此时的概率可能仅仅只是10%。

所以啊，要多学数学啊。

结论是：在第10秒爆炸概率最大，爆炸时间的期望是12.5331秒。

其实可以给出更为普适的答案：假设有很多个互不影响的炸弹，每个未炸炸弹在下一秒被引爆的概率随时间线性增长，m秒时达到1，那么炸弹爆炸概率（密度）最大的时间为 ，爆炸的时间期望是 。

这题在很多物理化学问题中都有实际应用：假设有某个化学反应，其反应速率正比于反应物剩余浓度C（相当于炸弹未爆炸的概率）和温度T（实际上随T指数增长，不过原理差不多）。那么在升温时，随时间的增长，一方面反应速率会先因温度的升高而升高；但另一方面，随着反应的进行，反应物浓度越来越低，逐渐限制反应进行，因此反应速率在某个温度达到峰值后，会开始随时间降低。

在我研究的氢脆/聚变堆材料领域当中，研究人员经常会测量升温时氢的脱附速率，通过确定热脱附谱（Thermal desorption spectroscopy）峰值的温度，来推算出反应的能量学信息。

典型的热脱附结果如下：

先来看离散情形吧，也就是把时间拆成一秒一秒来看。这种情形最简单，高中数学就能解决。

根据题设可知，某个炸弹在1~n-1秒都没爆炸的概率为

令炸弹在第n秒爆炸的概率为P(n)，有：

显然

该比值函数在研究范围内是个单调递减函数，所以当 时，P(n)取得最大值。对应的解为 。当m比较大时，可以近似认为 。

简单写了个程序（代码见答案底部），当m=100时，数值解出来的爆炸期望为12.21，概率随时间变化的函数长这样子，是不是和上面的热脱附谱长得很像？

再来看连续情形。假设炸弹在0~t秒内不爆炸的概率为C(t)，那么将炸弹在t~t+dt时间内爆炸的概率，记为 ，它等于

其实题中问的是单位时间内的爆炸概率，即概率密度，记为 ：

按照题中定义，它正比于时间t，以及未爆炸概率C(t)，

将这两个式子结合起来并积分，得到：

;

可以看出C(t)是个正态分布，P(t)是个瑞利分布（Rayleigh distribution）。

对P(t)求导，使其等于零，得爆炸概率密度的最大值

解得 ，与m较大时的离散解是一致的。

积分可得爆炸的时间期望：

当m=100时，时间期望E=12.5331，与离散情形稍有出入，应该是离散/连续模型中时间分割细致度不一样导致的。

相关收藏：认真写科普

上面数值计算的MATLAB代码为：

       clear;clc;close all;  days=100; P_explode(1:days)=0; P_explode(1)=1/days;  P_unexplode(1:days)=0; P_unexplode(1)=1-1/days;  E=P_explode(1);  for n=2:1:days;     P_unexplode(n)=P_unexplode(n-1)*(1-n/days);     %第n秒不爆炸的概率     P_explode(n)=P_unexplode(n-1)*(n/days);     %第n秒爆炸的概率     E=E+P_explode(n)*n; end  E=E  set(0,'defaultfigurecolor','w'); figure set(gcf,'Position',[300 300 800 600]); set(gca,'Position',[.10 .10 .85 .85]); plot(1:days,P_explode,'b-','linewidth',1.5); lx=xlabel('Time (s)','FontSize',16); ly=ylabel('Explosion rate','FontSize',16);