通过前面的学习,我们知道了线性回归模型的回归系数表达式是:
β
=
(
X
T
X
)
−
1
X
T
y
\beta =\left ( X^TX \right )^{-1}X^Ty
β=(XTX)−1XTy
现在的问题是,能保证回归系数一定有解吗?答案是不一定,这是有条件的,从该式可以看出来,必须确保矩阵
X
T
X
X^TX
XTX是满秩的,即
X
T
X
X^TX
XTX可逆。但在实际的数据当中,自变量之间可能存在高度自相关性,这样就会导致回归系数无解或结果无效。那么如何解决这个问题?
第一,可以根据业务知识,人工判断,将那些高自相关的变量进行删除;
第二,选用岭回归也能够避免
X
T
X
X^TX
XTX的不可逆。
一、岭回归
岭回归一般可以用来解决线性回归模型系数无解的两种情况,一方面是自变量间存在高度多重共线性,另一方面则是自变量个数大于等于观测个数。针对这两种情况,通过举例进行进行说明:
第一种情况(多重共线性):
X
=
[
1
2
2
2
5
4
2
3
4
]
X=\begin{bmatrix} 1 & 2 &2 \\ 2 & 5 &4 \\ 2 & 3 &4 \end{bmatrix}
X=⎣⎡122253244⎦⎤
X
T
X
=
[
1
2
2
2
5
3
2
4
4
]
[
1
2
2
2
5
4
2
3
4
]
=
[
9
18
18
18
38
36
18
36
36
]
X^TX=\begin{bmatrix} 1 & 2 &2 \\ 2 & 5 &3 \\ 2 & 4 &4 \end{bmatrix} \begin{bmatrix} 1 & 2 &2 \\ 2 & 5 &4 \\ 2 & 3 &4 \end{bmatrix} =\begin{bmatrix} 9 & 18 &18 \\ 18 & 38 &36 \\ 18 & 36 &36 \end{bmatrix}
XTX=⎣⎡122254234⎦⎤⎣⎡122253244⎦⎤=⎣⎡91818183836183636⎦⎤
计算
X
T
X
X^TX
XTX的行列式
∣
X
T
X
∣
=
0
\left | X^TX \right |=0
∣∣XTX∣∣=0
第二种情况(列数比行数多):
X
=
[
1
2
5
6
1
3
]
X=\begin{bmatrix} 1 & 2& 5\\ 6& 1& 3 \end{bmatrix}
X=[162153]
X
T
X
=
[
1
6
2
1
5
3
]
[
1
2
5
6
1
3
]
=
[
37
8
23
8
5
13
23
13
34
]
X^TX=\begin{bmatrix} 1 &6 \\ 2 &1 \\ 5 &3 \end{bmatrix}\begin{bmatrix} 1 & 2& 5\\ 6& 1& 3 \end{bmatrix} =\begin{bmatrix} 37 & 8 & 23\\ 8& 5 & 13\\ 23& 13 &34 \end{bmatrix}
XTX=⎣⎡125613⎦⎤[162153]=⎣⎡378238513231334⎦⎤
计算
X
T
X
X^TX
XTX的行列式
∣
X
T
X
∣
=
0
\left | X^TX \right |=0
∣∣XTX∣∣=0
所以,不管是高度多重共线性的矩阵还是列数多于观测数(PS:就是特征数量很多,但是样本数量很少)的矩阵,最终算出来的行列式都等于0或者是近似为0,类似于这样的矩阵,都会导致线性回归模型的回归系数无解或解无意义(因为矩阵行列式近似为0时,其逆将偏于无穷大,导致回归系数也被放大)。那如何来解决这个问题呢?1970年Heer提出了岭回归方法,巧妙化解了这一难题,即在 X T X X^TX XTX的基础上加上一个较小的 λ \lambda λ扰动 ,从而使得行列式不再为0。
将岭回归系数的求解表达式写成如下这个式子:
β
(
λ
)
=
(
X
T
X
+
λ
I
)
−
1
X
T
y
\beta \left ( \lambda \right )=\left ( X^TX+\lambda I \right )^{-1}X^Ty
β(λ)=(XTX+λI)−1XTy
回归系数
β
\beta
β的值将随着
λ
\lambda
λ的变化而变化,当
λ
=
0
\lambda=0
λ=0时,就退化为线性回归模型的系数值。
实际上,岭回归系数的解依赖于下面这个最优化问题:
β
^
r
i
d
g
e
=
a
r
g
m
i
n
β
{
∑
i
=
1
N
(
y
i
−
β
0
−
∑
j
=
1
p
x
i
j
β
j
)
2
+
λ
∑
i
=
1
p
β
j
2
}
\hat{\beta }^{ridge}=\underset{\beta}{argmin}\left \{ \sum_{i=1}^{N}\left ( y_i-\beta_0-\sum_{j=1}^{p}x_{ij}\beta_j \right )^2 +{\color{Red} \lambda \sum_{i=1}^{p}\beta_j^2}\right \}
β^ridge=βargmin⎩⎨⎧i=1∑N(yi−β0−j=1∑pxijβj)2+λi=1∑pβj2⎭⎬⎫
是不是有些眼熟呢?对!标红色的被称作目标函数的惩罚函数,它是一个
L
2
L_2
L2范数!它可以确保岭回归系数
β
\beta
β值不会变的很大,起到收缩的作用,这个收缩力度就可以通过
λ
\lambda
λ来平衡。见下图:
这幅图的横坐标是模型的复杂度,纵坐标是模型的预测误差,绿色曲线代表的是模型在训练集上的效果,蓝色曲线代表的是模型在测试集的效果。从预测效果的角度来看,随着模型复杂度的提升,在训练集上的预测效果会越来越好,呈现在绿色曲线上就是预测误差越来越低,但是模型运用到测试集的话,预测误差就会呈现蓝色曲线的变化,先降低后上升(过拟合);从模型方差角度(即回归系数的方差)来看,模型方差会随着复杂度的提升而提升。针对上面这个图而言,我们是希望平衡方差和偏差来选择一个比较理想的模型,对于岭回归来说,随着
λ
\lambda
λ的增大,模型方差会减小(因为矩阵
X
T
X
X^TX
XTX行列式在增大,其逆就是减小,从而使得岭回归系数在减小)而偏差会增大。故通过
λ
\lambda
λ来平衡模型的方差和偏差,最终得到比较理想的岭回归系数。
上面讲解了关于岭回归模型的参数求解,那么岭回归的几何意义是什么呢?
上式
β
^
r
i
d
g
e
\hat{\beta }^{ridge}
β^ridge还可以表示成如下形式:
β
^
r
i
d
g
e
=
a
r
g
m
i
n
β
{
∑
i
=
1
N
(
y
i
−
β
0
−
∑
j
=
1
p
x
i
j
β
j
)
2
}
\hat{\beta }^{ridge}=\underset{\beta}{argmin}\left \{ \sum_{i=1}^{N}\left ( y_i-\beta_0-\sum_{j=1}^{p}x_{ij}\beta_j \right )^2 \right \}
β^ridge=βargmin⎩⎨⎧i=1∑N(yi−β0−j=1∑pxijβj)2⎭⎬⎫
s
u
b
j
e
c
t
t
o
:
λ
∑
i
=
1
p
β
j
2
≤
C
subject\quad to:\lambda \sum_{i=1}^{p}\beta_j^2\leq C
subjectto:λi=1∑pβj2≤C
为什么要添加这个岭回归系数平方和的约束呢?我们知道,岭回归模型可以解决多重共线性的麻烦,正是因为多重共线性的原因,才需要添加这个约束,也就是说变量的回归系数之间定会存在相互抵消的作用。
如果把上面的等价目标函数展示到几何图形中的话,将会是如下的形式(这里以两个变量的回归系数为例):
其中,半椭圆体表现的是
a
r
g
m
i
n
β
{
∑
i
=
1
N
(
y
i
−
β
0
−
∑
j
=
1
p
x
i
j
β
j
)
2
}
\underset{\beta}{argmin}\left \{ \sum_{i=1}^{N}\left ( y_i-\beta_0-\sum_{j=1}^{p}x_{ij}\beta_j \right )^2 \right \}
βargmin⎩⎨⎧i=1∑N(yi−β0−j=1∑pxijβj)2⎭⎬⎫
圆柱体表现的是
λ
∑
i
=
1
p
β
j
2
≤
C
\lambda \sum_{i=1}^{p}\beta_j^2\leq C
λi=1∑pβj2≤C
重点来了:黄色的交点就是满足目标函数下的岭回归系数值。
进一步,可以将这个三维的立体图映射到二维平面中,表现的就更加直观了:
是不是又很面熟呢?
那么接下来就要讨论如何选取岭参数
λ
\lambda
λ呢?
由岭回归的目标函数可以看出,惩罚函数的系数也就是
λ
\lambda
λ越大,目标函数中惩罚函数所占的重要性越高,从而估计参数
β
\beta
β也就越小了。可以看出岭参数不是唯一的,所以我们得到的岭回归估计
β
^
r
i
d
g
e
\hat{\beta }^{ridge}
β^ridge实际是回归参数
β
\beta
β的一个估计族。
当矩阵存在奇异性时,由岭回归的参数估计结果可以看出来,刚开始岭参数不够大时,奇异性并没有得到太大的改变,所以随着岭参数的变化,回归的估计参数震动很大,当岭参数足够大时,奇异性的影响逐渐减少,从而估计参数的值变的逐渐稳定。
岭参数选择的一般原则:
1、各回归系数的岭估计基本稳定
2、不存在有明显不符合常理的回归参数,其岭估计的符号应当要变得合理
3、回归系数没有不合实际意义的绝对值
4、残差平方和增大不多
岭参数的选取方法有很多,目前我在论文里看到过的选取方法有岭迹法、方差扩大因子法、双 H H H公式法、 H o r e l − K e n n a r d Horel-Kennard Horel−Kennard公式法。
岭回归存在的问题:
1、岭参数计算方法太多,差异较大;
2、根据岭迹图进行变量筛选,随意性太大;
3、岭回归返回的模型(如果没有经过变量筛选)包含了所有的变量。
二、LASSO回归
LASSO回归其实与岭回归很相似,不同的是求解回归系数的目标函数中使用的惩罚函数是
L
1
L_1
L1范数,即
β
^
l
a
s
s
o
=
a
r
g
m
i
n
β
{
∑
i
=
1
N
(
y
i
−
β
0
−
∑
j
=
1
p
x
i
j
β
j
)
2
+
λ
∑
i
=
1
p
∣
β
j
∣
}
\hat{\beta }^{lasso}=\underset{\beta}{argmin}\left \{ \sum_{i=1}^{N}\left ( y_i-\beta_0-\sum_{j=1}^{p}x_{ij}\beta_j \right )^2 +{\color{Red} \lambda \sum_{i=1}^{p}\left |\beta_j \right |}\right \}
β^lasso=βargmin{∑i=1N(yi−β0−∑j=1pxijβj)2+λ∑i=1p∣βj∣}
同理,该目标函数也可以转化成约束的形式:
β
^
l
a
s
s
o
=
a
r
g
m
i
n
β
{
∑
i
=
1
N
(
y
i
−
β
0
−
∑
j
=
1
p
x
i
j
β
j
)
2
}
\hat{\beta }^{lasso}=\underset{\beta}{argmin}\left \{ \sum_{i=1}^{N}\left ( y_i-\beta_0-\sum_{j=1}^{p}x_{ij}\beta_j \right )^2 \right \}
β^lasso=βargmin⎩⎨⎧i=1∑N(yi−β0−j=1∑pxijβj)2⎭⎬⎫
s
u
b
j
e
c
t
t
o
:
λ
∑
i
=
1
p
∣
β
j
∣
≤
C
subject\quad to:\lambda \sum_{i=1}^{p}\left|\beta_j\right|\leq C
subjectto:λi=1∑p∣βj∣≤C
将上面的等价目标函数可以展示到二维几何图形中,效果如下图所示(这里以两个变量的回归系数为例):
其实,LASSO回归对于岭回归只是在惩罚函数部分有所不同,但这个不同却让LASSO明显占了很多优势,例如在变量选择上就比岭回归强悍的多。就以直观的图形为例,LASSO回归的惩罚函数映射到二维空间的话,就会形成“角”,一旦“角”与抛物面相交,就会导致
β
1
\beta_1
β1为
0
0
0(如上图所示),这样
β
1
\beta_1
β1对应的变量就是一个可抛弃的变量。但是在岭回归过程中,没有“角”的圆形与抛物面相交,出现岭回归系数为
0
0
0的概率还是非常小的。
(PS:这其实就是前面一章所讲的
L
1
L_1
L1范数稀疏的特性,用于对特征进行降维!)
那么如何对LASSO回归算法进行求解呢?
网上最多的是用LAR算法(最小角回归)来近似求解LASSO。
这篇文章很长,我就不记录了,发个链接吧。
使用LAR算法进行求解Lasso算法
未完待续。。。