从LASSO回归到结构性稀疏：线性回归的正则项都带来了什么？

会员服务 ·

从LASSO回归到结构性稀疏：线性回归的正则项都带来了什么？

2022 年 3 月 25 日 PaperWeekly

©作者 | 黄秋实

单位 | 香港中文大学（深圳）

研究方向 | 智能电网

本文我们主要关注于以下问题：

1. LASSO Regression 是什么？

2. 稀疏性的意义是什么？（从数学上证明）

3. 为什么 LASSO Regression 可以带来稀疏性？

4. 如何求解 LASSO Regression（附代码）？

5. 线性回归中不同的正则项都带来了什么？

6. 为什么需要结构性稀疏？

7. 如何得到结构性的稀疏？

8. 稀疏性的拓展一些常用的工具箱

LASSO是什么？

LASSO Regression 是一种特殊的线性回归模型。与常见的最小二乘（Least Square Regression）回归相比，LASSO 回归仅仅在损失函数上增加了一个对变量的 -norm（L1-范数）的惩罚：

其中是输入的样本，是一个给定的矩阵（一组给定的基），是我们需要求解的变量，是这个模型的惩罚系数。与最小二乘回归相比，LASSO 回归不仅可以帮助我们得到一个预测模型，同时也会对特征进行一定的选择。对特征的选择反应在变量上时就会呈现出向量具有稀疏性。这听起来很晦涩，所以我们在下面举一个例子来说明。

稀疏性的意义是什么？

我们来考虑一个商业分析场景。假设有一系列销售火爆的同类产品，因为我们也想做一个这样的产品来卖，所以我们很关注它们的销量。当前，主要有 A，B，C，... , Z 这 26 家公司在卖这个产品。我们令公司 A 的产品销量为，其余各家的产品销量以此类推，我们就可以得到一个产品销量的向量。

我们都知道要做一个好的产品，并让它有很好的销量，需要有很多方面的投入。例如：外观设计，功能研发，广告宣传等很多方面都需要进行投入。这样便引出了一个很自然的分析问题：在那些方面进行投入，我们可以获得更好的产品销量呢？这些公司已有的做法可能会给我们带来一定的启示。

我们假设得到了之前这 26 家公司在做产品时每个方面的投入大小，并且每个方面的投入与产品的销量之间呈线性关系，因此我们就可以用 LASSO 回归等线性模型来进行分析。从数学上来说，对于每一个公司的产品销量，由于线性关系的存在，我们可以使用如下线性模型来近似（假设有 100 个可投入方面）：

其中代表在第 i 方面收益/投入比，代表公司 A 在第 i 方面的投入大小。从上面的线性模型中，我们可以发现，我们想要分析的问题其实就是去求解每个可投入方面的收益/投入比。从这些我们知道的产品销量以及每个公司的投入细节中求解收益投入比的方法就是线性回归。

如果我们使用 LASSO 回归来求解上面的这个回归问题，那么我们其实是想要得到：

其中是记录 26 个公司产品销量的向量，是记录 26 个公司在 100 个方面投入大小的数据矩阵（矩阵的任意一行代表对应公司对给个方面的投入大小，矩阵 X 的每一列代表各个公司对对应方面的投入大小），是我们想要了解的每个可投入方面的收益/投入比向量。此时，我们可以看到如果我们最后得到的收益/投入比向量是稀疏的，那么相当于 LASSO 回归模型帮助我们对这 100 个可投入方面进行了选择，仅留下了一些值得投入，对产品销量有较大影响的方面。这便是稀疏性的意义：对特征进行选择，仅留下对分析目标相关性更强的特征。

为什么LASSO Regression可以带来稀疏性？

相比于其他的相对定性的诠释，在这里我们给出一个更加严谨的解释。为了从数学上证明 LASSO Regression 确实可以使优化的目标稀疏，我们首先要回顾一些数学知识。

3.1 次梯度（Subgradient）

对于一个凸函数，存在一个向量是在处的次梯度，那么它需要满足：

我们从函数图像的角度来解释次梯度的意义。我们知道是一条经过点的直线，其中这个向量就是这条直线的斜率（梯度）。所以对于一个次梯度来说，它需要满足它所形成的这条直线要一直处于原函数的下方，即。同时，与梯度不同的，对于函数上的一个点，可能会有很多满足以上条件的向量。如下图所示，在标出的点，可以有很多满足条件的直线，因此这个点的次梯度就会是一个集合。

▲ 一个次梯度的例子

我们知道对于可导的凸函数，有凸函数的一阶性质存在即：

由于凸函数的一阶性质保证了在 w 处能够满足如上性质的向量，只能是对应点的梯度即。因此对于凸函数中的可导部分，它的梯度就等于它的次梯度。而对于函数的不可导部分，函数的次梯度就会是一个集合。举一个最简单的例子，对于函数而言，它的次梯度就是：

对于任意凸函数，由于凸函数的性质，我们可以知道以下关于次梯度的性质：

1. 对于上每一点 w 都有次梯度存在；

2. 在可导的点上，它的次梯度就等于它的梯度；

3. 对于中不可导的点，它的次梯度就会是一个凸集；

4.当且仅当时，就是整个函数的全局最小点（global minimum）。

有了以上的数学知识，我们就可以来解释为什么 LASSO 回归可以带来稀疏性，亦即为什么 -norm 的正则项可以带来稀疏性。

3.2 LASSO回归与稀疏性

我们知道 LASSO 回归的损失函数是：。由于并不是全局可导，所以我们只能使用次梯度来对它进行分析。对于是否变量是否稀疏，我们一一分析变量的每一个分量会在什么情况下为 0 。对于损失函数，任意一个分量上在处有次梯度：

其中是矩阵的第 j 列。我们知道当时，为全局最小点。此时对于的每一个分量都有。基于这一点，我们知道当且仅当。由于是线性拟合时的残差（residual），我们令，于是。因此我们可以看到：

从之前的案例分析中，我们可以知道矩阵，其中 m 是数据（sample）的数量，而 n 是特征（feature）的数量。因此矩阵的第 j 列就是第 j 个特征的数据向量。结合上式我们便可分析得到，当数据特征向量与最终的残差更加正交的话，那么与它对应的变量就更容易等于 0。并且我们还可以看到，随着惩罚参数的增加，变量会愈加稀疏，因此它的每一个分量等于 0 的可能性都越大。

如何求解LASSO Regression？

基于之前对 LASSO 回归稀疏性的分析，最自然地求解 LASSO 回归的方法就是次梯度下降（subgradient descent）。但是次梯度下降的收敛速率很慢，并且需要计算次梯度，计算过程较为复杂。因此常用坐标轴下降法（coordinate descent）来解决 LASSO 回归。由于完整地讲述坐标轴下降法需要的篇幅过长，我们仅给出其对于 LASSO 回归问题的解法。

其中是软门限操作符（soft-thresholdin g），它的具体表示如下：

因此对于 LASSO 回归我们可以用以下代码来求解：

import numpy as np

def Lasso(X, y, Lambda):

    n = X.shape[1]

    w = np.random.randn(n)  # 随机初始化 w

    for i in range(10):  #迭代10次
        for k in range(n):
            y_res = y - np.dot(X, w) + X[:,k] * w[k]

            inner = np.dot(X[:,k], y_res)

            if inner > Lambda:
                w[k] = (inner - Lambda)/(np.linalg.norm(X[:,k])**2)
            elif inner < -Lambda:
                w[k] = (inner + Lambda)/(np.linalg.norm(X[:,k])**2)
            else:
                w[k] = 0

        print('第'+str(i+1)+'步损失函数的值:' + str(np.linalg.norm(y-np.dot(X,w))**2+Lambda*np.linalg.norm(w,ord=1))) # 输出损失函数大小

    return w

X = np.random.randn(10,20)
y = np.random.randn(10)
Lasso(X,y,0.2)

正则项与回归系数轨迹

岭回归，LASSO 回归都是在最小二乘回归的基础上对损失函数添加正则项得到的回归方法。其中岭回归在最初的最小二乘回归的损失函数上增加了一项对回归系数的惩罚，得到：

我们可以看到增加的惩罚项后，损失函数的 Hessian 矩阵从之前的一个半正定矩阵变为了这样一个正定矩阵，这让损失函数具有了强凸的性质（strongly convex）。因为一定可逆而不一定可逆，这很好地简化了问题的求解复杂度。同时，因为强凸性质的存在，岭回归有且仅有一个全局最优解，而与之对比的最小二乘回归，在很多情况下它可能会存在很多的全局最优解。

岭回归的这一系列优势都是由正则项的增加所带而来。同时我们之前分析过的 LASSO 回归，也是因为增加了一项针对回归系数的惩罚才具有变量稀疏与特征选择的功能。所以我们想知道，这些正则项（不仅限于这两种）都做了什么，带来了什么？

我们分别绘制与观察添加不同正则项时，回归系数结果（的每一个分量的大小）与正则项大小（）的关系图。通过观察不同的回归系数的轨迹图，我们可以看到不同正则项所带来的影响。我们分别展示正则项为，，，时，回归系数的轨迹图，这基本包括了所有常用的正则项。

我们可以看到，当添加 范数的平方 进行正则化的时候，随着增大，回归系数整体逐渐减小，并逐渐趋于 0，但是始终不会等于 0。我们从岭回归的解析解中也可以验证这一点，岭回归的解析解为：。当时，的大小完全决定了解析解的大小，因此会产生如上的系数回归轨迹。

而使用与之极为相似的 不平方的 范数进行正则化时，所得的结果却完全不同：

我们可以看到此时，当正则项系数增大到一定程度时，回归系数会变为。范数正则化与岭回归结果不同是因为当时，。因此我们如下分析当添加范数项做正则化时，在满足什么条件时会导致回归系数。我们假设当前为当前最优的选择，则：

由于我们只需要考虑这一部分，对于其他的，其必然满足以上条件。因此，进而我们可以得到当

时，为当前最优的选择。

对于添加范数做正则化的 LASSO 回归，它的回归系数轨迹图如下：

与我们预想中的一样，LASSO 回归会让系数一一分别趋于 0，进而得到一个系数的回归系数向量。我们可以运用次梯度在理论层面证明得到同样的结论，具体细节可以参见上一节。

而对于使用范数进行正则化时，为了凸显它与范数正则化的不同，我们设计让回归系数向量内的每两项受到正则化，即分别为如下优化问题：

在实验中，我们可以分别得到如下结果：

我们可以看到：使用范数进行正则化会让被正则化的系数先同时趋向于一个非零值，再一起趋向于 0（期间一直保持相等）；而使用范数进行正则化会让正则化的系数在同一时，一起变为 0。我们可以看到由于因此只有缩减正则项中最大的参量才会使得正则项大小减小。于是随着的增大，中较大的会被首先缩减，直到中所有参量相等。

为什么要从LASSO到结构性稀疏？

我们已经看到 LASSO 回归可以带来稀疏性，并对特征进行选择。但是，我们对 LASSO 回归所带来的稀疏性进行观察时，我们可以看到 LASSO 回归所带来的稀疏是全局稀疏，也就是说：在没有其他先验知识的情况下，所有变量的被置为 0 的概率相同。

如果我们仅有变量全局稀疏的先验知识，那么 LASSO 回归确实会很好地满足我们的需求。但是，如果我们有了比全局稀疏更多的先验知识，LASSO 回归却会出现一定的问题。例如，LASSO 回归并不会考虑变量之间的相关性（捆绑关系，即有几个变量要么同时为 0，要么同时不为 0）；同时，LASSO 回归也无法准确地进行多标签回归时的变量选择。这些都在一定程度上限制了 LASSO 回归被更为广泛与更为精准的应用。

▲ 结构性稀疏形成的一个凸的非零区域

如何得到结构性的稀疏？

结构性稀疏主要可以分为两种：一种是正则化范数不相互重合的小组稀疏（group sparsity）；一种是正则化范数相互重合的结构稀疏（structural sparsity）[1]。

带有小组稀疏的线性回归问题可以被表示成如下：

其中可以根据先验知识的不同分别选取等。我们可以看到小组稀疏虽然可以容纳能多的先验知识，但是对于回归变量 w 中的每一个分量，都只能允许单一先验知识的存在（正则化范数不互相重合）。而对于一些结构性的知识，小组稀疏就不能完美地满足相应的需求了。因此，我们需要允许正则化范数相互重合，进而得到了结构稀疏（structural sparsity）。结构性稀疏可以很好的使用一族并封闭的支集来表示（union-closed families of supports）即：

给定一个小组集，有：

则带有结构性稀疏的线性回归问题等价于求解如下优化问题：

通常，带有结构性稀疏的线性回归问题是一个凸优化问题，可以使用次梯度下降的方法进行求解。我们也可以使用 CVX，CVXPY 等凸优化问题的求解工具来对设计好的优化问题进行求解。

7.1 一个结构性稀疏的例子

我们简单比较一个带有结构性稀疏的线性回归和一个 LASSO 回归在先验知识上的区别。对于如下优化问题：

它优化得到结果只有；；这三种情况。与传统的 LASSO 回归相比，它相当于增加了一个先验知识即且的情况不存在。

稀疏性的扩展

其实稀疏性并不仅仅局限于参数选择这一个场景中，很多现实中的场景中都存在稀疏性或者稀疏性的变种（可以被转化为稀疏特征）。例如很多惯性比较大的系统（如：电网等）中的系统特征（如：城市区域的用电量）大多具有波动较小（相对稳定）的属性。这些在时间上波动较小，且相对稳定的属性，我们也可以将其转化为稀疏特征来处理。其中常见的一种是 Fused Lasso Signal Approximator：

其中：

这种方法可以将数据在时间上的相关特征作为先验知识加入到模型中去，因此 Fused Lasso的方法常被用来进行时间序列等序列信息的降噪处理。

同时，压缩感知（compressed sensing）中的矩阵补全，矩阵分解里常用的核范数（nuclear norm）的正则化方法也是一个变形的稀疏正则项。例如 Robust Principal Component Analysis [2]：

其中为核范数，它等于对应矩阵奇异值的绝对值的和。我们知道矩阵的非零奇异值个数是矩阵的秩。因此矩阵的低秩性可以使用核范数来正则化。

一些常用的工具包

一般而言，带有稀疏特性（包括小组稀疏，结构性稀疏）的线性回归会是一个凸优化问题，因此我们可以使用 CVX [3]，CVXPY [4] 等凸优化求解工具包来求解。同时一些使用较多的方法如：Fused Lasso Signal Approximator，Isotonic Signal Approximator 等方法可以使用 RegReg 工具包 [5] 来求解。