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

发布时间:2020-12-07 阅读 2267

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 回归的参数可以由如下最优化问题求解:

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

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

具体而言,相对于 OLS 方法, Lasso 方法添加了 λ(j=1p|βj|) 作为惩罚项(也称正则化函数),即

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

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

其中, J(θ) 为目标函数,L(yi,f(xi;θ)) 称作损失函数,衡量单个样本和拟合值之间的误差,1ni=1nL(yi,f(xi;θ)) 为代价函数,衡量整个样本集的平均误差, λΦ(θ) 为正则化函数或惩罚项,当 λ0 时,上述目标函数则从 OLS 回归转向惩罚回归。

另外,很多教材,包括 Stata16 的手册中,在代价函数前面都会乘以 12n,原因有二:

  • 乘以 1n 是为了针对不同的样本容量,λ 的可比性;
  • 乘以 12 纯属是为了求解最优化问题的方便,因为损失函数的指数为 2。

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

可能有人想问,为何是 L1 范数而不是 L2 范数或者 Lq 范数呢?原因主要是 L1 范数有如下特殊性:

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

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

不同范数惩罚项的约束区域
不同范数惩罚项的约束区域

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

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

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

Lasso回归(左)与岭回归(右)图示
Lasso回归(左)与岭回归(右)图示

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

综上,综合了解的稀疏性和计算的便捷性两个主要优点,Lasso 方法的惩罚项选择了 L1 范数。

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

2.2 从 Lasso 到弹性网回归

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

其中,α 为弹性网惩罚系数,设定 α=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

相关推文

相关课程

连享会-直播课 上线了!
http://lianxh.duanshu.com

免费公开课:


课程一览

支持回看

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

Note: 部分课程的资料,PPT 等可以前往 连享会-直播课 主页查看,下载。


关于我们

  • Stata 连享会 由中山大学连玉君老师团队创办,定期分享实证分析经验。直播间 有很多视频课程,可以随时观看。
  • 连享会-主页知乎专栏,300+ 推文,实证分析不再抓狂。
  • 公众号推文分类: 计量专题 | 分类推文 | 资源工具。推文分成 内生性 | 空间计量 | 时序面板 | 结果输出 | 交乘调节 五类,主流方法介绍一目了然:DID, RDD, IV, GMM, FE, Probit 等。
  • 公众号关键词搜索/回复 功能已经上线。大家可以在公众号左下角点击键盘图标,输入简要关键词,以便快速呈现历史推文,获取工具软件和数据下载。常见关键词:课程, 直播, 视频, 客服, 模型设定, 研究设计, stata, plus, 绘图, 编程, 面板, 论文重现, 可视化, RDD, DID, PSM, 合成控制法

连享会主页  lianxh.cn
连享会主页 lianxh.cn

连享会小程序:扫一扫,看推文,看视频……

扫码加入连享会微信群,提问交流更方便

✏ 连享会学习群-常见问题解答汇总:
https://gitee.com/arlionn/WD

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