Lasso算法

假设一个样本包括N种事件，每个事件包括p个协变量和一个输出值。让${\displaystyle y_{i}}$为输出值，并且${\displaystyle x_{i}:=(x_{1},x_{2},\ldots ,x_{p})^{T}}$为第i种情况的协变量向量，那么Lasso要计算的目标方程就是：

对所有 ${\displaystyle \sum _{j=1}^{p}|\beta _{j}|\leq t}$，计算 ${\displaystyle \min _{\beta _{0},\beta }\left\{{\frac {1}{N}}\sum _{i=1}^{N}(y_{i}-\beta _{0}-x_{i}^{T}\beta )^{2}\right\}}$[1]

这里 ${\displaystyle t}$ 是一个决定规则化程度的预定的自由参数。 设${\displaystyle X}$为协变量矩阵，那么 ${\displaystyle X_{ij}=(x_{i})_{j}}$，其中 ${\displaystyle x_{i}^{T}}$ 是 ${\displaystyle X}$的第 i 行，那么上式可以写成更紧凑的形式：

对所有 ${\displaystyle \|\beta \|_{1}\leq t}$，计算 ${\displaystyle \min _{\beta _{0},\beta }\left\{{\frac {1}{N}}\left\|y-\beta _{0}-X\beta \right\|_{2}^{2}\right\}}$

这里 ${\displaystyle \|\beta \|_{p}=\left(\sum _{i=1}^{N}|\beta _{i}|^{p}\right)^{1/p}}$ 是标准 ${\displaystyle \ell ^{p}}$ 范数，${\displaystyle 1_{N}}$是${\displaystyle N\times 1}$维的1的向量。

因为 ${\displaystyle {\hat {\beta }}_{0}={\bar {y}}-{\bar {x}}^{T}\beta }$，所以有

${\displaystyle y_{i}-{\hat {\beta }}_{0}-x_{i}^{T}\beta =y_{i}-({\bar {y}}-{\bar {x}}^{T}\beta )-x_{i}^{T}\beta =(y_{i}-{\bar {y}})-(x_{i}-{\bar {x}})^{T}\beta ,}$

对变量进行中心化是常用的数据处理方法。并且协方差一般规范化为 ${\displaystyle \textstyle \left(\sum _{i=1}^{N}x_{ij}^{2}=1\right)}$ ，这样得到的解就不会依赖测量的规模。

它的目标方程还可以写为：

${\displaystyle \min _{\beta \in \mathbb {R} ^{p}}\left\{{\frac {1}{N}}\left\|y-X\beta \right\|_{2}^{2}\right\}{\text{ subject to }}\|\beta \|_{1}\leq t.}$

其拉格朗日形式为：

${\displaystyle \min _{\beta \in \mathbb {R} ^{p}}\left\{{\frac {1}{N}}\left\|y-X\beta \right\|_{2}^{2}+\lambda \|\beta \|_{1}\right\}}$

其中 ${\displaystyle t}$ 和 ${\displaystyle \lambda }$ 的关系取决于数据特征。

现在考虑一些Lasso回归估计的基本性质。

首先假定所有的协变量都是正交的，即 ${\displaystyle (x_{i}\mid x_{j})=\delta _{ij}}$，其中${\displaystyle \delta _{ij}}$ 为 克罗内克δ函数。等价的矩阵写法为 ${\displaystyle X^{T}X=I}$，使用次梯度法可有如下的表达形式

${\displaystyle {\begin{aligned}{\hat {\beta }}_{j}={}&S_{N\lambda }({\hat {\beta }}_{j}^{\text{OLS}})={\hat {\beta }}_{j}^{\text{OLS}}\max \left(0,1-{\frac {N\lambda }{|{\hat {\beta }}_{j}^{\text{OLS}}|}}\right)\\&{\text{ 其中 }}{\hat {\beta }}^{\text{OLS}}=(X^{T}X)^{-1}X^{T}y\end{aligned}}}$ [1]

${\displaystyle S_{\alpha }}$ 用于表示软阈值算子，当这个值非常小的时候为0。一个与之相近的记号${\displaystyle H_{\alpha }}$用来表示硬阈值算子，将较小的数值记为0的同时保留原有的较大数值。

与岭回归相比较，其中岭回归的目标在于最小化

${\displaystyle \min _{\beta \in \mathbb {R} ^{p}}\left\{{\frac {1}{N}}\|y-X\beta \|_{2}^{2}+\lambda \|\beta \|_{2}^{2}\right\}}$

即有

${\displaystyle {\hat {\beta }}_{j}=(1+N\lambda )^{-1}{\hat {\beta }}_{j}^{\text{OLS}}.}$

因此岭回归是对OLS回归中所有的系数以一致的系数${\displaystyle (1+N\lambda )^{-1}}$缩放，并不会进行变量选择。

同样也可以对best subset selection算法进行比较，其目标在于最小化

${\displaystyle \min _{\beta \in \mathbb {R} ^{p}}\left\{{\frac {1}{N}}\left\|y-X\beta \right\|_{2}^{2}+\lambda \|\beta \|_{0}\right\}}$

其中 ${\displaystyle \|\cdot \|_{0}}$ 表示 "${\displaystyle \ell ^{0}}$ norm"，即0范数，被定义为该向量中非零元的个数。在这个例子中，可以得到

${\displaystyle {\hat {\beta }}_{j}=H_{\sqrt {N\lambda }}\left({\hat {\beta }}_{j}^{\text{OLS}}\right)={\hat {\beta }}_{j}^{\text{OLS}}\mathrm {I} \left(\left|{\hat {\beta }}_{j}^{\text{OLS}}\right|\geq {\sqrt {N\lambda }}\right)}$

其中 ${\displaystyle H_{\alpha }}$ 被称为软阈值算子，${\displaystyle \mathrm {I} }$ 为示性函数。

总的来说，Lasso估计量展现出了岭回归和最佳子划分算法的系数收缩的优点，使得部分系数为0。此外，在岭回归全部使用一个常数系数缩放的时候，Lasso回归会将一个接近0的系数变为0。

对于一般的情况中，不同的协变量之间可能并不是独立的，其中一种特例即为变量存在重复，例如变量j和变量k，有${\displaystyle x_{(j)}=x_{(k)}}$，${\displaystyle x_{(j),i}=x_{ij}}$。在这种情况下参数 ${\displaystyle \beta _{j}}$ 和 ${\displaystyle \beta _{k}}$的Lasso回归的估计量不是唯一确定的。

事实上，如果有一些 ${\displaystyle {\hat {\beta }}}$ 中存在 ${\displaystyle {\hat {\beta }}_{j}{\hat {\beta }}_{k}\geq 0}$，寻找一个${\displaystyle s\in [0,1]}$进行变换，将${\displaystyle {\hat {\beta }}_{j}}$变换为${\displaystyle s({\hat {\beta }}_{j}+{\hat {\beta }}_{k})}$ 的同时有 ${\displaystyle {\hat {\beta }}_{k}}$ 变换为 ${\displaystyle (1-s)({\hat {\beta }}_{j}+{\hat {\beta }}_{k})}$，并保留其他参数不变，此时Lasso回归具有有效的连续性质。一些基于Lasso回归的改进，例如弹性网络正则化，旨在解决这个缺点。

Forms of the constraint regions for lasso and ridge regression.

Lasso回归可以使得某些项系数为0，从几何上来看，不同约束边界形状的岭回归则不能。他们都可以解释为最小化相同的目标函数

${\displaystyle \min _{\beta _{0},\beta }\left\{{\frac {1}{N}}\left\|y-\beta _{0}-X\beta \right\|_{2}^{2}\right\}}$

但是有不同的约束条件：在Lasso回归中为 ${\displaystyle \|\beta \|_{1}\leq t}$ 而在岭回归中为 ${\displaystyle \|\beta \|_{2}^{2}\leq t}$。1-范数

The figure shows that the constraint region defined by the ${\displaystyle \ell ^{1}}$ norm is a square rotated so that its corners lie on the axes (in general a cross-polytope), while the region defined by the ${\displaystyle \ell ^{2}}$ norm is a circle (in general an n-sphere), which is rotationally invariant and, therefore, has no corners. As seen in the figure, a convex object that lies tangent to the boundary, such as the line shown, is likely to encounter a corner (or a higher-dimensional equivalent) of a hypercube, for which some components of ${\displaystyle \beta }$ are identically zero, while in the case of an n-sphere, the points on the boundary for which some of the components of ${\displaystyle \beta }$ are zero are not distinguished from the others and the convex object is no more likely to contact a point at which some components of ${\displaystyle \beta }$ are zero than one for which none of them are.

The lasso can be rescaled so that it becomes easy to anticipate and influence the degree of shrinkage associated with a given value of ${\displaystyle \lambda }$.[6] It is assumed that ${\displaystyle X}$ is standardized with z-scores and that ${\displaystyle y}$ is centered (zero mean). Let ${\displaystyle \beta _{0}}$ represent the hypothesized regression coefficients and let ${\displaystyle b_{OLS}}$ refer to the data-optimized ordinary least squares solutions. We can then define the Lagrangian as a tradeoff between the in-sample accuracy of the data-optimized solutions and the simplicity of sticking to the hypothesized values.[7] This results in

${\displaystyle \min _{\beta \in \mathbb {R} ^{p}}\left\{{\frac {(y-X\beta )'(y-X\beta )}{(y-X\beta _{0})'(y-X\beta _{0})}}+2\lambda \sum _{i=1}^{p}{\frac {|\beta _{i}-\beta _{0,i}|}{q_{i}}}\right\}}$

where ${\displaystyle q_{i}}$ is specified below. The first fraction represents relative accuracy, the second fraction relative simplicity, and ${\displaystyle \lambda }$ balances between the two.

Solution paths for the ${\displaystyle \ell _{1}}$ norm and ${\displaystyle \ell _{2}}$ norm when ${\displaystyle b_{OLS}=2}$ and ${\displaystyle \beta _{0}=0}$

Given a single regressor, relative simplicity can be defined by specifying ${\displaystyle q_{i}}$ as ${\displaystyle |b_{OLS}-\beta _{0}|}$, which is the maximum amount of deviation from ${\displaystyle \beta _{0}}$ when ${\displaystyle \lambda =0}$. Assuming that ${\displaystyle \beta _{0}=0}$, the solution path can be defined in terms of ${\displaystyle R^{2}}$:

${\displaystyle b_{\ell _{1}}={\begin{cases}(1-\lambda /R^{2})b_{OLS}&{\mbox{if }}\lambda \leq R^{2},\\0&{\mbox{if }}\lambda >R^{2}.\end{cases}}}$

If ${\displaystyle \lambda =0}$, the ordinary least squares solution (OLS) is used. The hypothesized value of ${\displaystyle \beta _{0}=0}$ is selected if ${\displaystyle \lambda }$ is bigger than ${\displaystyle R^{2}}$. Furthermore, if ${\displaystyle R^{2}=1}$, then ${\displaystyle \lambda }$ represents the proportional influence of ${\displaystyle \beta _{0}=0}$. In other words, ${\displaystyle \lambda \times 100\%}$ measures in percentage terms the minimal amount of influence of the hypothesized value relative to the data-optimized OLS solution.

If an ${\displaystyle \ell _{2}}$-norm is used to penalize deviations from zero given a single regressor, the solution path is given by

${\displaystyle b_{\ell _{2}}={\bigg (}1+{\frac {\lambda }{R^{2}(1-\lambda )}}{\bigg )}^{-1}b_{OLS}}$. Like ${\displaystyle b_{\ell _{1}}}$, ${\displaystyle b_{\ell _{2}}}$ moves in the direction of the point ${\displaystyle (\lambda =R^{2},b=0)}$ when ${\displaystyle \lambda }$ is close to zero; but unlike ${\displaystyle b_{\ell _{1}}}$, the influence of ${\displaystyle R^{2}}$ diminishes in ${\displaystyle b_{\ell _{2}}}$ if ${\displaystyle \lambda }$ increases (see figure).
Given multiple regressors, the moment that a parameter is activated (i.e. allowed to deviate from ${\displaystyle \beta _{0}}$) is also determined by a regressor's contribution to ${\displaystyle R^{2}}$ accuracy. First,

${\displaystyle R^{2}=1-{\frac {(y-Xb)'(y-Xb)}{(y-X\beta _{0})'(y-X\beta _{0})}}.}$

An ${\displaystyle R^{2}}$ of 75% means that in-sample accuracy improves by 75% if the unrestricted OLS solutions are used instead of the hypothesized ${\displaystyle \beta _{0}}$ values. The individual contribution of deviating from each hypothesis can be computed with the ${\displaystyle p}$ x ${\displaystyle p}$ matrix

${\displaystyle R^{\otimes }=(X'{\tilde {y}}_{0})(X'{\tilde {y}}_{0})'(X'X)^{-1}({\tilde {y}}_{0}'{\tilde {y}}_{0})^{-1},}$

where ${\displaystyle {\tilde {y}}_{0}=y-X\beta _{0}}$. If ${\displaystyle b=b_{OLS}}$ when ${\displaystyle R^{2}}$ is computed, then the diagonal elements of ${\displaystyle R^{\otimes }}$ sum to ${\displaystyle R^{2}}$. The diagonal ${\displaystyle R^{\otimes }}$ values may be smaller than 0 or, less often, larger than 1. If regressors are uncorrelated, then the ${\displaystyle i^{th}}$ diagonal element of ${\displaystyle R^{\otimes }}$ simply corresponds to the ${\displaystyle r^{2}}$ value between ${\displaystyle x_{i}}$ and ${\displaystyle y}$.

A rescaled version of the adaptive lasso of can be obtained by setting ${\displaystyle q_{{\mbox{adaptive lasso}},i}=|b_{OLS,i}-\beta _{0,i}|}$.[8] If regressors are uncorrelated, the moment that the ${\displaystyle i^{th}}$ parameter is activated is given by the ${\displaystyle i^{th}}$ diagonal element of ${\displaystyle R^{\otimes }}$. Assuming for convenience that ${\displaystyle \beta _{0}}$ is a vector of zeros,

${\displaystyle b_{i}={\begin{cases}(1-\lambda /R_{ii}^{\otimes })b_{OLS,i}&{\mbox{if }}\lambda \leq R_{ii}^{\otimes },\\0&{\mbox{if }}\lambda >R_{ii}^{\otimes }.\end{cases}}}$

That is, if regressors are uncorrelated, ${\displaystyle \lambda }$ again specifies the minimal influence of ${\displaystyle \beta _{0}}$. Even when regressors are correlated, the first time that a regression parameter is activated occurs when ${\displaystyle \lambda }$ is equal to the highest diagonal element of ${\displaystyle R^{\otimes }}$.

These results can be compared to a rescaled version of the lasso by defining ${\displaystyle q_{{\mbox{lasso}},i}={\frac {1}{p}}\sum _{l}|b_{OLS,l}-\beta _{0,l}|}$, which is the average absolute deviation of ${\displaystyle b_{OLS}}$ from ${\displaystyle \beta _{0}}$. Assuming that regressors are uncorrelated, then the moment of activation of the ${\displaystyle i^{th}}$ regressor is given by

${\displaystyle {\tilde {\lambda }}_{{\text{lasso}},i}={\frac {1}{p}}{\sqrt {R_{i}^{\otimes }}}\sum _{l=1}^{p}{\sqrt {R_{l}^{\otimes }}}.}$

For ${\displaystyle p=1}$, the moment of activation is again given by ${\displaystyle {\tilde {\lambda }}_{{\text{lasso}},i}=R^{2}}$. If ${\displaystyle \beta _{0}}$ is a vector of zeros and a subset of ${\displaystyle p_{B}}$ relevant parameters are equally responsible for a perfect fit of ${\displaystyle R^{2}=1}$, then this subset is activated at a ${\displaystyle \lambda }$ value of ${\displaystyle {\frac {1}{p}}}$. The moment of activation of a relevant regressor then equals ${\displaystyle {\frac {1}{p}}{\frac {1}{\sqrt {p_{B}}}}p_{B}{\frac {1}{\sqrt {p_{B}}}}={\frac {1}{p}}}$. In other words, the inclusion of irrelevant regressors delays the moment that relevant regressors are activated by this rescaled lasso. The adaptive lasso and the lasso are special cases of a '1ASTc' estimator. The latter only groups parameters together if the absolute correlation among regressors is larger than a user-specified value.[6]

Laplace distributions are sharply peaked at their mean with more probability density concentrated there compared to a normal distribution.

Just as ridge regression can be interpreted as linear regression for which the coefficients have been assigned normal prior distributions, lasso can be interpreted as linear regression for which the coefficients have Laplace prior distributions. The Laplace distribution is sharply peaked at zero (its first derivative is discontinuous at zero) and it concentrates its probability mass closer to zero than does the normal distribution. This provides an alternative explanation of why lasso tends to set some coefficients to zero, while ridge regression does not.[1]

Lasso can also be viewed as a convex relaxation of the best subset selection regression problem, which is to find the subset of ${\displaystyle \leq k}$ covariates that results in the smallest value of the objective function for some fixed ${\displaystyle k\leq n}$, where n is the total number of covariates. The "${\displaystyle \ell ^{0}}$ norm", ${\displaystyle \|\cdot \|_{0}}$, (the number of nonzero entries of a vector), is the limiting case of "${\displaystyle \ell ^{p}}$ norms", of the form ${\displaystyle \textstyle \|x\|_{p}=\left(\sum _{i=1}^{n}|x_{j}|^{p}\right)^{1/p}}$ (where the quotation marks signify that these are not really norms for ${\displaystyle p<1}$ since ${\displaystyle \|\cdot \|_{p}}$ is not convex for ${\displaystyle p<1}$, so the triangle inequality does not hold). Therefore, since p = 1 is the smallest value for which the "${\displaystyle \ell ^{p}}$ norm" is convex (and therefore actually a norm), lasso is, in some sense, the best convex approximation to the best subset selection problem, since the region defined by ${\displaystyle \|x\|_{1}\leq t}$ is the convex hull of the region defined by ${\displaystyle \|x\|_{p}\leq t}$ for ${\displaystyle p<1}$.