简单回顾一下多元回归分析里的参数估计问题：

对于多元线性回归模型 Y=Xm{eta}+varepsilon 的参数，我们有 least-square 估计量 m{hat{eta}}=(X^TX)^{-1}X^TY ，从而得到模型： hat{Y}=Xm{hat{eta}} = X(X^TX)^{-1}X^TY

需要注意的是建立以上模型的假设条件。其中最重要的是，对于随机误差向量 varepsilon ， 满足Evarepsilon=0 且有 Gauss-Markov 假设： [Cov(varepsilon_i,varepsilon_j)=left{egin{array}{ll} sigma^2& ext{$i=j$},\ 0& ext{$i≠j$,} end{array} ight.] 。在进行线性模型有效性检验和因变数与自变数的线性关系检验以及预测值区间估计时，还需要假设 varepsilon_i 服从正态分布。

回归诊断要研究的第一个问题就是考察我们的数据是否符合这些假设，如果假设不成立，探讨如何对数据进行修正以使它们（近似）满足这些假设。另一个重要问题是检测并处理对统计推断造成较大影响的数据点，也称为强影响点（influence case）。

对于后一个问题，可以在网上找到很多资料，包括基于统计推断的，基于距离的，基于聚类演算法的检测方法。这里只做简单的讨论。回归分析中，我们希望每组数据对未知参数的估计或其他推断有一定影响，但影响不要过大，否则得到的估计的稳定性较差。但是当我们检测出强影响点后，也不能一概的认为含强影响点的回归分析结果是不好的。如果强影响点是由系统误差造成的，这样的异常点（outliers）可以直接筛除掉。而如果强影响点由系统本身的运动变化产生，则它们很可能暗含著一些重要信息，如特殊环境下的变化趋势，解释新环境的出现等，这时可以采用各种稳健估计的方法缩小强影响点的作用，或者学习一个新的预测模型捕捉它们的信息，或者增加数据量。

这篇文章将从以残差为诊断统计量，以残差分析的角度讨论前面提到的第一个问题。

· 检测高杠杆点

若GM假设成立，且 varepsilon_1,...,varepsilon_n 有联合正态分布，则 m{ varepsilon }sim N(0,sigma^2m{I}) 。残差定义为 e=Y-hat{Y} ，服从 N(0,sigma^2(1-H)) ，并且 cov(hat{m{eta}}, e)=0 。其中，帽子矩阵 H= X(X^TX)^{-1}X^T（因为它对 Y 进行映射后变成 hat{Y} ）。

如果假设不成立，随机误差向量 m{varepsilon} 的分布就会发生变化，从而导致残差 e 的分布变化。残差分析的主要思想就是通过对 e 的表现进行研究，从而推断出关于 m{varepsilon} 的假设有哪些不满足的地方。

由 e 的分布可以知道，帽子矩阵 H 完全决定了残差的方差和协方差，所以可以猜测， H 将在回归诊断中起到重要作用。事实确实如此：

  • H 的对角线元素 h_{ii} 为第 i 个样本 x_i 到样本中心的 Mahalanobis 距离（加上常数1/n）；
  • 第 i 个样本的预测残差的方差是 h_{ii}(0leq h_{ii}leq1) 的单调递减函数，即Var(e_i)=sigma^2(1-h_{ii}) ，而且当h_{ii}=1时， Var(e_i)=0 。

这些事实说明，当 h_{ii} 很大，即当数据点距离样本中心很远的时候，残差趋于0。从几何上反映为，距离很远的点把回归直线拉向它自己。这种数据点称为高杠杆点（high leverage case），它们对参数的估计的影响很大。由这些结论知道，可以基于手头上的数据计算出矩阵 H ，然后根据对角元的取值大小检测高杠杆点是否存在。而如果 h_{ii} 的取值很大，但其对应的样本的残差也很大，说明当前数据不满足假设。

究竟 h_{ii} 大到多少算是取值较大呢？这个没有普遍适用的标准。一种做法是，把样本看做服从正态分布的的随即向量 X=(X_1,...X_p) 的简单随机样本（独立同分布），则有枢轴量 F=frac{n-p-1}{p}frac{h_{ii}-frac{1}{n}}{1-h_{ii}} sim F(p, n-p-1) ，h_{ii}很大等价于 F 很大。对给定置信水平 alpha ，从已知的 h_{ii} 算出的 F_{1-alpha}(p, n-p-1) 则认为h_{ii}在相应置信水平下显著的为较大值，对应的样本点 (x_i,y_i) 为高杠杆点。

· 残差的分布

我们直观地希望直接根据残差的大小来判定异常点，如果知道了残差的分布，我们就可以根据分位数等方式判别某个样本计算的残差是否统计意义上的异常。由于残差的方差与因变数的度量单位以及 h_{ii} 有关，所以将它们标准化： frac{e_i}{sigmasqrt{1-h_{ii}}} . 但是其中误差的方差 sigma 未知，所以用估计值： hat{sigma}=frac{||Y-Xm{hat{eta}}||^2}{sqrt{n-p-1}} ，得到「学生化残差」 r_i = frac{e_i}{hat{sigma}sqrt{1-h_{ii}}} 。

可惜的是 r_i 并不服从学生 t 分布，诸 r_i 也不彼此独立。不过，当假设 m{ varepsilon }sim N(0,sigma^2m{I}) 成立时， Er_i=0 ， Var(r_i)=1 ，在应用时将诸 r_i 看做来自 N(0,1) 的简单随机样本，然后用残差图进行分析诊断。

还有另一个学生化残差定义 {r_i^*} = frac{e_i}{hat{sigma(i)}sqrt{1-h_{ii}}} ，这里的 hat{sigma(i)} 是从剔除掉第 i 个样本后计算得到的回归模型导出来的。这样处理的意思是，对残差标准化时，排除自身在误差方差估计中的份额。在许多应用场合， r_i^* 近似服从 t(n-p-2) 分布；而在满足关于随机误差的假设时，则是严格服从 t(n-p-2)。 用于检测异常值时， r_i^* 相比 r_i 更有效。

· 预测残差

以上提到的残差都是从「拟合角度」提出的。因为回归模型最重要的应用在于预测，有必要从「预测角度」定义在 x_i 处的预测残差 delta_i=Y_i-x_i^That{eta(i)} ，这里 hat{eta(i)} 同样是在去掉第 i 个样本后，回归模型的估计参数。也就是说，前面的拟合残差是机器学习中通常所说的「训练误差」，预测残差则是「预测误差」。

预测残差和拟合残差具有关系： delta_i=frac{e_i}{1-h_{ii}} . 所以对于 h_{ii} 较大的样本点，其预测残差也较大。所以，根据前面的讨论，使用预测残差检测高杠杆点或检验假设合理性时，会更看重 h_{ii} 较大的数据，也就是远离样本中心的数据。

在假设 m{ varepsilon }sim N(0,sigma^2m{I}) 成立时，有 {delta_i}sim N(0, sigma^2/(1-h_{ii})) ，所以对预测残差进行标准化处理，就得到前面讨论的学生化残差 r_i 和 r_i^* .

预测残差还可以用于自变数选择。一种 PRESS (prediction error sum of square)准则认为，好的回归模型应该具有较小的 Sigma_idelta_i^2 。

（有空后面再补充一点残差图的内容。）


推荐阅读：
相关文章