Stata连享会 主页 || 视频 || 推文 || 知乎 || Bilibili 站
温馨提示: 定期 清理浏览器缓存,可以获得最佳浏览体验。
New!
lianxh
命令发布了:
随时搜索推文、Stata 资源。安装:
. ssc install lianxh
详情参见帮助文件 (有惊喜):
. help lianxh
连享会新命令:cnssc
,ihelp
,rdbalance
,gitee
,installpkg
⛳ Stata 系列推文:
作者:陈卓然 (中山大学岭南学院)
邮箱:chenzhr25@mail2.sysu.edu.cn
编者按:本文主要摘译自下文,特此致谢!
Source: Cattaneo M D, Jansson M, Ma X. Local regression distribution estimators[J]. Journal of Econometrics, 2021. -PDF-
目录
自从 Whitney Newey 开创性地将核平滑应用到非参数或半参数估计中以后,分布函数和密度函数的非参数核估计及其高阶导数在计量经济学中就一直扮演着重要的角色。本文在此基础上进一步研究局部回归分布估计量的大样本性质,以便更好地将其应用于非参/半参数的估计和推断当中。
局部回归分布估计量是通过对一个随机变量
相比于其他类似的估计量,局部回归分布估计量不需要进行预先数据的平滑,因此避开了调整参数选择问题,并且很容易实施和解释。
假定
局部回归分布估计量应用于当
其中,
其中,
局部回归分布估计量是通过经验分布函数的一个最小二乘映射而得到的,这种映射赋予每一个观测值相等的权重。因此将前述的局部回归估计量应用
本文的案例演示都以 Stata 为例,关于 Python 和 R 语言的应用,详见复制文件「CJM_2021_JOE」。
*命令安装
lxhget lpdensity.pkg, install replace
// 或者 net install lpdensity, from(https://raw.githubusercontent.com/nppackages/lpdensity/master/stata) replace
// 关于 Python 和 R 的命令安装包请参考:https://nppackages.github.io/lpdensity
*命令语法
lpdensity Var [if] [in] [,
grid(Var) bw(Var or #) p(#) q(#) v(#) kernel(KernelFn)
scale(#) nomasspoints
bwselect(BwMethod) nlocalmin(#) nuniquemin(#) noregularize nostdvar
cweights(Var) pweights(Var)
genvars(NewVarName)
rgrid(Var) rindex(Var) level(#) ciuniform cisimul(#) separator(#)
plot
estype(ESOpts)
esline_options(ESLineOpts) espoint_options(ESPointOpts)
citype(CIOpts)
ciregion_options(CIRegionOpts)
ciline_options(CILineOpts) ciebar_options(CIEbarOpts)
histogram hiplot_options(HistOpts)
graph_options(GraphOpts)
]
其中,
grid(var)
:指明密度被估计的网格,默认是选择数据集的 5% 到 95%,步长为 0.05;bw(var or #)
:指明用作估计的带宽 (可以是包含带宽的变量,也可以是一个数字),当被省略时,带宽将被通过 bwselect(BwMethod)
来计算;p(#)
:指明用作构建点估计的局部多项式的阶数,默认为 2 阶;q(#)
:指明用来构建置信区间的局部多项式的阶数,默认为 p(#)+1
,通常情况下它应该严格大于 p(#)
;v(#)
:指明待估计的分布函数的导数,v(0)
表示分布函数,v(1)
表示密度函数 (默认);kernel(KernelFn)
:指明用来构建局部多项式估计量的核函数:
triangular
:epanechnikov
:uniform
:scale(#)
:控制估计量的成倍增加或减少;nomassponts
:将不会调整点估计或者标准误即使数据集中有质点;bwselect(BwMethod)
:指明带宽选择的方法:
mse-dpi
:每一个网格点的均方误最优带宽;imse-dpi
:加总均方误最优带宽;mse-rot
:基于高斯推断模型的经验带宽;imse-rot
:加总的基于高斯推断模型的经验带宽。nlocalmin
:指明在每一个局部邻域内不同观测值的最小个数;noregularize
:将每一个局部邻域内不同观测值的最小个数设置为 0;nostdvar
:将不会根据带宽选择来标准化数据;cweights(var)
:指明用于反事实分布构建的权重;pweights(var)
:指明用于抽样的权重,应该是非负的;genvars(NewVarName)
:生成一个用来存储估计结果的变量;rgrid(var)
:指明用来展示结果的网格点;rindex(var)
:指明用来呈现结果的指数;level(#)
:控制置信区间的水平;ciuniform
:计算统一的置信区间,而非逐点置信区间;cisimul(#)
:指明用来构建临界值的模拟次数,默认是 2000 次;separator(#)
:在每 # 个变量之后画一条分割线;plot
:如果被指明的话表示点估计和区间估计都将被绘制;estype(ESOpts)
:指明点估计的绘制形状:
line
:一条曲线 (默认);points
:单独的点;both
:点和曲线;none
:将不会绘制点估计。esline_options(ESLineOpts)
:指明twoway line
的相关选项;espoint_options(ESPointOpts)
:指明twoway scatter
的相关选项;citype(CIOpts)
:指明置信区间的绘制形状:
region
:绘制阴影区域;line
:绘制上下边界;ebar
:误差柱 (error bars);all
:以上全部;none
:将不会绘制置信区间。ciregion_options(CILineOpts)
:指明 twoway rarea
的有关选项;ciline_options(CILineOpts)
:指明 twoway rline
的有关选项;ciebr_options(VIEbarOpts)
:指明 twoway rcap
的有关选项;histgram
:如被指明,将绘制一个直方图;hiplot_options(HistOpts)
:指明额外的 twoway histogram
的选项;graph_options(GraphOpts)
:指明额外绘制选项,例如图例和标签等。反事实分析中的一个关键特征是使用被估计出来的权重体系去估计在全部支撑集 (包括边界和边界邻域) 上的密度函数。为了构建一个反事实密度估计量,我们需要合理设定权重体系
结果变量
那么一个更具挑战性的问题出现了:如果处理组中的个体具有和控制组中的个体相同的协变量分布时,结果变量分布将会是怎样的呢?这时的处理组中结果变量的密度函数被称为反事实密度,记为
此外应该注意到的是,反事实分布还有一个更为有用的解释:假定结果变量是从潜在产出
在实践当中,这种权重体系是未知的,因为倾向得分
在这一节中我们主要展示边际密度,条件密度以及反事实密度是如何通过使用局部分布回归的方法去估计的。具体地,我们采用 Abadie 等 (2002) 研究教育对收入影响数据集中的一个子数据集。这份数据中包括没有登记在 Job Training Partnership Act (JTPA) 的全部个体。
主要的结果变量是在 30 个月的期限内收入的加总,这些个体根据其教育程度被分为两组:
样本量共有 5447,其中有 3927 个高中毕业生。该子样本的描述性统计我们放在下一个应用案例 (用到的是全样本) 中一起展示。众所周知的是教育对于劳动收入具有显著的影响,我们首先分别画出具有高中学历和不具有高中学历两个子样本的收入分布情况,即画出
// 调用 dta 格式数据
lxhuse jtpa.dta, clear
// 下载 csv 格式数据
lxhget jtpa.csv, replace
* Figure 2
capture drop grid // set up grid points
qui gen grid = 1.97 + 0.03 * _n if _n <= 101
capture drop eduFull_*
capture drop edu1_*
capture drop edu0_*
// full sample
qui lpdensity logincome if treatment == 0, bwselect(imse-dpi) grid(grid) ///
ciuniform genvars(eduFull)
// with HS or GED
qui lpdensity logincome if treatment == 0 & hsorged > 0, bwselect(imse-dpi) grid(grid) ///
cweights(hsorged) ciuniform genvars(edu1)
// no HS or GED
qui lpdensity logincome if treatment == 0 & hsorged == 0, bwselect(imse-dpi) grid(grid) ///
ciuniform genvars(edu0)
// Figure 2(a) * NewVarName_f_p 表示用 p(#) 多项式得到的点估计
twoway (line eduFull_f_p eduFull_grid, col(red)) ///
(line edu1_f_p edu1_grid, col(blue)) ///
(line edu0_f_p edu0_grid, col(green)) ///
(rarea eduFull_CI_l eduFull_CI_r eduFull_grid, col(red%20) lcolor(red%0)) ///
(rarea edu1_CI_l edu1_CI_r edu1_grid, col(blue%20) lcolor(blue%0)) ///
(rarea edu0_CI_l edu0_CI_r edu0_grid, col(green%20) lcolor(green%0)), ///
legend(cols(3) order(1 "Full" 2 "HS or GED" 3 "No HS or GED")) ///
xtitle("log(income)") ytitle("")
从上图可以看出高中毕业生和没有取得高中学历的毕业生之间的收入分布是非常不同的。
// counterfactual weights estimation
capture drop male2nonwhite
capture drop male2wkless13
capture drop nonwhite2wkless13
capture drop afdc2male
capture drop cweights
qui gen male2nonwhite = male * nonwhite
qui gen male2wkless13 = male * wkless13
qui gen nonwhite2wkless13 = nonwhite * wkless13
qui gen afdc2male = afdc * male
qui reg hsorged male nonwhite wkless13 married male2nonwhite male2wkless13 nonwhite2wkless13 afdc2male ///
age2225 age2629 age3035 age3644 age4554 if treatment == 0
qui predict cweights
qui replace cweights = (1 - cweights) / cweights
qui replace cweights = cweights * hsorged if hsorged > 0
// estimate counterfactual density
capture drop educ_*
qui lpdensity logincome if treatment == 0 & hsorged > 0, bwselect(imse-dpi) grid(grid) ///
cweights(cweights) ciuniform genvars(educ)
// Figure 2(b)
twoway (line educ_f_p educ_grid, col(red)) ///
(line edu1_f_p edu1_grid, col(blue)) ///
(line edu0_f_p edu0_grid, col(green)) ///
(rarea educ_CI_l educ_CI_r educ_grid, col(red%20) lcolor(red%0)) ///
(rarea edu1_CI_l edu1_CI_r edu1_grid, col(blue%20) lcolor(blue%0)) ///
(rarea edu0_CI_l edu0_CI_r edu0_grid, col(green%20) lcolor(green%0)), ///
legend(cols(3) order(1 "Counterfactual" 2 "HS or GED" 3 "No HS or GED")) ///
xtitle("log(income)") ytitle("")
从上图我们可以看到
自选择和处理效应异质性是因果推断和社会经济学中需要重点考虑的方面。众所周知的是,包括平均处理效应或处理组的处理效应在内的古典处理效应,由于存在不完全服从个体从而导致即便在充分随机分组的条件下也不能被识别。基于此,我们首先应用局部回归分布估计量的密度函数来进行工具变量的识别检验,然后估计一下服从者潜在产出的密度函数。
我们采用和 3.2 节中相同的数据集。JIPA 是一项大型的对经济困难或存在就业障碍的个体提供再就业培训的公共福利计划。个体被随机提供 JIPA 培训,但是在被提供 JIPA 培训的样本当中只有 67% 的个体接受。因此 JTPA 项目是一个有效的工具变量用以研究就业培训对收入的影响。
我们首先通过局部分布估计量
lxhuse jtpa.dta, clear
qui gen grid = 1.97 + 0.03 * _n if _n <= 101
// full sample
capture drop full_*
qui lpdensity logincome, bwselect(imse-dpi) grid(grid) genvars(full) ciuniform
// JTPA offer
capture drop offer1_*
capture drop offer0_*
// JTPA offer
qui lpdensity logincome if instrument == 1, bwselect(imse-dpi) grid(grid) ///
ciuniform genvars(offer1)
// no JTPA offer
qui lpdensity logincome if instrument == 0, bwselect(imse-dpi) grid(grid) ///
ciuniform genvars(offer0)
// JTPA Enroll
capture drop treat1_*
capture drop treat0_*
// Enrolled in JTPA
qui lpdensity logincome if treatment == 1, bwselect(imse-dpi) grid(grid) ///
ciuniform genvars(treat1)
// not Enrolled in JTPA
qui lpdensity logincome if treatment == 0, bwselect(imse-dpi) grid(grid) ///
ciuniform genvars(treat0)
// weights for compliers, Abadie (2003)
capture drop cweights
capture drop cweights0
capture drop cweights1
capture drop hsorged2male
capture drop hsorged2nonwhite
capture drop hsorged2wkless13
capture drop male2nonwhite
capture drop male2wkless13
capture drop nonwhite2wkless13
qui gen hsorged2male = hsorged * male
qui gen hsorged2nonwhite = hsorged * nonwhite
qui gen hsorged2wkless13 = hsorged * wkless13
qui gen male2nonwhite = male * nonwhite
qui gen male2wkless13 = male * wkless13
qui gen nonwhite2wkless13 = nonwhite * wkless13
qui reg instrument male hsorged nonwhite wkless13 married hsorged2male hsorged2nonwhite ///
hsorged2wkless13 male2nonwhite male2wkless13 nonwhite2wkless13
qui predict cweights
qui gen cweights0 = cweights
qui gen cweights1 = cweights
qui replace cweights = 1 - instrument * (1 - treatment) / cweights - ///
(1 - instrument) * treatment / (1 - cweights)
qui replace cweights1 = treatment * (instrument - cweights1) / (cweights1 * (1 - cweights1))
qui replace cweights0 = (1 - treatment) * (1 - instrument - 1 + ///
cweights0) / (cweights0 * (1 - cweights0))
// density for compliers
capture drop complier_*
qui lpdensity logincome, bwselect(imse-dpi) grid(grid) cweights(cweights) ///
ciuniform genvars(complier)
capture drop complier0_*
qui lpdensity logincome, bwselect(imse-dpi) grid(grid) cweights(cweights0) ///
ciuniform genvars(complier0)
capture drop complier1_*
qui lpdensity logincome, bwselect(imse-dpi) grid(grid) cweights(cweights1) ///
ciuniform genvars(complier1)
// Figure 4(a)
twoway (line full_f_p full_grid, col(red)) ///
(line complier_f_p complier_grid, col(blue)) ///
(rarea full_CI_l full_CI_r full_grid, col(red%20) lcolor(red%0)) ///
(rarea complier_CI_l complier_CI_r complier_grid, col(blue%20) lcolor(blue%0)), ///
legend(cols(3) order(1 "Full" 2 "Compliers")) ///
xtitle("log(income)") ytitle("") title("a") saving(a)
// Figure 4(b)
twoway (line offer0_f_p offer0_grid, col(red)) ///
(line offer1_f_p offer1_grid, col(blue)) ///
(rarea offer0_CI_l offer0_CI_r offer0_grid, col(red%20) lcolor(red%0)) ///
(rarea offer1_CI_l offer1_CI_r offer1_grid, col(blue%20) lcolor(blue%0)), ///
legend(cols(3) order(1 "JTPA Offer = 0" 2 "1")) ///
xtitle("log(income)") ytitle("") title("b") saving(b)
// Figure 4(c)
twoway (line treat0_f_p treat0_grid, col(red)) ///
(line treat1_f_p treat1_grid, col(blue)) ///
(rarea treat0_CI_l treat0_CI_r treat0_grid, col(red%20) lcolor(red%0)) ///
(rarea treat1_CI_l treat1_CI_r treat1_grid, col(blue%20) lcolor(blue%0)), ///
legend(cols(3) order(1 "JTPA Enrollment = 0" 2 "1")) ///
xtitle("log(income)") ytitle("") title("c") saving(c)
// Figure 4(d)
twoway (line complier0_f_p complier0_grid, col(red)) ///
(line complier1_f_p complier1_grid, col(blue)) ///
(rarea complier0_CI_l complier0_CI_r complier0_grid, col(red%20) lcolor(red%0)) ///
(rarea complier1_CI_l complier1_CI_r complier1_grid, col(blue%20) lcolor(blue%0)), ///
legend(cols(3) order(1 "x(0) for compliers" 2 "x(1) for compliers")) ///
xtitle("log(income)") ytitle("") title("d") saving(d)
graph combine a.gph b.gph c.gph d.gph
在图 (a) 中,我们绘制了全样本和只针对服从者的收入分布,两者看起来很相似,尽管后者有着更高的均值并且左尾更瘦。从图 (b) 中可以看到被提供 JIPA 和没有被提供 JIPA 的收入分布是显著不同的。在图 (c) 中,我们也绘制出是否接受 JIPA 对个体收入分布的影响,不出意外,两者的收入分布差距更大。
在图 (d) 中,我们构建出对于服从者而言的潜在收入分布,并估计了在这样的工具变量设定背景下被识别的分布处理效应。可以看出,JIPA 确实对收入有着正的显著的影响,而且它证实了自选择的存在:那些参与 JIPA 的个体平均而言会受益最多,其次是那些被认为在 “在无差异边缘” 的服从者。
Note:产生如下推文列表的 Stata 命令为:
lianxh 估计, m
安装最新版lianxh
命令:
ssc install lianxh, replace
免费公开课
最新课程-直播课
专题 | 嘉宾 | 直播/回看视频 |
---|---|---|
⭐ 最新专题 | 文本分析、机器学习、效率专题、生存分析等 | |
研究设计 | 连玉君 | 我的特斯拉-实证研究设计,-幻灯片- |
面板模型 | 连玉君 | 动态面板模型,-幻灯片- |
面板模型 | 连玉君 | 直击面板数据模型 [免费公开课,2小时] |
⛳ 课程主页
⛳ 课程主页
关于我们
课程, 直播, 视频, 客服, 模型设定, 研究设计, stata, plus, 绘图, 编程, 面板, 论文重现, 可视化, RDD, DID, PSM, 合成控制法
等
连享会小程序:扫一扫,看推文,看视频……
扫码加入连享会微信群,提问交流更方便
✏ 连享会-常见问题解答:
✨ https://gitee.com/lianxh/Course/wikis
New!
lianxh
和songbl
命令发布了:
随时搜索连享会推文、Stata 资源,安装命令如下:
. ssc install lianxh
使用详情参见帮助文件 (有惊喜):
. help lianxh