概率和统计间的旋转门 – 极限理论

P() 可以看作一个函数，但它是作用在事件而不是数字上的，事件是一个数学名词，定义为一个带有逻辑真假的命题（陈述句），对于逻辑为真的情况，可能会分配到一个非 0 值， 比如 P(骰子结果 1)=1/6, P(明天下雨)=0.4 。但对于逻辑为假的情况，概率均为 0, 比如在标准六面骰子下， P(骰子结果 7)=0; P(0>1)=0 。

对于逻辑值为真的事件，其具体的概率值来自于随机实验的属性，一般可以用概率密度函数来描述，也可以用累计分布函数，或者部分事件的初始概率赋值加上数学结构来隐式地编码整个概率模型，比如 Bayes 网络。

P(X>1) 表示随机实验结果大于 1 的概率，等于对概率密度函数某一侧的积分 \( \int_{x=1}^{\infty}p(x)dx \) 或求和, 如果 X 不可能取到大于 1 的值, 那么 P(X>1) 为 0 。

回到 P(|X-a|>=e), 从具体到抽象分情况理解：

• a=0 且 e 为具体数 1 的情况： P(|X|>=1)

  根据"绝对值"的逻辑定义： \(P(|X|>=1) = P(X>=1 \cup X\leq -1) \) 因此它是概率函数双侧的积分或求和： \( \int_{x=-\infty}^{-1}p(x)dx + \int_{x=1}^{\infty}p(x)dx \)

• a 是不为 0 的实数的情况： P(|X-a|>=1)

  得到是平移后的双侧积分或求和：

  \[ \int_{x=-\infty}^{a-1}p(x)dx + \int_{x=a+1}^{\infty}p(x)dx \]

  这个值可以通过 X 的概率密度函数而确定。

• e 是任意正数且不考虑极限的情况： \( \forall e>0, P(|X-a|>=e)=0 \)

  记 |X-a|>=e 为事件 A, 前文说到， P(A)=0 等同于 A 这个事件对应的命题在当前实验下是一个假命题，因此 P(|X-a|>=e)=0 意味着 |X-a|>=e 不可能成立。如果 e 是一个固定值，那么它说明 X 只能取到 (a-e, a+e) 之间的值，而由于 e 可以取任意小的正数，那么 (a-e, a+e) 区间宽度可以无限缩小但又不为 0 。

  这表明对 X 采样时几乎只能取到值 a，X 的概率密度函数中几乎所有的概率密度都集中到了 a 点上。X 从一个随机变量"退化" 到一个确定的实数 a 上，这种分布有时候也称为 delta 分布（Dirac delta distribution）。如下图（1）所示。

尽管后文的大数定律并不要求依概率收敛公式中的目标变量 Y 是一个真正的随机变量，但本节还是想要理清当 Y 是一般的诸如遵循二项分布或正态分布的随机变量时，足够大的 n 对应的 \( X_{n} \) 和 Y 是什么关系？

如何理解 P(|X-Y|>=e)

首先 X-Y 是两个随机变量的减法，如果要严格计算它们的分布，即便它们独立也需要用到卷积，但我们在讨论一般的变量 X 和 Y ，而没有给定具体分布，因此这里只是定性地去理解。

P(|X-Y|>=e) 实际是一个联合分布，如果 e 是一个比较小的固定值，比如 0.01, 同时将概率也限制在一个很小的值内，如 \( P(|X-Y|>=0.01) \le 0.001 \), 那么其意思是，从 X Y 中分别采样出值 x,y, 它们之间差值大于 0.01 的概率很低，其关系大致如上图（2）所示，它们的概率密度大多集中在 y=x 这条线的附近的狭缝中。

如果随机序列 \( \{X_{n}\} \) 依概率收敛于 Y, 图中两条直线之间的距离就会无限小，但它们不一定是线性关系，比如图（3）中，X,Y 的概率密度集中在狭缝的一个单位圆里，那么它们的相关系数实际为 0 。 所以尽管依概率收敛有一种 X Y 距离逐渐缩小的感觉，但是：

• 它不意味着 X 和 Y 相关系数接近 1, 也就是越来越线性相关

  线性相关体现的是两个随机变量"值"的依赖关系，图 (3) 是一个反例。

• 它不意味着 X 和 Y 所遵循的分布 P 和 Q 的 KL 散度会接近 0 。

  KL 散度的标准公式是 \[ D_{KL}(P||Q) =\int_z P(z)log \frac{P(z)}{Q(z)} dz \]

  首先它衡量的是两个概率分布的关系，这体现在，该公式完全不涉及随机变量 X 和 Y 的值（比如掷骰子中骰子的点数），而只涉及它们对应的概率值 P(x) 和 Q(y) 。 它还要求 Q(z) 不能为 0, 以及随机变量 X,Y 的取值范围是相同的，这里都用 z 表示。

  KL 散度为 0 意味着 P 和 Q 分布相等，它实际可以作为“独立同分布”中的同分布的一个判断依据，而 KL 越接近 0 则认为是两个分布形状上越相似，最终完全相等。但即便 X 和 Y 在依概率收敛语境下很接近（e 很小），那只说明了它们的联合分布概率密度只能在图 (2) 所示的狭缝里，而不能说明联合分布的在 X 和 Y 上的边缘分布 Px, Py 形状是相似的。

  比如把图（3）里的圆替换成一条接近竖线的矩形，所有概率密度集中在这个矩形中，该矩形可以看作是带宽度的直线 x=a 在 y=x+0.01 和 y=x 之间的线段，假设宽度是高度的 1/10, Q(Y) 是在 [a,a+0.01] 区间的均匀分布，每个点的概率密度是 1/0.01=100 ，P(x) 则是 [a,a+0.001] 区间的均匀分布，每个点的概率密度是 1000 。

  P, Q 的 KL 散度为：

  \[ \int_{z}P(z)\log P(z)dz - \int_{z}P(z)\log Q(z)dz \]

  第一项结果是 log100，第二项中，z 只有在 a 到 a+0.001 中 Q(z) 才是非 0 常数 1000, 结果为 0.1log1000, 二者相减结果与 0 相差甚远。

  import math
  math.log(100)-0.2*math.log(1000)
  

  3.2236191301916644
  

在依概率收敛描述中，出现了随机变量序列 \( \{X_{n}\} \) 的说法，这很突兀，因为在概率论基础知识里很少讨论随机变量的序列（更多时候是讨论某个实验里事件序列，比如抛硬币）。该序列还被参数 n 所控制，随着 n 变大它会更接近另外一个随机变量或实数。

合理的数学公式千千万，之所以写出以下公式，是为了捕捉从不确定性到确定性的“感觉”，上一节已经有了直观解释。

\[ \lim_{n \to \infty} P(|X_n - a| \geq \epsilon) = 0 \]

有了这种词汇，其真正的意图是把 \( X_{n} \) 替换成一个非常有实际意义的对象 – n 个随机变量的和 \( S_{n}=X_{1}+X_{2}+\dots +X_{n} \)或者平均 \( Z_{n}=\frac{S_{n}}{n} \) 。

这种情况在现实中非常常见，比如抛了 n 次硬币（且 n 可以任意多）；某种粒子随机游走了 n 步之后；信号传输中每个时刻都有微小的扰动，那么经过 t 时间后误差叠加了多少；这些问题都可以写成多个随机变量相加的场景，这是一种问题分解再组合的还原论或分治思想。

而 a 则对应现实中均值这种概念（比如平均身高，平均正面次数，平均回报，这些是人所想要知道的确定性的值）

因此，在定义中引入序列是在为了表达某种非常具体观念前的词汇准备，先用抽象的数学符号写出形式，比如 forall e 和 |x-a|<e 表示距离上的任意接近，lim P()=0 表示不确定性的逐渐消除。然后代入真正有意义的 x 和 a 。

要证明依概率收敛到一个具体实数，应该先找到一个形如 \( P(|X_n - a| \geq \epsilon) \le f(n) \) 的式子，然后两边取极限：

\[ \lim_{n\to \infty}P(|X_n - a| \geq \epsilon) \le \lim_{n\to \infty}f(n) \]

对于右边，只要能找一个极限为 0 的关于 n 的函数，比如最简单的 \( \frac{1}{n} \), 那么就可以得到 \( \lim P() \le 0 \), 由于概率 P 一定大于等于 0, 意味着 \( \lim P()=0 \) ，这就是依概率收敛的形式了。

对于不等式左边呢？ 由于它涉及到计算区间概率的问题，需要引入概率中的一些著名不等式， 在下一节证明大数定律中具体说明。

从概率论基础可以知道，对于任何 X, Y 有 E[X+Y] = E[X] + E[Y] 以及 Var(X+Y) = Var(X)+ Var(Y)+2Cov(X,Y), 如果扩展到更多的变量，那么对于期望来说，可以很自然地扩展到：

\[ E[X+Y+Z] = E[X]+E[Y]+E[Z] \]

这意味着，如果 X, Y, Z 都有相同的均值 \( \mu \)，那么：

\[ E[\frac{X+Y+Z}{3}] = \mu \]

继续扩展到 n 个变量：

\[ E[\frac{1}{n}\sum_{i}^{n}X_{i}] = \mu \]

这里核心是对 \( \frac{1}{n}\sum_{i}^{n}X_{i} \) 进行概念上的理解，它实际是一个随机变量，将其记做 \( \bar{X} \), 它的值可能是不稳定的，比如抛掷三个硬币，有的时候出现 2 个正面，有点时候出现 0 个正面，如果把正面映射到数字 1, 反面映射到数字 0, 那么 \( \bar{X} \) 在这两次的值分别是 \( \frac{2}{3} \) 和 0 。不过以上式子表明，不论 n 多大，该变量的期望 \( E[\bar{X}] \) 都等于单个随机变量的期望 \( E[X_{i}] \) 。

尽管从期望的视角来看， \( \bar{X} \) 随机变量没有随着 n 的变化而变化，但期望只是代表随机变量的某个性质，从另外一些视角上看， \( \bar{X} \) 是在随着 n 变化的，比如方差。方差本身表示各个值偏离均值的平均程度，或者说概率密度在均值点附近的集中程度，方差为 0 意味着随机变量是某个确定值，比如对于随机变量 X, 它在 X=1 的时候概率为 1, 其他为 0, 那么其方差为 0, 而 \( P(|X-1|\ge \epsilon ) \) 的概率都是 0, 这很符合依概率收敛的描述。因此方差越小，本身就体现了分布依概率“越”收敛于某个数值。

由于目前对 \( X_{i} \) 之间没有添加 "独立性" 约束，这使得对 n 个变量求平均后的方差很难计算，因此退后一步，加上"独立性"和同方差假设，可以得到：

\[ Var(\sum_{i}^{n}X_{i}) = nVar(X_{i}) \]

\[ Var(\frac{1}{n}\sum_{i}^{n}X_{i}) = \frac{1}{n^{2}} nVar(X_{i}) = \frac{Var(X_{i})}{n} \]

两边加上极限： \[ \lim_{n\to \infty}Var(\frac{1}{n}\sum_{i}^{n}X_{i}) = \lim_{n\to \infty} \frac{Var(X_{i})}{n} = 0 \]

以上等式为 0 的前提是单个变量的方差 Var(Xi) 不是无限大（当然如果从微积分角度， n/n 型是收敛的，而只要 Var(Xi) 比 n 的阶更低，比如 \( sqrt(n) \) 那么即便它是无限大，也是收敛到 0 的，但实践中更多情况方差都是固定的，不随着 n 变化）。在该前提下，同方差的假设还可以继续放宽，即便不是同方差，只要每个 \( X_{i} \) 方差是有限就可以从中选出一个方差的最大值记为 \( V_{\max } \) ，得到：

\[ \lim_{n\to \infty}Var(\frac{1}{n}\sum_{i}^{n}X_{i}) \le \lim_{n\to \infty} \frac{V_{\max} }{n} = 0 \]

如何把以上情况描述成依概率收敛的语言呢？核心在于要把 P(|X-a|>=e) 和方差建立一个等式或者不等式关系，而 Chebyshev 不等式补上了这层关系，最终得到了大数定律。

大数定律提供了一个经典的依概率收敛的例子（或者说依概率收敛的语法就是为大数定律创造的）。考虑一个随机变量序列 \(\{X_i\}\)，这些随机变量之间是独立的，且有相同的期望 \(\mu\) 和有限的方差 \(\sigma^2\) （并不一定同分布，或者说 Xi 之间 KL 散度不一定都要是 0）。定义样本均值 \(\overline{X}_n = \frac{1}{n}\sum_{i=1}^{n}X_i\)。

我们想证明 \(\overline{X}_n\) 依概率收敛到 \(\mu\)，即 \(\overline{X}_n \xrightarrow{P} \mu\) ，就要证明对于任何整数 \( \epsilon \) 有 \[ \lim_{n\to \infty}P(|\bar{X}_n - \mu | \geq \epsilon) = 0 \]

这需要用到切比雪夫不等式（Chebyshev's inequality）：

\[ P(|\overline{X}_n - \mu| \geq \epsilon) \leq \frac{\sigma^2}{n\epsilon^2} \]

先不用考虑该公式是如何的到的，只是观察这个公式（就像你从来不知道自己是怎么进化成当前这个样子的，但你至少有照镜子观察过自己），给定任何一个 \( \epsilon \), 只要 \(n\) 足够大，\(\frac{\sigma^2}{n\epsilon^2}\) 都会趋于 0 ，因此这个公式可以用来说明，样本均值 \(\overline{X}_n\) 依概率收敛到总体均值（或者说单个样本的均值） \(\mu\),这就是弱大数定律。

“定律”是早期该公式还没有形式化证明前的称呼（就像现在我们还没有证明它），表示这是一种观察到现象后归纳出来的规律，比如抛大量硬币后发现正面朝上的频率接近 1/2, 而这和人触摸硬币后觉得它是 "均匀" 的感觉是比较契合的，实际上在严格证明出大数定律之前，这条定律已经被人使用了很久，后文的中心极限定理也是如此。这些定律甚至可以作为公理来使用，但接下来会看到，用概率论和其他数学领域中更基础的概念和方法，是可以把它们证明出来，从而变成成为定理的。

尽管我们可能永远无法搞清楚自己是怎么来的（这似乎不妨碍人类的生活），但再阅读两节就可以知道大数定律是如何通过“证明”变成大数定理了。

对于非负的随机变量，假设期望 E(X) 存在（和无限方差类似，有些分布的 E(X) 也是无限的），那么 markov 不等式为：

\( P(X>=t) <= \frac{E(x)}{t} \)

可以用 P(X>t) 的定义来证明，由于 X 的非负性质，因此可以舍弃部分值来构造出不等式：

因此有 \[ P(x>=t) \leq \frac{E(X)}{t} \]

这个不等式算是比较弱的，因为经过了两次缩放。“弱”意味着从数学层面上，它能提供的信息量不是很大，比如说“事物是普遍联系的”并没有带来多少信息，它甚至像是同义反复：事物是普遍联系的，因为它们都是在同一个世界。当然从心理学或哲学上可能会有一些启发，或者在这个大前提下继续去说明事物之间到底有什么样的 "普遍" 联系。

尽管弱，但 Markov 不等式还是提供了对概率分布右侧尾部区间的概率估计的上界，还直观上表明：如果某个概率分布的期望很小，那么从中选出一个很大的数的概率会很低。

可以以它作为起点，引入更多外部信息以获得更强的不等式：

把 Markov 不等式中的随机变量写成 \( Y=(X-\mu)^{2} \) ，且 \(t = \epsilon ^{2} \)后得到：

\[ P((X-\mu)^{2}>=\epsilon^{2}) <= \frac{E[(X-\mu)^{2}]}{\epsilon ^{2}} = \frac{\sigma^{2}}{\epsilon ^{2}} \]

\( (X-\mu)^{2}>=\epsilon^{2} \) 等价于 \( |X-\mu |\ge \epsilon \)

因此：

\(\mu ,\sigma^{2} \) 是 X 变量的期望和方差。

可以看到这个变换把 Markov 不等式里关于 P(x>=t) 单侧概率的下界估计变成了 x 距离均值中心概率集中程度的估计。

另外，以上整个变换看上去没有增加什么信息，比如 \( Y=(X-\mu)^2 \) 里涉及到的均值 \( \mu \) 在 Markov 不等式里也存在 ，但关键是把它塞入到 E[] 后，方差 \( E[(X-\mu)^{2}] \) 并不能从 Markov 不等式中已有的符号中组合获得，尽管方差依赖于均值，但不能只通过均值来计算出方差，所以在应用这个不等式的时候，你需要额外地去测量或者查阅这个随机变量的方差，相比于 Markov 不等式，这是从外部世界新获得的新的观测信息，因此该不等式实际比 Markov 不等式更强。

可以用一个具体的例子来说明：

对于均值为 0 方差为 1 的任意分布，根据 Markov 不等式， P(X>=a)<=1/a ，但：

\[ P(X\ge a)=P(X-0 \ge a)\le P(|X-0|\ge a) \]

根据 Chebyshev 不等式，\(P(|X-0|\ge a) \le \frac{1}{a^{2}} \) ，随着 a 增大，它比 1/a 更快地趋近于 0 （类似算法中 O(n^2) 和 O(n) 效率的差别）。

当 x 取 1 的时候，两个不等式都给出了一个“无信息量的”估计，即 P(X>=1)<=1, 因为概率公理已经给出了该信息，但当 x 取 3 的时候， Markov 给出的是 P(X>=3)<=1/3, 而 Chebyshev 给出的是 P(X>=3)<=1/9, 如果这对应了押注获胜的概率，知道 Chebyshev 的估计可能会让你放弃这个决策。

以上 Chebyshev 不等式给出了一个框架，它在 Markov 不等式的基础上，通过引入外部世界的方差，进一步告诉我们可以用额外的测量信息来获得围绕均值的区间概率的估计。

这种思路可以不断应用，比如不等式中的 \( \epsilon \) 目前似乎是任意的，可以给它赋予人类实践中更容易理解的意义 – 标准差的倍数：

尽管这没有引入额外的信息，但它使得公式更加接近某些应用场景，即人类更喜欢用整数而不是实数去描述对象，比如谈论偏离 2 个标准差的概率，比如企业管理中有一种叫作 6 sigma 的理念，尽管我不知道具体是什么，但它的来源就是以上形式的不等式，我们可以知道，超过 6 个标准差的概率是 1/36, 这是一个比较小的概率（ <3%, 即 97% 的样本都在 6 个标准差中），体现的是要把质量提升到一个很高的程度，当然随着不等式越来越强，能估出的质量会越来越高。

总的来说，这种变换的目的不是做数学上优化，而是对数学表达的润色，使得它更契合人类在特定实践场景下的偏好。

还可以继续对 X 进行变换，将其替换成前文提到的更具有 "统计" 意义的对象 – \(\bar{X}_n=\frac{1}{n}\sum_{1}^{n}X_i \) ，这里 \( X_{i} \) 要有相同的均值，并且：

• 如果有相同的有限的方差 \( \sigma^{2} \), 由前文知道 \( Var(\bar{X}_{n}) = \frac{\sigma^{2} }{n} \), 那么

  \[ P(|\bar{X}_n-\mu|\geq \epsilon) \leq \frac{\sigma^{2} }{n\epsilon ^2} \]

• 如果不同的有限的方差，并且最大值仍然记为 \( \sigma^{2} \), 由前文知道 \( Var(\bar{X}_{n}) \le \frac{\sigma^{2} }{n} \), 还是可以得到以上不等式。

如果两边加上极限就得到了辛钦大数定理（弱大数定律被证明后的称呼）

\[ \lim_{n\to \infty} P(|\bar{X}_n-\mu|\geq \epsilon) = 0 \]

它表明 \(\overline{X}_n\) 依概率收敛到 \(\mu\) 。

这个式子引入了额外信息吗？ 没有，但也有：

• 没有是因为它和刚才把 \( \epsilon \) 替换成 \( k\sigma \) 有点像，引入了一个整数 n, 和 k 类似， 这种整数体现了人类能方便控制的部分，比如公司有能力要求令 k=6 去做严格的质量把关。同样的，人可以很容易重复 n 次抛硬币实验。
• 有是因为它又告诉了我们，整个上下文都变化了，我们不是关注单个随机变量，而是可以控制的 n 个随机变量，在不需要有额外测量信息的情况下就可以估计到多个重复实验的平均值偏离方差的概率，而且这个概率随着 n 会越来越小。更重要的是，这个定律的使用场景更多不是在知道单个 \( X_{i} \) 的均值和方差情况下去对多个 Xi 的平均进行区间概率估计，而是反过来， 在 \( X_{i} \) 有相同的均值和有限的方差、但均值未知的情况下，通过足够多的独立地对 \( X_{i} \) 的采样再求平均可以得到一个很接近 Xi 真实期望的值（依概率保证，即很大概率这个值不会超过某个范围），不关心方差具体值是因为它们在 n 无限大的情况下无关紧要。 这个问题会在蒙特卡洛采样一节中更详细展开。

我们继续对 chebyshev 不等式的“大数”版本进行润色：

对于任何重复发生的事情都可以用一个指示（indicator） 函数来把这个事件转成 0,1 的 Bernoulli 随机变量。 例如，当掷骰子时，一般我们会将该随机变量 X 的取值定义到 {1,2,3,4,5,6}, 但如果只关注其中某个事件，比如是否出现 6, 就可以用指示函数把它变成一个随机变量 Y ，取值为 {1, 0}, 1 表示 6 出现 – 概率为 5/6，0 表示不出现 – 概率为 1/6 。

这样的话，n 次随机事件的平均值 \( \bar{X}_{n} \) 就可以写成 F/n, 其中 F 是事件发生的次数（这是随机变量）， F/n 就是事件发生的频率。而 \( \bar{X}_{n} \) 的均值为 p ，即单次事件发生的理想概率（总体概率）。

弱大数定律（辛钦定理）去掉极限后的公式：

\[ P(|\bar{X}_n-\mu|\geq \epsilon) \leq \frac{\sigma^{2} }{n\epsilon ^2} \]

可以写成以下形式（取极限之后称为伯努利大数定律。）：

\[ P(|\frac{F}{n}-p|\geq \epsilon) \leq \frac{\sigma^{2} }{n\epsilon ^2} \]

由于有了特定的意义，我们可以重新解释它：

当 n 很大时，事件发生的频率较大偏离总体概率的可能性很小，或者说“概率即长期频率”，这实际是频率主义对“概率”的定义， 这是经验主义的，因为频率（empirical frequency）本身是对事件发生的经验观察或测量，这是哲学观点获得数学支撑的一个经典例子。

再回顾弱大数定律的基本形式：

\[ \forall \epsilon \lim_{n\to \infty}P(|\bar{X}_n - \mu | \geq \epsilon) = 0 \]

它表明，当 n 足够大时， \( \bar{X}_{n} \) 随机变量的概率密度会聚集到 \( \mu \) 点上，因此如果 n 足够大，从 \( \bar{X}_{n} \) 中采样出一个值之后，它有很大的概率接近 \( \mu \), 而根据 \( \bar{X}_{n} \) 的定义，它就是对对每个 X[i] 进行采样，然后求和并除以 n 的结果。 如果每个 X[i] 都相同，那么就是对 X 重复采样 n 次求平均，而 \( \mu \) 就是 X 的期望。这里尽管我们进行 n 次实验，但它只是对 \( \bar{X}_{n} \) 随机变量的一次采样。

蒙特卡洛（Monte Carlo method）算法的核心思想是利用随机抽样来估计一个我们很难直接计算的量（通常是期望值、积分或概率），它的名字来源于摩纳哥的蒙特卡罗赌场，表示它是基于概率和随机性的算法思想（或一类算法，类似分治算法，贪心算法），而不是具体某个算法。

由于弱大数定律直接就是关于期望估计准确性的描述，因此用 \( \bar{X}_{n} \) 作为估计器去采样就可以得到对期望的估计， 这里“估计器”（estimator）就是一个随机变量（往往是由更基础的随机变量组合而来），通过对它进行采样就可以间接地去得到想要估计的值，如 E(X) ，从估计器 \( \bar{X}_{n} \) 中进行一次采样在“物理”上实际是对 X 进行 n 次采样然后求平均，因此只要能够用计算机程序模拟出从 X 的分布进行采样的过程，那么数量越多，精度就越高。

比如估计 6 面骰子的期望：

import random
for n in [10, 100, 1000, 10000, 20000]:
    print(sum(random.randint(1, 6) for _ in range(n))/n)

3.4
3.41
3.507
3.5048
3.4684

随着 n 增大，这越来越接近 3.5 。

但为何估计积分也能用这种方法呢？积分估计和期望有什么关系？以单变量为例，积分一般是某个区间中函数曲线下 x 轴之上的面积（如果曲线在 x 轴下方，则计算负面积），但期望是给集合里元素加上一个概率后的加权平均，这似乎没有直觉上的关系， 然而它们的在数学描述上都是一种形式，即都可以写成定积分形式：

假设计算目标是 \( \int_{a}^{b}f(x)dx \) ，如果改写成 \( (b-a)\int_{a}^{b}f(x)\frac{1}{b-a}dx \) ，由于 1/(b-a) 是 a 到 b 区间均匀分布的概率，那么 \( \int_{a}^{b}f(x)\frac{1}{b-a}dx \) 实际就是在 [a,b] 区间上遵循均匀分布的随机变量 X 的函数 f(X) 的期望 E[f(X)] ，于是定积分 \( \int_{a}^{b}f(x)dx=(b-a)E[f(X)]=E[(b-a)f(X)] \)

可以看到，由于概率论中的期望（现实中可能是平均回报）和我们要计算的目标（比如 f(x) 曲线下的面积）都共享了数学中定积分这个结构，于是一个确定性的问题（面积）转换成了不确定性的问题 – 概率分布的期望。

而根据大数定律，概率分布的期望可以用模拟的方式来估计，那么 \( \frac{1}{N}\sum_{i=1}^{n}(b-a)f(X_{i}) \) 就是定积分 \( \int_{a}^{b}f(x)dx \) 的一个蒙特卡洛估计器。

比如估计 cos(x) 在 0 到 2 范围内的积分，先用微积分的方法可以知道 cos 的原函数是 sin(x), 那么定积分就是 sin(2), 再用以上蒙特卡洛估计器来估计：

import math
def cos_i02(n):
    return 2*sum(math.cos(random.random()*2) for _ in range(n))/n

for n in [100, 1000, 10000, 100000, 1000000]:
    print(cos_i02(n))

print(f"{math.sin(2)} <-- sin(2)")

0.8715420390618961
0.8780918181000097
0.8980362348055912
0.9123663503481906
0.908918270874923
0.9092974268256817 <-- sin(2)

可以将以上估计和基础的确定性的数值积分方法对比，即把 [a,b] 区间分成 N 等分，每等分长度是 (b-a)/N, 高度用区间右侧端点的值表示，那么数值积分公式为：

\[ \frac{(b-a)}{N}\sum_{i=1}^{n}f(X_{i}) \]

这和以上蒙特卡洛估计器完全一样，区别仅仅在于这里 Xi 不是随机变量，而是确定的值 \( X_{i} = a+i*\frac{(b-a)}{N} \)

import math

for n in [100, 1000, 10000, 100000, 1000000]:
    print(2/n*sum(math.cos(0+i*2/n) for i in range(1, n+1)))

print(f"{math.sin(2)} <-- sin(2)")

0.8951056483439147
0.9078809768899722
0.9091558091110356
0.9092832653270064
0.909296010678542
0.9092974268256817 <-- sin(2)

这种情况下，确定性的数值积分算法明显比随机的蒙特卡洛算法更快地接近理想值，但为什么要蒙特卡洛算法呢？

这主要是因为确定性数值积分在高维多重积分中会遇到维度爆炸的问题，在一维中，误差大小一般是和 n^{2} 成反比，但到了高维，它是是 \( n^{2/d} \), 因此如果 d 很大，收敛速度会很慢。

而对于蒙特卡洛估计器，它的误差的期望实际就是 \( \frac{(b-a)}{N}\sum f(X) \) 的标准差，而由于：

\[ Var(\frac{1}{n}\sum_{i}^{n}X_{i}) = \frac{1}{n^{2}} nVar(X_{i}) = \frac{Var(X_{i})}{n} \]

因此标准差和 \( n^{\frac{1}{2}} \) 成反比，它不随着维度的变化而变化，因此计算高维积分一般会选择蒙特卡洛算法。

\( \frac{1}{N}\sum_{i=1}^{n}(b-a)f(X_{i}) \) 被称为定积分 \( \int_{a}^{b}f(x)dx \) 的基础蒙特卡洛估计器，之所以是基础，是因为它只是无限多估计器中最简单的那一个，在模拟中使用到的也是最基础的均匀分布采样。

为了更一般地方式在定积分中引入概率，可以写成： \[ \int_{a}^{b}f(x)dx = \int_{a}^{b}\frac{f(x)}{p(x)}p(x)dx = E[\frac{f(x)}{p(x)}] \]

这里引入了一个在 [a,b] 区间上的概率分布 p(x), 前文中 b-a 只是 p(x) 的均匀分布特例, 这时候定积分的蒙特卡洛估计器就是：

\[ F_{N} = \frac{1}{N}\sum_{i=1}^{N}\frac{f(X_{i})}{p(X_{i})} \]

由于 \( \int_{a}^{b}f(x)dx = E[\frac{f(x)}{p(x)}] \), 根据期望线性性质，多个同分布的随机变量的平均值的期望和单个变量的期望是一样的，意味着该蒙特卡洛估计器是无偏的。

理论上，任何满足在 f(x) 不等于 0 时 p(x) 也不等于 0 的概率分布都是可行的，但实际上 p(x) 的选择不单要考虑计算机采样效率，还要追求更小的方差，而如果无偏性导致方差或效率达不到预期，还可以适当放弃无偏的追求，从而引入偏差-方差权衡，根据不同的限制和偏好引出非常多的采样策略和方法，可以搜索 "采样方法" 了解。

这可能不是一个好的问题，因为在中心极限定理被严格证明之前（就像大数定律被证明之前，赌徒早就知道可以抛 1000 次硬币后统计正面的频率来估计这个硬币是不是均匀的），费希尔关于似然函数的理论就默认建立在该定理之上，拉普拉斯在 19 世纪早期用中心极限定理假设证明了最小二乘法（参见《女士品茶》第 9 章）。

因此统计学本身可能不需要这些理论证明来支撑，它本身是关于现实数据的实践类学科，至少可以用经验手段来维持其有效性，然而极限理论的证明使得统计学获得了更加强稳固的数学支撑，反过来也表明概率公理这种形式系统建立得很合理，在极限情况下还能与现实中的观察保持一致。

大数定律建立了纯数学化的概率论和带有"经验" 色彩的统计学之间的关系，但它核心描述的还只是经验均值和理想期望在样本很大的情况下的渐进关系。

同样作为数学推导的结果、一种数学发现，中心极限定理则表明大量样本情况下样本均值随机变量会越来越接近正态分布。这意味着，我们不单单能够在大数定律的支持下，用大量抽样样本去估计均值这种参数，还能够谈论这个均值偏离理想期望的概率（或者说置信度）是多少，因为可以用近似的正态分布精确地计算出抽样结果偏离理想期望的概率，这是经典假设检验里的重要方法，工程方面，我们可以用高斯分布去对系统中大量却微小的背景扰动进行建模。

期望和方差可以看作是对某个随机变量遵循的分布的特定“切片”或特征，期望的重要性可能一方面来自人的直觉，比如用平均成绩衡量两个班级的水平高低，在掷骰子游戏里计算期望以判断自己是否能够长期获利。另一方面来自于数学符号表达的简洁性，许多常见的初等概率分布仅需要期望和方差就可以描述，比如正态分布、二项分布。

假设两个班平均成绩完全相同，这个时候可能就会去对比另外的指标，比如方差越小说明成绩分布越集中，学生的成绩更加均匀，这可能体现了老师对每个学生的关照程度，但如果方差也相同，还能用哪些指标呢？

矩（moment）提供了统一的形式来描述随机变量的分布特征。它分为原始矩和中心矩：

• 第 \(n\) 个原始矩定义为随机变量 \(X\) 的 \(n\) 次幂的期望值 \(E[X^n]\)
• 第 \(n\) 个中心矩定义为随机变量 \(X\) 与均值之差的 \(n\) 次幂的期望值 \(E[(X-E[X])^n]\)

可以看到，均值属于一阶原始矩，而方差则是二阶中心矩（一阶中心矩是多少？），它们都看作是矩的一种特殊形式

偏度（Skewness）是三阶中心矩（用方差常数进行了缩放），的计算公式为：

\[ \text{Skewness} = E\left[\left(\frac{X - \mu}{\sigma}\right)^3\right] = \frac{\sum_{i=1}^{n} (X_i - \overline{X})^3}{(n-1) \cdot \sigma^3} \]

由于这是 3 次方，因此它不像方差只能大于 0, 还可以小于 0 。大于 0 的偏度称为正偏度，表示分布的右尾（正方向）较长，负偏度则表示左尾（负方向）较长。

峰度（Kurtosis）定义如下：

\[ \text{Kurtosis} = E\left[\left(\frac{X - \mu}{\sigma}\right)^4\right] - 3 = \frac{\sum_{i=1}^{n} (X_i - \overline{X})^4}{(n-1) \cdot \sigma^4} - 3 \]

它是（经过方差缩放后的）四阶中心距减去 3 ，减 3 是为了让正态分布的峰度为 0 ，这样其他分布的峰度都是相对正态分布的值。峰度是衡量概率分布尖锐程度的统计量，也可以解释为分布尾部的厚度。

如果不减去 3, 那么超过 3 的峰度表示分布就比正态分布更尖锐，称为“尖峰态”；小于 3 的峰度表示分布比正态分布更平缓，称为“扁平态”。

对这两个度量的直观理解可以参考： 39_偏度和峰度 峰度.zh_en_哔哩哔哩_bilibili, 这里只是想表明它们是分布的矩的特例，不会具体讨论如何利用这两种度量进行分布分析。

对于连续分布的原始矩，计算方式为：

一种形式上和矩很相似的抽象是向量的 p-范数：

\[ \|v\|_p = \left( \sum_{i=1}^{n} |x_i|^p \right)^{\frac{1}{p}} \]

不同的 p-范数对应了向量不同视角下的一个度量。

另外一种内涵上和矩接近的类比是高阶导数：

考虑一个实数值随机变量的分布，其一阶矩（即期望）提供了分布中心的位置信息，就像一阶导数在某点处提供了函数曲线在该点附近的最佳线性近似的斜率。二阶矩（方差）描述了数据围绕中心值的离散程度，类似于二阶导数描述的函数曲线的凹凸性（曲率）。

泰勒展开（也称为麦克劳林级数）可以用函数某点的无穷阶导数的组合来"还原"出真实的函数点的所有信息，同样，后文会展示用分布的无穷阶矩的组合来还原或保存真实的分布函数（一些经典分布只需要少量矩来表示，例如正态分布只需一阶和二阶矩，它的高阶矩都是 0, 就像 f(x)=x 除了一阶导数外的高阶导数都是 0 一样）。

还是用上一节中的例子来说明，如果两个班级的平均成绩相等，可以比较成绩的方差，再相等则比较偏度、峰度，以此类推，也许可以一直对比下去，如果所有阶的 moment 都相等，那么"基本"可以认为这两个班级的成绩分布是一样的。

这类似于对比两个函数泰勒展开的每一项，如果二者都相等，那么可以认为这两个函数（在该点附近）基本相等。

那么我们就需要一个数学对象，它可以把随机变量的所有阶的 moment 都存储在该对象中，然后经过一定操作就可以获得任意阶的 moment 。当然，一个疑问是，我直接给出概率密度函数不就可以计算出任意阶的 moment 吗？为什么要捣鼓出一个新的数学对象。

这里有两种类比性的解释：

• 第一种来自编程：比如你已经有一个数列保存了所有需要的数据，但为了查询和插入的总体效率， 需要用平衡二叉树（比如红黑树）来重新保存这些数据。因此尽管保存的数据本身没有变， 但新的数据结构在实际场景中有很大的效率提升。
• 第二种来自数学本身：当解一个多元一次线性方程组的时候，可以用高斯消元法去解，但如果将其写成矩阵和向量形式 Ax=b, 那么我们可以用线性空间里的各种手段来探测、分解、重组 A 和 b 从而使得计算变得更加方便。因此尽管方程组和矩阵或向量都拥有一样的数据，但将方程组变换到一个具有丰富数学性质的空间里，我们有更多方便的工具来展现方程的特点和性质，也可以获得更灵活简便的计算方法。

对于分布的 moment, 我们将其保存在所谓的生成函数 (Generating Function) 中，有时候也翻译为母函数，它都表示一种展开、繁衍的意思。该函数的形式一般是多项式进行累加，就像一个数据结构中的链表。

例如，矩生成函数（矩母函数）如下：

\[ M_{X}(t) =\sum_{0}^{\infty}m_{n}\frac{t^{n}}{n!} \]

其中的 \( m_{n} \) 是 X 随机变量的 n 阶矩，对于原始矩就是 \( E[X^n] \) ， 注意这个形式实际可以看作函数 M 在 0 点附近的泰勒展开式，其中 \( m_{n} \) 是函数在 0 点的 n 阶导数，反过来这意味着，只要对 \( M_x(t) \) 求第 n 阶导数再取 t=0，就可以得到 X 的 n 阶矩，因此该函数把随机变量 X 的所有矩信息都保存住了。（所以这里实际把 n 阶导数和 n 阶矩性质绑定了）

由于期望的线性性质，可以将上式子继续写成：

\[ M_{X}(t) = \sum_{0}^{\infty}E[X^{n}]\frac{t^{n}}{n!} = E[\sum_{0}^{\infty}X^{n}\frac{t^{n}}{n!}] = E[\sum_{0}^{\infty}\frac{(tX)^{n}}{n!}] \]

最后一个式子内无穷级数实际等于 \( e^{tX} \) （X 是一个随机变量），因此有：

\[ M_{X}(t) = E[e^{tX}] \]

注意这是一种"高层"的符号表达，矩母函数的意义还是来自以上展开的多项式形式。但由于指数形式有比较好的运算性质（比如乘法），因此这种"高层" 的对随机变量求指数之后再求期望的操作也称为变换 transform ，有许多这种变换的例子，比如

比如针对离散随机变量的 z-变换, 针对连续随机变量的 拉普拉斯变换，信号处理领域的傅里叶变换等，本文不对此展开。

这里给出另外一种对分布的变换 – 特征函数：

\[ C_{X}(t) = E[e^{itX}] \]

它和矩生成函数唯一区别是在指数位置加入了虚数 i, 这是出于数学性质上的考虑，因为 \( e^{itX} \) 中指数是纯虚部的，模长 \( \|e^{itX}\| \) 就等于 1 ，因此对它求期望永远不会出现无穷大的情况，而如果是矩母函数，如果有些 X 随机变量取值可以无穷大，再取指数后期望可能不存在。因此后文给出的证明是基于特征函数，然而基于矩母函数也是可以的（不严格考虑极端情况或者加上一点假设限制），参考： 正题3：中心极限定理证明_哔哩哔哩_bilibili

无论是特征函数或者矩母函数，经过转换后的随机变量有很好的加法/乘法性质，假设随机变量 \( X_{i} \) 之间独立的话，会得到：

\[ M_{X_1+\dots +X_n}(t)=E[e^{t(X_1 +X_2+\dots)}] = E[e^{tX_{1}}]E[e^{tX_{2}}]\dots = \prod_{i=1}^{n}M_{X_i}(t) \]

注意是独立性使得 E[XY]=E[X][Y] 的，这个性质使得多个随机变量之间的加法在转换后变成了特征函数之间的乘法，而指数形式的乘法又有很好的性质，这是证明中心定理的关键技巧。

特征函数中涉及到复指数的性质，而后文还涉及到欧拉方程，完全严格的证明会需要用到数学分析、复分析等知识，相当于为了写一段 println("hello") 去重写操作系统和编译器，因此这里不打破抽象层去列出它们底层的合理性， 而只是做一点"高层"的解释

直接求出多个随机变量的和的分布公式是非常困难的，即便多个随机变量是独立的，那么也需要用到卷积。而两个变量的卷积就涉及积分或累加，若是 \(n\) 个变量卷积就更加困难了。根据上节矩生成函数（MGF）或特征函数（CF）的性质， 多个随机变量的和可以变成多个特征函数的乘积，那么只要这个结果和正态分布的特征函数一样，就可以认为多个随机变量和服从正态分布。

这种思路可以对应到两个类比：

• 数的乘法在对数空间里变成了加法
• 当我们检查某个软件下载后是否损坏，一般可以先计算该软件的签名，然后和官方提供的签名进行对比，如果一致，那么很大概率说明软件没有问题。但密码学是通过哈希函数几乎不会碰撞的性质来保证"签名" 的有效性，特征函数同样要有类似的性质，实际上特征函数有比哈希函数更好的性质，可以通过特征函数反向推导出唯一的“加密前的”概率分布函数，而不需要用概率来保证，该性质的证明参考 mathematical statistics - Proof that moment generating functions uniquely determine probability distributions - Cross Validated。

随机变量 \(X_1, X_2, ..., X_n\) 具有相同的期望 \(\mu\) 和方差 \(\sigma^2\)，我们首先关注它们的和 \(S_n = \sum_{i=1}^n X_i\) 的分布。

通过期望和方差的性质就可以知道 \( E[S_{n}] = n\mu \) 以及 \(Var(S_{n}) = n\sigma^{2} \), 因此，随着 n 的增大， \( S_{n} \) 的方差会越来越大，概率分布会越来越平坦，或者说越来越发散，因此 \( S_{n} \) 本身没有收敛性质，即便它能接近一个高斯分布，也是非常扁平的高斯分布， 而这时候用高斯分布去估计概率或者置信区间误差会很大。所以更好的考虑对象是 \( Z_{n} = \frac{S_{n}}{n} \), 根据期望和方差在线性变换下的性质就可以知道 \( Z_{n} \) 的期望是 \( \mu \), 方差是 \( \frac{\sigma^{2}}{n} \) 。

以上只是对均值（一阶矩）和方差（二阶矩）在 n 个分布累加情况下性质的论述，中心极限定理还给 \( Z_{n} \) 的分布加上了更精细的限定，即当 n 逐渐增大， \( Z_{n} \) 的分布会越来越接近以 \( \mu \) 为均值，以 \( \frac{\sigma^{2}}{n} \) 为方差的正态分布 \(N(\mu, \frac{\sigma^{2}}{n}) \) 。

如果对 \( Z_{n} \) 进行标准化：

\[ Z' = \frac{Z_{n }-\mu}{\sigma /\sqrt{n}} = \frac{S_n - n\mu}{\sigma\sqrt{n}} = \sum \frac{X_{i}-\mu }{\sigma \sqrt{n}} \]

那么中心极限定理表明 \( Z' \) 会随着 n 的增大而越来越接近标准正态分布 \( N(0, 1) \)

由于 \( Z' \) 可以写成多个标准化的变量 \( Y_{i}= \frac{X_{i}-\mu }{\sigma \sqrt{n}} \) 求和的结果，所以接下来要证明多个 \( Y_{i} \) 求和后的分布接近标准正态分布，也就是多个 \( Y_{i} \) 的特征函数的乘积接近标准正态分布的特征方程 \( \) \( e^{-t^{2} /{2}} \) 。

先准备 \( Y_{i} \) 的几个矩信息，证明中会使用到：

• \( E[Y_{i}] = 0 \), 因为 \( Y_{i} \) 本身就是对 \( X_{i} \) 进行标准化基础上再除以 \( \sqrt{n} \) 的结果。
• \( Var[Y_{i}] = Var(\frac{1}{\sigma \sqrt{n}}X_{i}) = \frac{1}{\sigma ^{2}n} \sigma ^{2} = \frac{1}{n}\)
• \( E[Y_{i}^{2}] = Var[Y_{i}] + E[Y_{i}]^{2} = \frac{1}{n} + 0 = \frac{1}{n}\)

接着是中心极限定理的证明：

1. 写出单个变量 \( Y_{i} \) 的特征方程：

   \[ C_{Y_i}(t) = E[e^{tY_i}] \]

   继续对 e 展开成级数形式就得到： \[ C_{Y_i}(t) = E[e^{tY_i}] = E[1+itY_{i}-\frac{1}{2}t^{2}Y_{i}^{2}+O(t^{3})] \]

   这里 \( O(t^{3}) \) 表示更高阶乘的余项，由于我们只关注 \( t^{2} \), 因此只展开到 2 阶，根据期望的线性性质，可以继续展开：

   \[ C_{Y_i}(t) = E[1+itY_{i}-\frac{1}{2}t^{2}Y_{i}^{2}+O(t^{3})] = 1+itE[Y_{i}]-\frac{1}{2}t^{2}E[Y_{i}^{2}]+O(t^{3}) = 1+0-\frac{1}{2}\frac{t^{2}}{n} + O(t^{3}) \]

   这些步骤实际上是回到矩母函数的"底层"定义，即用多项式的系数来保存分布的各阶矩，前文把它写成 e 的指数形式是为了之后能够更好地处理"高层"的乘法。

   整理以上过程得到： \[ C_{Y_i}(t) = 1-\frac{1}{2}\frac{t^{2}}{n^{2}} + O(t^{3}) \]

   这里 \( O(t^{3}) \) 表示更高阶乘的余项的期望了，同样不关心它的具体值。

2. 对 n 个 \( Y_{i} \) 求和的结果 \( Z' \) 计算特征方程：

   \[ C_{Z'}(t) = E[e^{tZ'}] = E[\prod e^{tY_{i}}] \]

   由于 \( Y_{i} \) 之间的完全独立性，因此根据 E[XYZ] = E[X]E[Y]E[Z], 最后一个式子求期望可以写到累积内部：

   \[ C_{Z'}(t) = E[\prod e^{tY_{i}}] = \prod E[e^{tY_{i}}] = \prod C_{Yi}(t) \]

   完全独立和两两独立

   E[XYZ] = E[X]E[Y]E[Z] 实际上是在说三个随机变量的行为互不影响，即它们之间不存在任何形式的相关性或依赖性。 因为 E[XYZ] 可以写成 E[UZ], 这里 U = XY, 那么 E[XYZ]=E[UZ]=E[U]E[Z] = E[X]E[Y]E[Z]

   这意味着 Z 和 U 也相互独立，根据乘法交换律，X 要和 YZ 独立， Y 和 XZ 独立，因此它们完全独立。

   再根据 \( Y_{i} \) 是同分布假设，所有的 \( C_{Y_{i}}(t) \) 都相同，因此：

   \[ C_{Z'}(t) = C_{Yi}(t)^{n} = (1-\frac{1}{2}\frac{t^{2}}{n} + O(t^{3}))^{n} \]

   对比 e 的定义 \( e=\lim\limits_{n\to \infty} (1+\frac{1}{n})^{n} \), 它可以扩展到：

   \( e^{x}=\lim\limits_{n\to \infty} (1+\frac{x}{n})^{n} \)

   那么 \( (1-\frac{1}{2}\frac{t^{2}}{n})^{n} \) 就应该是 \( e^{-\frac{1}{2}t^{2}} \), 这就是标准高斯分布的特征函数， 而 \( 1-\frac{1}{2}\frac{t^{2}}{n} + O(t^{3}) \) 中的 \( O(t^{3}) \) 里由于分母中都是 n 的更高的乘方，因此它不会对极限产生影响，最终

   \[ \lim\limits_{n \to \infty} C_{Z'}(t) = \lim\limits_{n \to \infty} (1-\frac{1}{2}\frac{t^{2}}{n} + O(t^{3}))^{n} = e^{-\frac{1}{2}t^{2}} \]

   这就说明了 \( Z' \) 随着 n 越大越趋近标准正态分布。

总结"独立" 和 "同分布" 在证明中的作用：

• 独立性： 这是使得特征函数（的指数形式）可以进行累乘的关键，中心极限定理有其他更宽松的变种，比如不需要同分布， 但独立性假设仍然无法去掉。
• 同分布： 这使得随机变量卷积在特征函数空间变成了某个函数的 n 次方，这和 e 的极限定义匹配上了，最终导出了最终 \( e^{t^{2}} \) 形式，而这是正态分布的标准"签名" 。

前文通过中心极限定理得到了 \( e^{-\frac{x^{2}}{2}} \) 形式的公式，但它并不是一个合格的概率分布，所谓合格，要保证其积分为 1 ， 也就是 \( \int_{-\infty}^{\infty} e^{-\frac{x^{2}}{2}} = 1 \) 是否成立，如果不成立，比如等于某个其他的值 A, 那么应该用 \( \frac{1}{A}e^{-\frac{x^{2}}{2}} \) 来表示最终的分布。当然更坏的情况是这个积分不存在，比如是无穷大。

求解积分 \(\int_{-\infty}^{\infty} e^{-\frac{x^2}{2}} \, dx\) 需要用到平方的技巧，即不直接计算该积分，而是计算积分的平方：

\[ I^2 = \left(\int_{-\infty}^{\infty} e^{-\frac{x^2}{2}} \, dx\right)^2 \]

1. 把式子变成二维平面上关于 \( x \) 和 \( y \) 的高斯函数的积分。： \[ I^2 = \left(\int_{-\infty}^{\infty} e^{-\frac{x^2}{2}} \, dx\right) \left(\int_{-\infty}^{\infty} e^{-\frac{y^2}{2}} \, dy\right) \]
2. 将 \( x \) 和 \( y \) 的平方和重写： \[ I^2 = \int_{-\infty}^{\infty} \int_{-\infty}^{\infty} e^{-\frac{x^2 + y^2}{2}} \, dx \, dy \]
3. 转换到极坐标： 在极坐标中，\( x = r \cos \theta \) 和 \( y = r \sin \theta \)，而 \( dx \, dy \) 转换为 \( r \, dr \, d\theta \)。因此，积分变为： \[ I^2 = \int_{0}^{2\pi} \int_{0}^{\infty} e^{-\frac{r^2}{2}} r \, dr \, d\theta \]
4. 计算内层积分： \[ \int_{0}^{\infty} e^{-\frac{r^2}{2}} r \, dr \] 使用变量替换 \( u = \frac{r^2}{2} \)，则 \( du = r \, dr \)。因此， \[ \int e^{-\frac{r^2}{2}} r \, dr = \int e^{-u} \, du = -e^{-u}\Big|_{0}^{\infty} = 1 \]
5. 计算外层积分： \[ I^2 = \int_{0}^{2\pi} 1 \, d\theta = 2\pi \]

因此有 \( I^2 = 2\pi \), 从而得出 \( I = \sqrt{2\pi} \)。

因此，标准正态分布的公式为： \[ \frac{1}{\sqrt{2 \pi}}e^{-\frac{x^2}{2}} \]

正态分布公式的推导最早来源于高斯对天文观测的误差的分析，而不是概率论。在天文观测中记录了许多数据， 不同观测数据往往存在误差，这些误差可能由观测人员操作的精细程度、 观测工具的精度、天气等共同作用。比如一年前出现了一个新的星体，已经跟踪了一年， 每个月各个天文台里不同设备都记录了数据，现在要预测下个月星体出现的位置，那么这就是根据过去的一组 带有误差的数据预测新的数据，高斯在研究中发明了最小二乘法和误差函数：

他假设了一个函数 \( \phi(x) \) 来描述误差的分布，从直觉和数据统计的结果上，该函数应该满足以下三个性质：

• 小的误差比大的误差更多
• 正误差出现的次数和负误差次数基本相同，所以 \( \phi(x) \) 应该是偶函数。
• 预测出的结果应该和下个月从不同天文台使用不同工具获得的数据的均值很接近

为了找出满足这些性质的函数，高斯假设 \( \phi(x) \) 的导数和 \( x \) 的关系是线性的：

\[ \frac{\phi'(x)}{\phi(x)} = kx \]

这体现了如果误差越大，那么它更强的衰减倾向，即更不可能出现更大误差

对两边进行积分：

\[ \int \frac{1}{\phi(x)} d\phi(x) = k \int x \, dx \]

解得：

\[ \ln |\phi(x)| = \frac{kx^2}{2} + C \]

两边取指数：

\[ \phi(x) = A \exp\left( \frac{kx^2}{2} \right) \]

参考 The Origins of the Normal Distribution | by William Sundstrom | Medium

假设 n 非常大， p 非常小，并且 \( np= \lambda \) 是一个常数，它可以表示单位时间里发生的次数是固定的，一个例子是，cpu 每秒 10^9 时钟周期，其中只有 1 次或 0 次会出现周期偏移。

n 可以认为是趋于无穷大，而 k 一般是一个远小于 n 的值，那么最后一行公式的后半部分取极限情况下 \[ (1-\frac{\lambda}{n})^{n-k} = (1-\frac{\lambda}{n})^n (1-\frac{\lambda}{n})^{-k} \simeq e^{-\lambda}\cdot 1 = e^{-\lambda}\],

前半部分如下：

因此在词极限下，二项分布 \( {n\choose k} p^{k}(1-p)^{n-k} \\ \) 就等价于泊松分布 \( (\frac{\lambda^{k}}{k!})e^{-\lambda} \)

Poisson 分布是假设大量的 p 非常小的伯努利分布叠加的结果。 我们从 \( (1+\frac{1}{n})^{n} \) 的极限角度计算出了 Possion 分布公式，但中心极限定理表明大量的独立同分布累加起来会收敛到正态分布形式，其期望也是 \( np=\lambda \), 这样我们就可以 得到一个以 \( \lambda \) 为均值，为方差的正态分布 \[ \frac{1}{\sqrt{2\pi \lambda }}e^{\frac{-(x-\lambda)^{2}}{2\lambda }} \]

这里把 Possion 分布里的 k 也用连续值插值了 这个分布应该和 Possion 分布在极限情况下近似，于是有：

设 \( k=\lambda \) ，那么有

这就是 Stirling 对阶乘的估计公式， \( x! \approx e^{-x}x^{x}\sqrt{2\pi x } \)

取 ln 形式为： \( \ln x! \approx x\ln x -x +\frac{1}{2}\ln 2\pi x \)