斯坦福-CS109-计算机科学家的概率论-三-

斯坦福 CS109：计算机科学家的概率论（三）

原文：chrispiech.github.io/probabilityForComputerScientists/en/index.html

译者：飞龙

协议：CC BY-NC-SA 4.0

信息论

原文：chrispiech.github.io/probabilityForComputerScientists/en/part4/information_theory/

---

信息论是一个极其强大的视角，在众多算法中扮演着核心角色，包括决策树、WordleBot、自适应测试、最佳扑克玩法甚至数据压缩（如霍夫曼编码或 JPEG 文件）！本章的目标是在展示信息论的强大功能的同时，尽可能保持内容的直接性。为此，一个很好的起点是思考如何编写一个能够玩“想一个动物”问答游戏的机器人。

想一个动物！

“想一个动物”游戏是这样的：人类将要想一个动物。我们假设他们选择动物的概率分布是已知的（基于四岁孩子对动物的喜爱程度）：

你的算法任务是选择下一个要问的问题。假设你被提供了一个包含经典问题（如）的“是”或“否”问题库：

• 它是宠物吗？

• 它住在水里吗？

• 你在想着狗吗？

选择最佳问题

我们如何选择最佳问题来提问呢？考虑一个简化版的游戏，有五种动物（狗、猫、大象、熊和猴子），并且只有两个问题可供选择：“它是宠物吗？”和“它是狗吗？”对于每个问题，我们可以使用全概率公式来思考回答“是”或“否”的概率。甚至更好！对于每个可能的问题和每个可能的答案，如果我们看到那个问题的答案，我们可以思考由此产生的动物概率质量函数！

我们已经非常接近了！如果我们能够量化四个结果概率质量函数的不确定性，我们就可以简单地选择一个最小化我们对人们所想动物预期不确定性的问题。这种不确定性的量化形式上称为“熵”，它是信息论中的关键概念。

从高层次衡量不确定性

让 \(X\) 是任何随机变量。量化 \(X\) 的概率质量函数所代表的多少不确定性的一个非常优雅的方法是思考 \(X\) 可能取的所有值，并对每个值计算如果 \(X\) 确实取了那个值你会多么惊讶。如果你对这些惊讶值的加权求和，你会得到随机变量的预期惊讶值

\[\begin{align*} &\text{不确定性}(X) \\ &= E[\text{惊讶}(X)] \\ &= \sum_x \text{惊讶}(x) \cdot P(X=x) \end{align*}$$ 这就是核心思想！主要剩余的任务是定义我们所说的“事件惊讶”是什么。 ### 事件惊讶度衡量 如果我们被告知 $X$ 取值 $x$，我们应该如何量化我们会感到的惊讶程度？有几种方法可以量化惊讶，所有这些方法都是基于事件 $\P(X=x)$ 的概率。以下是三种合理的惊讶函数：  所有这些函数都符合以下要求： + 低概率事件令人惊讶 + 高概率事件并不令人惊讶 + 惊讶是概率的单调递减函数 由于许多原因，信息理论使用中间方程的变体来定义惊讶：$$\text{Surprise}(E) = \log_2 \frac{1}{P(E)}\]

观察事件 \(P(E)\) 和 \(\text{Surprise}(E)\) 的值之间的关系，我们可以观察到一些直观的关系：

事件概率 \(P(E)\)	事件惊讶 \(\text{Surprise}(E) = \log_2 \frac{1}{P(E)}\)
1/2	1
1/4	2
1/8	3
1/16	4
1/32	5
1/64	6
1/128	7

如果一个概率为 \(P(E) = 1/16\) 的事件发生，我们将比一个概率为 \(P(E) = 1/2\) 的事件发生时感到四倍的“惊讶”。这感觉很好！

定义： 事件的惊讶（信息含量）

事件 E 的信息含量，也称为惊讶或自信息，是一个随着事件概率的降低而增加的函数。当概率接近 1 时，事件的惊讶值较低，但如果概率接近 0，事件发生的惊讶值较高。这种关系由以下函数描述：

\[\begin{align*} \text{Surprise}(E) &= \log_2 \Big(\frac{1}{\p(E)} \Big) \end{align*}$$ 这可以写成（以一种更令人困惑的方式）：$$\begin{align*} \text{Surprise}(E) &= \log_2 \Big(\frac{1}{\p(E)} \Big) \\ &= \log_2 \p(E)^{-1} \\ &= -\log_2 \p(E) \end{align*}$$ 在事件惊讶的初始定义中，函数名 $I$ 被用作惊讶的简称。$I$ 代表“信息含量”或“自信息”，这是惊讶的另一种名称。 我们还可以讲述许多其他的故事，来说明为什么 $\log_2 \frac{1}{\p(E)}$ 是我们衡量惊讶程度的绝佳选择。信息理论的奠基人克劳德·香农选择了以 2 为底的对数，因为这可以使你用比特（如计算机中使用的 0 和 1 值）来表达你的惊讶程度。信息理论有许多应用，但最初是在他试图找到一种方法来最优地压缩基于文本的数据时发明的！ ### 随机变量的不确定性（熵） 现在我们已经有了惊讶的正式定义，我们可以重新审视随机变量不确定性的计算。 ***定义：*** 随机变量的不确定性（熵） 让 $H$ 成为衡量我们对随机变量 $X$ 不确定性的度量。定义 $H$ 为观察 $X$ 分配的预期惊喜。$H(X) = E[\text{Surprise}(X)]$。使用无意识统计学家定律和事件惊喜的定义，我们可以展开 $H$ 的公式如下： $$ H(X) = \sum_{x \in X} \log_2 \Big(\frac{1}{\p(X=x)}\Big) \cdot \p(X=x) $$ 不确定度 ($H$) 也可以（以更令人困惑的方式）重写为：$$\begin{align*} \text{H}(X) &= \sum_{x \in X} \text{Surprise}(X=x) \cdot P(X=x) \\ &= \sum_{x \in X} \log_2 \frac{1}{P(X=x)} \cdot P(X=x) \\ &= \sum_{x \in X} \log_2 P(X=x) ^{-1}\cdot P(X=x) \\ &= \sum_{x \in X} - \log_2 P(X=x)\cdot P(X=x) \\ &= - \sum_{x \in X} \log_2 P(X=x)\cdot P(X=x) \\ &= - \sum_{x \in X} \log_2 P(X=x)\cdot P(X=x) \end{align*}\]

这里是计算随机变量不确定度 \((H)\) 的代码，基于其概率质量函数：

import numpy as np

def calc_uncertainty(pmf):
  """
  Calculate how much uncertainty is represented by this 
  probability mass function. Also known as the Entropy of a
  random variable, H(X).
  """
  uncertainty = 0
  for x in pmf:
      p_x = pmf[x]
      # skip zero probabilities
      if p_x == 0: continue
      suprise_x = np.log2(1/p_x)
      uncertainty += suprise_x * p_x
  return uncertainty 

现在我们已经拥有了在“想一个动物”游戏中选择最佳问题的所有理论工具。对于每个可能的概率质量函数，我们都可以计算该 PMF 的不确定度 \((H)\)：

问题“它是宠物吗？”的预期不确定度为 1.3。问题“它是狗吗？”的预期不确定度为 1.7。因此，如果我们问“它是宠物吗？”这个问题，我们预期对朋友所想的动物会更有信心。

这只是信息论众多应用之一！

分布之间的距离

原文：chrispiech.github.io/probabilityForComputerScientists/en/part4/divergence/

---

这里有一个重要的问题！你如何量化两个分布之间的差异？换句话说，如果我有两个概率质量函数，我能计算出一个数字来表示它们彼此之间有多大的差异？这里有三种合理的方式来量化分布之间的距离：

总变差距离

遍历所有可能的值并计算概率的绝对差异。

定义： 总变差距离

设 \(X\) 和 \(Y\) 为离散随机变量

\[\begin{align*} \text{TV}(X, Y) = \frac{1}{2} \sum_{i} |\P(X = i) - \P(Y = i)| \end{align*} \]

地球迁移距离

想象一个分布是一堆泥土。要让它看起来和另一个分布一样，需要做多少工作？也称为 Wasserstein 度量。

这个值非常有意义，但它没有封闭形式的方程。相反，它通过求解线性规划来计算。如果随机变量可以取的值的范围大小为 \(n\)，线性规划的计算时间为 \(O(n³ \cdot \log n)\)，这非常慢。

库尔巴克-莱布勒散度

设 \(X\) 和 \(Y\) 为离散随机变量。计算当实际概率质量函数为 \(X\) 时，使用 \(Y\) 作为概率质量函数带来的“额外惊喜”的期望值。

定义： 库尔巴克-莱布勒散度

设 \(X\) 和 \(Y\) 为离散随机变量

\[\begin{align*} \text{KL}(X, Y) &= \sum_{x \in X} \text{ExcessSurprise}(x) \cdot P(X=x) \\ &= \sum_{x \in X} \Big[\text{Surprise}_X(x) - \text{Surprise}_Y(x)\Big] \cdot P(X=x) \\ &= \sum_{x \in X} \Big[\log_2\frac{1}{P(Y=x)} - \log_2 \frac{1}{P(X=x)}\Big] \cdot P(X=x) \\ &= \sum_{x \in X} - \log_2 P(Y=x) + \log_2 P(X=x) \cdot P(X=x) \\ &= \sum_{x \in X} \log_2 \frac{P(X=x)}{P(Y=x) } \cdot P(X=x) \\ \end{align*} \]

这里是计算每年观察到的飓风概率质量函数与预测的泊松分布差异的示例代码

from scipy import stats
import math

def kl_divergence(predicted_lambda, observed_pmf):
    """
    We predicted that the number of hurricanes would be 
    X ~ Poisson(predicted_lambda) and observed a real world
    number of hurricanes Y ~ observed_pmf
    """
    X = stats.poisson(predicted_lambda)
    divergence = 0
    # loop over all the values of hurricanes
    for i in range(0, 40):
      pr_X_i = X.pmf(i)
      pr_Y_i = observed_pmf[i]
      excess_surprise_i = math.log(pr_X_i / pr_Y_i)
      divergence += excess_surprise_i * pr_X_i
    return divergence 

维度诅咒

原文：chrispiech.github.io/probabilityForComputerScientists/en/examples/curse_of_dimensionality/#part4_examples/

---

在机器学习中，就像计算机科学的许多领域一样，经常涉及高维点，而高维空间具有一些令人惊讶的概率性质。

一个随机值 \(X_i\) 是一个均匀分布的值：Uni(0, 1)。

一个 \(d\) 维的随机点是一个 \(d\) 个随机值的列表：\([X_1 \dots X_d]\)。

一个随机值 \(X_i\) 如果 \(X_i\) 小于 0.01 或者 \(X_i\) 大于 0.99，则该值接近边缘。一个随机值接近边缘的概率是多少？

设 \(E\) 为随机值接近边缘的事件。\(P(E) = P(X_i < 0.01) + P(X_i > 0.99) = 0.02\)

一个三维的随机点 \([X_1, X_2, X_3]\) 如果其任何一个值接近边缘，则该点接近边缘。一个三维点接近边缘的概率是多少？

该事件等价于点的任何维度都不接近边缘的补事件，即：\(1 - (1 - P(E))³ = 1 - 0.98³ \approx 0.058\)

一个 \(100\) 维的随机点 \([X_1, \dots X_{100}]\) 如果其任何一个值接近边缘，则该点接近边缘。一个 \(100\) 维点接近边缘的概率是多少？

同样，它也是：\(1 - (1 - P(E))^{100} = 1 - 0.98^{100} \approx 0.867\) 高维点的许多其他现象：例如，点之间的欧几里得距离开始收敛。

Thompson 采样

原文：chrispiech.github.io/probabilityForComputerScientists/en/examples/thompson/

---

想象一下你必须做出以下一系列决定。你有两种药物可以给药，药物 1 或药物 2。最初你不知道哪种药物更好。你想要知道哪种药物最有效，但与此同时，探索的成本很高。

这里有一个例子：

Welcome to the drug simulator. 
There are two drugs: 1 and 2.

Next patient. Which drug? (1 or 2): 1
Failure. Yikes!

Next patient. Which drug? (1 or 2): 2
Success. Patient lives!

Next patient. Which drug? (1 or 2): 2
Failure. Yikes!

Next patient. Which drug? (1 or 2): 1
Failure. Yikes!

Next patient. Which drug? (1 or 2): 1
Success. Patient lives!

Next patient. Which drug? (1 or 2): 1
Failure. Yikes!

Next patient. Which drug? (1 or 2): 2
Success. Patient lives!

Next patient. Which drug? (1 or 2): 2
Failure. Yikes!

Next patient. Which drug? (1 or 2): 

这个问题出奇地复杂。有时它被称为“多臂老丨虎丨机问题！”事实上，对这个问题的完美答案可能非常难以计算。有许多近似解，并且这是一个活跃的研究领域。

一种解决方案已经变得相当受欢迎：Thompson 采样。它易于实现，易于理解，有可证明的保证 [1]，并且在实践中表现良好 [2]。

你对选择的了解

Thompson 采样的第一步是表达你对选择（以及你不知道的）的了解。让我们回顾一下上一节中提到的两种药物。到结束时，我们测试了药物 1 四次（1 次成功）和药物 2 四次（2 次成功）。一种复杂的方式来表示我们对药物 1 和药物 2 后面隐藏的概率的信念是使用 Beta 分布。让 \(X_1\) 表示对药物 1 概率的信念，让 \(X_2\) 表示对药物 2 概率的信念。$$ $$\begin{align} X_1 \sim \Beta(a = 2, b = 4)\ X_2 \sim \Beta(a = 3, b = 3) \end{align}$$ $$

回想一下，在具有均匀先验的 Beta 分布中，第一个参数 \(a\) 是观察到的成功次数 + 1。第二个参数 \(b\) 是观察到的失败次数 + 1。从图形上观察这两个分布是有帮助的：

如果我们不得不猜测，药物 2 看起来更好，但仍然有很多不确定性，这体现在这些信念的高方差中。这是一个有用的表示。但我们如何利用这些信息来做出下一个药物的好决定。

做出选择

确定正确的选择很难！如果你只有一个更多的病人，那么你该做什么就很明确了。你应该计算 \(X_2 > X_1\) 的概率，如果这个概率超过 0.5，那么你应该选择 \(a\)。然而，如果你需要持续给药，那么正确的选择就不那么明确了。如果你选择 1，你将错过了解 2 的机会。我们该怎么办？我们需要平衡“探索”的需求和利用已知信息的需要。

汤普森采样的简单思想是随机地根据其成为最优选择的概率来做出选择。在这种情况下，我们应该选择药物 1，其选择的概率是 1 大于 2。人们在实践中是如何做到这一点的呢？他们有一个非常简单的公式。从每个 Beta 分布中随机采样。选择其采样值较大的选项。

sample_a = sample_beta(2, 4)
sample_b = sample_beta(3, 3)
if sample_a > sample_b:
    choose choice a 
else:
    choose choice b 

采样意味着什么？这意味着根据概率密度函数（或概率质量函数）选择一个值。所以在我们上面的例子中，我们可能会为药物 1 采样 0.4，为药物 2 采样 0.35。在这种情况下，我们会选择药物 1。

在开始时，汤普森采样“探索”了相当长的时间。随着它越来越确信一种药物比另一种药物更好，它将开始大多数时候选择那种药物。最终，它将收敛到知道哪种药物最好，并且它将始终选择那种药物。

夜视功能

原文：chrispiech.github.io/probabilityForComputerScientists/en/examples/night_sight/

---

在这个问题中，我们探讨如何使用概率论在黑暗中拍照。数码相机有一个传感器，在照片拍摄的持续时间中捕获光子以生成图片。然而，这些传感器容易受到“噪声”的影响，即撞击镜头的光子数量的随机波动。在问题的范围内，我们只考虑单个像素。噪声光子到达表面的到达率是恒定的。

左：使用标准照片捕获的照片。右：使用快门连拍的同张照片 [1].

对于噪声，标准差是关键！为什么？因为如果相机可以计算预期的噪声量，它可以直接减去它。但是围绕平均值的波动（以标准差衡量）会导致相机无法简单减去的测量变化。

第一部分：标准照片

首先，让我们计算以标准方式拍照时的噪声量。如果照片拍摄的持续时间是 1000 𝜇s，那么像素在单次拍照中捕获的光子数量的标准差是多少？请注意，噪声光子以每微秒 10 光子的速率（𝜇s）落在特定的像素上。

标准照片中的噪声：正如你可能猜到的，因为照片以恒定的速率撞击相机，并且彼此独立，所以撞击任何像素的噪声照片数量被建模为泊松分布！对于给定的噪声率，设\(X\)为撞击像素的噪声照片数量：$$X \sim \Poi(\lambda = 10,000).$$ 注意，10,000 是 1000𝜇s 内撞击的平均光子数（持续时间以微秒为单位乘以每微秒的光子数）。泊松的标准差简单地等于其参数的平方根，\(\sqrt{\lambda}\)。因此，捕获的噪声光子的标准差为 100（相当高）。

第二部分：快门拍摄

为了减轻噪声，斯坦福大学的毕业生意识到你可以进行快门拍摄（连续快速拍摄多张照片）并累加捕获的光子数量。由于手机相机的限制，相机在 1000μs 内最多可以拍摄 15 张照片，每张照片的持续时间为 66μs。如果我们对 15 张照片的快门拍摄进行平均，光子数量的标准差是多少？

快门拍摄中的噪声：

设 \(Y\) 为单个像素捕获的 15 张照片中光子散粒噪声的平均数量。我们想要计算 \(\text{Var}(Y)\)。具体来说，\(Y = \frac{1}{15}\sum_{i=1}^{15} X_i\)，其中 \(X_i\) 是第 i 张照片中的散粒噪声光子数量。类似于前一部分：$$X_i \sim \Poi(\lambda = 66 \cdot 10)$$ 由于 \(X_i\) 是泊松分布，\(\E[X_i] = 660\) 和 \(\Var(X_i) = 660\)。

由于 \(Y\) 是独立同分布随机变量的平均值，中心极限定理将会发挥作用。此外，根据中心极限定理的规则，\(Y\) 的方差将等于 \({1}/{n}\cdot \Var(X_i)\)。 $$\begin{align} \Var(Y) &= {1}/{n}\cdot \Var(X_i) \ &= 1/15 \cdot 660 = 44 \end{align}$$ 因此，标准差将是这个方差的平方根 \(\Std(Y) = \sqrt{44}\)，这大约是 6.6。这是一个巨大的散粒噪声减少！

---

问题由 Will Song 和 Chris Piech 提出。由谷歌的 Night Sight 提供。

P-Hacking

原文：chrispiech.github.io/probabilityForComputerScientists/en/examples/p_hacking/

---

结果表明，科学有一个缺陷！如果你测试了许多假设，但只报告了具有最低 p 值的那个，那么你更有可能得到一个虚假的结果（一个由机会产生而不是真实模式的结果）。

回顾 p 值：p 值原本是用来表示虚假结果概率的。它是如果两个总体实际上完全相同，则在数据集中观察到至少与观察到的均值差异（或你测量的任何统计量）一样大的差异的概率。p 值 < 0.05 被认为是“统计上显著的”。在课堂上，我们比较了两个总体的样本均值并计算了 p 值。如果我们有 5 个总体并寻找具有显著 p 值的成对比较，这被称为 p-hacking！

为了探索这个想法，我们将寻找一个完全随机的数据集中的模式——每个值都是均匀分布（0,1），并且与其他任何值独立。在这个玩具数据集中，任何均值差异都没有明显的意义。然而，我们可能只是偶然地发现一个看似统计上显著的结果。以下是一个包含 5 个随机总体、每个总体有 20 个样本的模拟数据集的例子：

上表中的数字只是为了演示目的。你不应该基于这些数字得出结论。我们称每个总体为随机总体，以强调其中没有模式。

有许多比较

你有多少种方法可以从一组五个总体中选择一对进行比较？总体内元素的价值以及成对的顺序都不重要。

\[\binom{5}{2} \]

理解从均匀分布（0,1）中抽取 20 个样本的均值分布的近似

均匀分布（0, 1）的方差是多少？

设 \(Z \sim \Uni(0, 1)\) $$\begin{align} \Var(Z) &= \frac{1}{12}(\beta - \alpha) \ &=\frac{1}{12} (1-0)\ &= \frac{1}{12} \end{align}$$

从均匀分布（0,1）中抽取 20 个样本的均值分布的近似是什么？

设 \(Z_1...Z_n\) 是独立同分布的 \(\Uni(0,1)\)。设 \(\bar{X} = \frac{1}{n} \sum_{i=1}^{n} Z_i\)。 [\E[X] = \frac{1}{n} \sum_{i=1}^{n} E[Z_i] = \frac{1}{n} \sum_{i=1}^{n} 0.5 = \frac{n}{n} 0.5 = 0.5] $$\begin{align} \Var(X) &= \Var\left(\frac{1}{n} \sum_{i=1}^{n} Z_i\right) \ &= \frac{1}{n²} \Var\left( \sum_{i=1}^{n} Z_i\right)\ &= \frac{1}{n²} \sum_{i=1}^{n} \Var\left(Z_i\right) \ &= \frac{1}{n²} \sum_{i=1}^{n} v \ &= \frac{n}{n²} v = \frac{v}{n} = \frac{v}{20} = \frac{1}{240} \end{align}$$ 使用中心极限定理，\(\bar{X}\sim N\left(\mu = 0.5 , \sigma² = \frac{1}{240}\right)\)

从一个总体减去另一个总体的均值分布的近似是什么？注意：如果第一个总体的均值小于第二个总体，则此值可能为负。

设 \(X_1\) 和 \(X_2\) 为总体均值。

\(X_1\sim N(\mu = 0.5 , \sigma² = \frac{1}{240})\)

\(X_2\sim N(\mu = 0.5 , \sigma² = \frac{1}{240})\) 期望值计算简单，因为 $$E[X_1 - X_2] = E[X_1] - E[X_2] = 0$$ $$\begin{align} \Var(X_1 - X_2) &= \Var(X_1) + \Var(X_2) \ &= \frac{1}{120} \end{align}$$ 独立正态分布的和（或差）仍然是正态分布：\fbox{\(Y \sim N(\mu = 0, \sigma² = \frac{v}{10})\)}(8 分) 如果只有两个种群，那么最小的均值差异\(k\)，使得差异在统计上显著是多少？换句话说，观察到均值差异为\(k\)或更大的概率小于 0.05。这个问题的一个棘手之处在于要认识到距离的双向性。如果\(P(Y<-k)\)或\(P(Y > k)\)，我们会认为这是一个显著的距离。 $$\begin{align} P(Y < -k) + P(Y > k) &= 0.05 \ F_Y(-k) + (1 - F_Y(k)) &= 0.05 \ (1-F_Y(k)) + (1 - F_Y(k)) &= 0.05 \ 2 - 2F_Y(k) &= 0.05 \ F_Y(k) &= 0.975 \end{align}$$ 现在我们需要求逆\(\Phi\)来得到\(k\)的值。 $$\begin{align} 0.975 &= \Phi\Big(\frac{k - 0}{\sqrt{v/10}}\Big) \ \Phi^{-1}(0.975) &= \frac{k}{\sqrt{v/10}} \ k &= \Phi^{-1}(0.975)\sqrt{v/10} \end{align}$$(5 分) 给出一个表达式，表示 5 个随机种群中最小的样本均值小于 0.2 的概率。设\(X_i\)为种群\(i\)的样本均值。 $$\begin{align} P(min{X_1 ... X_n} < 0.2) &= P\left(\bigcup_{i=1}^{5} X_i < 0.2\right) \ &= 1 - P\left(\left(\bigcup_{i=1}^{5} X_i < 0.2\right)^{\complement}\right) \ &= 1 - P\left(\bigcap_{i=1}^{5} X_i \geq 0.2\right) \ &= 1 - \prod_{i=1}⁵ P(X_i \ge 0.2) \ &= 1 - \prod_{i=1}⁵ 1 - \Phi\Big(\frac{0.2 - 0.5}{\sqrt{v/20}}\Big) \ \end{align}$$(7 分) 使用以下函数编写代码，估计在 5 个种群中找到均值差异的概率，这种差异被认为是显著的（使用设计用于比较两个种群的 bootstrap 方法）。至少运行 10,000 次模拟来估计你的答案。你可以使用以下辅助函数。

# the smallest difference in means that would look statistically significant
k = calculate_k()

# create a matrix with n_rows by n_cols elements, each of which is Uni(0, 1)
matrix = random_matrix(n_rows, n_cols)

# from the matrix, return the column (as a list) which has the smallest mean
min_mean_col = get_min_mean_col(matrix)

# from the matrix, return the row (as a list) which has the largest mean
max_mean_col = get_max_mean_col(matrix)

# calculate the p-value between two lists using bootstrapping (like in pset5)
p_value = bootstrap(list1, list2)

编写伪代码：

n_significant = 0
k = calculate_k()
for i in range(N_TRIALS): 
    dataset = random_matrix(20, 5)
    col_max = get_max_mean_col(dataset)
    col_min = get_min_mean_col(dataset)}
    diff = np.mean(col_max) - np.mean(col_min)}
    if diff >= k: 
        n_significant += 1}

print(n_significant / N_TRIALS)

差分隐私

原文：chrispiech.github.io/probabilityForComputerScientists/en/examples/differential_privacy/

---

最近，许多组织已经发布了在大型数据集上训练的机器学习模型（GPT-3、YOLO 等...）。这对科学是一个巨大的贡献，并简化了现代人工智能研究。然而，公布这些模型可能会允许对模型进行潜在的“逆向工程”，以揭示模型的训练数据。具体来说，攻击者可以下载一个模型，查看参数值，然后尝试重建原始训练数据。这对于在敏感数据（如健康信息）上训练的模型来说尤其糟糕。在本节中，我们将使用随机性作为防御算法“逆向工程”的方法。

注入随机性

一种对抗算法逆向工程的方法是在已经存在的数据集中添加一些随机元素。设

\[\begin{align*} X_1, \dots, X_i \overset{\text{i.i.d}} \sim \text{Bern}(p) \end{align*} \]

代表一组真实的人类数据。考虑以下代码片段：

def calculateXi(Xi): 
	return Xi

简单来说，攻击者可以针对所有 100 个样本调用上述操作，揭示所有 100 个数据点。相反，我们可以注入一些随机性：

def calculateYi(Xi): 
	obfuscate = random() # Bern with parameter p=0.5
	if obfuscate:
		return indicator(random())
	else:
		return Xi

攻击者可以期望调用新函数 100 次，并得到其中 50 个的正确值（但他们不知道是哪 50 个）。

恢复 \(p\)

现在考虑如果我们发布 calculateYi 函数，一个对样本均值感兴趣的研究者如何获取有用的数据？他们可以查看：

\[\begin{align*} Z = \sum^{100}_{n=1} Y_i. \end{align*} \]

其期望值为：

\[\begin{align*} E[Z] = E\left[\sum^{100}_{n=1} Y_i\right] = \sum^{100}_{n=1} E[Y_i] = \sum^{100}_{n=1}\left(\frac{p}{2} + \frac{1}{4}\right) = 50p + 25 \end{align*} \]

然后为了揭示一个估计值，科学家可以这样做，

\[\begin{align*} p \approx \frac{Z - 25}{50} \end{align*} \]

然后继续进行更多研究！

第五部分：机器学习

参数估计

原文：chrispiech.github.io/probabilityForComputerScientists/en/part5/parameter_estimation/

---

我们已经学习了多种随机变量的分布，所有这些分布都有参数：当你定义随机变量时作为输入提供的数字。到目前为止，当我们处理随机变量时，要么明确告知参数的值，要么通过理解生成随机变量的过程来推断这些值。

如果我们不知道参数的值，也无法从自己的专业知识中估计它们，那会怎样？如果我们知道的不是随机变量，而是大量使用相同潜在分布生成的数据示例，那又会怎样？在本章中，我们将学习从数据中估计参数的正式方法。

这些思想对于人工智能至关重要。几乎所有现代机器学习算法都是这样工作的：(1)指定一个具有参数的概率模型。(2)从数据中学习这些参数的值。

参数

在我们深入参数估计之前，首先让我们回顾一下参数的概念。给定一个模型，参数是产生实际分布的数字。在伯努利随机变量的情况下，单个参数是值 \(p\)。在均匀随机变量的情况下，参数是定义最小值和最大值的 \(a\) 和 \(b\) 值。以下是随机变量及其对应参数的列表。从现在开始，我们将使用符号 \(\theta\) 来表示所有参数的向量：

Distribution	Parameters
Bernoulli(\(p\))	\(\theta = p\)
Poisson(\(\lambda\))	\(\theta = \lambda\)
Uniform(\(a, b\))	\(\theta = [a, b]\)
Normal(\(\mu, \sigma²\))	\(\theta = [\mu, \sigma²]\)

在现实世界中，通常你不知道“真实”的参数，但你能够观察到数据。接下来，我们将探讨如何使用数据来估计模型参数。

估计参数值的方法不止一种。主要有两种思想流派：最大似然估计（MLE）和最大后验估计（MAP）。这两种思想流派都假设你的数据是独立同分布（IID）的样本：\(X_1, X_2, \dots X_n\)。

最大似然估计

原文：chrispiech.github.io/probabilityForComputerScientists/en/part5/mle/

---

我们用于估计参数的第一个算法被称为最大似然估计（MLE）。MLE 背后的核心思想是选择那些使观察到的数据最可能的参数（\(\theta\)）。

我们将要用来估计参数的数据将是 \(n\) 个独立同分布（IID）的样本：\(X_1, X_2, \dots X_n\)。

似然

我们假设我们的数据是同分布的。这意味着它们必须有相同的概率质量函数（如果数据是离散的）或相同的概率密度函数（如果数据是连续的）。为了简化我们关于参数估计的讨论，我们将使用符号 \(f(X=x|\Theta = \theta)\) 来指代这个共享的 PMF 或 PDF。我们的新符号有两个有趣之处。首先，我们现在已经包含了关于 \(\theta\) 的条件，这是我们表示不同值 \(X\) 的似然依赖于我们的参数值的方式。其次，我们将使用相同的符号 \(f\) 来表示离散和连续分布。

似然是什么意思，它与“概率”有何不同？在离散分布的情况下，似然是数据概率质量或联合概率质量的同义词。在连续分布的情况下，似然指的是数据的概率密度。

由于我们假设每个数据点是独立的，因此我们所有数据的似然是每个数据点似然的乘积。从数学上讲，我们的数据给出参数 \(\theta\) 的似然是：$$\begin{align} L(\theta) = \prod_{i=1}^n f(X_i = x_i|\Theta = \theta) \end{align}$$

对于不同的参数值，我们数据的似然将不同。如果我们有正确的参数，我们的数据将比我们有错误的参数时更有可能。因此，我们将似然写成参数（\(\theta\)）的函数。

最大化

在最大似然估计（MLE）中，我们的目标是选择参数（\(\theta\)）的值，以最大化上一节中的似然函数。我们将使用符号 \(\hat{\theta}\) 来表示参数的最佳选择值。形式上，MLE 假设：$$\begin{align} \hat{\theta} = \underset{\theta}{\operatorname{argmax }} \text{ }L(\theta) \end{align}$$

Argmax 是“最大值之论证”的缩写。函数的 argmax 是函数在该域上最大化的值。这适用于任何维度的域。

argmax 的一个酷特性是，由于对数是一个单调函数，函数的 argmax 与其对数的 argmax 是相同的！这很好，因为对数使数学更简单。如果我们找到似然的对数的 argmax，它将等于似然的 argmax。因此，对于 MLE，我们首先写出对数似然函数（\(LL\)） $$\begin{align} LL(\theta) &= \log L(\theta) \ &= \log \prod_{i=1}^n f(X_i=x_i|\Theta = \theta) \ &= \sum_{i=1}^n \log f(X_i = x_i|\Theta = \theta) \end{align}$$

要使用最大似然估计器，首先写出给定参数的数据的对数似然。然后选择使对数似然函数最大化的参数值。argmax 可以通过多种方式计算。我们在这个课程中涵盖的所有方法都需要计算函数的第一导数。

二项式 MLE 估计

在我们的第一个例子中，我们将使用 MLE 来估计伯努利分布的 \(p\) 参数。我们将基于 \(n\) 个数据点进行估计，我们将这些数据点称为独立同分布（IID）随机变量 \(X_1, X_2, \dots X_n\)。这些随机变量中的每一个都被假定为来自相同的伯努利分布，具有相同的 \(p\)，\(X_i \sim \text{Ber}(p)\)。我们想找出那个 \(p\) 是多少。

MLE 的第一步是将伯努利分布的似然写成我们可以最大化的函数。由于伯努利是一个离散分布，似然是概率质量函数。

伯努利 \(X\) 的概率质量函数可以写成 \(f(x) = p^{x}(1-p)^{1-x}\)。哇！这是怎么回事？这是一个方程，它允许我们说 \(X = 1\) 的概率是 \(p\)，而 \(X = 0\) 的概率是 \(1- p\)。请你自己证明当 \(X_i=0\) 和 \(X_i=1\) 时，PMF 返回正确的概率。我们这样写 PMF 是因为它可导。

现在我们来进行一些 MLE 估计： $$\begin{align} L(\theta) &= \prod_{i=1}^n p^{x_i}(1-p) && \text{首先写出似然函数} \ LL(\theta) &= \sum_{i=1}^n \log p^{x_i}(1-p) && \text{然后写出对数似然函数}\ &= \sum_{i=1}^n x_i (\log p) + (1 - x_i) log(1-p) \ &= Y \log p + (n - Y) \log(1-p) && \text{其中 } Y = \sum_{i=1}^n x_i \end{align}$$

天哪！我们得到了对数似然方程。现在我们只需要选择使我们的对数似然最大化的 \(p\) 值。正如你的微积分老师可能教你的那样，找到使一个函数最大化的值的一种方法就是找到该函数的第一导数并将其设为 0。 $$\begin{align} \frac{\delta LL(p)}{\delta p} &= Y \frac{1}{p} + (n - Y) \frac{-1}{1-p} = 0 \ \hat{p} &= \frac{Y}{n} = \frac{\sum_{i=1}^n x_i}{n} \end{align}$$

所有这些工作，结果发现最大似然估计（MLE）仅仅是样本均值...

正态 MLE 估计

实践是关键。接下来，我们将尝试估计正态分布的最佳参数值。我们所能访问的只有来自正态分布的 \(n\) 个样本，我们将其称为独立同分布的随机变量 \(X_1, X_2, \dots X_n\)。我们假设对于所有 \(i\)，\(X_i \sim N(\mu = \theta_0, \sigma² = \theta_1)\)。这个例子似乎更复杂，因为正态分布有两个参数需要估计。在这种情况下，\(\theta\) 是一个包含两个值的向量，第一个是均值（\(\mu\)）参数。第二个是方差（\(\sigma²\)）参数。 $$\begin{align} L(\theta) &= \prod_{i=1}^n f(X_i|\theta) \ &=\prod_{i=1}^n \frac{1}{\sqrt{2\pi\theta_1}} e^{-\frac{(x_i - \theta_0)²}{2\theta_1}} && \text{连续变量的似然是概率密度函数}\ LL(\theta) &= \sum_{i=1}^n \log \frac{1}{\sqrt{2\pi\theta_1}} e^{-\frac{(x_i - \theta_0)²}{2\theta_1}} && \text{我们想要计算对数似然} \ &= \sum_{i=1}^n\left[ - \log(\sqrt{2\pi\theta_1}) - \frac{1}{2\theta_1}(x_i - \theta_0)² \right] \end{align}$$

再次，MLE 的最后一步是选择最大化对数似然函数的 \(\theta\) 的值。在这种情况下，我们可以计算 \(LL\) 函数相对于 \(\theta_0\) 和 \(\theta_1\) 的偏导数，将两个方程都设为等于 0，然后求解 \(\theta\) 的值。这样做会导致最大化似然的 \(\hat{\mu} = \hat{\theta_0}\) 和 \(\hat{\sigma²} = \hat{\theta_1}\) 的值。结果是：\(\hat{\mu} = \frac{1}{n}\sum_{i=1}^n x_i\) 和 \(\hat{\sigma²} = \frac{1}{n}\sum_{i=1}^n (x_i - \hat{\mu})²\)。

线性变换加噪声

最大似然估计（MLE）是一种可以用于任何具有可导似然函数的概率模型的算法。以一个例子来说，让我们估计一个模型中的参数 \(\theta\)，其中存在一个随机变量 \(Y\)，使得 \(Y = \theta X + Z\)，\(Z \sim N(0, \sigma²)\)，而 \(X\) 是一个未知的分布。

在你被告知 \(X\) 的值的情况下，\(\theta X\) 是一个数字，而 \(\theta X + Z\) 是高斯分布和一个数字的和。这意味着 \(Y|X \sim N(\theta X, \sigma²)\)。我们的目标是选择一个 \(\theta\) 的值，以最大化独立同分布的概率：\((X_1, Y_1), (X_2, Y_2), \dots (X_n, Y_n)\)。

我们通过首先找到给定 \(\theta\) 的数据的对数似然函数来解决这个问题。然后我们找到最大化对数似然函数的 \(\theta\) 的值。首先，使用正态分布的概率密度函数来表达 \(Y|X,\theta\) 的概率：$$\begin{align} f(Y_i | X_i , \theta) = \frac{1}{\sqrt{2\pi} \sigma} e^{-\frac{(Y_i - \theta X_i)²}{2\sigma²}} \end{align}$$

现在我们已经准备好编写似然函数，然后取其对数得到对数似然函数：$$\begin{align} L(\theta) &= \prod_{i=1}^n f(Y_i, X_i | \theta) && \text{让我们分解这个联合函数}\ &= \prod_{i=1}^n f(Y_i | X_i, \theta)f(X_i ) && f(X_i) \text{ 与 } \theta \text{ 独立}\ &= \prod_{i=1}^n\frac{1}{\sqrt{2\pi} \sigma} e^{-\frac{(Y_i - \theta X_i)²}{2\sigma²}}f(X_i) &&\text{代入 }f(Y_i | X_i) \text{ 的定义} \end{align}$$ $$\begin{align} LL(\theta) &= \log L(\theta) \ &= \log \prod_{i=1}^n \frac{1}{\sqrt{2\pi} \sigma} e^{-\frac{(Y_i - \theta X_i)²}{2\sigma²}} f(X_i) &&\text{代入 }L(\theta)\ &= \sum_{i=1}^n \log \frac{1}{\sqrt{2\pi} \sigma} e^{-\frac{(Y_i - \theta X_i)²}{2\sigma²}} + \sum_{i=1}^n \log f(X_i) &&\text{乘积的对数是各对数之和}\ &=n \log \frac{1}{\sqrt{2\pi}\sigma} - \frac{1}{2\sigma²} \sum_{i=1}^n (Y_i - \theta X_i)² + \sum_{i=1}^n \log f(X_i ) \end{align}$$

移除常数乘数和不含 \(\theta\) 的项。我们剩下的是尝试找到一个 \(\theta\) 的值，使得最大化：$$\begin{align} \hat{\theta} &= \underset{\theta}{\operatorname{argmax}} - \sum_{i=1}^m (Y_i - \theta X_i)²\ &= \underset{\theta}{\operatorname{argmin}} \sum_{i=1}^m (Y_i - \theta X_i)² \end{align}$$

这个结果说明，使数据最可能的 \(\theta\) 值是使 \(Y\) 的预测平方误差最小化的值。我们将在几天后看到，这是线性回归的基础。

最大后验概率

chrispiech.github.io/probabilityForComputerScientists/en/part5/map/

---

MLE（最大似然估计）很棒，但它并不是估计参数的唯一方法！本节介绍另一种算法，最大后验概率（Maximum A Posteriori，MAP）。MAP 的范式是我们应该选择最有可能的数据对应的参数值。乍一看，这似乎与 MLE 相同，然而请注意，MLE 选择的是使数据最有可能的参数值。形式上，对于独立同分布（IID）的随机变量 \(X_1, \dots, X_n\)：$$\begin{align} \theta_{\text{MAP}} =& \underset{\theta}{\operatorname{argmax }} \text{ } f(\Theta = \theta | X_1=x_1, X_2=x_2, \dots X_n=x_n) \end{align}$$在上面的方程中，我们试图计算给定观测随机变量的未观测随机变量的条件概率。当这种情况发生时，想想贝叶斯定理！使用贝叶斯定理的连续版本展开函数 \(f\)。$$\begin{align} \theta_{\text{MAP}} =& \underset{\theta}{\operatorname{argmax }} \text{ } f(\Theta = \theta | X_1=x_1, X_2=x_2, \dots X_n=x_n) \ =& \underset{\theta}{\operatorname{argmax }} \text{ } \frac{f(X_1=x_1, X_2=x_2, \dots, X_n=x_n |\Theta = \theta) f(\Theta = \theta)}{f(X_1=x_1, X_2=x_2, \dots X_n=x_n)} && \text{贝叶斯定理} \end{align}$$请注意，\(f\) 都是概率密度函数或概率质量函数。现在我们将利用两个观察结果。首先，假设数据是 IID 的，因此我们可以分解给定 \(\theta\) 的数据密度。其次，分母相对于 \(\theta\) 是常数。因此，它的值不影响 argmax，我们可以省略这一项。数学上：$$\begin{align} \theta_{\text{MAP}} =& \underset{\theta}{\operatorname{argmax }} \text{ } \frac{f(\Theta = \theta) \cdot \prod_{i=1}^n f(X_i =x_i| \Theta = \theta) }{f(X_1 = x_1, X_2 = x_2, \dots X_n = x_n)} && \text{由于样本是 IID 的}\ =& \underset{\theta}{\operatorname{argmax }} \text{ } f(\Theta = \theta) \cdot\prod_{i=1}^n f(X_i =x_i |\Theta = \theta) && \text{分母相对于 \(\theta\) 是常数} \end{align}$$与之前一样，找到 MAP 函数的对数的 argmax 将更方便，这给出了参数 MAP 估计的最终形式。 $$\begin{align} \theta_{\text{MAP}} =& \underset{\theta}{\operatorname{argmax }} \text{ } \left( \log f(\Theta = \theta) + \sum_{i=1}^n \log(f(X_i =x_i|\Theta = \theta)) \right) \end{align}$$使用贝叶斯术语，MAP 估计是 \(\theta\) 的“后验”分布的众数。如果你将这个方程与 MLE 方程并排放置，你会注意到 MAP 是相同函数的 argmax 加上一个对数先验项。

参数先验

为了准备好 MAP 估计的世界，我们需要复习我们的分布。我们需要为我们的不同参数中的每一个找到一个合理的分布。例如，如果你正在预测泊松分布，\(\lambda\) 的先验分布应该选择哪种合适的随机变量类型？

对于先验分布的一个期望是，得到的后验分布具有相同的函数形式。我们称这些为“共轭”先验。在你多次更新你的信念的情况下，共轭先验使得在数学方程中进行编程变得容易得多。

这里是一个不同参数及其通常用于其先验分布的分布列表：$$\begin{align} &\text{参数 }&& \text{分布}\ &\text{伯努利 } p && \text{Beta}\ &\text{二项式 } p && \text{Beta}\ &\text{泊松 } \lambda && \text{Gamma}\ &\text{指数 } \lambda && \text{Gamma}\ &\text{多项式 } p_i && \text{Dirichlet}\ &\text{正态 } \mu && \text{Normal}\ &\text{正态 } \sigma² && \text{Inverse Gamma} \end{align}$$ 你只需要对新的分布有一个高层次的理解。你不需要了解逆伽马分布。我包括它是为了完整性。

用于表示你对随机变量“先验”信念的分布通常会具有自己的参数。例如，Beta 分布使用两个参数 \((a, b)\) 定义。我们是否必须使用参数估计来评估 \(a\) 和 \(b\)？不。这些参数被称为“超参数”。这是一个我们保留用于在运行参数估计之前固定的模型参数的术语。在你运行 MAP 之前，你将决定 \((a, b)\) 的值。

Dirichlet

Dirichlet 分布以与多项式分布推广伯努利分布相同的方式推广了 Beta 分布。一个服从 Dirichlet 分布的随机变量 \(X\) 被参数化为 \(X \sim \text{Dirichlet}(a_1, a_2, \dots, a_m)\)。该分布的 PDF 为：$$\begin{align} f(X_1 = x_1, X_2 = x_2, \dots, X_m = x_m) = K \prod_{i=1}^m x_i^{a_i - 1} \end{align}$$ 其中 \(K\) 是归一化常数。

你可以直观地理解 Dirichlet 分布的超参数：想象你已经看到了 \(\sum_{i=1}^m a_i - m\) 个虚拟试验。在这些试验中，你有 \((a_i - 1)\) 个值为 \(i\) 的结果。作为一个例子，考虑估计一个六面偏斜骰子（每个面是不同形状）得到不同数字的概率。我们将通过反复掷骰子 \(n\) 次来估计掷出每个面的概率。这将产生 \(n\) 个独立同分布的样本。对于 MAP 模式，我们需要为每个参数 \(p_1 \dots p_6\) 的信念设定一个先验。我们希望表达出我们轻信每次掷骰子的概率是相等的。

在你掷骰子之前，让我们想象你已经掷了六次骰子，并且得到了每种可能值的一个。因此，“先验”分布将是 Dirichlet(\(2, 2, 2, 2, 2, 2\))。在观察到 \(n_1 + n_2 + \dots + n_6\) 个新的试验结果，其中 \(n_i\) 是结果 \(i\) 的次数后，“后验”分布是 Dirichlet(\(2 + n_1, \dots 2 + n_6\))。使用代表每个结果的一个想象观察的先验称为“拉普拉斯平滑”，并且它保证了你的概率不会是 0 或 1。

Gamma

Gamma(\(k, \theta\)) 分布是泊松分布 \(\lambda\) 参数的共轭先验分布（它也是指数分布的共轭分布，但我们将不会深入探讨这一点）。

超参数可以被解释为：你在 \(\theta\) 个想象的时间段内观察到了 \(k\) 个总想象事件。在接下来的 \(t\) 个时间段内观察到 \(n\) 个事件后，后验分布是 Gamma(\(k + n, \theta + t\))。

例如，Gamma(10, 5) 将代表在 5 个时间段内观察到 10 个想象事件。这就像想象一个有某种置信度的 2 的速率。如果我们以那个 Gamma 作为先验，然后在接下来的 2 个时间段内观察到 11 个事件，我们的后验是 Gamma(21,7)，这相当于一个更新的速率为 3。

机器学习

原文：chrispiech.github.io/probabilityForComputerScientists/en/part5/machine_learning/

---

机器学习是计算机科学的一个子领域，它使计算机能够在没有明确编程的情况下执行任务。机器学习领域包含多个不同的任务，以及多种不同的“学习”算法。在本章中，我们将重点关注分类以及两种经典的分类算法：朴素贝叶斯和逻辑回归。

分类

在分类任务中，你的任务是使用带有特征/标签对 (\(\mathbf{x}\), y) 的训练数据来估计一个函数 \(\hat{y} = g(\mathbf{x})\)。这个函数可以用来进行预测。在分类中，\(y\) 的值取自一个 \textbf{离散} 的数值集合。因此，我们通常选择 \(g(\mathbf{x}) = \underset{y}{\operatorname{argmax }}\text{ }\hat{P}(Y=y|\mathbf{X} = \mathbf{x})\)。

在分类任务中，你将得到 \(N\) 个训练对：\((\mathbf{x}^{(1)},y^{(1)}), (\mathbf{x}^{(2)},y^{(2)}), \dots , (\mathbf{x}^{(N)},y^{(N)})\) 其中 \(\mathbf{x}^{(i)}\) 是第 \(i\) 个训练示例的 \(m\) 个离散特征向量，而 \(y^{(i)}\) 是第 \(i\) 个训练示例的离散标签。

在我们介绍机器学习时，我们将假设训练数据集中的所有值都是二元的。虽然这不是一个必要的假设，但它使得学习核心概念变得容易得多。具体来说，我们假设所有标签都是二元的 \(y^{(i)} \in \{0, 1\} \text{ }\forall i\)，并且所有特征都是二元的 \(x^{(i)}_j \in \{0, 1\} \text{ }\forall i, j\)。

Naïve Bayes

原文：chrispiech.github.io/probabilityForComputerScientists/en/part5/naive_bayes/

---

Naive Bayes 是一种用于“分类任务”的机器学习算法。它做出了一个实质性的假设（称为朴素贝叶斯假设），即所有特征在给定分类标签的情况下彼此独立。这个假设是错误的，但它允许算法快速且高效，通常很有用。为了实现朴素贝叶斯，你需要学习如何训练你的模型以及如何使用它进行预测，一旦模型训练完成。

训练（又称参数估计）

训练的目标是估计所有 \(0 < i \leq m\) 个特征的概率 \(P(Y)\) 和 \(P(X_i | Y)\)。我们使用符号 \(\hat{p}\) 来明确表示概率是一个估计值。

使用 MLE 估计：$$\begin{align} \hat{p}(X_i = x_i | Y = y) = \frac{ (\text{训练样本中 \(X_i = x_i\) 且 \(Y = y\) 的数量})}{(\text{训练样本中 \(Y = y\) 的数量})} \end{align}$$

使用拉普拉斯 MAP 估计：$$\begin{align} \hat{p}(X_i = x_i | Y = y) = \frac{ (\text{训练样本中 \(X_i = x_i\) 且 \(Y = y\) 的数量}) + 1 }{(\text{训练样本中 \(Y = y\) 的数量}) + 2} \end{align}$$

使用 MLE 估计训练 \(Y\) 的先验概率：$$\begin{align} \hat{p}(Y = y) = \frac{ (\text{训练样本中 \(Y = y\) 的数量})}{(\text{训练样本数量})} \end{align}$$

预测

对于 \(\mathbf{x} = [x_1, x_2, \dots , x_m]\) 的例子，估计 \(y\) 的值如下：$$\begin{align} \hat{y} &= \argmax_{y = {0, 1}} \text{ } \log \hat{p}(Y = y) + \sum_{i=1}^m \log \hat{p}(X_i = x_i | Y = y) \end{align}$$ 注意，对于足够小的数据集，你可能不需要使用 argmax 的对数版本。

理论

在分类的世界里，当我们做出预测时，我们希望选择 \(y\) 的值，使其最大化 \(P(Y=y|\mathbf{X})\)。 $$\begin{align} \hat{y} &= \argmax_{y = {0, 1}} P(Y = y|\mathbf{X} = \mathbf{X}) && \text{我们的目标}\ &= \argmax_{y = {0, 1}} \frac{P(Y=y)P(\mathbf{X} =\mathbf{x}| Y = y)}{P(\mathbf{X} =\mathbf{x})} && \text{根据贝叶斯定理}\ &= \argmax_{y = {0, 1}} P(Y=y)P(\mathbf{X} =\mathbf{x}| Y = y)) && \text{因为 \(P(\mathbf{X} =\mathbf{x})\) 与 \(Y\) 无关} \end{align}$$

使用我们的训练数据，我们可以将 \(\mathbf{X}\) 和 \(Y\) 的联合分布解释为一个巨大的多项式分布，每个 \(\mathbf{X}=\mathbf{x}\) 和 \(Y=y\) 的组合都有一个不同的参数。例如，如果输入向量只有长度为 1。换句话说 \(|\mathbf{x}| = 1\)，\(x\) 和 \(y\) 可以取的值数量很少，比如二进制，这是一个完全合理的做法。我们可以使用 MLE 或 MAP 估计器来估计多项式分布，然后在我们的表中进行几次查找来计算 argmax。

当特征数量变得很大时，糟糕的情况就会出现。回想一下，我们的多项式分布需要为向量 \(\mathbf{x}\) 和值 \(y\) 的每个唯一组合估计一个参数。如果存在 \(|\mathbf{x}| = n\) 个二进制特征，那么这种策略将需要 \(\mathcal{O}(2^n)\) 的空间，并且可能有许多参数在没有与相应赋值匹配的训练数据的情况下被估计。

朴素贝叶斯假设

朴素贝叶斯假设是，给定 \(y\)，\(\mathbf{x}\) 的每个特征都是相互独立的。朴素贝叶斯假设是错误的，但很有用。这个假设允许我们使用与特征大小成线性关系的空间和数据来做出预测：如果 \(|\mathbf{x}| = n\)，则为 \(\mathcal{O}(n)\)。这使我们能够训练和预测具有巨大特征空间的数据，例如具有互联网上每个单词指示符的特征空间。使用这个假设，预测算法可以被简化。 $$\begin{align} \hat{y} &= \argmax\limits_{y = {0, 1}} \text{ }P(Y = y)P(\mathbf{X} = \mathbf{x}| Y = y) && \text{正如我们上次留下的}\ &= \argmax\limits_{y = {0, 1}} \text{ } P(Y= y)\prod_i P(X_i = x_i| Y = y) &&\text{朴素贝叶斯假设}\ &= \argmax\limits_{y = {0, 1}} \text{ } \log P(Y = y) + \sum_i \log P(X_i = x_i| Y = y) && \text{为了数值稳定性}\ \end{align}$$

在最后一步，我们利用了函数的 argmax 等于函数的 log 的 argmax 的性质。这个算法在训练和预测时都既快又稳定。

逻辑回归

原文：chrispiech.github.io/probabilityForComputerScientists/en/part5/log_regression/

---

逻辑回归是一种分类算法（我知道，这个名字很糟糕。或许逻辑分类会更好一些），它通过尝试学习一个近似 \(\P(y|x)\) 的函数来工作。它假设 \(\P(y|x)\) 可以近似为应用在输入特征线性组合上的 sigmoid 函数。学习逻辑回归非常重要，因为逻辑回归是人工神经网络的基本构建块。

从数学上讲，对于一个单个训练数据点（\(\mathbf{x}, y\)），逻辑回归假设：$$\begin{align} P(Y=1|\mathbf{X}=\mathbf{x}) &= \sigma(z) \text{ 其中 } z = \theta_0 + \sum_{i=1}^m \theta_i x_i \end{align}$$ 这个假设通常以等价的形式写出：$$\begin{align} P(Y=1|\mathbf{X}=\mathbf{x}) &=\sigma(\mathbf{\theta}^T\mathbf{x}) &&\text{ 其中我们总是将 \(x_0\) 设置为 1}\ P(Y=0|\mathbf{X}=\mathbf{x}) &=1-\sigma(\mathbf{\theta}^T\mathbf{x}) &&\text{ 根据概率的总法则} \end{align}$$ 使用这些方程来表示 \(Y|X\) 的概率，我们可以创建一个算法，该算法选择 \(\theta\) 的值，以最大化所有数据的概率。我首先将陈述对数概率函数和关于 \(\theta\) 的偏导数。然后我们将在（a）展示一个可以选择 \(\theta\) 的最优值的算法，以及（b）展示这些方程是如何推导出来的。

一个重要的事实是：给定参数（\(\theta\)）的最佳值，逻辑回归通常可以很好地估计不同类别标签的概率。然而，给定不良的，甚至随机的 \(\theta\) 值，它的工作表现会很差。逻辑回归机器学习算法的“智能”程度取决于 \(\theta\) 的良好值。

符号

在我们开始之前，我想确保我们对于符号的使用是一致的。在逻辑回归中，\(\theta\) 是一个长度为 \(m\) 的参数向量，我们将根据 \(n\) 个训练样本来学习这些参数的值。参数的数量应该等于每个数据点的特征数量。

在逻辑回归中，我们经常使用的一些符号你可能不太熟悉：$$\begin{align} \mathbf{\theta}^T\mathbf{x} &= \sum_{i=1}^m \theta_i x_i = \theta_1 x_1 + \theta_2 x_2 + \dots + \theta_m x_m && \text{点积，也称为加权求和}\ \sigma(z) &= \frac{1}{1+ e^{-z}} && \text{sigmoid 函数} \end{align}$$

对数似然

为了选择逻辑回归参数的值，我们使用最大似然估计（MLE）。因此，我们将有两个步骤：（1）写出对数似然函数；（2）找到最大化对数似然函数的 \(\theta\) 值。

我们预测的标签是二元的，我们的逻辑回归函数的输出应该是标签为 1 的概率。这意味着我们可以（并且应该）将每个标签解释为一个伯努利随机变量：\(Y \sim \text{Bern}(p)\) 其中 \(p = \sigma(\theta^T \textbf{x})\)。

首先，这里有写出一个数据点的概率的超级简洁方法（回想这是伯努利概率质量函数的方程形式）：$$\begin{align} P(Y=y | X = \mathbf{x}) = \sigma({\mathbf{\theta}^{T\mathbf{x}})}y \cdot \left[1 - \sigma({\mathbf{\theta}^{T\mathbf{x}})\right]} \end{align}$$

现在我们知道了概率质量函数，我们可以写出所有数据的似然：$$\begin{align} L(\theta) =& \prod_{i=1}^n P(Y=y^{(i)} | X = \mathbf{x}^{(i)}) && \text{独立训练标签的似然}\ =& \prod_{i=1}^n \sigma({\mathbf{\theta}^T\mathbf{x}{(i)}})^{y{(i)}} \cdot \left[1 - \sigma({\mathbf{\theta}^T\mathbf{x}{(i)}})\right]^{(1-y)} && \text{代入伯努利的似然} \end{align}$$ 如果你取这个函数的对数，你得到逻辑回归报告的对数似然。对数似然方程是：$$\begin{align} LL(\theta) = \sum_{i=1}^n y^{(i)} \log \sigma(\mathbf{\theta}^T\mathbf{x}) + (1-y^{(i)}) \log [1 - \sigma(\mathbf{\theta}^T\mathbf{x})] \end{align}$$

回想一下，在最大似然估计（MLE）中，唯一剩下的步骤就是选择参数（\(\theta\)）以最大化对数似然。

对数似然的梯度

现在我们有了对数似然函数，我们只需要选择最大化它的\(\theta\)值。我们可以通过使用优化算法找到\(\theta\)的最佳值。然而，为了使用优化算法，我们首先需要知道对数似然相对于每个参数的偏导数。首先，我将给你偏导数（这样你可以看到它是如何被使用的）。然后，我将向你展示如何推导它：$$\begin{align} \frac{\partial LL(\theta)}{\partial \theta_j} = \sum_{i=1}^n \left[ y^{(i)} - \sigma(\mathbf{\theta}^T\mathbf{x}) \right] x_j^{(i)} \end{align}$$

梯度下降优化

我们的目标是选择参数（\(\theta\)）以最大化似然，并且我们知道对数似然相对于每个参数的偏导数。我们已经准备好使用我们的优化算法了。

在逻辑回归的情况下，我们不能从数学上求解\(\theta\)。相反，我们使用计算机来选择\(\theta\)。为此，我们采用一个称为梯度下降（优化理论中的经典算法）的算法。梯度下降背后的思想是，如果你持续地向下走小步（在负梯度的方向上），你最终会到达局部最小值。在我们的情况下，我们想要最大化我们的似然。正如你所想象的，最小化我们似然的负值将等同于最大化我们的似然。

导致每个小步骤的参数更新可以计算如下：$$\begin{align} \theta_j^{\text{ new}} &= \theta_j^{\text{ old}} + \eta \cdot \frac{\partial LL(\theta^{\text{ old}})}{\partial \theta_j^{\text{ old}}} \ &= \theta_j^{\text{ old}} + \eta \cdot \sum_{i=1}^n \left[ y^{(i)} - \sigma(\mathbf{\theta}^T\mathbf{x}) \right] x_j^{(i)} \end{align}$$其中\(\eta\)是我们采取的步长的大小。如果你继续使用上述方程更新\(\theta\)，你将收敛到\(\theta\)的最佳值。你现在有一个智能模型。以下是逻辑回归的梯度上升算法的伪代码：

小贴士：别忘了，为了学习\(\theta_0\)的值，你可以简单地定义\(\textbf{x}_0\)始终为 1。

导数

在本节中，我们提供了对数似然梯度的数学推导。这些推导值得了解，因为这些思想在人工神经网络中得到了广泛的应用。

我们的目标是计算对数似然相对于每个 theta 的导数。首先，这是 sigmoid 函数相对于其输入的导数的定义：$$\begin{align} \frac{\partial}{\partial z} \sigma(z) = \sigma(z)[1 - \sigma(z)] && \text{要得到相对于\(\theta\)的导数，使用链式法则} \end{align}$$花点时间来欣赏 sigmoid 函数导数的美丽。sigmoid 具有如此简单的导数的原因在于 sigmoid 分母中的自然指数。

由于似然函数是对所有数据的求和，而在微积分中，和的导数是各个导数的和，因此我们可以专注于计算一个示例的导数。theta 的梯度仅仅是每个训练数据点的这个项的总和。

首先，我将向您展示如何通过困难的方式计算导数。然后，我们将探讨一种更简单的方法。对于一个数据点 \((\mathbf{x}, y)\) 的梯度导数：$$\begin{align} \frac{\partial LL(\theta)}{\partial \theta_j} &= \frac{\partial }{\partial \theta_j} y \log \sigma(\mathbf{\theta}^T\mathbf{x}) + \frac{\partial }{\partial \theta_j} (1-y) \log [1 - \sigma(\mathbf{\theta}^T\mathbf{x}] && \text{项之和的导数}\ &=\left[\frac{y}{\sigma(\theta^T\mathbf{x})} - \frac{1-y}{1-\sigma(\theta^T\mathbf{x})} \right] \frac{\partial}{\partial \theta_j} \sigma(\theta^T \mathbf{x}) &&\text{log \(f(x)\) 的导数}\ &=\left[\frac{y}{\sigma(\theta^T\mathbf{x})} - \frac{1-y}{1-\sigma(\theta^T\mathbf{x})} \right] \sigma(\theta^T \mathbf{x}) [1 - \sigma(\theta^T \mathbf{x})]\mathbf{x}_j && \text{链式法则 + sigma 的导数}\ &=\left[ \frac{y - \sigma(\theta^{T\mathbf{x})}{\sigma(\theta}T \mathbf{x}) [1 - \sigma(\theta^T \mathbf{x})]} \right] \sigma(\theta^T \mathbf{x}) [1 - \sigma(\theta^T \mathbf{x})]\mathbf{x}_j && \text{代数操作}\ &= \left[y - \sigma(\theta^T\mathbf{x}) \right] \mathbf{x}_j && \text{消去项} \end{align}$$

《微分不再流泪》

那是艰难的方法。逻辑回归是人工神经网络的基石。如果我们想要扩展规模，我们就必须习惯于一种更简单的方式来计算导数。为此，我们必须欢迎我们老朋友链式法则的回归。根据链式法则：$$\begin{align} \frac{\partial LL(\theta)}{\partial \theta_j} &= \frac{\partial LL(\theta)}{\partial p} \cdot \frac{\partial p}{\partial \theta_j} && \text{其中 } p = \sigma(\theta^T\textbf{x})\ &= \frac{\partial LL(\theta)}{\partial p} \cdot \frac{\partial p}{\partial z} \cdot \frac{\partial z}{\partial \theta_j} && \text{其中 } z = \theta^T\textbf{x} \end{align}$$ 链式法则是微积分的分解机制。它允许我们通过将其分解成更小的部分来计算复杂的偏导数（\(\frac{\partial LL(\theta)}{\partial \theta_j}\)）。$$\begin{align} LL(\theta) &= y \log p + (1-y) \log (1 - p) && \text{其中 } p = \sigma(\theta^T\textbf{x}) \ \frac{\partial LL(\theta)}{\partial p} &= \frac{y}{p} - \frac{1-y}{1-p} && \text{通过求导得到} \end{align}$$ $$\begin{align} p &= \sigma(z) && \text{其中 }z = \theta^T\textbf{x}\ \frac{\partial p}{\partial z} &= \sigma(z)[1- \sigma(z)] && \text{通过求导得到 sigmoid 函数} \end{align}$$ $$\begin{align} z &= \theta^T\textbf{x} && \text{如前所述}\ \frac{\partial z}{\partial \theta_j} &= \textbf{x}_j && \text{只有 \(\textbf{x}_j\) 与 \(\theta_j\) 交互} \end{align}$$ 这些导数都更容易计算。现在我们只需将它们相乘。$$\begin{align} \frac{\partial LL(\theta)}{\partial \theta_j} &= \frac{\partial LL(\theta)}{\partial p} \cdot \frac{\partial p}{\partial z} \cdot \frac{\partial z}{\partial \theta_j} \ &= \Big[\frac{y}{p} - \frac{1-y}{1-p}\Big] \cdot \sigma(z)[1- \sigma(z)] \cdot \textbf{x}_j && \text{通过代入每个项} \ &= \Big[\frac{y}{p} - \frac{1-y}{1-p}\Big] \cdot p[1- p] \cdot \textbf{x}_j && \text{因为 }p = \sigma(z)\ &= [y(1-p) - p(1-y)] \cdot \textbf{x}_j && \text{相乘得到} \ &= [y - p]\textbf{x}_j && \text{展开} \ &= [y - \sigma(\theta^T\textbf{x})]\textbf{x}_j && \text{因为 } p = \sigma(\theta^T\textbf{x}) \end{align}$$

扩散

原文：chrispiech.github.io/probabilityForComputerScientists/en/part5/diffusion/

---

扩散任务

目标： 创建一个可以从“树照片分布”生成树木图片的模型

数据： 许多树木图片：

整体思路

通过迭代地向像素添加高斯噪声来补充你的数据集，然后训练一个深度学习模型来去除噪声。

合理的步骤数量是每次添加 10% 的噪声，这样在 10 个时间步之后，每个像素都是完全噪声。

关键任务是训练一个深度神经网络来预测像素的“去噪”值：

损失：预测像素和真实颜色之间的均方误差。将你的神经网络的参数设置为最小化损失。然后你就有了一个可以一次去除 10% 噪声的模型。从随机噪声开始，然后运行它通过你的去噪神经网络 10 次。

扩散模型的原理

扩散模型的魔力在于高斯噪声过程。让我们分解一下：

1. 添加噪声：正向过程

在每个步骤 ( t ) 中，我们对像素值添加高斯噪声：[ x_{t+1} = x_t + n_t \quad \text{where } n_t \sim N(0, \sigma²). ] 这逐渐将原始图像转换为纯噪声。

2. 去除噪声：逆向过程

要逆转这个过程，我们需要条件分布 ( x_{t-1} | x_t )。这里令人惊讶的部分是：

关键事实 #1： ( x_{t-1} | x_t ) 是具有已知方差的 Gaussian

如果噪声方差 ( \sigma² ) 足够小，( x_{t-1} | x_t ) 的分布可以近似为：[ x_{t-1} | x_t \sim N(\mu_{t-1}(x_t), \sigma²), ] 其中：

• ( \mu_{t-1}(x_t) )：Gaussian 的均值，它依赖于 ( x_t )。

• ( \sigma² )：噪声的已知方差。

这是一个好消息！这意味着我们只需要估计均值 ( \mu_{t-1}(x_t) ) 来完全描述 ( x_{t-1} | x_t )。

3. 训练神经网络

关键事实 #2： 你只需要标准回归！

为了训练神经网络，我们需要它预测高斯均值 ( \mu_{t-1}(x_t) )。我们如何衡量预测的质量？

预测的高斯 ( q_\theta(x_{t-1} | x_t) )（来自神经网络）和真实的高斯 ( p(x_{t-1} | x_t) ) 之间的差异可以通过 KL 散度 来衡量。幸运的是，在这种情况下：[ \text{最小化 KL 散度 } \Leftrightarrow \text{最小化均方误差 (MSE)}. ] 这只因为分布是 Gaussian。因此，我们可以简单地训练神经网络来最小化其预测的像素值（( \mu_{t-1}(x_t) )）和真实像素值之间的 MSE。

一旦训练完成，神经网络可以迭代地去除图像噪声，从随机噪声开始，直到生成清晰、逼真的图像。

4. 完整的扩散算法

这里是扩散模型的完整工作流程：

1. 正向过程： 向图像添加高斯噪声，将其转换为纯噪声。

2. 反向过程： 训练一个神经网络来预测均值 ( \mu_{t-1}(x_t) ) 并逐步去除噪声。

3. 图像生成： 从随机噪声开始，反向运行神经网络 ( T ) 次以生成逼真的图像。

这种优雅的方法结合了简单的高斯噪声和深度学习的力量，以生成惊人的结果！

---

关键思想 #1 的证明草图

主张： 当噪声方差 ( \sigma² ) 足够小时，( x_{t-1} | x_t ) 大约是高斯分布。

我们从贝叶斯定理开始，表示条件概率：[ \p(x_{t-1} | x_t) = \frac{ \p(x_t | x_{t-1}) \p(x_{t-1}) }{ \p(x_t) } ] 我们将考虑这个表达式的对数。这是因为高斯的对数是二次函数，这将使我们的数学更容易。我们可以写出：[ \log \p(x_{t-1} | x_t) = \log \p(x_t | x_{t-1}) + \log \p(x_{t-1}) - \log \p(x_t) ] 让我们分解这个表达式中的各项。

正向过程似然：

在正向过程中，给定 ( x_{t-1} ) 的 ( x_t ) 是高斯分布：[ \p(x_t | x_{t-1}) = \frac{1}{\sqrt{2 \pi \sigma²}} \exp\left(-\frac{(x_t - x_{t-1})²}{2 \sigma²}\right). ] 取对数，我们得到：[ \log \p(x_t | x_{t-1}) = -\frac{(x_t - x_{t-1})²}{2 \sigma²} + \text{constant} ] 关于 ( x_{t-1} ) 的先验：

先验 ( p(x_{t-1}) ) 是 ( x_{t-1} ) 在前一步的概率。这个很难知道！一棵树的像素的先验分布是什么？然而，我们采用一个非常巧妙的技巧。其对数密度可以在 ( x_t ) 附近对 ( x_{t-1} ) 进行泰勒展开，假设 ( x_{t-1} ) 接近 ( x_t )：[ \log \p(x_{t-1}) \approx \log \p(x_t) + \Big[ \frac{\partial}{\partial x} \log \p(x_t) \Big] (x_{t-1} - x_t) ] 完成平方：

上面的表达式中涉及两个涉及 ( x_{t-1} ) 和 ( x_t ) 差别的项。作为一个有用的步骤，我们将对这些项进行平方完成。

对这些项的总和进行平方完成：[ -\frac{(x_t - x_{t-1})²}{2 \sigma²} + \Big[ \frac{\partial}{\partial x} \log \p(x_t) \Big] (x_{t-1} - x_t) ] 首先，重写 ( (x_t - x_{t-1})² )：[ (x_t - x_{t-1})² = (x_{t-1} - x_t)² ] 这使我们能够将总和重写为：[ -\frac{1}{2 \sigma²} (x_{t-1} - x_t)² + \Big[ \frac{\partial}{\partial x} \log \p(x_t) \Big] (x_{t-1} - x_t) ]

提取 ( -\frac{1}{2 \sigma²} ) 以使二次项更加明确：[ -\frac{1}{2 \sigma²} \left[ (x_{t-1} - x_t)² - 2 \sigma² \Big[ \frac{\partial}{\partial x} \log \p(x_t) \Big] (x_{t-1} - x_t) \right] ]

括号内的表达式是 ( (x_{t-1} - x_t) ) 的二次表达式。让我们对以下表达式完成平方：[ (x_{t-1} - x_t)² - 2 \sigma² \Big[ \frac{\partial}{\partial x} \log \p(x_t) \Big] (x_{t-1} - x_t) ]

这使我们能够将二次表达式重写为：[ \left[ (x_{t-1} - x_t) - \sigma² \Big[ \frac{\partial}{\partial x} \log \p(x_t) \Big] \right]² - \left( \sigma² \Big[ \frac{\partial}{\partial x} \log \p(x_t) \Big] \right)² ]

将其代入。让 ( K ) 代表常数：$$\begin{align} \log \p(x_{t-1} | x_t) &= \log \p(x_t | x_{t-1}) + \log \p(x_{t-1}) - \log \p(x_t) \ &= -\frac{(x_t - x_{t-1})²}{2 \sigma²} + \Big(\log \p(x_t) + \Big[ \frac{\partial}{\partial x} \log \p(x_t) \Big] (x_{t-1} - x_t)\Big) - \log \p(x_t) + K \ &= -\frac{(x_t - x_{t-1})²}{2 \sigma²} + \Big[ \frac{\partial}{\partial x} \log \p(x_t) \Big] (x_{t-1} - x_t)+ K \ &= -\frac{1}{2 \sigma²} \left[ \left( x_{t-1} - x_t - \sigma² \Big[ \frac{\partial}{\partial x} \log \p(x_t) \Big] \right)² \right] + K \end{align}$$

最终结果：回想一下高斯概率密度函数的对数看起来是这样的：

设 ( X \sim N(\mu, \sigma²) )。( X ) 的概率密度函数的对数是什么？$$\log \P(X=x) = -\frac{1}{2 \sigma²} (x - \mu)² + K$$从上面的公式中，我们可以看到 ( x_{t-1} | x_t ) 是高斯分布。我们如何知道这一点？分布与正态分布的对数密度相同，除了加性因子。 [ x_{t-1} | x_t \sim N\left( \mu_{t-1}, \sigma² \right), ] 其中：[ \mu_{t-1} = x_t + \sigma²\Big[ \frac{\partial}{\partial x} \log \p(x_t) \Big] ] 方差保持为 ( \sigma² )，这是从正向过程中固定的。

维度诅咒

原文：chrispiech.github.io/probabilityForComputerScientists/en/examples/curse_of_dimensionality/#part5_examples/

---

在机器学习中，就像计算机科学的许多领域一样，经常涉及高维点，而高维空间具有一些令人惊讶的概率性质。

一个随机值 \(X_i\) 是一个均匀分布（Uni）的（0, 1）。

一个维度为 \(d\) 的随机点是一个 \(d\) 个随机值的列表：\([X_1 \dots X_d]\)。

一个随机值 \(X_i\) 如果 \(X_i\) 小于 0.01 或者 \(X_i\) 大于 0.99，则认为它接近边缘。随机值接近边缘的概率是多少？

设 \(E\) 为随机值接近边缘的事件。\(P(E) = P(X_i < 0.01) + P(X_i > 0.99) = 0.02\)

一个维度为 \(3\) 的随机点 \([X_1, X_2, X_3]\) 如果其任何值接近边缘，则认为该点接近边缘。一个三维点接近边缘的概率是多少？

事件等价于点的所有维度都不接近边缘的补集，即：\(1 - (1 - P(E))³ = 1 - 0.98³ \approx 0.058\)

一个随机点 \([X_1, \dots X_{100}]\) 的维度为 \(100\)，如果其任何值接近边缘，则认为该点接近边缘。100 维点接近边缘的概率是多少？

同样，它是：\(1 - (1 - P(E))^{100} = 1 - 0.98^{100} \approx 0.867\)高维点的许多其他现象：例如，点之间的欧几里得距离开始收敛。

MLE 正态分布演示

原文：chrispiech.github.io/probabilityForComputerScientists/en/examples/mle_demo/

---

让我们手动执行最大似然估计。你的任务是选择参数值，使得数据看起来尽可能可能。这里有 20 个数据点，我们假设它们来自正态分布

数据 = [6.3 , 5.5 , 5.4, 7.1, 4.6, 6.7, 5.3 , 4.8, 5.6, 3.4, 5.4, 3.4, 4.8, 7.9, 4.6, 7.0, 2.9, 6.4, 6.0 , 4.3]

选择你的参数估计值

参数 \(\mu\): 参数 \(\sigma\):

给定你的参数的数据似然

似然：对数似然：最佳观测：

你的高斯分布

帕累托分布的最大似然估计

原文：chrispiech.github.io/probabilityForComputerScientists/en/examples/mle_pareto/

---

你正在创建不同大小的圆形艺术作品，这些圆形遵循帕累托分布：

\(X \sim \text{Pareto}(\alpha)\)

帕累托分布由单个参数 \(\alpha\) 定义，并具有概率密度函数

\[\begin{align*} f(x) = \frac{\alpha}{x^{\alpha + 1}} \end{align*} \]

你希望你的艺术作品中的 \(\alpha\) 与你当地海滩的沙子中的 \(\alpha\) 相匹配。你去了海滩，收集了 100 粒沙子并测量了它们的大小。将测量的半径称为 \(x_1 , \dots, x_{100}\)：

 observations = [1.677, 3.812, 1.463, 2.641, 1.256, 1.678, 1.157, 
1.146, 1.323, 1.029, 1.238, 1.018, 1.171, 1.123, 1.074, 1.652, 
1.873, 1.314, 1.309, 3.325, 1.045, 2.271, 1.305, 1.277, 1.114, 
1.391, 3.728, 1.405, 1.054, 2.789, 1.019, 1.218, 1.033, 1.362, 
1.058, 2.037, 1.171, 1.457, 1.518, 1.117, 1.153, 2.257, 1.022, 
1.839, 1.706, 1.139, 1.501, 1.238, 2.53 , 1.414, 1.064, 1.097, 
1.261, 1.784, 1.196, 1.169, 2.101, 1.132, 1.193, 1.239, 1.518, 
2.764, 1.053, 1.267, 1.015, 1.789, 1.099, 1.25 , 1.253, 1.418, 
1.494, 1.015, 1.459, 2.175, 2.044, 1.551, 4.095, 1.396, 1.262, 
1.351, 1.121, 1.196, 1.391, 1.305, 1.141, 1.157, 1.155, 1.103, 
1.048, 1.918, 1.889, 1.068, 1.811, 1.198, 1.361, 1.261, 4.093, 
2.925, 1.133, 1.573]

根据你收集的数据，推导出 \(\alpha\) 的最大似然估计公式。

编写对数似然函数

在最大似然估计（MLE）中，第一个主要目标是为我们数据构造一个对数似然表达式。为此，我们首先写出如果我们被告知 \(\alpha\) 的值，我们的数据集看起来有多可能：\begin{aligned} L(\alpha) = f(x_1\dots x_n) = \prod_{i=1}^{n\frac{\alpha}{x_i}{\alpha+1}} \end{aligned}

如果我们尝试优化对数似然，优化将变得容易得多：\begin{aligned} LL(\alpha) &= \log L(\alpha) = \log \prod_{i=1}^{n\frac{\alpha}{x_i}{\alpha+1}} \ &= \sum_{i=1}^n \log \frac{\alpha}{x_i^{\alpha+1}} \ &= \sum_{i=1}^n \log \alpha - (\alpha +1)\log x_i \ &= n \log \alpha - (\alpha +1) \sum_{i=1}^n\log x_i \ \end{aligned}

选择 \(\alpha\)

我们将选择 \(\alpha\) 以最大化对数似然。为此，我们需要 \(\alpha\) 对 LL 的导数 \begin{aligned} \frac{\partial LL(\alpha)}{\partial \alpha} &= \frac{\partial LL(\alpha)}{\partial \alpha} \Big( n \log \alpha - (\alpha +1) \sum_{i=1}^n\log x_i \Big) \ &= \frac{n}{\alpha} - \sum_{i=1}^n\log x_i \end{aligned}

一种优化方法是求导数并将其设为零：\begin{aligned} 0=\frac{n}{\alpha} - \sum_{i=1}^n\log x_i\ \sum_{i=1}^n\log x_i = \frac{n}{\alpha} \ \alpha = \frac{n}{ \sum\limits_{i=1}^n\log x_i} \end{aligned}

到目前为止，我们有一个可以用来计算 \(\alpha\) 的公式！哇哦

将其放入代码

import math

def estimate_alpha(observations):
    # This code computes the MLE estimate of alpha
    log_sum = 0
    for x_i in observations:
        log_sum += math.log(x_i)
    n = len(observations)
    return n / log_sum

def main():
    observations = [1.677, 3.812, 1.463, 2.641, 1.256, 1.678, 1.157, 1.146, 
    1.323, 1.029, 1.238, 1.018, 1.171, 1.123, 1.074, 1.652, 1.873, 1.314, 
    1.309, 3.325, 1.045, 2.271, 1.305, 1.277, 1.114, 1.391, 3.728, 1.405, 
    1.054, 2.789, 1.019, 1.218, 1.033, 1.362, 1.058, 2.037, 1.171, 1.457, 
    1.518, 1.117, 1.153, 2.257, 1.022, 1.839, 1.706, 1.139, 1.501, 1.238, 
    2.53 , 1.414, 1.064, 1.097, 1.261, 1.784, 1.196, 1.169, 2.101, 1.132, 
    1.193, 1.239, 1.518, 2.764, 1.053, 1.267, 1.015, 1.789, 1.099, 1.25 , 
    1.253, 1.418, 1.494, 1.015, 1.459, 2.175, 2.044, 1.551, 4.095, 1.396, 
    1.262, 1.351, 1.121, 1.196, 1.391, 1.305, 1.141, 1.157, 1.155, 1.103, 
    1.048, 1.918, 1.889, 1.068, 1.811, 1.198, 1.361, 1.261, 4.093, 2.925, 
    1.133, 1.573]
    alpha = estimate_alpha(observations)
    print(alpha)

if __name__ == '__main__':
    main()

高斯混合

原文：chrispiech.github.io/probabilityForComputerScientists/en/examples/mixture_models/

---

数据 = [6.47, 5.82, 8.7, 4.76, 7.62, 6.95, 7.44, 6.73, 3.38, 5.89, 7.81, 6.93, 7.23, 6.25, 5.31, 7.71, 7.42, 5.81, 4.03, 7.09, 7.1, 7.62, 7.74, 6.19, 7.3, 7.37, 6.99, 2.97, 3.3, 7.08, 6.23, 3.67, 3.05, 6.67, 6.5, 6.08, 3.7, 6.76, 6.56, 3.61, 7.25, 7.34, 6.27, 6.54, 5.83, 6.44, 5.34, 7.7, 4.19, 7.34]参数 \(t\): 参数 \(\mu_a\): 参数 \(\sigma_a\): 参数 \(\mu_b\): 参数 \(\sigma_b\):

似然：对数似然：最佳观察：

什么是高斯混合？

高斯混合描述了一个随机变量，其概率密度函数可能来自两个高斯分布（或更多，但在这个演示中我们只使用两个）。样本来自第一个高斯分布的概率是确定的，否则来自第二个。它有五个参数：四个用于描述两个高斯分布，一个用于描述两个高斯分布的相对权重。

生成代码

from scipy import stats
def sample():
   # choose group membership
   membership = stats.bernoulli.rvs(0.2)
   if membership == 1:
      # sample from gaussian 1
      return stats.norm.rvs(3.5,0.7)
   else:
      # sample from gaussian 2
      return stats.norm.rvs(6.8,0.7)

概率密度函数

\[\begin{align*} f(X=x) = t \cdot f(A=x) + (1-t) \cdot f(B=x) \end{align*}$$ $\text{st}$ $A \sim N(\mu_a, \sigma_a²)$ $B \sim N(\mu_b, \sigma_b²)$ 将所有这些放在一起，高斯混合的概率密度函数（PDF）为： $$\begin{align*} f(x) = t \cdot \Big(\frac{1}{\sqrt{2\pi}\sigma_a} e^{-\frac{1}{2}(\frac{x-\mu_a}{\sigma_a})²}\Big) + (1-t) \cdot \Big(\frac{1}{\sqrt{2\pi}\sigma_b} e^{-\frac{1}{2}(\frac{x-\mu_b}{\sigma_b})²}\Big) \end{align*}\]

高斯混合的最大似然估计

特别说明：尽管生成故事有一个伯努利（组成员）属性，但它从未被观察到。MLE 最大化的是观察数据的似然。

令 \(\vec{\theta} = [t, \mu_a,\mu_b,\sigma_a, \sigma_b]\) 为参数。因为数学表达式会变得很长，我将使用 \(\theta\) 作为符号来代替 \(\vec{\theta}\)。只需记住它是一个向量。

最大似然估计（MLE）的想法是选择使对数似然最大的 \(\theta\) 值。所有优化方法都需要我们计算我们想要优化的东西（对数似然）相对于我们可以改变的值（我们的参数）的偏导数。

似然函数

\[\begin{align*} L(\theta) &= \prod_i^n f(x_i | \theta) \\ &= \prod_i^n \Big[t \cdot \Big(\frac{1}{\sqrt{2\pi}\sigma_a} e^{-\frac{1}{2}(\frac{x_i-\mu_a}{\sigma_a})²}\Big) + (1-t) \cdot \Big(\frac{1}{\sqrt{2\pi}\sigma_b} e^{-\frac{1}{2}(\frac{x_i-\mu_b}{\sigma_b})²}\Big) \Big] \end{align*} \]

对数似然函数

\[\begin{align*} LL(\theta) &= \log L(\theta) \\ &= \log \prod_i^n f(x_i | \theta) \\ &= \sum_i^n \log f(x_i | \theta) \\ \end{align*} \]

现在已经足够了，但如果你想要扩展这个术语，你会得到：$$\begin{align} LL(\theta) &= \sum_i^n \log \Big[t \cdot \Big(\frac{1}{\sqrt{2\pi}\sigma_a} e^{-\frac{1}{2}(\frac{x_i-\mu_a}{\sigma_a})²}\Big) + (1-t) \cdot \Big(\frac{1}{\sqrt{2\pi}\sigma_b} e^{-\frac{1}{2}(\frac{x_i-\mu_b}{\sigma_b})²}\Big) \Big] \end{align}$$

关于 \(\theta\) 的 LL 的导数

这里是一个计算关于参数之一 \(\mu_a\) 的偏导数的例子。你需要为所有参数都找到这样的导数。

注意：当我最初写这个演示时，我以为这只是一个简单的导数。它并不简单，因为对数中有一个和。因此，对数项不会减少。对数仍然用于将外层的 \(\prod\) 转换为 \(\sum\)。因此，\(LL\) 的偏导数是可解的，但证明使用了相当多的链式法则。要点：本节的主要要点（如果你想跳过导数证明）是，得到的导数足够复杂，我们希望有一种方法来计算 argmax，而无需将那个导数设为零并求解 \(\mu_a\)。梯度下降法登场！

在对对数似然函数进行大量求导时，一个好的第一步是考虑单个数据点的似然对数的导数。这是对数似然表达式中的内层求和：$ \frac{\d}{\d \mu_a} \log f(x_i|\theta) $

在开始之前：请注意，\(\mu_a\) 在 $f(x_i|\theta) $ 这个项中不出现：$(1-t) \cdot \Big(\frac{1}{\sqrt{2\pi}\sigma_b} e^{-\frac{1}{2}(\frac{x_i-\mu_b}{\sigma_b})²}\Big) =K $ 在证明中，当我们遇到这个项时，我们将把它视为一个常数，我们称之为 \(K\)。好的，让我们开始吧！ $$\begin{align} \frac{\d}{\d \mu_a} & \log f(x_i|\theta) \ &= \frac{1}{f(x_i|\theta)} \frac{\d}{\d \mu_a} f(x_i|\theta) &&\text{链式法则在 }\log \ &= \frac{1}{f(x_i|\theta)} \frac{\d}{\d \mu_a} \Big[t \cdot \Big(\frac{1}{\sqrt{2\pi}\sigma_a} e^{-\frac{1}{2}(\frac{x_i-\mu_a}{\sigma_a})²}\Big) + K\Big] &&\text{代入 }f(x_i|\theta)\ &= \frac{1}{f(x_i|\theta)} \frac{\d}{\d \mu_a} \Big[t \cdot \Big(\frac{1}{\sqrt{2\pi}\sigma_a} e^{-\frac{1}{2}(\frac{x_i-\mu_a}{\sigma_a})²}\Big)\Big] &&\frac{\d}{\d \mu_a} K = 0\ &= \frac{t}{f(x_i|\theta) \sqrt{2\pi}\sigma_a} \cdot \frac{\d}{\d \mu_a} e^{-\frac{1}{2}(\frac{x_i-\mu_a}{\sigma_a})²} &&\text{提取常数}\ &= \frac{t}{f(x_i|\theta) \sqrt{2\pi}\sigma_a} \cdot e^{-\frac{1}{2}(\frac{x_i-\mu_a}{\sigma_a})²} \cdot \frac{\d}{\d \mu_a} -\frac{1}{2}(\frac{x_i-\mu_a}{\sigma_a})²&&\text{对 }e^x \text{ 求链式}\ &= \frac{t}{f(x_i|\theta) \sqrt{2\pi}\sigma_a} \cdot e^{-\frac{1}{2}(\frac{x_i-\mu_a}{\sigma_a})²} \cdot \Big[-(\frac{x_i-\mu_a}{\sigma_a}) \frac{\d}{\d \mu_a} (\frac{x_i-\mu_a}{\sigma_a})\Big]&&\text{对 }x² \text{ 求链式}\ &= \frac{t}{f(x_i|\theta) \sqrt{2\pi}\sigma_a} \cdot e^{-\frac{1}{2}(\frac{x_i-\mu_a}{\sigma_a})²} \cdot \Big[ -(\frac{x_i-\mu_a}{\sigma_a}) \cdot \frac{-1}{\sigma_a}\Big] &&\text{最终导数} \ &= \frac{t}{f(x_i|\theta) \sqrt{2\pi}\sigma_a³} \cdot e^{-\frac{1}{2}(\frac{x_i-\mu_a}{\sigma_a})²} \cdot (x_i-\mu_a) &&\text{简化} \ \end{align}$$

那是针对单个数据点的。对于完整的数据集： $$\begin{align} \frac{\d LL(\theta)}{\d \mu_a} &= \sum_i^n \frac{\d}{\d \mu_a} \log f(x_i|\theta)\ &= \sum_i^n \frac{t}{f(x_i|\theta) \sqrt{2\pi}\sigma_a³} \cdot e^{-\frac{1}{2}(\frac{x_i-\mu_a}{\sigma_a})²} \cdot (x_i-\mu_a) \end{align}$$

这个过程应该对所有五个参数重复进行！现在，我们应该如何找到一个 \(\mu_a\) 的值，在存在其他参数设置和数据的情况下，使得这个导数为零？将导数设为 0 并求解 \(\mu_a\) 是不会奏效的。

使用优化器估计参数

一旦我们有了 \(LL\) 函数以及 \(LL\) 对每个参数的导数，我们就可以使用优化器来计算 argmax。在这种情况下，最佳选择可能是梯度上升（或负对数似然梯度下降）。 $$\begin{align} \nabla_{\theta} LL(\theta) &=\begin{bmatrix} \frac{\d LL(\theta)}{\d t}\ \frac{\d LL(\theta)}{\d \mu_a}\ \frac{\d LL(\theta)}{\d \mu_b}\ \frac{\d LL(\theta)}{\d \sigma_a}\ \frac{\d LL(\theta)}{\d \sigma_b}\ \end{bmatrix} \end{align}$$