1. 命令安装与方法介绍
2. 模拟生成非线性相关数据
3. 进行非参数估计
4. 总结
5. 参考资料
6. 相关推文
资源共享
关于我们
RDD：断点回归的非参数估计及Stata实现
2020-05-27
崔颖
10990

连享会   主页 || 推文 || 知乎 || B 站 || 在线课堂

New！ 搜推文，找资料，用 lianxh 命令：
安装： ssc install lianxh, replace
使用： lianxh 合成控制
       lianxh DID + 多期, w

作者：崔颖（中央财经大学）
邮箱：cuiying@email.cufe.edu.cn

Source: Non-Parametric Regression Discontinuity (Francis, 2013)

 

  • Title: RDD：断点回归的非参数估计及Stata实现
  • Keywords: 非参数, 断点回归, RDD, rd, rdrobust, Cattaneo, rcspline, 最优带宽

 

断点回归 (RDD) 的应用日益广泛。在 「连享会-RDD 专辑1」, 曹昊煜 (2021), 亢延锟 (2020), 张子楠 (2020), 张子楠 (2021), 李鑫 (2020), 谭睿鹏 (2021), 谭睿鹏 (2020) 中，我们对 RDD 的原理、平滑性检验方法等主题进行了介绍。

本文介绍 RDD 的非参数估计在 Stata 中的实现方法，主要采用 rdrobust 命令。该命令是由哥伦比亚大学 Sebastian Calonico 教授、普林斯顿大学 Matias D. Cattaneo 教授及众合作者共同开发。Google 的网页 RD software package 提供了丰富的学习资料，包括许多相关论文的原始数据及复制结果代码。

 

1. 命令安装与方法介绍

net install rdrobust, from(http://www-personal.umich.edu/~cattaneo/rdrobust)

RD 可以用来识别自然实验或结构性政策变化附近的局部处理效应。

例如，如果你关心政府奖金对大学入学情况有怎样的影响，你可能会想要将那些获得政府奖金的学生和未获得政府奖金的学生进行比较。但这种方法是存在问题的，因为获得政府奖金的低收入家庭学生与未获得政府奖金的学生可能在多方面均存在差异。

应用 RD 方法的前提条件是个人不能通过合理低报收入水平而获得政府奖金，那些在断点附近的人自报收入分布情况应该和非断点附近的人基本上保持一致。

如果政府奖金资格确定的收入线是未知的，那么，此前提条件可能是合理的。即使学生会系统性地低报他们的收入，但因他们并不知道实际确认资格的收入分界线，可以认为那些收入在断点上下的学生随机抽取自相同的池子，仅是否收到政府奖金这一项差异。

缺乏实验数据的计量经济学识别方法往往需要建立在外生性假定基础之上。也就是说，x 对 y 的影响与误差项 u 不相关。在外生变量直接导致被解释变量变化的情况下，回归识别因果效应才是充分有效的。

在上述例子中，显然，不能简单地将 y (GPA、出勤率、毕业率等)的变化归结为政府奖金的功劳，因为那些收到奖金和未收到奖金的学生存在多方面差异。然而，由于确认资格的收入分界线是未知的，在断点两侧小邻域内的个体可以被视为是相同的。因此，我们有理由认为未知的收入线外生随机地将断点附近的个体分成了两组，一组收到了政府奖金，一组未收到。

 

2. 模拟生成非线性相关数据

这里，我们假设被解释变量与收入的关系是非线性的 (线性相关性的举例和分析可以参见 Sharp Regression Discontinuity Example and Limitations)。Stata 随机生成一些非线性相关的自变量收入 income 和因变量学习表现 perfo 并绘制二者散点关系图。

clear
set obs 10000

gen income=3^((runiform()-0.75)*4)
label var income "Reported Income"

sum income
    Variable |        Obs        Mean    Std. Dev.       Min        Max
-------------+---------------------------------------------------------
      income |     10,000    .6789349    .7606786   .0370671   2.999232

gen perfo=ln(income)+sin((income-r(min))/r(max)*4*_pi)/3+3
label var perfo "Performance Index - Base"

scatter perfo income

下图展示了自变量与因变量的非线性关系：

现在，让我们加入一些随机扰动。

gen perf1=perfo+rnormal()*0.5
label var perf1 "Performance Index - with noise"

scatter perf1 income

接着，使用命令rcspline，可以将局部平均的学习表现视作收入的三次样条 (Cubic Spline) 函数。

ssc install rcspline

rcspline perf1 income, nknots(7) showknots title(Cubic Spline)

此时，样条曲线是平滑的。接下来，让我们在 0.5 处设置一个断点。

gen grant=income<0.5
sum grant

    Variable |        Obs        Mean    Std. Dev.       Min        Max
-------------+---------------------------------------------------------
       grant |     10,000       .5921     .491469          0          1

*样本中大约有59%低收入学生是具备获得政府奖金资格的
*现在加入政府奖金对学生表现的正向效应
*首先生成以政府奖金资格确认收入线为中心的收入变量
gen income_center=income-0.5
gen perf2=perf1+0.5*grant-0.1*income_center*grant
*这样政府奖金对低收入学生将更加有效
label var perf2 "Observed Performance"

 

3. 进行非参数估计

首先，让我们尝试使用传统的普通最小二乘回归。

reg perf2 income grant

      Source |       SS           df       MS      Number of obs   =    10,000
-------------+----------------------------------   F(2, 9997)      =   5734.77
       Model |   7041.8845         2  3520.94225   Prob > F        =    0.0000
    Residual |  6137.80314     9,997  .613964504   R-squared       =    0.5343
-------------+----------------------------------   Adj R-squared   =    0.5342
       Total |  13179.6876     9,999  1.31810057   Root MSE        =    .78356

------------------------------------------------------------------------------
       perf2 |      Coef.   Std. Err.      t    P>|t|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
      income |   .9003367   .0168801    53.34   0.000     .8672482    .9334251
       grant |  -.3767682   .0261265   -14.42   0.000    -.4279814    -.325555
       _cons |   1.924569   .0267009    72.08   0.000      1.87223    1.976909
------------------------------------------------------------------------------

显然，估计结果是错误的。政府奖金 grant 的估计系数为负，表明政府奖金会阻碍学习表现，这显然与常理相违背。

接下来，我们尝试使用 RD 命令 rdrobust。

*-默认断点在 0 点处，因此我们使用中心化后的变量 income_centered
rdrobust perf2 income_center

Sharp RD estimates using local polynomial regression.

      Cutoff c = 0 | Left of c  Right of c            Number of obs =      10000
-------------------+----------------------            BW type       =      mserd
     Number of obs |      5921        4079            Kernel        = Triangular
Eff. Number of obs |       683         530            VCE method    =         NN
    Order est. (p) |         1           1
    Order bias (q) |         2           2
       BW est. (h) |     0.129       0.129
       BW bias (b) |     0.197       0.197
         rho (h/b) |     0.652       0.652

Outcome: perf2. Running variable: income_center.
--------------------------------------------------------------------------------
            Method |   Coef.    Std. Err.    z     P>|z|    [95% Conf. Interval]
-------------------+------------------------------------------------------------
      Conventional | -.48486     .06467   -7.4971  0.000   -.611614     -.358102
            Robust |     -          -     -6.1641  0.000   -.633493     -.327828
--------------------------------------------------------------------------------

回归结果中估计系数为负，原因是 RD 方法通常默认断点处处置变量由 0 变为 1，与本案例中「收入高于收入线，获得政府奖金的资格取消」正好相反。因此，我们需要改变 RD 估计量的符号方向。这可以通过将收入取相反数来实现。

gen nincome_center=income_center*(-1)
rdrobust perf2 nincome_center

Sharp RD estimates using local polynomial regression.

      Cutoff c = 0 | Left of c  Right of c            Number of obs =      10000
-------------------+----------------------            BW type       =      mserd
     Number of obs |      4079        5921            Kernel        = Triangular
Eff. Number of obs |       530         683            VCE method    =         NN
    Order est. (p) |         1           1
    Order bias (q) |         2           2
       BW est. (h) |     0.129       0.129
       BW bias (b) |     0.197       0.197
         rho (h/b) |     0.652       0.652

Outcome: perf2. Running variable: nincome_center.
--------------------------------------------------------------------------------
            Method |   Coef.    Std. Err.    z     P>|z|    [95% Conf. Interval]
-------------------+------------------------------------------------------------
      Conventional |  .48486     .06467   7.4971   0.000    .358102      .611614
            Robust |     -          -     6.1641   0.000    .327828      .633493
--------------------------------------------------------------------------------

同时，我们可以将 rdrobust 新命令的回归结果与 Stata 传统 RD 回归命令 rd 的结果相比较。

rd perf2 nincome_center

Two variables specified; treatment is
assumed to jump from zero to one at Z=0.

 Assignment variable Z is nincome_center
 Treatment variable X_T unspecified
 Outcome variable y is perf2

Estimating for bandwidth .1832282339354582
Estimating for bandwidth .0916141169677291
Estimating for bandwidth .3664564678709164
------------------------------------------------------------------------------
       perf2 |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
       lwald |   .4860547   .0539727     9.01   0.000     .3802703    .5918392
     lwald50 |   .4929885   .0744274     6.62   0.000     .3471136    .6388635
    lwald200 |   .5737225   .0377565    15.20   0.000     .4997212    .6477238
------------------------------------------------------------------------------

rd 命令的好处在于它可以同时汇报出不同带宽下的估计结果。此处，带宽 (100 50 200) 分别代表最小化 MSE (mean squared error) 的最优带宽 (100%) 及其一半与两倍的带宽 (50%, 200%)。

我们可以将不同带宽的估计结果绘制在一张图上。

gen effect_est = .
label var effect_est "Estimated Effect"

gen band_scale = .
label var band_scale "Bandwidth as a Scale Factor of Bandwidth that Minimizes MSE"

forv i = 1/16 {
  rd perf2 nincome_center, mbw(100 `=`i'*25')
    if `i' ~= 4 replace effect_est = _b[lwald`=`i'*25'] if _n==`i'
    if `i' == 4 replace effect_est = _b[lwald] if _n==`i'
    replace band_scale = `=`i'*25'     if _n==`i'
}

gen true_effect = .5
label var true_effect "True effect"

two (scatter effect_est band_scale) (line true_effect band_scale)

从图上可以看出，最小化 MSE (即 100%) 带宽估计出的结果最为准确，且估计系数在其附近区间内也相对稳定。

 

4. 总结

本文介绍了 Stata 实现断点回归的最新命令 rdrobust。选取政府奖学金影响学生学业表现的案例，通过数值模拟随机生成观测数据，分别运用 rdrobust 和 rd 命令对模拟数据进行回归分析并比较回归结果。

 

5. 参考资料

  • RD Software Packages https://sites.google.com/site/rdpackages/.
  • Calonico, S., Cattaneo, M. D., Farrell, M. H., & Titiunik, R. (2017). Rdrobust: Software for Regression-discontinuity Designs. The Stata Journal: Promoting Communications on Statistics and Stata, 17(2), 372–404. -Link-, -PDF-.

 

6. 相关推文

Note：产生如下推文列表的 Stata 命令为：
  lianxh
安装最新版 lianxh 命令：
  ssc install lianxh, replace

资源共享

尊敬的老师 / 亲爱的同学们：
连享会致力于不断优化和丰富课程内容，以确保每位学员都能获得最有价值的学习体验。为了更精准地满足您的学习需求，我们诚挚地邀请您参与到我们的课程规划中来。
请您在下面的问卷中，分享您 感兴趣的学习主题或您希望深入了解的知识领域 。您的每一条建议都是我们宝贵的资源，将直接影响到我们课程的改进和创新。
我们期待您的反馈，因为您的参与和支持是我们不断前进的动力。感谢您抽出宝贵时间，与我们共同塑造更加精彩的学习旅程！https://www.wjx.cn/vm/YgPfdsJ.aspx# 再次感谢大家宝贵的意见！

关于我们

  • Stata连享会 由中山大学连玉君老师团队创办，定期分享实证分析经验。more……
  • 扫码加入连享会微信群，提问交流更方便

⚡ 最新课程推荐：
 

课程一：结构模型 · 2025 · 现场班

课程二：连享会 · 面板数据因果推断

FAQ：往期课程板书和答疑文档

选题征集，点击 此处 或者扫码填写问卷留下您的宝贵意见，感谢！

NEW！视频课堂：lianxh-class.cn

⚡ 新版 lianxh 命令：使用说明
 . ssc install lianxh, replace
 . lianxh 多期DID

⚡ ihelp：帮助文件 PDF 版
 . ssc install ihelp, replace
 . ihelp xtreg

lianxh_cn_saoma

⭕ 连玉君 · Stata 33 讲 - 免费
  o 点击观看，dofile 下载, --Book--

推荐专题