离散随机变量

回到抛硬币的例子。假设硬币是特殊铸造的，以便抛出正面的概率为\(p\)。在这种情况下，我们抛三次。使用正面数量作为随机变量\(X\)，它包含四个可能的结果：0、1、2或3。

我们可以计算每个输出结果的可能性。由于每次投掷都是一个独立的试验，三个正面的可能性 P(X=2) 是 \(p^3\)。两个正面包括三种情况：HHT、HTH、THH，每种情况的概率为 \(p^2(1-p)\)，总共 \(P(X=2) = 3p^2(1-p)\)。类似地，\(P(X=1)=3p(1-p)^2\)，\(P(X=0)=(1-p)^3\)。

形式上，考虑一系列\(n\)个独立的试验，每个试验包含两个可能的结果，并且感兴趣的结果发生的概率为\(p\)，则随机变量\(X\)的概率分布是（\(X\)是感兴趣结果的数量）：

\[P(X=k) = {N\choose k} p^k(1-p)^{n-k}.\] {#eq:stats:binomial_pdf}

这种分布称为二项概率分布。我们可以使用Stats.binomial_rvs函数模拟抛硬币的过程。假设抛硬币的概率为0.4，抛10次。

let _ =
  let toss = Array.make 10 0 in
  Array.map (fun _ -> Stats.binomial_rvs 0.3 1) toss
>- : int array = [|0; 0; 0; 0; 0; 0; 0; 0; 1; 0|]

方程[@eq:stats:binomial_pdf]称为此二项分布的概率密度函数（PDF）。形式上，随机变量X的PDF表示为\(p_X(k)\)，定义为：

\[p_X(k)=P({s \in S | X(s) = k}),\]

其中\(S\)是样本空间。这也可以用代码表示：

let x = [|0; 1; 2; 3|]
>val x : int array = [|0; 1; 2; 3|]
let p = Array.map (Stats.binomial_pdf ~p:0.3 ~n:3) x
>val p : float array =
>  [|0.342999999999999916; 0.440999999999999837; 0.188999999999999918;
>    0.0269999999999999823|]
Array.fold_left (+.) 0. p
>- : float = 0.999999999999999778

除了PDF之外，另一个相关且经常使用的想法是查看随机变量\(X\)落在某个范围内的概率：\(P(a \leq X \leq b)\)。它可以重写为\(P(X \leq b) - P(X \leq a - 1)\)。这里术语\(P(X \leq t)\)称为随机变量\(X\)的累积分布函数。对于二项分布，它的CDF是：

\[p(X\leq~k)=\sum_{i=0}^k{N\choose i} p^k(1-p)^{n-i}.\]

我们可以再次使用代码计算3次抛硬币问题的CDF。

let x = [|0; 1; 2; 3|]
>val x : int array = [|0; 1; 2; 3|]
let p = Array.map (Stats.binomial_cdf ~p:0.3 ~n:3) x
>val p : float array = [|0.342999999999999972; 0.784; 0.973; 1.|]

连续随机变量

与离散随机变量不同，连续随机变量有无限数量的可能结果。例如，在均匀分布中，我们可以选择0到1之间的一个随机实数。显然可以有无限数量的输出。

最广泛使用的连续分布之一无疑是高斯分布。它的概率函数是一个连续的函数：\[p(x) = \frac{1}{\sqrt{2\pi~\delta}}e^{-\frac{1}{2}\left(\frac{t - \mu}{\sigma}\right)^2}\] {#eq:stats:gaussian_pdf}

这里的\(\mu\)和\(\sigma\)是参数。根据它们，\(p(x)\)可以采取不同的形状。让我们看一个例子。

在这个例子中，我们生成了两个数据集，都包含了从不同的高斯分布\(\mathcal{N} (\mu, \sigma^{2})\)中抽取的999个点。对于第一个数据集，配置为\((\mu = 1, \sigma = 1)\)；而对于第二个数据集，配置为\((\mu = 12, \sigma = 3)\)。

let noise sigma = Stats.gaussian_rvs ~mu:0. ~sigma;;
let x = Array.init 999 (fun _ -> Stats.gaussian_rvs ~mu:1. ~sigma:1.);;
let y = Array.init 999 (fun _ -> Stats.gaussian_rvs ~mu:12. ~sigma:3.);;

我们可以使用直方图绘制数据集，如下所示。在调用histogram时，我们还明确指定了30个bin。您还可以使用spec命名参数微调图形，以指定颜色、x范围、y范围等。我们将在单独的章节中详细讨论如何使用Owl绘图。

(* 将数组转换为矩阵 *)

let x' = Mat.of_array x 1 999;;
let y' = Mat.of_array y 1 999;;

(* 绘制图形 *)

let h = Plot.create ~m:1 ~n:2 "plot_02.png" in

Plot.subplot h 0 0;
Plot.set_ylabel h "频率";
Plot.histogram ~bin:30 ~h x';
Plot.histogram ~bin:30 ~h y';

Plot.subplot h 0 1;
Plot.set_ylabel h "PDF p(x)";
Plot.plot_fun ~h (fun x -> Stats.gaussian_pdf ~mu:1. ~sigma:1. x) (-2.) 6.;
Plot.plot_fun ~h (fun x -> Stats.gaussian_pdf ~mu:12. ~sigma:3. x) 0. 25.;

Plot.output h;;

在子图1中，我们可以看到第二个数据集的分布范围更广。在子图2中，我们还绘制了两个数据集的概率密度函数。

高斯的CDF可以通过无限求和，即积分来计算：

\[p(x\leq~k)=\frac{1}{\sqrt{2\pi}}\int_{-\infty}^k~e^{-t^2/2}dt.\]

我们可以使用gaussian_cdf观察这个函数。

let h = Plot.create "plot_gaussian_cdf.png" in
Plot.set_ylabel h "CDF";
Plot.plot_fun ~h ~spec:[ RGB (66,133,244); LineStyle 1; LineWidth 2.; Marker "*" ] (fun x -> Stats.gaussian_cdf ~mu:1. ~sigma:1. x) (-2.) 6.;
Plot.plot_fun ~h ~spec:[ RGB (219,68,55);  LineStyle 2; LineWidth 2.; Marker "+" ] (fun x -> Stats.gaussian_cdf ~mu:12. ~sigma:3. x) 0. 25.;
Plot.(legend_on h ~position:SouthEast [|"mu=1,sigma=1"; "mu=12, sigma=3"|]);
Plot.output h

描述性统计

一个随机变量描述一个单独的事件。一整组对某些感兴趣的个体的集合称为总体。总体可以用多个描述性统计来描述。其中最常用的两个是均值和方差。总体的均值\(X\)，具有\(n\)个元素，定义为：

\[E(X) = \frac{1}{n}\sum_{i}x_i,\] {#eq:stats:mean} 其中\(x_i\)是总体中的第\(i\)个元素。方差的定义类似：

\[Var(X) = \frac{1}{n}\sum_{i}(x_i - E(X))^2.\] {#eq:stats:variance}

一个类似且常用的概念是标准差，它是方差的平方根。均值（或期望值）和方差的含义是很明显的，前者是总体的代表性中心值，后者是值围绕中心期望的分散程度。

这些定义是针对离散随机变量的，但它们可以很容易地扩展到连续情况。为了使其更一般，我们将实变量关于值\(X\)的n阶矩定义为：

\[M_n(X) = \int_x~(x_i - c)^2~f(x_i)dx,\] {#eq:stats:moment}

其中\(f(x)\)是变量\(X\)的连续函数，\(c\)是某个常数。您可以看到，均值实际上是一阶矩，而方差是二阶矩。
三阶矩称为偏度，表示实随机变量的概率分布的不对称性。四阶矩称为峰度，显示概率分布的“尾巴”有多长。

让我们看一个简单的例子。我们首先生成一百个在0到10之间均匀分布的随机数。在这里，我们使用Stats.uniform_rvs函数生成遵循均匀分布的数字。

let data = Array.init 100 (fun _ -> Stats.uniform_rvs 0. 10.);;

然后我们使用mean函数计算样本平均值。如预期的那样，它大约是5。我们还可以使用相应的函数轻松计算其他更高的矩。我们可以对这些结果进行一个非常粗糙而快速的解释。它具有广泛的分布（左右约3），根据非常小的偏度数字，分布不对称。最后，小的峰度表明分布没有明显的尾部。

Stats.mean data
>- : float = 5.18160409659184573
Stats.std data
>- : float = 2.92844832850280135
Stats.var data
>- : float = 8.57580961271085229
Stats.skew data
>- : float = -0.109699186612116223
Stats.kurtosis data
>- : float = 1.75165078829330856

以下代码计算分布的不同中心矩。中心矩是关于随机变量均值的概率分布的矩。零阶中心矩始终为1，第一阶接近零，第二阶接近方差。

Stats.central_moment 0 data
>- : float = 1.
Stats.central_moment 1 data
>- : float = -3.13082892944294137e-15
Stats.central_moment 2 data
>- : float = 8.49005151658374224
Stats.central_moment 3 data
>- : float = -2.75496511397836663

除了矩之外，我们经常使用顺序统计量来了解数据。顺序统计量和秩统计量是非参数统计和推断中最基本的工具之一。统计样本的第\(k^{th}\)顺序统计量等于它的第k小值。例如：

在Owl的Stat模块中有许多有序统计函数供您探索。其中一些最常用的如下所示：

Stats.min;;
Stats.max;;
Stats.median;;
Stats.quantile;;
Stats.first_quartile;;
Stats.third_quartile;;
Stats.percentile;;

min和max很容易使用。median是整个样本中排序列表中的中间数字。它有时比mean更能描述数据，因为后者更容易受到异常值的影响。

类似的概念是quartile：样本中有75%的测量值大于第一四分位数，有25%大于第三四分位数。median也是第二四分位数。更一般的概念是percentile，即总值中低于该度量的百分比。例如，第一四分位数也是第25百分位数。

无偏估计

有多种抽样方法。随机抽样是一个常见选择。类似的方法是“分层随机抽样”，它首先将总体分成几个组，然后在每个组内随机选择。例如，在设计问卷调查时，你希望各个年龄组的人都能得到平等代表，那么分层随机抽样就是一种更合适的方法。当然，只要样本具有代表性，也有其他合理的抽样方法，这意味着总体中的任何成员被选入样本的可能性是相等的。

在选择适当的样本之后，下一步是用样本描述总体。诸如均值和方差等统计数据仍然非常有用，但我们能否直接使用样本的统计数据并声明它们也可以用来代表整个总体呢？事实上，这取决于统计数据是否是一个无偏估计，即其值的期望值是相应的总体参数。

例如，让我们取一个包含\(n\)个元素的样本，其均值为\(m\)：

\[m = \frac{1}{n}\sum_{i=1}^n~x_i,\]

其中\(x_i\)是样本中的一个元素。将总体表示为\(\mu\)，可以进一步证明：\(E(m) = \mu\)。因此，样本均值是总体的无偏估计。

方差就不能这样说了。样本方差为：

\[v = \frac{1}{n}\sum_{i=1}^n(x_i - m)^2.\]

假设总体的方差为\(\sigma^2\)，则可以证明\(E(v) = \frac{n - 1}{n}\sigma^2\)。因此，总体方差的无偏估计不是样本的\(v\)，而是\(\frac{n}{n-1}v\)。

推断总体参数

在前一节中，我们已经展示了如何在给定来自总体的样本的情况下得到均值和方差的期望值。但我们也许需要了解的不仅仅是期望值。例如，我们能否确定一个区间，我们可以相当确信总体均值位于其中？本节探讨了这个问题。

首先，我们需要解释中心极限定理。它指出，如果你有一个总体，并从总体中取足够大的随机样本（带有替换），则样本均值的分布将近似正态分布。如果样本大小足够大（如\(n \lt 20\)），无论总体分布如何，该定理都成立。

具体而言，假设我们反复从相同大小的子集中抽样，然后我们可以定义随机变量\(X\)表示每个抽样子集的均值。根据中心极限定理，可以推导出，假设总体的均值为\(\mu\)，方差为\(\sigma^2\)，均未知，那么\(X\)遵循均值为\(\mu\)，方差为\(\frac{\sigma^2}{n}\)的正态分布。

由于总体的均值和方差都是未知的，显然我们不能同时解决这两个未知数。为了更精确地估计总体均值\(\mu\)，让我们首先假设可以直接用样本方差计算总体方差：\(\sigma^2 = \frac{1}{n-1}\sum_{i=1}^n(x_i - m)^2\)。这个假设在实践中是合理的，当\(n\)足够大时。

既然我们知道\(X\)遵循正态分布，我们可以利用其中的一些好的性质。例如，我们知道概率质量的95%位于均值的1.96标准差之内。我们可以使用正态分布的CDF函数来验证这一点：

let f = Stats.gaussian_cdf ~mu:0. ~sigma:1. in
f 1.96 -. f (-1.96)
>- : float = 0.950004209703559

因此，对于\(X\)中的任何值\(x\)，我们知道：\[P(\mu - 1.96~\frac{\sigma}{\sqrt{n}} \le x \le \mu + 1.96~\frac{\sigma}{\sqrt{n}}).\]

稍作变化，它变成了：

\[P( x - 1.96~\frac{\sigma}{\sqrt{n}} \le \mu \le x + 1.96~\frac{\sigma}{\sqrt{n}}).\]

这意味着在给定样本均值 \(m\) 的情况下，总体均值 \(\mu\) 位于这个范围内 [ \(m - 1.96~\frac{\sigma}{\sqrt{n}}\), \(m + 1.96~\frac{\sigma}{\sqrt{n}}\) ] 的概率为 95%。这被称为其置信区间。同样，总体方差 \(\sigma^2\) 直接使用样本的无偏估计。

让我们回到 1.96 这个数字。我们使用这个范围是因为假设 \(X\) 遵循正态分布。变量 \(\frac{x - \mu}{\sigma/\sqrt{n}}\) 遵循标准正态分布。它被称为标准 Z 变量。我们可以查看标准正态分布表，找到对应于 95% 置信水平的范围。然而，正如我们所解释的，当 \(n\) 较小时，这并不成立，因为我们实际上使用的是 \(\frac{x-\mu}{\sqrt{\frac{\sum_{i}(x - m)^2}{n(n-1)}}}\) 而不是真正的 \(z\) 变量。后者被称为标准 t 变量，它遵循自由度为 \(n-1\) 的 t 分布。当 \(n\) 很大时，t 分布几乎与正态分布相同。因此，如果 \(n\) 较小，我们需要查阅 t 表。例如，如果 \(n=17\)，则范围参数约为 2.12，可以验证为：

let f x = Stats.t_cdf x ~df:16. ~loc:0. ~scale:1. in 
f 2.12 -. f (-2.12)
>- : float = 0.950009071286895823

这就是总体均值的全部内容。总体方差范围的估计使用 \(\chi\)-square 分布，但在实践中很少使用。因此，在本节中我们省略了它。

理论

虽然描述性统计只关注观察数据的属性，统计推断着重于研究数据集是否是从一个更大的总体中抽样的。换句话说，统计推断对总体提出了假设。假设检验是推断性统计分析中的一种重要方法。

让我们考虑经典的抛硬币实验。假设我们有一个基本假设/假设，即大多数硬币在投掷时可以以50/50的概率出现正面或反面。现在，如果你抛一枚硬币3次，得到3个正面，你能否声称这枚硬币不是一枚正常的硬币？以多大的概率？另一个例子是，假设你声称你提出的算法提高了最先进的运行速度，你有两个关于使用两种不同算法的执行时间的样本，然后你如何确信你的声明是合理的。这就是我们需要假设检验的地方。

关于数据集之间的统计关系，提出了两个假设。

• 零假设 \(H_0\)：两个数据集之间没有关系。
• 备择假设 \(H_1\)：两个数据集之间存在统计显著关系。

在假设成立的情况下，出现某个结果的概率称为其 p 值。在实践中，p 值被设定为 5%（1% 也经常使用）。例如，如果我们相信一枚硬币是正常的，但实验结果只有很低的概率，比如 0.03，我们可以在 5% 的置信水平下拒绝这个假设（“这枚硬币是正常的”）。

请注意，如果我们不拒绝一个假设，并不意味着它被接受。如果我们抛硬币三次得到三个正面。在假设这枚硬币是正常的情况下，这个结果的概率是 12.5%，因此我们不能拒绝这个假设。
不拒绝并不意味着我们非常确定硬币完全没有偏倚。因此，在选择假设时需要非常小心。\(H_0\) 应该是我们认为足够坚实以解释数据的东西，除非被观察到的数据强烈挑战。此外，将零假设尽可能明确地表述也有助于使其尽可能不可否认。

在假设检验中的正态分布

最常见的测试之一是查看观察数据是否来自某个正态分布。这被称为“z-测试”。现在让我们看看如何在 Owl 中执行 z-测试。我们首先生成两个数据集，两者都是从正态分布中抽取的，但参数不同。第一个 data_0 从 \(\mathcal{N}(0, 1)\) 抽取，而第二个 data_1 从 \(\mathcal{N}(3, 1)\) 抽取。

let data_0 = Array.init 10 (fun _ -> Stats.gaussian_rvs ~mu:0. ~sigma:1.);;
let data_1 = Array.init 10 (fun _ -> Stats.gaussian_rvs ~mu:3. ~sigma:1.);;

我们的假设是数据集是从正态分布 \(\mathcal{N}(0, 1)\) 中抽取的。从我们生成合成数据的方式来看，显然 data_0 将通过测试，但让我们看看 Owl 是否会使用其 Stats.z_test 函数测试我们。

Stats.z_test ~mu:0. ~sigma:1. data_0
>- : Owl_stats.hypothesis =
>{Owl.Stats.reject = false; p_value = 0.289340080583773251;
> score = -1.05957041132113083}

返回的结果是以下类型定义的记录。字段是自说明的：reject 字段告诉我们是否拒绝零假设，以及使用给定数据集计算的 p 值和分数。

type hypothesis = {
  reject : bool;
  p_value : float;
  score : float;
}

从前面的结果中，我们可以看到 reject = false，表示拒绝零假设，因此数据集 data_0 是从 \(\mathcal{N}(0, 1)\) 中抽取的。那么第二个数据集呢？

Stats.z_test ~mu:0. ~sigma:1. data_1
>- : Owl_stats.hypothesis =
>{Owl.Stats.reject = true; p_value = 5.06534675819424548e-23;
> score = 9.88035435799393547}

正如我们预期的那样，零假设在一个非常小的p值下被接受。这表明data_1不是从假定的\(\mathcal{N}(0, 1)\)分布中抽取的。

在前一节中，我们介绍了z变量和t变量。除了z检验之外，另一种经常使用的测试是检查在零假设下给定的数据是否遵循学生t分布。在Stats模块中，t_test ~mu ~alpha ~side x函数返回t检验的测试决策，这是当总体标准差未知时位置参数的参数测试。这里mu是总体均值，alpha是显著性水平。

双样本推断

另一种常见的测试类型是测试两个样本是否来自同一总体。例如，我们需要测试对算法的改进是否真的有效。默认的零假设是，两个样本是从同一总体中抽取的。测试策略取决于样本大小。

如果两个样本的大小相同，那么可以通过逐个相减两个集合中的元素来简化测试。测试然后变成检查结果集是否来自均值为0的总体。这个问题可以使用前一节介绍的方法通过使用t检验来解决。具体而言，我们提供了配对样本t检验函数。t_test_paired ~alpha ~side x y返回一个测试决策，用于零假设，即x – y中的数据来自均值为零且方差未知的正态分布。

当这两个子集的大小不相同时，问题变得更加复杂。然后，我们需要讨论关于它们的置信区间的两种情况。如果两个区间都不重叠，显然我们可以相当肯定两个样本来自不同的总体，因此拒绝零假设。如果区间重叠，我们需要确保其他一些条件成立。例如，如果我们可以假设两个总体的方差相同，那么可以证明变量：

\[\frac{\bar{x} - \bar{y}}{\sqrt{\frac{\sum_{i=1}^a~(x_i - \bar{x})^2 + \sum_{i=1}^b~(y_i - \bar{y})^2}{a + b - 2}(\frac{1}{a} + \frac{1}{b})}},\]

是一个具有\(a + b - 2\)自由度的标准t变量，其中\(a\)和\(b\)是样本\(X\)和\(Y\)的长度。

这个想法被实现为不配对样本t检验。函数t_test_unpaired ~alpha ~side ~equal_var x y返回一个测试决策，用于零假设，即向量x和y中的数据来自具有相等均值和相等但未知方差的独立正态分布的随机样本。这里equal_var指示两个样本是否具有相同的方差。如果两个方差不同，我们需要使用非参数测试，如Welche的t检验。

其他测试类型

除了我们到目前为止提到的内容之外，Owl中的Stats模块还支持许多其他不同类型的假设检验。本节将对它们进行简要介绍。

• Kolmogorov-Smirnov测试：ks_test ~alpha x f返回一个测试决策，用于零假设，即向量x中的数据来自具有CDF f的分布的独立随机样本。备择假设是x中的数据来自不同的分布。结果(h,p,d)：d是Kolmogorov-Smirnov测试统计量。

• 双样本Kolmogorov-Smirnov测试：ks2_test ~alpha x y返回一个测试决策，用于零假设，即向量x和y中的数据来自相同分布的独立随机样本。

• 卡方方差测试var_test ~alpha ~side ~variance x返回一个测试决策，用于零假设，即x中的数据来自具有输入variance的正态分布，使用卡方方差测试。备择假设是x来自具有不同方差的正态分布。

• Jarque-Bera测试jb_test ~alpha x返回一个测试决策，用于零假设，即数据x来自具有未知均值和方差的正态分布，使用Jarque-Bera测试。

• Wald–Wolfowitz Runs测试runs_test ~alpha ~v x返回一个测试决策，用于零假设，即数据x是随机排序的，反对它们不是，通过运行Wald–Wolfowitz测试。该测试基于连续值在均值上方或下方的运行的数量。~v是参考值，默认值是x的中位数。

• Mann-Whitney秩测试mannwhitneyu ~alpha ~side x y对样本x和y执行Mann-Whitney秩测试。如果每个样本的长度小于10且没有结合，则使用精确测试，否则使用渐近正态分布。