【数理统计基础】 03 - 参数估计

合集 - 随机数学(11)

1.【初等概率论】 01 - 不确定中的确定性2017-02-072.【初等概率论】 02 - 条件概率和随机变量2017-02-073.【初等概率论】 03 - 常见分布和问题举例2017-02-074.【初等概率论】 04 - 数字特征2017-02-075.【初等概率论】 05 - 极限定理和正态分布2017-02-076.【数理统计基础】 01 - 关于预测的方法学2017-05-227.【数理统计基础】 02 - 统计量和三大分布2017-05-22

8.【数理统计基础】 03 - 参数估计2017-05-22

9.【数理统计基础】 04 - 假设检验2017-05-2210.【数理统计基础】 05 - 回归分析2017-05-2211.【数理统计基础】 06 - 相关分析和方差分析2017-05-24

收起

　　数理统计要解决的问题是，根据样本的信息猜测随机变量的信息。随机变量的分布可能完全未知，也可能已经判定为某类分布$f(x,{\theta }_1,\cdots ,{\theta }_k)$，但有未知参数$\overline{\theta }=({\theta }_1,\cdots ,{\theta }_k)$，这是数理统计中最常研究的情景。

1. 点估计

　　一类最简单的问题是，要求给出参数函数$g(\overline{\theta })$的一个统计量估计$\hat{g}(X_1,\cdots ,X_n)$。因为对于某次试验，估计量$\hat{g}$是确定的值，它被称为$g$的点估计。需要估计的函数$g$往往是参数$\overline{\theta }$本身，这里只讨论常用分布的参数估计。

1.1 矩估计

　　点估计就是值的近似，而分布的矩和样本矩正是一对现成的近似关系，而且大数定理也保证了它们的极限关系。被估函数$g$和估计量$\hat{g}$分别取矩和样本矩的方法就叫矩估计法，只需要有$k$个不同的联立方程（1）便能得到参数的估算值${\hat{\theta }}_1,\cdots ,{\hat{\theta }}_k$（但不一定所有时候都要求估计参数）。而常见的分布只有一个或两个参数（正态分布），往往用样本均值$\overline{X}$或样本方差$S^2$就能得到很好的估计。

$$\begin{array}{}\text{(1)}& {\hat{g}}_i(X_1,\cdots ,X_n)=g_i({\theta }_1,\cdots ,{\theta }_k),(i=1,\cdots ,k)\end{array}$$

　　对于常用分布的矩估计，这里就不一一列举了，计算并无本质困难。这里仅给出均匀分布$[{\theta }_1,{\theta }_2]$的矩估计（均值和方差）结果（式（2）），请注意和其它方法的比较。

$$\begin{array}{}\text{(2)}& {\hat{\theta }}_1=\overline{X}-\sqrt{3}S;{\hat{\theta }}_2=\overline{X}+\sqrt{3}S\end{array}$$

1.2 最大似然估计

　　矩估计法是一个非常符合直观的方法，但并不是唯一的方法，换一个思维，也能找到别的“合乎情理”的方法。站在概率论的角度，我们希望能找到参数“可能性最大”的值，但关于可能性的度量还比较含糊。教材上直接告诉我们，如果分布的密度函数是$f(x,\overline{\theta })$，则整个样本的密度函数为式（3），对于确定的采样$(x_1,\cdots ,x_n)$，合适的${\hat{\theta }}_i$应当使得$L(\overline{\theta })$达到最大值。用这样的${\hat{\theta }}_i$作为$\overline{\theta }$的点估计的方法叫做最大似然估计法，$L(\overline{\theta })$则称为似然函数。

$$\begin{array}{}\text{(3)}& L(x_1,\cdots ,x_n,{\theta }_1,\cdots ,{\theta }_k)=\prod _{i=1}^{n}f(x_i,{\theta }_1,\cdots ,{\theta }_k)\end{array}$$

　　这个方法看起来很合理，但和矩估计法差别太大，甚至让人担心会得出相矛盾的结论。并且如果仔细推敲，这样的似然函数其实并不合理，因为$L(\overline{\theta })$最大和最可能的$\overline{\theta }$还是有很大差别的，$L(\overline{\theta })$并不能轻易看成$\overline{\theta }$的密度函数。讲到这里，其实我们已经触及到数理统计中的一个争论点了，就是$\overline{\theta }$究竟应该看成确定值还是随机变量？

　　矩估计法就是把$\overline{\theta }$看做确定值，这也是符合直觉的。但如果非要讨论最大可能性的$\overline{\theta }$，就不得不把它当做随机变量看待，这就是贝叶斯思想。似然函数本质上应当是个条件概率$P(A|B)$，条件$B$就是观察值$x_1,\cdots ,x_n$，但初始概率$P(A)$是什么？这就是问题的关键，$\overline{\theta }$应当有个初始分布，答案很简单，初始的$\overline{\theta }$默认是均匀分布的。这正是最大似然法适用的场合与实际意义，使用时一定要确保这个假设是成立的，详细的贝叶斯法在下一段展开讨论。

　　为了计算方便，一般是求解$\ln L(\overline{\theta })$的极值，即求联立方程（4）的解，得到的是函数不动点，有时还要论证是否为最值。常用分布的最大似然估计大多与矩估计法的结果一样，这是一种巧合，但也说明了最大似然估计法的有效性。不过也有与矩估计法不同的，比如正态分布的方差，得到的估计是$m_2$，而非修正的样本方差。再比如均匀分布$[{\theta }_1,{\theta }_2]$，由于密度函数非零的部分是$({\theta }_2-{\theta }_1{)}^{-n}$，显然在${\theta }_2-{\theta }_1$最小时取得最大值，故有式（5）的估计。

$$\begin{array}{}\text{(4)}& \frac{\partial [\ln L({\theta }_1,\cdots ,{\theta }_k)]}{\partial {\theta }_i}=0,(i=1,\cdots ,k)\end{array}$$

$$\begin{array}{}\text{(5)}& {\hat{\theta }}_1=min\{X_i\};{\hat{\theta }}_2=max\{X_i\}\end{array}$$

1.3 优良性准则

　　对于同一个参数，可以有不同的点估计，在具体的场景下应当如何选择，是很重要的问题。在制定准则时，有两点需要注意：一是判定准则是根据具体需求制定，好坏并不是绝对的；二是判定准则往往是针对全体样本的，某个具体样本的好坏不足以说明问题。

　　最简单的准则就是无偏性，它要求$E(\hat{g}(\theta ))=g(\theta )$。无偏性适用于多次误差可以补偿的情况，比如买东西的重量，误差造成的双方损失可以互补。但对于精度要求高的场景，还希望$\hat{g}(\theta )$尽量聚拢在$g(\theta )$周围，也就是说它的方差还要尽量小。在所有无偏估计中，方差最小的称为最小方差无偏估计，简称MVU估计。

　　MVU估计比较难找，甚至根本不存在。有一个朴素的思想是，如果能得到$D(\hat{g}(\theta ))$的一个下界，并且正好找到了这样的$\hat{g}(\theta )$，那么就是找到了MVU估计。这个看似异想天开的方法，居然还真有比较好的结论，下面来看看。结论的灵感来自于不等式${\text{Cov}}^2(\xi ,\eta )⩽D(\xi )D(\eta )$，其中等式成立的充要条件是$\xi ,\eta$（中心化后）有简单的线性关系。

　　构造的思路是这样的，需要选择一个统计量$G$，它使得$\text{Cov}(\hat{g},G)$是与$\hat{g}$无关的常量，然后所有$\hat{g}$中还存在与$G$有简单线性关系的统计量。对$\hat{g}$的唯一限制条件是其期望为$g$，它也是唯一可利用的等式，因此$G$要取与样本分布密度$p=\prod f(X_i,\theta )$有关的函数。这个问题正面求解似乎很难，我们不妨从最简单的场景入手，就拿正态分布的均值估计$\overline{X}$为例，与$p$有关且与$\overline{X}$有线性关系的量是式（6）中的$G$。

$$\begin{array}{}\text{(6)}& G(X_1,\cdots ,X_n)=\sum _{i=1}^n[\ln f(X_i,\theta ){]}_{\theta }^{\prime }=\sum _{i=1}^n{\displaystyle \frac{f_{\theta }^{\prime }(X_i,\theta )}{f(X_i,\theta )}}\end{array}$$

　　可以计算得到，$\text{Cov}(\hat{g},G)=g^{\prime }(\theta )$，然后还能得到$E(G)=0$，接下来由$X_i$的独立性可以把$D(G)$记作$nI(\theta )$。最终得到式（7）的克拉美-劳不等式，其中$I(\theta )$被称为费歇尔信息量（如果$X_i$是离散的也有类似的表达式），也许在信息论中会有更好的阐述，这里就不探究了。证明中还需满足一些一致性的要求，这里也省略了，请自行参考教材。

$$\begin{array}{}\text{(7)}& D[\hat{g}(X_1,\cdots ,X_n)]⩾{\displaystyle \frac{[g^{\prime }(\theta ){]}^2}{nI(\theta )}};I(\theta )={\int }_{-\mathrm{\infty }}^{\mathrm{\infty }}{\displaystyle \frac{[f_{\theta }^{\prime }(x,\theta ){]}^2}{f(x,\theta )}}\text{d}x\end{array}$$

　　式（6）对任何统计量都成立，我们更应当关注等号成立的条件，即$\hat{g}$和$G$有线性关系。对于均值估计$\overline{X}$，要想它是MVU估计，只需$[\ln f(x,\theta ){]}_{\theta }^{\prime }$是$x$的线性函数（当然还要验证一致性，这里略去），容易验证常见分布一般都满足这个条件。另外还可以证明，方差$S^2$也是${\sigma }^2$的MVU估计。

2. 区间估计

　　参数估计的目的是对参数更多的了解，点估计的结果虽然直接易用，但却丢失了太多参数的信息，使用上也没有灵活性。为了包含参数的更多信息，我们希望找到两个统计量${\hat{g}}_1,{\hat{g}}_2$，以区间形式估计参数，并达到一定的概率要求。一般是对给定足够小的$\alpha >0$，要找到尽量小的区间$[{\hat{g}}_1,{\hat{g}}_2]$，使它能以$1-\alpha$的概率包含$g(\theta )$（式（8））。

$$\begin{array}{}\text{(8)}& P[{\hat{g}}_1(X_1,\cdots ,X_n)⩽g(\theta )⩽{\hat{g}}_2(X_1,\cdots ,X_n)]=1-\alpha \end{array}$$

　　这样的估计方法叫区间估计，其中$[{\hat{g}}_1,{\hat{g}}_2]$叫置信区间，而$1-\alpha$是区间的置信系数。有两点需要强调：一个是这里仍然是把参数当做确定值，把样本当做随机变量，所以置信系数的意义是“区间能包含参数”的概率，而非“参数落在区间里”的概率；另一个是区间长度越小越好，但不做强求，因为区间本身就是随机变量，对它最小值的讨论比较困难。

　　为了构造统计量${\hat{g}}_1,{\hat{g}}_2$，观察式（8），其中只包含待估参数和样本值，以及它们之间的概率不等式。一种比较方便的构造方法是这样的，找一个变量$G(g(\theta ),X_1,\cdots ,X_n)$，它服从一个比较简单的分布$F$。为了生成置信区间，一般把变量值限定在$E(G)$的两侧，每测的概率分别取$(1-\alpha )/2$。如果用$f(\alpha )$表示$F$的$\alpha$分位点，则建立不等式（9），整理后便能得到式（8）置信区间。

$$\begin{array}{}\text{(9)}& E(G)-f(1-\frac{\alpha }{2})⩽G(g(\theta ),X_1,\cdots ,X_n)⩽E(G)+f(\frac{\alpha }{2})\end{array}$$

　　这里的$G$就是前面提到过的枢轴变量，因此该方法也叫枢轴变量法。很多分布的枢轴变量比较难构造或者计算量大，甚至有时对分布完全未知，这时如果样本足够大，可以利用中心极限定理，以标准正态分布作为枢轴变量。这里我们只讨论正态分布，它的常用枢轴变量还有上篇介绍的三大变量，请先回顾相关性质。以下讨论仅给出枢轴变量，具体置信区间请自行计算，并没有本质困难。

　　先讨论单样本的正态分布$X\sim N(\mu ,{\sigma }^2)$。估计均值$\mu$时，$\sigma$可能已知也可能未知，上篇的公式（8）和（15）便是对这两种情况的枢轴变量。再估计方差${\sigma }^2$，同样分为$\mu$已知和未知两种情况，$\mu$已知的情况比较简单，未知时上篇的公式（17）便是我们要的枢轴变量。

　　再讨论两样本的正态分布，一般是根据两个随机变量的观察值，比较它们的参数。设两个随机变量为$X\sim N({\mu }_1,{\sigma }_1^2)$和$X\sim N({\mu }_2,{\sigma }_2^2)$，样本分别是$X_1,\cdots ,X_m$和$Y_1,\cdots ,Y_n$。一种是要考察${\mu }_1-{\mu }_2$的大小，一般的做法当然是用$\overline{X}-\overline{Y}$去估计它。当${\sigma }_i$都已知时，$\overline{X}-\overline{Y}$的方差为${\sigma }^2={\displaystyle \frac{{\sigma }_1^2}{m}}+{\displaystyle \frac{{\sigma }_2^2}{n}}$，枢轴变量比较显然。

　　当${\sigma }_i$都未知时，暂时没办法把$\sigma$消除（即使用$S_1^2+S_2^2$也不行），这里只讨论${\sigma }_1={\sigma }_2$的场景。为了能使用$t$分布，直接使用式（10）中的$S^2$来近似方差，容易得到枢轴变量（11）。当${\sigma }_1\ne {\sigma }_2$时，暂时没有完美的解决方法，该问题称为贝伦斯-费歇尔问题。

$$\begin{array}{}\text{(10)}& S^2={\displaystyle \frac{\sum _{i=1}^m(X_i-\overline{X}{)}^2+\sum _{i=1}^n(Y_i-\overline{Y}{)}^2}{m+n-2}}\end{array}$$

$$\begin{array}{}\text{(11)}& {\displaystyle \frac{(\overline{X}-\overline{Y})-({\mu }_1-{\mu }_2)}{\sqrt{1/m+1/n}\cdot S}}\sim t_{m+n-2}\end{array}$$

　　最后来比较$X,Y$的方差，直接作差比较难处理，而且意义也不明显。一般是估计${\sigma }_1^2/{\sigma }_2^2$的大小，它可以直接使用上一篇的式（19）作为枢轴变量。

3. 贝叶斯估计

　　现在来正式讨论贝叶斯估计，它的模型直接从事件的条件概率扩展而来，只不过由事件概率扩展为分布密度（同样适用于离散分布）。贝叶斯法的最大特点就是把参数$\theta$看做一个随机变量，如何理解这一点非常关键。现实中参数$\theta$一定是确定的，只不过我们不知道它的信息。但根据过去的认识或合理的假设，对$\theta$的所有可能值会有个评估，这样的评估就使得$\theta$有了随机变量的性质。需要着重强调的是，随机变量不是变量，它只是对不同可能性的一种描述。

　　脑洞再开大一点，则可以认为我们以一定概率处在不同的平行时空中，而$\theta$在每一个时空中都有一个确定的值。在得到观察$X_1,\cdots ,X_n$后，我们需要重新评估处在不同时空的概率。这是典型的条件概率问题，但要注意，这时讨论的样本空间是$\theta ,x_1,\cdots ,x_n$。假设随机变量$X$的密度函数为$f(x,\theta )$，参数$\theta$的先验分布的密度函数为$h(\theta )$，容易得到$\theta$的后验概率的密度函数（式（12））。

$$\begin{array}{}\text{(12)}& h(\theta |x_1,\cdots ,x_n)={\displaystyle \frac{p(\theta )}{\int p(\theta )\text{d}\theta }},p(\theta )=h(\theta )\prod _{i=1}^{n}f(x_i,\theta )\end{array}$$

　　前面说过，最大似然估计本质上也是贝叶斯思想，只是先验分布采用的是均匀分布。这里有个很现实的问题，如何在无限区间（比如整个实数域、所有正数等）上定义均匀分布？这个对我们还是太困难，也许测度论中会有完美解释？这就不得而知了。但还是可以换个思路，$h(\theta )$在条件分布中本质上起到的是“权重”的作用，也就是说它的根本意义在于表明概率之间的“比重”。比如对于均匀分布，只需取$h(\theta )=1$就足以说明“均匀”的性质，不必要求$h(\theta )$是一个严格的密度函数（但由式（12）易知后验概率一定是密度函数）。

　　但无论如何，在式（12）的使用过程中，必须先给出一个先验分布$h(\theta )$，这个分布的选择非常影响估计结果。很多时候先验分布难以确定，只能凭主观经验，这给贝叶斯方法带来了很多诟病。但由于贝叶斯思想的有效性和方便性，它在数理统计中仍然大行其道，甚至形成了所谓的贝叶斯学派，以区别于坚持频率方法的学者。一种和解的方法是承认两个模型本质的不同，并且互相补充、互相学习。但以个人粗浅的了解，我觉得贝叶斯思想是对传统模型的扩充，它是用先验概率把传统模型补充完整而已。这个补充就如同虚数于实数系统一样，是打破直觉却非常必要的抽象，是现代数学所具有的特征。

　　贝叶斯模型是完整的，且逻辑自洽的，方法本身不应该被诟病。既然问题出在先验概率的选择，那么在使用时挑选最合适的即可。这就是另外一个问题了，需要更多的理论分析和支持，不能把这部分工作的欠缺怪罪于贝叶斯模型本身。在大部分场合，$h(\theta )$一般遵循“同等无知”原则，这个原则的缺点也是显然的：如果对$\theta$是同等无知的，则对它的函数则基本不满足。这时选择对谁同等无知就很关键，在正态分布中一般取$h(\sigma )={\sigma }^{-1}$，在指数分布中一般取$h(\lambda )={\lambda }^{-1}$。

　　后验概率可以看做包含了参数的所有信息，它可以被用作点估计和区间估计。在点估计时，最合理的应当是取期望值，而非最大似然法中的最大值。伯努利分布（先验概率取均匀分布）在贝叶斯方法中的期望值是$\hat{p}={\displaystyle \frac{N+1}{n+2}}$，这在$n$很小时明显更合理（最大似然法得到${\displaystyle \frac{N}{n}}$）。

　　后验概率上的区间估计实现起来非常方便，它有统一的计算过程，而不依赖于具体分布。针对指定的置信系数，寻找最小的置信区间，是存粹的分析学计算问题。后验分布有对称轴时，最小置信区间一般也是对称的。其它情况下，可以先固定一个边界以确定另一个边界，然后变动第一个边界寻找最小区间。实在复杂的计算，也可以直接交给计算机完成。最后需要提醒一下，这里的置信区间和第二段中的置信区间有着本质的区别，一个是$\theta$自身（随机变量）的取值区间，一个是可以包含$\theta$（确定值）的区间，请仔细体会。