5.4 天文学中的概率统计

蒙特卡罗思想

蒙特卡罗（Monte Carlo）模拟名字源自摩纳哥的赌城-蒙特卡洛，20 世纪 40 年代，由 John von Neumann，Stanislaw Ulam 和 Nicholas Metropolis 在曼哈顿计划中，为模拟中子扩散发展出的一种统计方法。 字如其名，蒙特卡洛方法本质上是跟赌博一样利用随机特性获得具有统计意义的结果。通过一个简单的例子，就能让你理解蒙卡的思想精髓：

如果 (x,y) 是独立地从 0 到 1 之间均匀分布抽样出的一系列的数对，那么这些随机的位置坐标 (x,y) 落在 1 为半径圆弧内的概率应该是：四分之一圆的面积/整个正方形的面积。

而因为 (x,y) 是 0 到 1 的均匀分布，所以这个概率当抽样足够多的时候就等于圆内的点数除以总共点数。这样一来，只要采样足够多，就可以得到无限趋近于 π 的值。

这种模拟方法称作舍选法（取舍法）。

蒙特卡洛方法还有很多其他的形式，如直接抽样法（反函数法）：http://hepg.sdu.edu.cn/zhangxueyao/Education/MonteCarlo/2013/ch03/MonteCarlo-ch03-03.pdf

中心极限定理和大数定律

统计就是能够从少量数据中得出结论的方法。例如我们只需要对 1000 个美国人进行电话调查，就能预测美国总统大选的得票数。有如此自信得到可靠结论源自：中心极限定理和大数定律。

中心极限定理：样本的平均值约等于总体的平均值。即不管总体是什么分布，任意一个样本平均值都会围绕在整体平均值周围，并且呈正态分布。

大数定律：如果统计数据足够大，那么事物出现的频率就能无限接近期望值。

小数定律：如果统计数据很少，那么事件就表现为各种极端情况，而这些情况都是偶然事件，跟它的期望值一点关系都没有。

正态分布

正态分布（Normal distribution）是一个非常常见的连续概率分布。正态分布在统计学上十分重要，经常用在自然和社会科学来代表一个不明的随机变量。

若随机变量$X$服从一个平均数μ、标准差为σ的正态分布，则记为：

$X \sim N\left(\mu, \sigma^2\right)$

则其概率密度函数为

$$f(x)=\frac{1}{\sigma \sqrt{2 \pi}} e^{-\frac{(x-\mu)^2}{2 \sigma^2}}$$

正态分布的数学期望值或期望值μ，决定分布的位置。其方差$\sigma^2$的平方根或标准差$\sigma$决定了分布的幅度。

期望、标准差、标准误和方差

• 期望（E）

期望是离散的随机变量的期望值，是试验中每次可能的结果乘以其结果概率的总和。

例如掷骰子（甲）的期望$E=\frac{1}{6}\times1+\frac{1}{6}\times2+\frac{1}{6}\times3+\frac{1}{6}\times4+\frac{1}{6}\times5+\frac{1}{6}\times6$

• 标准差

标准差，又称标准偏差、均方差（standard deviation, SD, σ）。用作描述一组测量值的离散程度。标准差的数学表达如下。N 为测量次数，xbar 为期望。

$$\sigma=\sqrt{\frac{1}{N} \sum_{i=1}^N\left(x_i-\bar{x}\right)^2}$$

例如，还有一个骰子（乙），它的六面均为 3.5，那么投 100 次骰子，甲乙的期望显然一样，但是甲的投掷结果显然更加离散，反映在坐标上就是甲的分布函数比乙宽。

标准差可以当作不确定性的一种测量。在物理中，做重复性测量，测量值集合的标准差代表这些测量的精确度。如果测量平均值与预测值相差太远（同时与标准差数值做比较），则认为测量值与预测值互相矛盾（即测量值都落在一定数值范围之外）。

• 标准误差

均值分布的标准差又叫做标准误差，或者标准误（standard error）

• 方差

方差（variance, Var）是标准差的平方。即方差的正平方根称为该随机变量的标准差。观测值通常是从真实世界的系统中测量的。如果给出系统的所有可能的观测，则它们算出的方差称为总体方差；然而，一般情况下我们只使用总体的一个子集（样本），由此计算出的方差称为样本方差。用样本计算出的方差可认为是对整个总体的方差的估计量。方差的表达式和展开式写作：

$$\begin{aligned} \operatorname{Var}(X) & =\mathrm{E}\left[(X-\mathrm{E}[X])^2\right] \\ & =\mathrm{E}\left[X^2-2 X \mathrm{E}[X]+\mathrm{E}[X]^2\right] \\ & =\mathrm{E}\left[X^2\right]-2 \mathrm{E}[X] \mathrm{E}[X]+\mathrm{E}[X]^2 \\ & =\mathrm{E}\left[X^2\right]-\mathrm{E}[X]^2 \end{aligned}$$

即方差也等于“方均减均方”。

卡方分布

所有数据形成的分布叫总体分布，在总体分布中随机抽取一批样本，这些样本形成了一个均值为 xbar1 的分布，再次抽取一批，得到均值为 xbar2 的分布，以此类推，多次取样后得到多个 xbar，且 xbar 形成一个新的分布，这个分布称为均值分布。

均值分布和总体分布的均值都是 μ，因为中心极限定理。当总体分布的标准差是 σ 时，那么均值分布的标准差（标准误）记作：$\sigma_{\bar{x}}=\frac{\sigma}{\sqrt{n}}$

现在在均值分布中随机抽取一个值 xbar1，这个值与均值之间存在一个距离 xbar1-μ，将这个距离与标准差作比得：

$$Z_1=\frac{\bar{x}_1-\mu}{\sigma_{\bar{x}}}=\frac{\bar{x}_1-\mu}{\frac{\sigma}{\sqrt{n}}}$$

且每次抽取都会得到一个Z值，那么足够多次抽样后，便得到了均值为 0 标准差为 1 的 Z 分布。

在这个分布中，继续抽样，每次都会得到一个 Z 值，由于 Z 分布是一个正态分布（有正有负），所以用 Z 的平方来表示与期望的距离，并将这些平方值按频率直方图的形式绘制便得到了自由度为 1 的卡方分布：

假设随机抽取时，一次抽取三个值且求得这三个值的平方和，这样多次抽样后，便得到了服从自由度为 3 的卡方分布：

$${\sum_{i=1}^{3}}Z_i^2\sim \chi_3^2$$

以此类推，得到服从自由度为 n 的卡方分布：

$${\sum_{i=1}^{n}}Z_i^2\sim \chi_n^2$$

又因：

$$Z_1=\frac{\bar{x}_1-\mu}{\sigma_{\bar{x}}}$$

自由度为 n 的卡方分布有可以写作：

$${\sum_{i=1}^{n}}Z_i^2=\frac{\sum(X-\mu)^2}{\sigma^2}\sim \chi_n^2$$

但是总体均值 μ 是不易获得的，而样本均值 xbar 是容易获得的。所以考虑：

$$\begin{aligned} & \frac{\sum(X-\bar X)^2}{\sigma^2} \\ & s^2=\frac{\sum(X-\bar{X})^2}{n-1} \Rightarrow s^2(n-1)=\sum(X-\bar{X})^2 \end{aligned}$$

则有服从 n-1 的自由度的卡方分布 ：

$$\frac{\mathrm{s}^2(\mathrm{n}-1)}{\sigma^2} \sim X_{n-1}^2$$

假设检验和显著水平

小概率原理：一个事件的发生概率很小，那么它在一次试验中是几乎不可能发生的，但在多次重复试验中是必然发生的。

假设检验的基本思想是小概率原理中的反证法——认为小概率事件在一次实验中基本不会发生（发生的概率记为P，称作P值）。在一次统计中，依据实际情况提出假设（称作原假设，记作H0，与之相反的称作备选假设，记作H1），当不符合假设的“极端”事件出现时，通过P值检验，判断是否应该相信假设。

判断逻辑是P值越小，越有理由拒绝假设H1。（极端事件发生的概率越小，越说明极端事件不应该发生，而统计过程发生了小概率事件，说明统计假设不成立。极端事件发生概率越大，说明该事件并不罕见）

在一个抛硬币的真实案例中，假设硬币为材质均匀的钱币，那么每次抛硬币，正反面出现的概率应该各为50%。但在某次实验中，该硬币连续七次的结果均为正面。此刻，是否应该相信硬币的材质均匀？

思路：

事件：“连续七次抛掷，均正面朝上”

原假设：“假设硬币材质均匀”，记为H0

备选假设：“假设硬币材质不均匀”，记为H1

事件发生的概率P值：P=0.00781

由于P值极小，所以有理由怀疑和拒绝原假设H1转而相信H2，即拒绝硬币材质均匀的假设，转而相信硬币材质不均匀。且误判的概率是0.00781。（连续七次正面朝上显然是一个小概率事件，但居然发生了，所以我们有理由怀疑，这个硬币材质不均匀）

在不同的研究领域中，有不同的标准来判断P值是否足够小。例如在物理学中，发现新粒子是小概率事件，那么当观察到新粒子信号超出时，应该先假设“是误差，不存在新粒子”，记作H0。

根据李马公式计算得显著性只有1时，此刻P值(S>=1)约为0.1586552538，这是一个非常大的概率了，所以不能拒绝原假设，或原假设成立，即没有发现新粒子，可能是误差。

当显著性达到5倍超出时，此刻P值(S>=5)足够小，小于0.0000003，意味原假设不成立，应该相信真的发现了粒子。且信号非常强烈。所以5倍显著性便被视为发现新信号的最低显著性。

这里有一个P值在线计算公式，方便理解P值与显著性https://tools.city/statistics/p-value

李马公式

李马公式

李马公式是高能天体物理中关于信号置信度判断的公式，由中国科学院高能物理研究所李惕碚院士和马宇倩于1983年共同提出的统计判据。

$$S=\sqrt{2}\left\{N_{o n} \ln \left[\frac{1+\alpha}{\alpha}\left(\frac{N_{o n}}{N_{o n}+N_{o f f}}\right)\right]+N_{o f f} \ln \left[(1+\alpha)\left(\frac{N_{o f f}}{N_{o n}+N_{o f f}}\right)\right]\right\}^{1 / 2}$$

其中$N_{on}$是向源窗口观测时间内观测到的粒子数，$N_{off}$是背源窗口观测到的粒子数。$\alpha$是向源窗口观测时间与背源窗口观测时间之比。当$N_{on}\gg10$,$N_{off}\gg10$时，可以简化为：

$$S=\frac{N_{o n}-\alpha N_{o f f}}{\sqrt{N_{o n}+\alpha^2 N_{o f f}}}$$

P值检验的步骤是先给出显著性水平α，它代表在原假设为真时，错误地拒绝的概率。

置信度$ζ=1-\alpha$，显著水平和置信度同样是概率，注意区分！

标准正态分布表：https://www.shuxuele.com/data/standard-normal-distribution-table.html

误差传递

https://en.wikipedia.org/wiki/Propagation_of_uncertainty

---

参考

https://www.zhihu.com/tardis/bd/ans/250046834

https://space.bilibili.com/446313875/channel/collectiondetail?sid=288143