Stata：拉索开心读懂-Lasso入门

发布时间：2020-12-07 阅读 9430

Stata 连享会主页 || 视频 || 推文 || 知乎

温馨提示： 定期清理浏览器缓存，可以获得最佳浏览体验。

New！ lianxh 命令发布了：
随时搜索连享会推文、Stata 资源，安装命令如下：
. ssc install lianxh
使用详情参见帮助文件 (有惊喜)：
. help lianxh

✌ 课程详情： https://gitee.com/arlionn/Course | lianxh.cn

✨ 课程主页： https://gitee.com/arlionn/Course

作者：杨继超 (中山大学)
E-Mail: yangjch7@mail2.sysu.edu.cn

1. 引言

你是否总是为摸不准你女神或男神喜好只能感叹恋爱不易？

你是否觉得自己已经用尽全力，可换来的却仍是一句“你是个好人”？

好消息！好消息！好消息！

今日我们请来一位大师将手把手帮你解决你的烦恼。

在传授大师的秘籍之前，先设想一个例子：

假设有小 A 、小 B 和小 C 三人，各方面基本条件都差不多，也都谈过一些恋爱，而且此前谈过的对象类型和数量几乎无差别。

但这三人在失恋后重新找到新对象的速度和难易程度却完全不同，具体为小 B 最快，小 A 次之，小 C 最差。

小 A 和小 C 就如同正在看推文的你，百思不得其解。

原来，小 A 属于“钢铁直男”，虽然也交往过不少女友，但是从头到尾就学会了一招——送花加多喝热水。

而小 B 却不同，他细心且善总结，思维方式不单一，总结出女生大多都不仅喜欢玫瑰花，还都喜欢会使用 Stata 的男生等其它关键点（假设别当真^_^）。

小 C 属于“蠢直男”，小 C 和小 B 不同，思维方式比小 B 还丰富，小 C 摸清了每一位女友的偏好。例如，有的喜欢头发茂密的男生，有的喜欢学过连享会课程的男生，甚至有的只喜欢发过 AER 的男生等等，但是没有总结让女生喜欢的共通点。就好像某些搞题海战术的同学，做过的每道题都会了，但遇到新题又不会了。

这三位的情况就如下图所示的三种情况一般，左边的图和小 A 类似，不管什么数据直接上如玫瑰花和多喝热水一般的 OLS，结果通常为“欠拟合”。中间的图和小 B 类似，如小 B 一样改变了思维方式，即变换了回归函数的形式，结果通常为“好的拟合”。右边的图和小 C 类似，如小 C 一样摸清了每一个女友的偏好，数据拟合得很好，但是没有总结出任何规律，这便为我们所熟知的高维数据（自变量相对于数据量而言过大）模型，结果通常为“过拟合”。

图片来源：Overfitting and Underfitting.pdf

对小 A 而言，解决办法较为简单，就是稍微把这位“钢铁直男”的思维“掰弯”一点，转变一下思维方式，即调整设定的回归函数。

而对于小 C 而言，其思维方式已经够丰富了，解决办法有两个，第一种办法就是加大数据量的训练。毕竟小 C 极其细心，只要给他足够多的机会交女朋友，极限情况就是让他摸清全世界所有类型的女生偏好，他必将成为恋爱场上的王者（希望他不要成为渣男）。

显然，第一种方法在现实世界里很不现实，想得倒是很美！

于是就有了第二种办法，就是寻求恋爱大师的帮助。其实小 C 最致命的缺点并非过分细心，而是缺乏总结的能力，让大师教他如何总结，便可将缺点转为优点。

下面有请大师登场！

不好意思搞错了。。。

今天要跟大家隆重介绍的恋爱大师是 Lasso，全称为 Least absolute shrinkage and selection operator，他将教你抓住主要矛盾，即如何在纷繁复杂的男神和女神偏好中找到最为关键的偏好。

回到正题，Lasso 当然主要还是为了解决在高维数据下的过拟合问题。下面将具体讲解 Lasso 回归的理论和实操。

2. Lasso 究竟为何方神圣

2.1 从 OLS 到 Lasso

在介绍 Lasso (Tibshirani, 1996) 方法之前，仍然还是从大家最为熟悉的 OLS 方法说起。通过简单的推导，OLS 回归的参数可以由如下最优化问题求解：

{\hat{β}}_{o l s} = \arg min_{β} {\frac{1}{2 n} \sum_{i = 1}^{n} {(y_{i} - x_{i}^{'} β)}^{2}} (1)

一般而言，通过（1）解出的系数都不为 0。如果自变量个数 $p$ 相对于样本量 $n$ 而言过大，即 $p > n$ ，则 OLS 估计的参数没有唯一解，出现前文提到的过拟合问题。

因此，上述估计过程需要一定的约束（或称为正则化）作为惩罚，剔除部分自变量，以解决自变量相对于样本量过大而产生的过拟合问题。

具体而言，相对于 OLS 方法， Lasso 方法添加了 $λ (\sum_{j = 1}^{p} | β_{j} |)$ 作为惩罚项（也称正则化函数），即

\hat{β} = \arg min_{β} {\frac{1}{2 n} \sum_{i = 1}^{n} {(y_{i} - x_{i}^{'} β)}^{2} + λ (\sum_{j = 1}^{p} | β_{j} |)} (2)

其中， $\sum_{j = 1}^{p} | β_{j} |$ 为 $β$ 的 $L_{1}$ 范数（ $L_{q}$ 范数则为 $\sum_{j = 1}^{p} {| β_{j} |}^{q}$ ）。 $λ (\sum_{j = 1}^{p} β_{j})$ 为惩罚项， $λ$ 为调节参数 (tuning parameter)，控制惩罚的力度，其取值的确定通常采用交叉验证 (cross-validation)、 adaptive 和 plugin 三种方法来确定，即找寻使得样本外均方误差 MSE 最小的 $λ$ 。

其实，OLS 和 Lasso 都可以划归为一类问题，具体为：

min J (θ) = \frac{1}{n} \sum_{i = 1}^{n} L (y_{i}, f (x_{i}; θ)) + λ Φ (θ) (3)

其中， $J (θ)$ 为目标函数， $L (y_{i}, f (x_{i}; θ))$ 称作损失函数，衡量单个样本和拟合值之间的误差， $\frac{1}{n} \sum_{i = 1}^{n} L (y_{i}, f (x_{i}; θ))$ 为代价函数，衡量整个样本集的平均误差， $λ Φ (θ)$ 为正则化函数或惩罚项，当 $λ \neq 0$ 时，上述目标函数则从 OLS 回归转向惩罚回归。

另外，很多教材，包括 Stata16 的手册中，在代价函数前面都会乘以 $\frac{1}{2 n}$ ，原因有二：

乘以 $\frac{1}{n}$ 是为了针对不同的样本容量， $λ$ 的可比性；
乘以 $\frac{1}{2}$ 纯属是为了求解最优化问题的方便，因为损失函数的指数为 2。

公式(2)亦可以等价为：

min_{β} {\frac{1}{2 n} \sum_{i = 1}^{n} {(y_{i} - x_{i}^{'} β)}^{2}} (4)

s t . \sum_{j = 1}^{p} | β_{j} | \leq t (t \geq 0)

可能有人想问，为何是 $L_{1}$ 范数而不是 $L_{2}$ 范数或者 $L_{q}$ 范数呢？原因主要是 $L_{1}$ 范数有如下特殊性：

如果 $t$ 足够小， $L_{1}$ 范数可以产生稀疏解（某些自变量系数为 0 的解），而且仅仅只有部分自变量系数为 0，这一点在 $q > 1$ 的 $L_{q}$ 范数中是无法保证的；
在 $q < 1$ 的 $L_{q}$ 范数中，虽然能够得到稀疏解，但是上述最优化问题将不再是一个凸规划的最优化求解问题，这会在计算上造成非常大的困难。

如下图所示，在只有两个待估参数的模型中，不同范数惩罚项的约束区域可用如下图形表示：

图片来源：Statistical Learning with Sparsity: The Lasso and Generalizations.pdf

下面以 $L_{1}$ 范数和 $L_{2}$ 范数为例，以图形的方式直观地分析 Lasso 方法为何更容易得到稀疏解，其中 $L_{2}$ 范数作为惩罚项的方法被称作岭回归 (Ridge Regession) 方法 (Hoerl and Kennard, 1970)。

如下图所示， $\hat{β}$ 为无惩罚项的参数估计值，即 OLS 估计值。在有约束的情况下，最优化问题会类似于微观经济学里的餍足偏好。其中，岭回归的约束集为圆形，最优解通常不在坐标轴上，得不到稀疏解，即无法通过约束使得某些变量系数为 0 来简化模型。而 Lasso 约束集为菱形，恰好可以弥补这一缺点。

图片来源：Statistical Learning with Sparsity: The Lasso and Generalizations.pdf

综上，综合了解的稀疏性和计算的便捷性两个主要优点，Lasso 方法的惩罚项选择了 $L_{1}$ 范数。

另外，Lasso 还有另外一种用法，称为 postselection。简单来讲，先采用 Lasso 回归筛选模型自变量，再采用无惩罚项的模型进行回归，这种方法的结果会作为各种模型结果对比的参照，且此方法也可简单用于变量的筛选。

2.2 从 Lasso 到弹性网回归

Lasso 回归与岭回归回归相比虽然达到了变量选择的目的，但是也有一定的局限性。通俗来讲，Lasso 丢变量太狠，而岭回归又太温柔，为了综合二者，于是 Zou and Hastie (2005) 提出了弹性网回归，具体为

\hat{β} = \arg min_{β} {\frac{1}{2 n} \sum_{i = 1}^{n} {(y_{i} - x_{i}^{'} β)}^{2} + λ [α \sum_{j = 1}^{p} | β_{j} | + \frac{(1 - α)}{2} \sum_{j = 1}^{p} β_{j}^{2}]} (5)

其中， $α$ 为弹性网惩罚系数，设定 $α = 0$ 则为岭回归，设定 $α = 1$ 则为 Lasso 回归。

3. Stata 范例

本节将采用 Stata16 软件中的系统数据auto.dta进行演示，虽然这个数据不是采用 Lasso 回归的最优数据，但本节仅仅只是用此数据做演示之用，以期让大家能对 Lasso 回归的基本操作有一个入门的了解。

Step1: 将数据分为训练集与验证集

sysuse "auto.dta", clear//调用系统数据
splitsample , generate(sample) split(.75 .25) rseed(12345)//将数据随机按照75%和25%分成两部分，前者为训练集用于回归，后者为验证集用于验证，rseed(12345)是为了让结果可重复
label define slabel 1 "Training" 2 "Validation"
label values sample slabel//给训练集添加标签为Training，给验证集添加标签为Validation

拆分数据后用tabulate看下样本结构,

tabulate sample

结果如下:

  sample |      Freq.     Percent        Cum.
------------+-----------------------------------
   Training |         56       75.68       75.68
 Validation |         18       24.32      100.00
------------+-----------------------------------
      Total |         74      100.00

Step2: 用 OLS 做基准回归

quietly regress mpg i.foreign i.rep78 headroom weight turn gear_ratio price trunk length displacement if sample==1//先采用OLS回归
estimates store ols

跑完 OLS，用lassogof计算均方误差，

lassogof ols, over(sample)//计算均方误差，over(sample)表示同时计算训练集和验证集的MSE

结果为：

Penalized coefficients
-------------------------------------------------------------
Name             sample |         MSE    R-squared        Obs
------------------------+------------------------------------
ols                     |
               Training |    4.490814       0.8710         56
             Validation |     33.2273      -0.5475         18
-------------------------------------------------------------

不出所料，训练集的均方误差比验证集的要小很多，毕竟这些都是谈过的女朋友，脾气还是摸得比较准的。^_^

Step3: 基于交叉检验的 Lasso 回归

lasso linear mpg i.foreign i.rep78 headroom weight turn gear_ratio price trunk length displacement if sample==1, nolog rseed(12345)//调节参数默认采用cv交叉检验
estimates store cv

结果为：

Lasso linear model                          No. of obs        =         51
                                            No. of covariates =         14
Selection: Cross-validation                 No. of CV folds   =         10

--------------------------------------------------------------------------
         |                                No. of      Out-of-      CV mean
         |                               nonzero       sample   prediction
      ID |     Description      lambda     coef.    R-squared        error
---------+----------------------------------------------------------------
       1 |    first lambda    4.898601         0       0.0282     37.17473
      32 |   lambda before    .2738715         8       0.7521     8.962883
    * 33 | selected lambda    .2495415         8       0.7525     8.948511
      34 |    lambda after    .2273729         9       0.7517     8.977777
      36 |     last lambda     .188769         8       0.7502     9.031694
--------------------------------------------------------------------------
* lambda selected by cross-validation.

根据回归结果，交叉检验选择 ID=33 的调节参数，模型最终选择了 8 个自变量。

也可以用cvplot这个命令直观地看 $λ$ 取值不同对样本外均方误差的影响，

cvplot

结果如图所示：

可以看出，CV 交叉检验选择最优的 $λ$ 为 0.25.

不仅如此，我们可以运行coefpath命令直观地观察自变量系数随着不同 $λ$ 取值的系数路径图。

coefpath, lineopts(lwidth(thick)) legend(on position(3) cols(1))

结果如图所示：

此外，还可以进一步用lassoknots看看不同 $λ$ 取值的不同选择过程，

lassoknots

结果为：

-------------------------------------------------------------------------------------------------------------------
       |              No. of   CV mean |
       |             nonzero     pred. |                       Variables (A)dded, (R)emoved,
    ID |   lambda      coef.     error |                            or left (U)nchanged
-------+-------------------------------+---------------------------------------------------------------------------
     2 | 4.463423          1  34.71904 | A weight
     3 | 4.066905          2  31.20129 | A length
     4 | 3.705612          3  28.30943 | A 5.rep78
    12 | 1.760466          4  13.16168 | A turn
    15 | 1.331728          5  11.10733 | A displacement
    25 |  .525261          6  9.177736 | A price
    27 | .4360809          7  9.096081 | A gear_ratio
    28 | .3973407          8   9.08659 | A headroom
  * 33 | .2495415          8  8.948511 | U
    34 | .2273729          9  8.977777 | A 2.rep78
    35 | .2071737          8  9.007081 | R weight
    36 |  .188769          8  9.031694 | U
-------------------------------------------------------------------------------------------------------------------
* lambda selected by cross-validation.

最后，可以根据上述结果进一步做敏感性分析，敏感性分析的目的是考察 $λ$ 一个微小的变化是否会引起结果产生巨大的差异。

先采用lassoselect手动选择一个最优 $λ$ 值附近的调节参数，

lassoselect id = 34
estimates store hand

然后采用lassogf命令比较手动选择的结果与交叉检验的结果差异，

lassogof hand cv if sample==2

结果如下：

Penalized coefficients
-------------------------------------------------
       Name |         MSE    R-squared        Obs
------------+------------------------------------
       hand |    35.06322      -0.6330         18
         cv |    34.89029      -0.6249         18
-------------------------------------------------

-Step3: 选择调节参数

上一步中提到，lasso linear默认的是用 CV 交叉检验方法选择调节参数，下面将采用 adaptive 和 plugin 两种方法重新选择，最后对比三种方法所对应的样本外均方误差，再确定最优的调节参数。

lasso linear mpg i.foreign i.rep78 headroom weight turn gear_ratio price trunk length displacement if sample==1, nolog rseed(12345) selection(adaptive)//采用adaptive方法
estimates store adaptive
lasso linear mpg i.foreign i.rep78 headroom weight turn gear_ratio price trunk length displacement if sample==1, nolog rseed(12345) selection(plugin)//采用plugin方法
estimates store plugin

为节省篇幅，省去了上述命令的结果，此处只给出不同方法的对比结果，仍然采用lassogof命令，

lassogof ols cv adaptive plugin if sample==2

结果为：

Penalized coefficients
-------------------------------------------------
       Name |         MSE    R-squared        Obs
------------+------------------------------------
        ols |     33.2273      -0.5475         18
         cv |    34.89029      -0.6249         18
   adaptive |    35.85554      -0.6699         18
     plugin |    18.06601       0.1586         18
-------------------------------------------------

从对比结果可见，采用 plugin 方法的均方误差最小。

Step4: 弹性网回归

这步将采用elasticnet做弹性网回归，弹性网回归里面内嵌了 Lasso 和岭回归，所以本步骤中也会给出岭回归代码，具体如下：

elasticnet linear mpg i.foreign i.rep78 headroom weight turn gear_ratio price trunk length displacement if sample==1, alpha(.25 .5 .75) nolog rseed(12345)//alpha分别取0.25、0.5和0.75的弹性网回归
estimates store enet
elasticnet linear mpg i.foreign i.rep78 headroom weight turn gear_ratio price trunk length displacement if sample==1, alpha(0) nolog rseed(12345)//岭回归
estimates store ridge

其中选项alpha()设定 $α$ 的取值，取值为 0 则为岭回归，取值为 1 为 Lasso，其他取值为弹性网回归，为节约篇幅，大家可以自行运行上述命令。

Note: 调节参数的选取默认为 CV 交叉检验。

Step5: 模型对比

根据本文所介绍的内容，调节参数有 3 种选择方法，除去 ols,模型有 lasso、岭回归和弹性网回归 3 种，模型做法又可以分为 postselection 和非 postselection 两种做法，因此组合起来一共是 18 种组合，大家可以自行探索，其核心还是比较各种模型的均方误差大小，从而选择最优的模型结果，本文只演示部分对比结果。

首先仍然是采用lassogof进行对比，

lassogof cv adaptive enet ridge plugin if sample==2

结果为：

Penalized coefficients
-------------------------------------------------
       Name |         MSE    R-squared        Obs
------------+------------------------------------
         cv |    34.89029      -0.6249         18
   adaptive |    35.85554      -0.6699         18
       enet |    31.30704      -0.4580         18
      ridge |    29.47994      -0.3729         18
     plugin |    18.06601       0.1586         18
-------------------------------------------------

上述结果表明应该选择 plugin 方法的结果。

接下来，采用 postselection 模型进行对比，

lassogof cv adaptive enet plugin if sample==2, postselection

结果为：

Penalized coefficients
-------------------------------------------------
       Name |         MSE    R-squared        Obs
------------+------------------------------------
         cv |     36.9961      -0.7230         18
   adaptive |    36.90559      -0.7188         18
       enet |    38.19441      -0.7788         18
     plugin |    39.96924      -0.8614         18
-------------------------------------------------

上述结果表明应该选择 adaptive 方法的结果。

学到这里，有了 Lasso 这位爱情大师坐镇，有数据的朋友可以丢进 Stata 跑一跑了，都说女人心海底针，在 Lasso 大师面前都能准确拿捏了。当然，女生也可以借力 Lasso 大师成功 pick 你心中的欧巴。

那么问题来了，正在看推文的你有“数据”吗？

【预告】 下篇推文将讲解如何 Lasso 做因果推断，敬请期待！

4. 参考文献和资料

Tibshirani, R. Regression shrinkage and selection via the lasso. Journal of the Royal Statistical Society, Series B, 1996, 58(1): 267–288. -PDF-
Zou, H., and T. Hastie. Regularization and variable selection via the elastic net[J]. Journal of the Royal Statistical Society, Series B, 2005, 67(2): 301–320.-PDF-
Hoerl, A., and Kennard, R. Ridge Regression: Biased Estimation for Nonorthogonal Problems[J]. Technometrics, 42(1), 80-86. -PDF-
David Drukker, The Stata Blog, An introduction to the lasso in Stata
Hastie, T. et al. Statistical Learning with Sparsity: The Lasso and Generalizations[M]. CRC Press, 2015. -PDF-
陈强，计量经济学及 Stata 应用推文，Stata 16 新功能之 Lasso 系列（一）：Lasso Basics

专题	嘉宾	直播/回看视频
⭐ 最新专题		因果推断, 空间计量，寒暑假班等
⭕ 数据清洗系列	游万海	直播, 88 元，已上线
研究设计	连玉君	我的特斯拉-实证研究设计，-幻灯片-
面板模型	连玉君	动态面板模型，-幻灯片-
面板模型	连玉君	直击面板数据模型 [免费公开课，2 小时]

Stata：拉索开心读懂-Lasso入门

1. 引言

2. Lasso 究竟为何方神圣

2.1 从 OLS 到 Lasso

2.2 从 Lasso 到弹性网回归

3. Stata 范例

4. 参考文献和资料

相关推文

相关课程

课程一览

关于我们

热门资讯