天天看点

Lasso回归理论Python实现

理论

Lasso回归在最小二乘法的基础上加上了一个 l 1 l_1 l1​惩罚项

损失函数: J ( θ ) = 1 2 m ∑ i = 1 m ( h θ ( x ( i ) ) − y ( i ) ) 2 + λ ∑ j = 1 n ∣ θ j ∣ J(\theta)=\frac 1 {2m}\sum_{i=1}^m(h_{\theta}(x^{(i)})-y^{(i)})^2+\lambda \sum_{j=1}^n|\theta_j| J(θ)=2m1​i=1∑m​(hθ​(x(i))−y(i))2+λj=1∑n​∣θj​∣

相比岭回归可以直接通过矩阵运算得到回归系数相比,LASSO的计算变得相对复杂。由于惩罚项中含有绝对值,此函数的导数是连续不光滑的,所以无法进行求导并使用梯度下降优化。

两种求解Lasso回归的方法:

  • 坐标下降法
  • 最小角回归

坐标下降法

坐标下降法,是沿着坐标轴的方向去下降。

坐标下降法的数学依据是:

         一个可微的凸函数 J ( θ ) J(\theta) J(θ),其中 θ \theta θ是 n ∗ 1 n*1 n∗1的向量,即有 n n n个维度。如果在某一点 θ ‾ \overline \theta θ,使得 J ( θ ) J(\theta) J(θ)在每一个坐标轴 θ ‾ i ( i = 1 , 2 , . . . , n ) \overline \theta_i(i=1,2,...,n) θi​(i=1,2,...,n)上都是最小值,那么 J ( θ ‾ i ) J(\overline \theta_i) J(θi​)就是一个全局最小值

于是,我们的优化目标是:在 θ \theta θ的 n n n个坐标轴上,对损失函数做迭代的下降,当所有的坐标轴上的 θ i ( i = 1 , 2 , . . . , n ) \theta_i(i=1,2,...,n) θi​(i=1,2,...,n)都收敛,此时损失函数最小,此时的 θ \theta θ即为我们要求的结果。

具体算法流程:

1、首先,初始化 θ \theta θ向量,随机取值即可,即为 θ ( 0 ) \theta^{(0)} θ(0),上面的括号里的数字表示当前迭代的轮数。

2、对于第 k k k轮的迭代,我们从 θ 1 ( k ) \theta_1^{(k)} θ1(k)​开始,到 θ n ( k ) \theta_n^{(k)} θn(k)​为止,依次求 θ i ( k ) \theta_i^{(k)} θi(k)​。 θ i ( k ) \theta_i^{(k)} θi(k)​的表达式如下:

θ i ( k ) ∈ a r g m i n ⎵ θ i J ( θ 1 ( k ) , θ 2 ( k ) , . . . θ i − 1 ( k ) , θ i , θ i + 1 ( k − 1 ) , . . . , θ n ( k − 1 ) ) \theta_i^{(k)} \in \underbrace{argmin}_{\theta_i} J(\theta_1^{(k)}, \theta_2^{(k)}, ... \theta_{i-1}^{(k)}, \theta_i, \theta_{i+1}^{(k-1)}, ..., \theta_n^{(k-1)}) θi(k)​∈θi​

argmin​​J(θ1(k)​,θ2(k)​,...θi−1(k)​,θi​,θi+1(k−1)​,...,θn(k−1)​)

θ \theta θ向量的 n n n个维度的迭代式如下:

θ i ( k ) ∈ a r g m i n ⎵ θ i J ( θ 1 ( k ) , θ 2 ( k ) , . . . θ i − 1 ( k ) , θ i , θ i + 1 ( k − 1 ) , . . . , θ n ( k − 1 ) ) \theta_i^{(k)} \in \underbrace{argmin}_{\theta_i} J(\theta_1^{(k)}, \theta_2^{(k)}, ... \theta_{i-1}^{(k)}, \theta_i, \theta_{i+1}^{(k-1)}, ..., \theta_n^{(k-1)}) θi(k)​∈θi​

argmin​​J(θ1(k)​,θ2(k)​,...θi−1(k)​,θi​,θi+1(k−1)​,...,θn(k−1)​)

也就是说 θ i ( k ) \theta_i^{(k)} θi(k)​是使得 J ( θ 1 ( k ) , θ 2 ( k ) , . . . θ i − 1 ( k ) , θ i , θ i + 1 ( k − 1 ) , . . . , θ n ( k − 1 ) ) J(\theta_1^{(k)}, \theta_2^{(k)}, ... \theta_{i-1}^{(k)}, \theta_i, \theta_{i+1}^{(k-1)}, ..., \theta_n^{(k-1)}) J(θ1(k)​,θ2(k)​,...θi−1(k)​,θi​,θi+1(k−1)​,...,θn(k−1)​)最小化的 θ i \theta_i θi​的值。此时 J ( θ ) J(\theta) J(θ)只有 θ i ( k ) \theta_i^{(k)} θi(k)​是变量,其余都是常量,所以这就是一个 J ( θ ) J(\theta) J(θ)关于 θ i ( k ) \theta_i^{(k)} θi(k)​的一元函数,很容易通过求导求得最小值。

如果上面这个式子不好理解,我们具体一点,在第 k k k轮, θ \theta θ向量的 n n n个维度的迭代式如下:

θ 1 ( k ) ∈ a r g m i n ⎵ θ 1 J ( θ 1 , θ 2 ( k − 1 ) , . . . , θ n ( k − 1 ) ) \theta_1^{(k)} \in \underbrace{argmin}_{\theta_1} J(\theta_1, \theta_2^{(k-1)}, ... , \theta_n^{(k-1)}) θ1(k)​∈θ1​

argmin​​J(θ1​,θ2(k−1)​,...,θn(k−1)​)

θ 2 ( k ) ∈ a r g m i n ⎵ θ 2 J ( θ 1 ( k ) , θ 2 , θ 3 ( k − 1 ) . . . , θ n ( k − 1 ) ) \theta_2^{(k)} \in \underbrace{argmin}_{\theta_2} J(\theta_1^{(k)}, \theta_2, \theta_3^{(k-1)}... , \theta_n^{(k-1)}) θ2(k)​∈θ2​

argmin​​J(θ1(k)​,θ2​,θ3(k−1)​...,θn(k−1)​)

. . . ... ...

θ n ( k ) ∈ a r g m i n ⎵ θ n J ( θ 1 ( k ) , θ 2 ( k ) , . . . , θ n − 1 ( k ) , θ n ) \theta_n^{(k)} \in \underbrace{argmin}_{\theta_n} J(\theta_1^{(k)}, \theta_2^{(k)}, ... , \theta_{n-1}^{(k)}, \theta_n) θn(k)​∈θn​

argmin​​J(θ1(k)​,θ2(k)​,...,θn−1(k)​,θn​)

3、检查 θ ( k ) \theta^{(k)} θ(k)向量和 θ ( k − 1 ) \theta^{(k-1)} θ(k−1)向量在各个维度上的变化情况,如果在所有维度上变化都足够小,则 θ ( k ) \theta^{(k)} θ(k)为最终结果,否则转入2,继续迭代

最小角回归

在介绍最小角回归之前,需要先看看两个预备算法:

  • 前向选择算法
  • 前向梯度算法

前向选择算法

前向选择算法的原理是一种典型的贪心算法。要解决的问题是:

        对于 Y = X θ Y=X\theta Y=Xθ这样的线性关系,如何求解系数 θ \theta θ。其中 Y Y Y是 m ∗ 1 m*1 m∗1的向量, X X X是 m ∗ n m*n m∗n的矩阵, θ \theta θ为 n ∗ 1 n*1 n∗1的向量。 m m m为样本数量, n n n为特征维度。

把矩阵 X X X看成 n n n个 m ∗ 1 m*1 m∗1的向量 X i ( i = 1 , 2 , . . . , n ) X_i(i=1,2,...,n) Xi​(i=1,2,...,n)。在这 n n n个向量中选择一个与目标 Y Y Y的余弦距离最大的一个 X k X_k Xk​,用 X k X_k Xk​来逼近 Y Y Y,得到下式:

Y ‾ = X k θ k \overline Y=X_k\theta_k Y=Xk​θk​,其中 θ k = &lt; X k , Y &gt; ∣ ∣ X k ∣ ∣ 2 \theta_k=\frac {&lt;X_k,Y&gt;}{||X_k||_2} θk​=∣∣Xk​∣∣2​<Xk​,Y>​

即 Y ‾ \overline Y Y是 Y Y Y在 X k X_k Xk​上的投影。那么,可以定义残差: Y y e s = Y − Y ‾ Y_{yes}=Y-\overline Y Yyes​=Y−Y。由于是投影,可知 Y y e s Y_{yes} Yyes​和 X k X_k Xk​是正交的。再以 Y y e s Y_{yes} Yyes​作为新的因变量,去掉 X k X_k Xk​后的剩下的自变量的集合 X i ( i = 1 , 2 , . . . , k − 1 , k + 1 , . . . , n ) X_i(i=1,2,...,k-1,k+1,...,n) Xi​(i=1,2,...,k−1,k+1,...,n)作为新的自变量集合,重复刚才投影和残差的操作,直到残差为0,或者所有的自变量都用完了,才停止算法。

Lasso回归理论Python实现

当 X X X只有2维时,如上图所示,和 Y Y Y最接近的是 X 1 X_1 X1​,首先在 X 1 X_1 X1​上投影,残差如上图长虚线。此时 X 1 θ 1 X_1\theta_1 X1​θ1​模拟了 Y Y Y, θ 1 \theta_1 θ1​模拟了 θ \theta θ(仅仅模拟了一个维度)。接着发现最接近的是 X 2 X_2 X2​,此时用残差接着在 X 2 X_2 X2​投影,残差为图中短虚线。由于没有其他自变量了,此时 x 1 θ 1 + x 2 θ 2 x_1\theta_1+x_2\theta_2 x1​θ1​+x2​θ2​模拟了 Y Y Y,对应的模拟了两个维度的 θ \theta θ即为最终结果。

此算法对每个变量只需执行一次操作,效率高,运算快。但,当自变量不是正交的时候,每次都在做投影,所以算法只能给出一个局部近似解。这个简单的算法太粗糙,不能直接用于Lasso回归。

前向梯度算法

前向梯度算法和前向选择算法有类似的地方,也是在 n n n个 X i X_i Xi​中选择和目标 Y Y Y最为接近(余弦距离最大)的一个变量 X k X_k Xk​,用 X k X_k Xk​来逼近 Y Y Y。但前向梯度算法不是粗暴的用投影,而是每次在最为接近的自变量 X t X_t Xt​的方向移动一小步,然后再看残差 Y y e s Y_{yes} Yyes​和哪个 X i X_i Xi​最为接近。此时我们也不会把 X t X_t Xt​去除,因为我们只前进了一小步,有可能下面最接近的自变量还是 X t X_t Xt​。如此进行下去,直到残差 Y y e s Y_yes Yy​es减小到足够小,算法停止。

Lasso回归理论Python实现

当 X X X只有2维时,例子如上图,和 Y Y Y最接近的是 X 1 X_1 X1​,首先在 X 1 X_1 X1​上面走一小段距离,此处 ϵ \epsilon ϵ为一个较小的常量,发现此时的残差还是和 X 1 X_1 X1​最接近。那么接着沿 X 1 X_1 X1​走,一直走到发现残差不是和 X 1 X_1 X1​最接近,而是和 X 2 X_2 X2​最接近,此时残差如上图长虚线。接着沿着 X 2 X_2 X2​走一小步,发现残差此时又和 X 1 X_1 X1​最接近,那么开始沿着 X 1 X_1 X1​走,走完一步后发现残差为0,那么算法停止。此时 Y Y Y由刚才所有的所有步相加而模拟,对应的算出的系数 θ \theta θ即为最终结果。

最小角回归算法

最小角回归对前向梯度和前向选择做了这种,保留了前向梯度算法一定程度上的精确性,同时简化了前向梯度算法一步步迭代的过程:

首先,找到与因变量 Y Y Y最接近或相关度最高的自变量 X k X_k Xk​,使用类似于前向梯度算法中的残差计算方法,得到新的目标 Y y e s Y_yes Yy​es,此时不用和前向梯度算法一样小步小步的走,而是直接向前走直到出现一个 X t X_t Xt​,使得 X t X_t Xt​和 Y ( y e s ) Y_{(yes)} Y(yes)​的相关度和 X k X_k Xk​与 Y ( y e s ) Y_{(yes)} Y(yes)​的相关度是一样的,此时残差 Y y e s Y_{yes} Yyes​就在 X t X_t Xt​和 X k X_k Xk​的角平分线上,此时我们开始沿着这个残差角平分线走,直到出现第三个特征 X p X_p Xp​和 Y y e s Y_yes Yy​es的相关度等于 θ t , θ k \theta_t,\theta_k θt​,θk​与 Y y e s Y_{yes} Yyes​的一样。将其也加入到 Y Y Y的逼近特征集合中,并用 Y Y Y的逼近特征集合的共同角分线,作为新的逼近方向,循环直到 Y y e s Y_{yes} Yyes​足够小或者所有变量都取完位置。

Lasso回归理论Python实现

当 θ \theta θ只有2维时。例子如上图,和 Y Y Y最接近的是 X 1 X_1 X1​,首先在 X 1 X_1 X1​上走一段距离,直到残差在 X 1 X_1 X1​和 X 2 X_2 X2​的角平分线上,此时沿着角平分线走,直到残差足够小才停止。此时对应的系数 β \beta β即为最终的结果。

最小角回归法是一个适用于高维数据的回归算法,其主要的优点有:

1)特别适合于特征维度n 远高于样本数m的情况。

2)算法的最坏计算复杂度和最小二乘法类似,但是其计算速度几乎和前向选择算法一样

3)可以产生分段线性结果的完整路径,这在模型的交叉验证中极为有用

主要的缺点是:

由于LARS的迭代方向是根据目标的残差而定,所以该算法对样本的噪声极为敏感

Python实现

坐标下降法 TODO

import numpy as np

           

最小角回归法 TODO

import numpy as np

           

继续阅读