【概率论与数理统计】小结9-2 - 点估计
注：点估计是参数估计中的一种。点估计常用的方法有两种：矩估计和最大似然估计。之所以要做估计，最本质的问题是我们能获得的信息量（样本的数量）有限，因此只能在有限的信息中，用合理的方法、在可接受的精度或置信度下做近似计算，以便对总体有一个大概的认识，也就是将某种在有限样本中获得的规律，推广到更大的样本量上。
0. 矩估计
0.1 引例
某大学新生有4000人参加第一学期的《微积分》考试。现随机选出100名同学，计算得到他们的平均成绩为68.5分，标准差为12.2分，试估计全体同学的平均成绩。
如果我们使用这100名同学的平均成绩来估计全体同学的平均成绩，那么我们就完成了一次矩估计。
- 记总体（所有4000名学生）的平均成绩为$\mu$，此时的$\mu$就是总体的一阶原点矩$\alpha_1$；
- 100名学生的平均成绩就是样本的一阶原点矩$A_1$；
- 用样本的一阶原点矩估计总体的一阶原点矩，得到$\hat{\mu} = 68.5$.
根据小结5中的约定，我们使用以下符号表示不同的矩：
- $\alpha_k$: 总体k阶原点矩；
- $\beta_k$: 总体k阶中心矩；
- $A_k$: 样本k阶原点矩；
- $B_k$: 样本k阶中心矩.
0.2 统计思想
矩估计是英国统计学家卡尔·皮尔逊于1894年提出的。是基于一种简单的 “替换” 思想建立起来的一种估计方法。其基本思想是以样本矩估计总体矩，以样本矩的函数估计总体矩的函数。
0.3 理论依据
辛钦大数定律和依概率收敛的性质。
假设$\alpha_j = E(X^j)$存在，其中$j = 1, 2, ..., k$. 则
$$\hat{\alpha_j} = A_j = \frac{1}{n}\sum_{i=1}^{n}{X_i^j}, \ j = 1, 2, ..., k, \ \to \ \alpha_j \ with \ probability \ p , \ j = 1, 2, ..., k$$
或
$$h(\hat{\alpha_1}, \ \hat{\alpha_2}, \ ..., \ \hat{\alpha_k}) = h(A_1, \ A_2, \ ..., \ A_k) \ \to \ h(\alpha_1, \ \alpha_2, \ ..., \ \alpha_k) \ with \ probability \ p$$
上面用公式表示了用样本的原点矩估计总体的原点矩或用样本原点矩的函数估计总体原点矩的函数，将原点矩换成中心矩也同样成立。
0.4 建立矩估计的步骤
简单来说，只要可以将一个待估计的参数用总体矩（单个矩或不同矩的函数）表示出来，然后用对应的样本矩替换总体矩就可以了。
设总体有k个未知参数，$\theta_1, \theta_2, ..., \theta_k$, $X_1, X_2, ..., X_n$是来自总体$X$的样本，假设总体的前k阶矩存在。
矩估计步骤：
（1）建立$(\theta_1, ..., \theta_k)$与$(\alpha_1, .., \alpha_k)$的联系：求总体前$k$阶原点矩关于$k$个参数的函数(即，用未知参数表示总体矩)，
$$\alpha_i = E(X^{(i)}) = h_i(\theta_1, ..., \theta_k), i = 1, ..., k.$$
（2）求各参数关于$k$阶原点矩的反函数(即，用总体距表示未知参数)，
$$\theta_i = g_i(\alpha_1, ..., \alpha_k), i = 1, ..., k$$
（3）以样本各阶矩$A_1, ..., A_k$代替总体$X$各阶矩$\alpha_1, ..., \alpha_k$，得到各参数的矩估计：
$$\hat{\theta_i} = g_i(A_1, ..., A_k)$$
在实际应用时，为求解方便，也可用总体中心矩$\beta_i$替换总体原点矩$\alpha_i$，相应的，以样本中心矩$B_i$估计总体中心矩$\beta_i$.
Notice: 采用的矩不同，得出的矩估计也可能不同。
0.5 常见的矩估计量
参考小结7，part2部分，总结了常见的统计量（所有的样本矩都是统计量，有些统计量本身也是样本矩）以及使用样本矩估计总体矩的例子。
0.6 总体方差的估计
总体方差的估计是一个比较特殊的例子，值得单独拿出来说一说。这是因为总体方差与总体的二阶中心矩相同，但是样本方差的计算公式与样本的二阶中心矩并不相同。
设$X$是一个随机变量，且方差存在，总体和样本的均值分别为$\mu, \ \bar{X}$. 则下面是该随机变量4种不同的数字特征的计算方式：
- 总体方差用$\sigma^2$表示，按照方差的定义可得$\sigma^2 = E[X - E(X)]^2 = \frac{1}{n}\sum_{i=1}^{n}{(X_i - \mu)^2}$;
- 按照中心矩的定义可得，总体二阶中心矩$\beta_2 = E[X - E(X)]^2 = \frac{1}{n}\sum_{i=1}^{n}{(X_i - \mu)^2}$，同$\sigma^2$;
- 同样，样本二阶中心矩$B_2 = E[X - \bar{X}]^2 = \frac{1}{n}\sum_{i=1}^{n}{(X_i - \bar{X})^2}$;
- 样本方差$S^2 = \frac{1}{n - 1}\sum_{i=1}^{n}{(X_i - \bar{X})^2}$，参考小结5中，2.3 样本方差的计算.
此时我们有两种方式来对总体方差进行估计，样本的二阶中心矩（矩估计）或样本方差（非矩估计）。那么哪一个估计更好呢？此时可以使用小结9中"点估计的评价"来比较这两种方法。参考wiki中的证明过程，我们可以知道：
- 使用样本的二阶中心矩估计总体方差时，估计量的期望$E(B_2) = (1 - \frac{1}{n})\sigma^2 \leq \sigma^2$，由此可以看到该方法的估计始终是偏小的，存在系统性误差；
- 使用样本的方差估计总体方差时，估计量的期望$E(S^2) = \sigma^2 = \sigma^2$，因此样本的方差是总体方差的无偏估计，比样本的二阶中心矩更优。
样本二阶中心矩与样本方差之间，只有分母上有差别，当n比较小的时候，$n$与$(n-1)$之间的差别对结果的影响很大，但是当$n \to +\infty$时，这种差异就可以忽略了。以上就是我们为什么使用样本方差而不是样本二阶中心矩来估计总体方差的原因（此时并不是矩估计）。
1. 最大似然估计
1.1 引例（基本原理）
考察以下例子：假设一个罐中放着许多白球和黑球，并假定已经知道两种球的数目之比是$1 : 3$，但不知道哪种颜色的球多。如果采用有放回抽样方法从罐中取5个球，观察结果为：黑、白、黑、黑、黑，估计取到黑球的概率$p$.
在该题目中，总体服从参数为p的伯努利分布。
解：设$X = 1, 取到黑球; X = 0, 取到白球$, 则$X \sim B(1, p)$.
p为取到黑球的概率，未知，$p = 1/4 或 3/4$.
抽取容量为5的样本$X_1, ..., X_5$，观察值为$1, 0, 1, 1, 1$.
当$p = 1/4$时，出现本次观察结果的概率为$(\frac{1}{4})^4(\frac{3}{4}) = \frac{3}{1024}$;
当$p = 3/4$时，出现本次观察结果的概率为$(\frac{3}{4})^4(\frac{1}{4}) = \frac{81}{1024}$.
由于$3/1024 < 81/1024$, 因此认为$p = 3/4$比$p = 1/4$更有可能，于是$\hat{p} = 3/4$更为合理.
由上面的求解过程可以看到，我们尝试了所有可能的取值（$p = 1/4$或$p = 3/4$），然后选择了一个最有可能产生出观察到的样本（$1, 0, 1, 1, 1$）的取值（$p = 3/4$）。
1.2 定义
设离散型总体$X \sim p(x; \theta), \theta 未知$. $X_1, ..., X_n$为样本，其观察值为$x_1, ..., x_n$，则事件$\{X_1 = x_1, ..., X_n = x_n\}$发生的概率为：
$$p(X_1,X_2, ..., X_n) = L(\theta) = \prod_{i=1}^{n}p(x_i; \theta) \quad \cdots \ (1-1)$$
式1-1表示所有观察值$x_1, ..., x_n$同时发生的概率（所有单个观察值发生概率的乘积）. 这个乘积有一个专门的名字：似然函数.
参考1.1中的例子，我们希望参数$\theta$的估计值$\hat{\theta}$能够使得所有观察值同时发生的概率最大，也就是似然函数取到最大值：
$$L(\hat{\theta}(x_1, ..., x_n)) = max L(\theta) \quad \cdots \ (1-2)$$
式子1-2就是极大似然估计的基本原理：将似然函数取到最大值时参数的取值作为最优的估计值；反过来也可以尝试参数所有可能的取值，使得似然函数取到最大值的取值就是最优估计值。
上面的不同思路也反映了求解的两种方式：第一种是直接求似然函数的极大值；第二种是带入所有参数比较函数值的大小。
1.3 求最大似然估计的一般步骤（以下以离散总体为例）
1.3.1 写出似然函数
$L = L(\theta_1, \theta_2, ..., \theta_m) = L(x_1, x_2, ..., x_n; \theta_1, \theta_2, ..., \theta_m) = \prod_{i=1}^{n}p(X_i = x_i; \theta_1, \theta_2, ..., \theta_m) $
由于每个样本之间都是独立同分布的，因此这些样本的联合概率密度函数就是所有单个样本概率密度函数的乘积（即，这些样本所代表的事件同时发生）
1.3.2 对似然函数取对数
$ln L = \sum_{i=1}^{n}{ln P(X_i = x_i; \theta_1, \theta_2, ..., \theta_m)}$
取对数是为了计算方便（对多项式加法求导比对多项连乘求导简单）
1.3.3 建立似然方程
$$\frac{\partial lnL(\theta_1, ..., \theta_m)}{\partial \theta_j} = 0, (j = 1, ..., m)$$
对各个参数求偏导，令偏导数等于0.
1.3.4 解似然方程
解上面的方程，就可以求出参数$\theta_j$的最大似然估计。若似然函数不可微，需要用定义求。
1.4 几点说明
- 极大似然估计值是一个数值，极大似然估计量是一个随机变量（同时也是一个统计量）；
- 似然函数是能够反映试验结果发生的概率大小的函数；
- 从极大似然估计的定义可以看出来，使用该方法时总体的概率密度函数形式是明确的（包含未知参数）；
- 未知参数可能不止一个，设为$\theta = (\theta_1, \theta_2, ..., \theta_k)$；
- 若$L(\theta)$关于某个$\theta_i$是单调递增(减)函数，则$\theta_i$的极大似然估计为$\theta_i$的最大(小)值(与样本无关)。
1.5 经典例题
例1
设$X$的概率密度为$f(x, \theta) = \sqrt{\theta} x^{\sqrt{\theta} - 1}, 0 \leq x \leq 1; f(x, \theta) = 0, 其他$. $X_1, ..., X_n$是样本，求$\theta$的极大似然估计量.
若已获得$n = 10$的样本值如下：
0.43, 0.01, 0.30, 0.04, 0.54, 0.14, 0.99, 0.18, 0.98, 0.02
求$\theta$的极大似然估计值.
解：$L(\theta) = \prod_{i=1}^{n} f(x_i; \theta) = \prod_{i=1}^{n} \sqrt{\theta} x_i^{\sqrt{\theta} - 1} = \theta^{\frac{n}{2}} (\prod_{i=1}^{n} x_i)^{\sqrt{\theta} - 1}$
取对数后，$ln \ L(\theta) = \frac{n}{2} ln \ \theta + (\sqrt{\theta} - 1) \sum_{i=1}^{n}{ln \ x_i}$
建立似然方程，$\frac{d \ ln \ L(\theta)}{d \theta} = \frac{n}{2} \cdot \frac{1}{\theta} + \frac{1}{2 \sqrt{\theta}} \sum_{i=1}^{n}{ln \ x_i} = 0$
$\Rightarrow \frac{n}{\sqrt{\theta}} = - \sum_{i=1}^{n}{ln \ x_i} \Rightarrow \sqrt{\theta} = -n / \sum_{i=1}^{n}{ln \ x_i}$
最终可以得到参数$\theta$的极大似然估计量为：$\hat{\theta} = \frac{n^2}{(\sum_{i=1}^{n}{ln \ X_i})^2}$
将上面的样本值带入估计量，可以得到$\theta$的极大似然估计值为：$\hat{\theta} = 0.305$
例2
设总体$X \sim N(\mu, \sigma^2) $, $X_1, ..., X_n$是样本，$\mu, \sigma^2$均未知. 求$\mu, \sigma^2$的极大似然估计.
解：$L(\mu, \sigma^2) = (\frac{1}{\sqrt{2\pi\sigma^2}})^n \cdot e^{- \sum_{i=1}^{n}{\frac{(x_i - \mu)^2}{2\sigma^2}}}$
取对数得，$ln \ L(\mu, \sigma^2) = n \ ln\frac{1}{\sqrt{2\pi}} - \frac{n}{2}ln \ \sigma^2 - \sum_{i=1}^{n}{\frac{(x_i - \mu)^2}{2\sigma^2}}$
建立似然方程，$\frac{\partial }{\partial \mu}\ ln \ L(\mu, \sigma^2) = \frac{1}{\sigma^2} \sum_{i=1}^{n}(x_i - \mu) = 0$,
$\frac{\partial }{\partial \sigma^2}\ ln \ L(\mu, \sigma^2) = \frac{n}{2\sigma^2} + \frac{1}{2\sigma^4} \sum_{i=1}^{n}(x_i - \mu)^2 = 0$
解方程得，$\hat{\mu} = \bar{X}, \hat{\sigma^2} = B_2$
这里利用最大似然估计得到的结果与矩估计得到的结果相同。
2. 矩估计与最大似然估计的比较
矩估计和最大似然估计是两个不同的点估计方法。有时候这两种方法得到的结果相同，有时候不同。
使用矩估计时，不需要知道随机变量的概率密度函数（有些情况可以假设为正态分布），需要将未知参数表示为总体矩的函数；
使用最大似然估计时，需要知道随机变量概率密度函数的确定形式（可以包含参数），然后通过求似然函数的极值的方式得到最优解。
从它们的定义和求解过程来看，矩估计更简单一些但是适用范围也比较窄（无法表示为总体矩的函数就无法使用该方法）；最大似然估计需要的条件更多，适用范围也更广一些。
欢迎阅读“概率论与数理统计及Python实现”系列文章
Reference
https://en.wikipedia.org/wiki/Bias_of_an_estimator#Sample_variance
https://zh.wikipedia.org/wiki/%E7%9F%A9%E4%BC%B0%E8%AE%A1
中国大学MOOC：浙江大学&哈工大，概率论与数理统计
