LASSO的解法

2021-06-17 本文已影响0人 Boye0212

LASSO非常实用，但由于它的惩罚项不可以常规地进行求导，使得很多人以为它无法显式地求出解析解。但其实并不是这样的。

1 单变量情形：软阈值法

1.1 软阈值的分类讨论

将 $N$ 个样本的真实值记为 $N$ 维向量 $y$ ，将 $N$ 个样本的自变量记为 $z$ ，假设我们已经将自变量做过标准化，即 $z' \ell_n=0$ ， $z'z/N=1$ ，这也意味着在LASSO模型中截距项为 $0$ 。系数 $\beta$ 是要优化的参数，惩罚项参数为 $\lambda\gt 0$ 。

LASSO就是要求解
$\argmin_\beta \dfrac{1}{2N}(y-z\beta)'(y-z\beta)+\lambda |\beta| \tag{1}$

忽略常数项后，上式等价于
$\argmin_\beta \dfrac{1}{2}\beta^2 -\dfrac{y'z}{N}\beta+ \lambda |\beta|$

将损失函数写成分段函数形式：
$f(\beta)=\begin{cases} f_1(\beta) = \dfrac{1}{2}\beta^2 -\left(\dfrac{y'z}{N} + \lambda\right)\beta , \beta \lt 0\\ f_2(\beta) = \dfrac{1}{2}\beta^2 -\left(\dfrac{y'z}{N}- \lambda\right)\beta, \beta \geq 0 \end{cases}$

分类讨论：

若 $\dfrac{y'z}{N}\gt \lambda$ ，则 $f_1(\beta) \gt 0$ ， $f_2(\beta)$ 在 $\hat\beta=\dfrac{y'z}{N}- \lambda$ 处取到最小值 $f_2(\hat\beta)\lt 0$ ，因此解为 $\hat\beta=\dfrac{y'z}{N}- \lambda$ ；
若 $\left|\dfrac{y'z}{N}\right|\leq \lambda$ ，则 $f_1(\beta) \geq 0$ ， $f_2(\beta) \geq 0$ ，且在 $\hat\beta=0$ 处有 $f_1(\hat\beta)=f_2(\hat\beta)=0$ ，因此解为 $\hat\beta=0$ ；
若 $\dfrac{y'z}{N}\lt -\lambda$ ，则 $f_2(\beta) \gt 0$ ， $f_1(\beta)$ 在 $\hat\beta=\dfrac{y'z}{N}+\lambda$ 处取到最小值 $f_1(\hat\beta)\lt 0$ ，因此解为 $\hat\beta=\dfrac{y'z}{N}+\lambda$ 。

利用软阈值算子（soft-thresholding operator） $S_\lambda(x)=\text{sign}(x)(|x|-\lambda)_+$ ，可将以上三种解统一为
$\hat\beta=S_\lambda(y'z/N)$

其实在我们的设定下，OLS估计量为 $\tilde\beta=y'z/N$ ，因此，将OLS估计量通过一个软阈值算子的操作，就变成了LASSO估计量。

1.2 次梯度

如果引入次梯度（subgradient）的概念，可以更直接地求解 $(1)$ 式。设 $|\beta|$ 的次梯度为 $s\in \text{sign}(\beta)$ ，它的形式是，当 $\beta \neq 0$ 时有 $s= \text{sign}(\beta)$ ，当 $\beta = 0$ 时有 $s\in [-1,1]$ 。根据凸优化（convex optimization）理论，求解 $(1)$ 相当于求解
$-\dfrac{1}{N}z'(y-z\beta) +\lambda s =0$
的解 $(\hat\beta,\hat\lambda)$ 。化简后得到 $y'z/N = \beta+\lambda s\in\beta+\lambda \text{sign}(\beta)$ ，最终同样可以解出 $\hat\beta=S_\lambda(y'z/N)$ 。比如 $\beta=0$ 时，就意味着 $y'z/N \in[-\lambda,\lambda]$ 。

2 多变量情形：循环坐标下降法

我们来看多变量的完整版LASSO问题。将自变量排成 $N\times p$ 的矩阵 $X$ ，我们要求解的是
$\argmin_{\beta\in \mathbb{R}^p} \dfrac{1}{2N}\left\Vert y-X\beta\Vert_2^2+\lambda \right\Vert\beta\Vert_1$

在这里，我们使用循环坐标下降法（cyclic coordinate descent），它的思想是，按一定顺序依次对 $p$ 个参数进行优化，比如按 $j=1,\ldots,p$ 的顺序，在第 $j$ 次优化时，保持其他所有系数不变，变动 $\beta_j$ 使损失函数最小化。

根据以上思想，我们将第 $j$ 次的最优化目标写为
$\argmin_{\beta_j\in \mathbb{R}^p} \dfrac{1}{2N}\left\Vert y-\sum_{k\neq j}x_{\cdot k}\beta_k-x_{\cdot j}\beta_j\right\Vert^2_2+\lambda \sum_{k\neq j}|\beta_k| +\lambda |\beta_j|$
记 $r^{(j)}=y-\sum_{k\neq j}x_{\cdot k}\hat{\beta}_k$ ，这称为partial residual，那么根据第1.1节中的结果，我们可以得出
$\hat\beta_j = S_\lambda (x_{\cdot j}' r^{(j)}/N)$

记 $r = r^{(j)}-x_{\cdot j}\hat\beta_j$ ，上式相当于更新规则
$\hat\beta_j \leftarrow S_\lambda (\hat\beta_j +x_{\cdot j}' r/N)$

由于目标函数是凸的，没有局部的极小值，因此，每次更新一个坐标，最终可以收敛到全局最优解。

Pathwise coordinate descent（逐路径坐标下降）：可以先取一个使 $\hat\beta=0$ 的最小的 $\lambda$ ，然后，略微减小 $\lambda$ 的值，并以上一次得到的 $\hat\beta$ 作为“warm start”，用坐标下降法迭代直到收敛，不断重复这个过程，最终就可以得到在 $\lambda$ 的一系列变化范围内的解。

那么，怎样才能使 $\hat\beta=0$ ？利用次梯度，我们可以知道，对于 $\hat\beta_j=0$ ，必有 $x_{\cdot j}'y /N \in [-\lambda,\lambda]$ ，即要求 $\lambda \geq |x_{\cdot j}'y| /N$ ，若要使整个 $\hat\beta=0$ ，则可取 $\lambda =\max_j |x_{\cdot j}'y| /N$ ，这就是使 $\hat\beta=0$ 的最小的 $\lambda$ 。

3 其他解法

求解LASSO还有其他的解法，如homotopy method，它可以从 $0$ 开始，得到序列型的解的路径，路径是分段线性的。

还有LARS（least angle regression）算法，这是homotopy method之一，可以有效得到分段线性路径。

这里不作展开。

4 正交基

在上面的过程中，如果将自变量正交化，可以大大简化计算。如在第2节中，如果自变量之间是正交的，则 $x_{\cdot j}' r^{(j)}/N = x_{\cdot j}' y/N$ ，此时 $\hat\beta_j$ 就是将 $y$ 对 $x_{\cdot j}$ 做回归的OLS解，通过软阈值算子后的值。