下载掌阅APP，畅读海量书库
立即打开
线性回归模型（Linear Regression Model，以下简称LM模型）是最简单的机器学习算法，但也是很多机器学习算法的基础，很多模型都可以看作线性模型的拓展。因为是线性关系，统计检验理论在LM模型上做得也非常充分。深刻理解与掌握LM模型对于掌握其他机器学习理论也有很多帮助。R语言里面有很多算法包提供了与LM模型相关的函数，见表3-4。
表3-4 与LM模型相关的函数
在介绍线性模型算法前要清楚一个概念： 线性是什么？ 线性模型中的线性，并不指因变量与自变量是线性关系，而是指因变量与参数是线性关系。也就说对于自变量来说，完全可以先对其进行非线性变换，再进行线性组合。从这个角度来说，线性模型完全具有描述非线性的能力。例如， y = β 0 + β 1 x + β 2 x 2 + β 3 x 3 也是线性模型，因为这里的线性并不指 y 与 x 呈线性关系，而是指 y 与参数 β 0 、 β 1 、 β 2 、 β 3 是线性关系。由此也可以扩展出广义线性模型，每个组分可以非线性的，但组分间是线性关系，例如：
1）相加模型，不同组分间是线性关系。
Y t =Trend t +Seasonal t +Irregular t
2）相乘模型，经过Log运算，也可以转化为线性关系。
Y t =Trend t ×Seasonal t ×Irregular t
1.普通最小二乘算法
普通最小二乘法（Ordinary Least Squares，OLS）的目标是极小化残差（预测值与真实值的差）的平方和。也就是说，在一维因变量
y
与
p
个自变量
x
1
，
x
2
，…，
x
p
的线性函数结构约束下（可以支持自变量间多项式、交叉等，只要因变量与参数间是线性关系），寻找一组最佳的参数
β
0
，
β
1
，
β
2
，…，
β
p
（其中
β
0
称为截距，其他参数称为斜率），在目前数据集的
N
个样本上计算出预测值
，即
使得残差平方和（Reidual Sum of Square，RSS）最小：
在解析表达上，通常采用矩阵表示，在自变量矩阵中，增加一个恒定为1的截距项的首列，构成 N ×( p +1)的矩阵 X ，参数向量 β 是 p +1维：
y = Xβ
如果
X
是满秩，最优解
可以用矩阵形式表示为
对于线性算法的推导，很多算法图书都有详细的介绍 [10-12] 。对于lm()函数的使用，特别是其参数formula的使用方法，可以参阅《R语言实战》 [10] 、《应用预测建模》 [13] 。
2.参数的置信区间
参数的置信区间（Confidence Interval）指的是在一定的置信度（常用95%）样本均会落在该区间，在正态分布下，这意味着1.96倍标准误差范围。
下面借用《R语言实战》 [10] 一书中的例子，confint()函数的输出结果表明，谋杀率斜率参数的95%置信区间为［2.38，5.90］。另外，因为Frost的置信区间包含0且在0附近，可以得出结论，当其他变量不变时，严寒天数（Frost）与谋杀率无关。
3.参数与模型的显著度
显著性检验（significance test）用于判断实验处理组与对照组的差异是否显著的方法。常把一个要检验的假设记作H0，称为原假设（或零假设）（null hypothesis），与H0对立的假设记作H1，称为备择假设（alternative hypothesis）。有两类错误的可能，在原假设为真时，决定放弃原假设，称为第一类错误；在原假设不真时，决定不放弃原假设，称为第二类错误。显著性检验中通常只限定犯第一类错误的最大概率 α ，不考虑第二类错误概率。概率 α 称为显著性水平。0.05是常用的显著度水平，意味着有5%的可能犯第一类错误。
OLS模型的参数显著度采用t检验，模型整体采用F检验。以上面的模型为例：
可以看出，在 α =0.05的水平下，只有Population、Illiteracy两个变量参数是显著的，二者的第一类错误概率分别为0.017、0.000022，模型整体也是显著的。
除了线性拟合模型，两个变量的相关性也可以直接使用cor.test()函数检验（Pearson、Spearman和Kendall相关系数），使用格式为
4.模型的拟合优度
通常采用
R
2
来评价OLS模型的拟合优度，看残差平方和与总体变化的比例。考虑到统计自由度，有调整后的
，它比
R
2
要小一些，
考虑了模型的复杂度，在同样的RSS下（
R
2
相同），参数越多的模型，
越小。
为了和大部分统计书籍中的公式保持一致，这里的 p 表示待回归参数的个数（包括截距项）。
另外，需要注意
R
2
或
并不能表示模型的优劣，这个后面讨论。
除了lm()函数的标准输出的置信度等信息外，还可以对模型和变量做更深入的检验与分析 [14] 。除了交叉检验（后面将统一讨论）、特征变量选择（见第2章）等通用分析外，OLS模型还有特定的检验，这是本小节讨论的重点。
1.OLS的4个统计假设的检验
OLS算法基于如下4个统计假设，否则统计显著性检验、置信区间等结果可能不精确。R语言的gvlma包的gvlma()提供了线性模型的综合检验。
·正态性：对于给定的自变量值，因变量值呈现正态分布，或者说拟合误差符合正态分布；
·独立性： y i 之间相互独立；
·线性：因变量与自变量为线性相关；
·同方差：因变量的方差不随自变量的水平不同而变化。
1）正态性检验：Q-Q图通过把测试样本数据的分位数（quantile）与已知分布相比较，从而来检验数据的分布情况。Q-Q图是一种散点图，对应于正态分布的Q-Q图，就是由标准正态分布的分位数为横坐标，样本值为纵坐标的散点图。要利用Q-Q图鉴别样本数据是否近似于正态分布，只需看Q-Q图上的点是否近似地在一条直线附近，图形是直线说明是正态分布。OLS模型的Q-Q图检验一般针对残差数据，R里面提供了qqplot()函数。另外，可以对残差序列做统计分布（见3.2节），看残差是否符合正态分布。
2）独立性检验：残差独立性的检验依赖于数据收集方式。例如，时序数据通常存在时间自相关（相隔近的相关性大），常采用Durbin-Watson方法进行检验，R里面有durbinnWatsonTest()函数。
3）线性检验：成分残差图（component plus residual plot）或偏残差图（partial residual plot）展示各个变量与残差的关系，若存在显著的非线性关系，说明模型的自变量不够充分，需要引入一些非线性项。R里面有crPlots()函数。
4）同方差检验：同方差的反面是异方差。二者的检验，除了查看残差图（残差-单个自变量），还有Breusch-Pagan检验（R中的car包有ncvTest函数，lmtest包有bptest函数）、Goldfeld-Quandt检验（R中的lmtest包有gqtest函数）、Portmanteau Q检验等方法。
2.多重共线性的VIF检验
多重共线性指的是自变量存在着线性关系，会造成矩阵非满秩或接近病态，对矩阵求逆等数值计算不稳定。造成多重共线性有很多原因，例如，解释变量有共同的时间趋势，一个解释变量滞后于另一个变量，两者遵循一个趋势；数据集太小或者解释变量间近似为线性关系。
判断是否为多重共线性可以用cor()函数分析回归系数的正负号是否符合预期，观察解释变量是否 R 2 高而 t 值低，或者通过对解释变量的添加和删除来重新评估回归效果。Pearson相关系数的一个明显缺陷就是只能观测单个变量间的相关性。需要和多重共线性的检验结合使用。
对于多重共线性有一个常用的检验方法是VIF（Variance Inflation Factor，方差扩大因子）检验，它是容忍度（tolerance）的倒数。一般认为VIF值大于4时，存在严重的多重共线性。
自变量
x
j
用其他自动变量回归建模，假设其拟合优度为
，一个自变量
x
j
的VIF定义为
上面拟合模型中各个变量的VIF数值都不大，可以认为它们之间不存在多重共线性。
3.异常值
异常值包括离群（outlier）点、高杠杆（high leverage）值点、强影响（high influence）点。离群点指的是模型预测的 y 值与真实的 y 值相差非常大的点，所以也称为 y 方向的离群值（y-outlier）；高杠杆值点指的是 x 值比较异常，通常与响应变量值 y 没有关系，所以也常被称为 x 方向的离群值（x-outlier）；强影响点是指对统计推断有影响的点。在介绍完这3个量的概念后，辨析它们之间的关系。
（1）离群点
通常检测离群点的方法有：
1）用Q-Q图检测，落在置信区间外的点通常被认为是离群点；
2）标准差检测，通常认为残差数值超过2倍或3倍标准差是离群点；
3）IQR（Inter Quartile Range，四分位间距）检测，IQR是25%分位数和75%分位数的差，以中值为中心，超出1.5倍IQR的点通常认为是离群点。
outliers包提供了几个有用的函数来系统地检测出离群值。其中outliers()会返回和平均值相比较后最极端的观测；scores()函数不仅可以计算规范化得分（诸如 z 得分），它还可以基于得分值，返回那些得分在相应分布百分位数之外的观测值。car包也提供了一种离群点的统计检验方法。outlierTest()函数可以求得最大标准化残差绝对值Bonferroni调整后的 p 值。
（2）高杠杆值点
高杠杆值点可通过帽子统计量（Hat Statistics）判断。对于一个给定的数据集，帽子统计量的均值为 p/N ，其中 p 是模型估计的参数数目（包含截距项）， N 是样本量。一般来说，若观测点的帽子值大于帽子均值的2倍或3倍，即可以认定为高杠杆值点。
在假定 X 中各变量已经中心化的情形下，帽子矩阵的计算公式为
帽子统计值（Hat Value）就是矩阵 H 的对角线上的元素。在图3-14所示的例子中，Alaska与California是显著高杠杆值点，仔细查看具体数值，可以发现：Alaska收入比其他州高很多，人口和温度要低很多；California的人口、收入和温度都比其他州高很多。
图3-14 高杠杆值点的例子
（3）强影响点
对模型参数的估计产生的影响过大，非常不成比例，当移除模型的一个观测点时模型会发生巨大的改变，这就属于强影响点。强影响点可以通过库克距离（Cook Distance）、DFFITS、DFBETAS来鉴别。
库克距离描述了单个样本对整个回归模型的影响程度。库克距离越大，说明影响越大。库克距离也可以用来检测异常点。在理想的情况下，每个样本对模型的影响是相等的。如某个样本的库克距离非常大，这个点就是强影响点，默认阈值设置为4 / ( N-p -1)。
第 i 个样本的库克距离定义如下：
式中，
表示去掉第
i
个样本后训练的OLS对第
j
个样本的预测值；
是全体样本训练模型的预测残差的标准差。
从上面的概念看，高杠杆值点同时又是离群点的样本，肯定是强影响点。特别大的离群点很有可能是强影响点，高杠杆但非离群的点，可能影响很小。在单变量情形下，杠杆可以简单理解为 x 距离质心的距离。简单示意如图3-15所示。
R的car包中有一个influencePlot()函数，可以把离群点、高杠杆值点、强影响点都整合在一个图上。图横坐标为帽子统计值，可以判断哪些点是高杠杆值点；纵坐标为学生化残差，因此纵坐标超过+2或者-2的点被认为是离群点；图中越往右上角的点，越有可能是强影响点。
4.思考
这里简单探讨几个概念，目的是让读者对各种统计指标的适用前提有个客观认知，避免盲从。更多的讨论可查阅参考文献［15］。
图3-15 离群点、高杠杆值点、强影响点间的关系示意图
1. R 2 小的LM模型一定差吗
不， R 2 小的模型也可能是很好的模型。有些研究领域本身就有大量无法解释的变异，这时 R 2 一定会很低。例如，人类行为学解释模型的 R 2 通常都低于0.5。如果模型的 R 2 很低，但自变量在统计上是显著的，那么仍然可以得出关于变量之间关系的重要结论。在这里，研究的重点是相关系数是否能够表示因变量的单位变化。换句话说，如果主要目标是理解数据中关系的本质，低 R 2 值可能不是问题，但如果需要相对精确的预测，那么较低的 R 2 可能会造成问题。
2. R 2 大的LM模型一定好吗
不，一个 R 2 值高的回归模型可能会产生许多问题。下面举个例子，如图3-16和图3-17所示。
拟合线图中的数据遵循一个非常低的噪声关系， R 2 是98.5%。然而，回归线始终过低或过高预测沿曲线的数据，始终在偏离。残差与拟合图有更明显的体现。非随机残差模式表明，尽管 R 2 较高，但拟合效果并不好。
造成这一点的原因是，当线性模型未进行详细解释时，容易出现这种规律性偏差。换句话说，它缺少重要的独立变量、多项式项和交互项。为了产生随机残差，可以尝试向模型中添加项或拟合一个非线性模型。
3.LM模型统计显著就可信吗
显著性检验结论中的“差异显著”或“差异极显著”不应该误解为相差很大或非常大，也不能认为在实际应用上一定就有重要或很重要的价值。“显著”或“极显著”是指表面差异为试验误差可能性小于0.05或0.01，已达到了可以认为存在真实差异的显著水平。有些试验结果虽然表面差异大，但由于试验误差大，也许还不能得出“差异显著”的结论，而有些试验的结果虽然表面差异小，但由于试验误差小，反而可能推断为“差异显著”。
图3-16 拟合线图
图3-17 残差与拟合图
导致“差异不显著”有两种可能：一是无本质差异；二是有本质差异，但被试验误差所掩盖，表现不出差异的显著性来。如果减小试验误差或增大样本容量，则可能表现出差异显著性。显著性检验只是用来确定无效假设能否被否定（证伪），而不能证明无效假设是正确的。同时还应该结合实际问题背景来考虑。
OLS的目标函数是残差平方和最小，但容易受到离群点的影响，包括 y 方向的离群值（y-outlier）和 x 方向的离群值（x-outlier），鲁棒线性回归（Robust LM）是降低这些异常值影响的一种方法。
1.分位数回归
分位数回归（Quantile Regression）把分位数的概念融入到普通的线性回归。 τ （0< τ <1）分位回归意味着在回归曲线直线之下的误差距离和与总误差距离和的比例为 τ ，例如，所谓的0.9分位数回归，就是希望回归曲线之下的距离和是其上部的距离和的9倍。 τ =0.5就是L1误差指标（Least Absolute Error，LAE）。 τ 分位回归的目标函数（极小化）如下式所示：
这里借用quantreg包里面的例子，如下述代码和图3-18、图3-19所示，来演示rq()函数在不同 τ 参数下，线性回归系数的不同。
图3-18 不同参数tau下分位数回归的结果对比
图3-19 不同参数tau下回归模型的截距与斜率
2.M估计与MM估计
降低 y 方向离群点的一种方法就是改变误差的惩罚，适当降低对大偏差点的惩罚。M估计与MM估计是一种常用的方法，R中的MASS包里面有rlm函数对两种方法有很好的支撑。
M估计由Huber在1964年提出，对残差进行饱和函数运算，限制单个数据点对于误差函数的影响力。其思想类似于L1的稀疏性。使用的函数称为M估计子。优化目标是极小化M估计子函数 ρ (·)意义下全体样本误差惩罚，即
式中，
R
i
（
β
；
m
）表示数据集
m
第
i
个样本的残差；
是一个尺度系数，通常选取绝对中位差（Median Absolute Deviation，MAD）指标，
常用的M估计子包括Huber、Hampel、Tukey等函数，见表3-5。得分函数（Score Function） ψ 表示每个误差的实际惩罚，是目标函数（Score Function） ρ 的微分，权重函数（Weight Function） w 表示每个样本折合的权重（如果采用误差平方和的目标函数）。
表3-5 常用的M估计子函数
上面的优化问题，通常采用迭代权重最小二乘法（Iteratively Reweighted Least-Squares，IRLS）求解。可以采用普通最小二乘的结果作为初始解 [16] 。
从上面的公式可以看出，M估计可以降低
y
方向的离群点，但不能消除高杠杆值点（
x
方向的离群点）的影响，为此，Mallow提出了GM估计子（Generalize M-estimator），根据
x
i
的帽子统计值（Hat Value），在目标函数中对每个样本的偏差引入权重
w
（
x
i
），消除高杠杆值点的影响。另外还有S估计子等提法，详情请参阅参考文献［16］。在这些估计子的基础上，Yohai于1987年提出MM估计子，是目前最常用的鲁棒回归技术。它首先采用S估计子获得一个强壮的解
β
，计算其预测残差和尺度参数
，后面采用类似IRLS的迭代过程，只不过尺度参数
保持不变。
3.LMS与LTS算法
最小平方中位数（Least Median of Squares，LMS）由Rousseeuw在1984年提出，它的优化目标为残差平方的中值（而是所有残差平方的和），这样提高了对离群点的鲁棒性，但由于缺乏梯度信息，在数值计算的收敛性要慢很多，但LMS可以为MM估计提供一个很好的初始值。
截断最小二乘估计（Least Trimmed Squares，LTS）回归也是Rousseeuw在1984年提出的，它将目标函数变为前 N-q 个最小残差平方的和（其中 q 是配置参数），即将 N 个样本预测残差平方按从小到大排序，忽略最高的 q 个数值，这样避免离群点的影响。与LMS类似，LTS在计算效率上有一定损失，通常作为其他估计初始值的计算方法。
R中的MASS包里面的lqs()函数提供了LTS、LMS等算法支持。鲁棒性的度量指标包括影响函数（Influence Function，IF）、崩溃点（Breakdown Point，BP）等指标，不同方法在鲁棒性测度和计算效率的对比可以参阅参考文献［16］，LTS、LMS、S估计子的BP值（简单理解，可以容忍的极端异常样本的比例）比较高，也被称为有抵抗力的回归模型（Resistant Regression）。
正则化指的是在RSS之上添加一个模型结构复杂度（待估计的系数）的惩罚项，LASSO回归是加入一个L1范数，岭（Ridge）回归是加入L2范数，弹性网络（Elastic Net）是L1范数与L2范数的加权。如果数据变异较大，那么在做这些正则化回归前最好进行数据标准化处理，否则结构惩罚侧重于系数（绝对值）较大的项。
岭回归的损失函数为
LASSO回归的损失函数定义为
弹性网络中的结构惩罚是L1范数与L2范数的加权。R里面glmnet包有glmnet()函数可以支持三种算法，glmnet()函数的参数值alpha=0，为LASSO回归，alpha=1，为岭回归，在0与1之间为弹性网。glmnet包的cv.glmnet()函数利用交叉检验，分别用不同的lambda值来观察模型误差。
前面讲述的模型都是一次性对原始变量 x 进行的计算，可以采用多次对 x 变换的方法，包括主成分回归（Principle Components Regression，PCR）、偏最小二乘回归（Partial Least Square，PLS）。除了全局采用一套参数，也可以采用局部信息的非参数化方法。
1.PLS
当自变量的个数大于观测值的个数，OLS的矩阵求逆就不成立。在样本量偏少（特别在高维回归时）的时候，OLS的应用也比较受限。PLS可以部分解决这些问题。PLS利用的自变量与预测量的相关性关系构建组分（Component），用这些组分去预测 y 。这与PCR不同，PCR主要考虑自变量 x 间关系，寻找方差最大的方向，而PLS侧重于与 y 最相关的方向。详细过程见参考文献［11］。R中的pls包有plsr()函数。
2.LOESS回归
LOESS（Locally Weighted Scatterplot Smoothing，局部散点平滑，也称局部加权多项式）回归是一种非参数化拟合方法，主要思想为在数据集合的每一点用低维多项式拟合数据点的一个子集，并估计该点附近自变量数据点所对应的因变量值，该多项式是用加权最小二乘法来拟合；离该点越远，权重越小，该点的回归函数值就是由这个局部多项式来得到，而用于加权最小二乘回归的数据子集是由最近邻方法确定。
可以看到，损失函数加上权重之后，在最小化损失函数时，就会更多地考虑权重大的点，希望它们更优，这样拟合出来的结果，自然就更加偏向权重大的点了，也就是说，距离拟合点更近距离的散点，对拟合直线的影响更大。
这种回归方法的优点是不需要事先设定一个全局函数来对所有数据拟合一个模型。并且可以对同一数据进行多次不同的拟合，先对某个变量进行拟合，再对另一变量进行拟合，以探索数据中可能存在的某种关系，这是普通的回归拟合无法做到的，具体的例子见下述代码和图3-20。
以mtcars数据集为例，mtcars是32个车型的设计与性能指标数据（11个指标），想探索车的重量（字段名为wt）与单位燃油行驶距离（字段名为mpg）的关系，LOESS回归可以比线性回归更光滑拟合实际数据关系。
图3-20 不同车型的车重量与单位燃油行驶距离的关系图
LOESS算法早期版本是LOWESS，LOWESS仅支持单变量（在二维图形上常用），而LOESS支持多变量，在局部权重和算法实现上也有一些差别，可以简单认为LOESS是LOWESS的扩展 [17] 。R语言也有两个对应的函数lowess()和loess()。
3.Kernel方法
Kernel回归（Kernel Regression）可以看作k近邻算法的扩展。k近邻算法中，一个点的预测值是基于 k 个最相近的样本点的平均。而核函数对此进一步细化，距离近的权重更高一些，并且不限于前 k 个最近样本点，相当于做了平滑。核密度估计是一种常见的非参数化概率密度估计方法。
式中， K 为核密度函数； h 为设定的窗宽。
4.思考
建议读者思考一下LM与GLM、GAM、ANN等算法的关系 [18] ，详细的分析可以阅读《图解机器学习》 [19] 或EoSL [11] 。