断点回归和读者的提问解答
本文包括两部分：政策评估方法里的断点回归设计（regression discontinuity design），附加了部分倾向匹配分析方法，和读者3个提问的解答（文章后面）。
断点回归是一种准实验设计。如果政策在一个关于个人背景的连续的变量(例如考试成绩、家庭人均收入等)上设定一个临界值(Cutoff/Threshold)，使得在临界值一侧的个体接受政策干预，而在临界值另一侧的个体不接受干预，则在临界值附近就构成了一个准实验。我们把这个决定了是否接受干预的连续变量叫做强制变量(Forcing Variable)，由于强制变量是连续的，所以在临界值两侧的个体应该是类似的、可比的，则这两侧的个体在产出上的差异就应该是干预造成的差异。
当个体是否接受政策干预由强制变量值与临界值之间的关系决定时，我们可以用如下数学表达式：
现在假设设立了奖学金，且只有成绩高于一个临界点的学生才能获得，则获得这个奖学金对上大学概率的影响可以用公式(10)来表达：
应用断点回归的一个经典研究是Lemieux&Milligan(2008)(17)。他们研究社会救助会不会影响就业率。劳动力经济学家根据理论推测，增加社会救济会减少接受救济的人群工作的必要性，从而减少劳动力供给、降低就业率。
Lemieux&Milligan(2008)研究的这个社会救助项目规定30岁以下的人只能获得185美元，而一旦超过30岁，就可以获得507美元，这是一个巨大的差额。因此年龄就是这个政策的强制变量，临界点是30岁。图4展示了1986年人口普查时30岁以下和30以上的人群获得社会救助的额度。可以看到在临界点两侧，人们的救助收入有一个飞跃。因此，实际情况完全符合政策设计。
图5展示了人口普查当天在临界点30岁附近，就业率的情况。可以看到，在30岁附近的就业率确实有一个跳跃。这就是在30岁时大幅增加社会救助的干预效应——降低就业率。
断点回归的分析必须经受多方面的假设检验，才能保证其内部有效性。比如，我们应该检验临界值两侧的两群人是否在一些重要的背景变量上分布连续，强制变量本身在临界值附近是否连续。此外，由于断点回归是在临界值附近构建准实验，因此，我们要界定“附近”的区域到底有多宽。如果这个区域太窄，则样本量很小，可能影响估计的精确度和统计推断力。如果这个区域太宽，则我必须控制其他因素，以保证两边人群的可比性。顺便指出，因为断点回归是在一个临界值附近估算干预效应，而不是在整个定义域内估计平均的干预效应(Average Treatment Effect，ATE)，因此这个估计值是局部平均干预效应(Local Average Treatment Effect, LATE)。
断点回归包括参数回归和非参数回归。参数回归主要强调对临界点附近的函数类型做出正确的假设。如图6(18)所示，第一张图是线性模型，第二张图是非线性模型，在临界点附近都有一个跳跃，因此可以用断点回归来识别这个干预效应。但第三张图实际上是一个在临界点附近没有跳跃的非线性模型，如果把这个非线性方程错当做线性方程对待，而且使用断点回归模型去估计，也能估算出一个跳跃，但这个估计结果是一个严重的错误。因此，在断点回归中，观察散点图的形状非常重要，而且我们往往把方程(11)写成一个更笼统的函数形式：
非参数回归则不对临界点附近的函数形式做出假设。最简单的非参数回归就是选择一个临界点附近的区域，分别计算临界点之上和之下的观测点的平均产出，然后比较均值。这种方法不够精确，且忽略了强制变量本身对产出的影响。其改进方法是Kernel方法，即在选定的区域中对观测点进行加权平均，权重大小由不同的分布函数来表示。第三种非参数回归则是局部线性回归，即在临界点左侧的局部区域和临界点右侧的局部区域分别进行线性回归，根据回归系数计算两个方程在临界点的估算值，估算值之差即为干预效应。这种方法的基础是随着强制变量从两边逼近临界值，临界值两边的方程形式已经不重要，而且可以用线性函数来近似。
图6 清晰的断点回归模型设计
那么，如何定义和寻找模拟的对照组成员呢？一般来讲我们都希望模拟对照组与实验组在一组变量上的统计值都无显著差异。这是比较难做到的。因为，当两组人在第一个变量上无显著差异时，可能在第二个变量上显著不同；当调整模拟对照组成员使得两组人在第二变量上无显著差异时，可能又在第三个变量上显著不同。因此，我们主要通过分析样本中的每个人选择接受干预的概率来进行配对。个体选择接受还是不接受干预，可以通过Probit模型或者Logit模型来估算，我们所关心的主要特征变量都应该是影响接受干预概率的变量。
这个配对的理念虽然比较简单，但具体的配对方法却比较复杂。例如，用于计算倾向分数的选择模型(14)是否有效？在配对时，如果一个没有接受干预的个体已经被配给某个接受干预的个体作为其对照，那么这个个体能否再配给另一个接受干预的个体？如果一个接受干预的个体在未接受干预的群体中找到若干个倾向分数相近的个体，是否都要把他们作为对照、或者只是随机选一个？如果都选为对照，给每个人赋予的权重是多少？Heckman(23)、Dehejia(24)等对倾向分数配对法都做了详细而深入的讨论。Altonji等(25)应用倾向分数配对法评估了天主教学校是否比公立学校更能提高学生的学业表现。读者可以参考。
断点回归设计示例和运行命令
断点回归由Thistlewaite and Campbell(1960)首次使用，但直到1990年代末才引起经济学家的重视。Hahn et al(2001)提供了断点回归在计量经济学理论基础。目前，断点回归在教育经济学、劳动经济学、健康经济学、政治经济学以及区域经济学的应用仍方兴未艾。参见Imbens and Lemieux(2008)，Van Der Klaauw(2008)以及Lee and Lemieux(2010)的文献综述。
例1： Thistlewaite and Campbell(1960)使用断点回归研究奖学金对于未来学业成就的影响。由于奖学金由学习成绩决定，故成绩刚好达到获奖标准与差一点达到的学生具有可比性。
例2： Angrist and Lavy(1999)在研究班级规模对成绩的影响时，利用以色列教育系统的一项制度进行断点回归；该制度限定班级规模的上限为40名学生，一旦超过40名学生（比如41名学生），则该班级被一分为二。
断点回归可以分为两种类型：
精确断点回归(Sharp Regression Discontinuity，SRD)，其特征是断点x=c处，个体得到处理的概率从0跳跃到1。
模糊断点回归(Fuzzy Regression Discontinuity，FRD)，其特征是断点x=c处，个体得到处理的概率从a跳跃到b，其中0<a<b<1。
1.精确断点回归
考察上大学对工资收入的影响。假设上大学与否（Di）完全取决于由高考成绩xi是否超过500分：
由于此回归存在一个断点，故称为“断点回归”(Regression Discontinuity)或称“断点回归设计”(Regression Discontinuity Design)。由于在断点附近仿佛存在随机分组，故一般认为断点回归是内部有效性比较强的一种准实验。
断点回归可视为“局部随机试验”(Local randomized experiment)；可通过考察协变量在断点两侧的分布是否有差异来检验随机性。但断点回归仅推断在断点处的因果关系，不能推广到其他样本值，故外部有效性受局限。
精确断点回归:
考察带宽的选择:
由于断点回归在操作上存在不同的选择，实践中一般建议同时汇报以下各种情况，以保证稳健性。
Ø分别汇报三角核和矩形核的局部线性回归结果（后者等价于线性参数回归）。
Ø分别汇报使用不同带宽的结果（比如，最优带宽及其二分之一或两倍带宽）。
Ø分别汇报包含协变量与不包含协变量的情形。
Ø进行模型设定检验，包括检验分组变量与协变量的条件密度是否在断点处连续。
2.模糊断点回归
例3: 高考成绩上线并不能完全保证上大学，能否上大学还取决于填报志愿，甚至有些上线考生放弃上大学 的机会；而即使成绩未上线，但也可能因某种特长而得到加分，从而得到上大学的机会。上大学的概率确实在分数线的位置上有一个不连续的跳跃。
在模糊断点的情况下，处理变量D不完全由分组变量x所决定。
断点回归的Stata实例
断点回归可通过非官方Stata命令rd来实现，基本句式：
rd D x , z0(real) strineq mbw(numlist) graph bdep oxline kernel(rectangle) cov(varlist) /// x(varlist)
其中，“y”为结果变量，“D”为处理变量，而“x”为分组变量。
选择项“z0(real)”用来指定断点位置，默认值为“z0(0)”，即断点为原点。
选择项“mbw(numlist)”用来指定带宽，默认值为“mbw(50 100 200)”。
选择项“graph”表示根据所选的每一 带宽，画出其局部现行回归图。
选择项“bdep”表示通过画图来考察断点回归估计量对带宽 的依赖性。
选择项“oxline”表示在此图的默认带宽上画一条直线，以便识别。
选择项“Kernel(rectangle)”表示使用矩形核，默认使用三角核。
选择项“cov(varlist)”用来指定加入局部线性回归的协变量。
选择项“x(varlist)”表示检验这些协变量是否在断点处有跳跃（估计其跳跃值和显著性）:
例4：考察美国国会选区（如果有）一名民主党众议员对该选区联邦支出的影响。
传统上，民主党倾向于大政府，故一个选区如果有民主党众议员，则该议员可能为该选区争取更多的联邦支出。然而，直接对二者进行回归可能存在遗漏变量问题或双向因果问题。为此，使用民主党候选人的得票率作为分组变量，以0.5作为断点（在两党政治中，得票率大于或等于0.5则当选，反之落选），进行断点回归。数据集votex.dta，其中结果变量为lne（选区联邦开支的对数）、分组变量为d（民主党候选人得票率减去0.5）、处理变量win（民主党候选人当选），以及一系列协变量。
先来做一个一般的lne对win的OLS回归，可以与断点回归作比较。
Stata命令：reg lne win i votpop black blucllr famer fedwrkr forborn ///
manufunemplyd union urban verterans
首先，使用最优带宽及默认的三角核进行精确断点回归，并画图。
Stata命令：rd lne d, gr mbw(100)
Stata命令：rd lne d, mbw(100) cov(i votpop black blucllr famer fedwrkr forborn manuf unemplyd union urban verterans)
再次，去掉协变量，同时估计三种带宽，并画出估计值对带宽的依赖性。
Stata命令：rd lne d , gr bdep oxline
进行断点回归后，还需要对其设定进行检验。先检验协变量在断点处的条件密度是否存在跳跃。
Stata命令：rd lne d , mbw(100) x(i votpop black blucllr farmer fedwrkr forborn manuf unemplyd union urban veterans)
McCrary检验可以通过非官方Stata命令Dcdensity来实现，ado文件下载地址为http://emlab.berkeley.edu/~jmccrary/DCdensity/。将“Dcdensity.ado”文件下载到文件夹“\\ado\\plus”即可。
命令Dcdensity基本句式为:
Dcdensity assign_var,breakpoint(#) generate(Xj Yj r0 fhat se_fhat) graphname(filename)
其中“assign_var”为分组变量，必选项“breakpoint(#)”指定断点位置，必选项“generate(Xj Yj r0 fhat se_fhat)”用来指定输出变量名，而选择项“graphname(filename)”用来指定密度函数图的文件名。
运用Dcdensity命令检验数据集(votex.dta)的分组变量d是否连续。
Stata命令：DCdensityd, breakpoint(0) generate(Xj Yj r0 fhat se_fhat) /// graphname(rd.eps)
模糊断点回归:
根据模糊断点回归的定义，随机生成一个新的处理变量ranwin。变量ranwin不由分组变量d完全决定，但与原来的处理变量win高度相关(win完全由d决定)。
Stata命令如下：
Set seed 10101
gen ranwin = cond(uniform () < .1 , 1 – win , win)
其中，“uniform()”表示生成一个在[0 , 1)区间服从均匀分布的随机变量，命令“cond(uniform () < .1 , 1 – win , win)”表示，如果此随机变量小于0.1，则变量ranwin= 1 – win；否则，变量ranwin= win。因此，ranwin也是虚拟变量，且与win高度相关。
使用最有带宽与默认的三角核进行模糊断点回归。
Stata命令：rd lne ranwin d , mbw(100)
模拟断点回归：
假设断点c=140
结果变量 y = 100+80*T+2*x+rnormal(0,20)
其中，分组变量x是服从均值为100、方差为50的正态分布(小于0的用0替代，删除大于280的数)。处理变量T是x的函数。当x>140时，T=1，否则T=0。
没有处理变量T的情况，y对x是一条没有断点的直线:
有高次项的情况：
例5：
数据集：LMB Data.dta
分组变量：民主党的得票率(demvoteshare)。当得票率大于50%时，则民主党获胜。
结果变量：来自ADA的自由投票得分。
考察在得票率相近的选举中获胜的候选人是否较少了他们的国会得票。
画散点图：
scatter score demvoteshare, msize(small) xline(0.5)
xtitle("Democrat vote share") ytitle("ADA score")
scatter score demvoteshare, msize(small) xline(0.5) xtitle("Democrat vote share")
ytitle("ADA score") jitter(2)
在断点两端添加拟合曲线：
加入高次项再回归：
依次添加高次项的综合图：
添加高次项以及他们的交互项（democrat与demvoteshare）:
用cmogram命令来画断点为0.5的断点回归图：
下列第一个图：断点为0.5的直线回归（flit）
下列第二个图：断点会0.5的局部光滑线性回归(lowess)
利用核函数做局部多项回归：
在断点0.5附近选择不同的带宽：
做精确断点回归：
论文：Evidence on the impact of sustained exposure to air pollution on life expectancy from China’s HuaiRiver policy 这篇文章由清华大学的李宏彬、北京大学的陈玉宇和另外两位作者共同完成。发表在《美国国家科学院院刊》（PNAS）的一篇论文，向已经逐渐意识到空气污染危害的中国民众再次展现了残酷的现实：以淮河为界，烧煤供暖的中国北方地区空气污染水平高于中国南方，北方5亿居民因严重的空气污染，平均每人失去5年寿命。使用的方法就是断点回归（RD）
而空气污染造成的结果，淮河两岸的人均寿命差异则是惊人的5年。
同时，不考虑空气污染，用其他变量预测出来的人均寿命则几乎没有差别。这说明了连续性假设是成立的，淮河两岸的两个十分接近的地区几乎拥有相同的特征。作者还发现，淮河以北的超额死亡率，正是由和空气污染相关的心肺疾病相关，并进而与空气总悬浮物相关，而与谋杀、自杀等因素无关。
作者给出结论：烧煤将使空气总悬浮物大幅度上升，并使中国北方地区的居民人均寿命降低5年。该文亦引起一些争论：
报告“科学性”的争论
亦有学者对该研究报告的科学性提出了质疑。
北京大学公共卫生学院教授潘小川指出，只有PM10以下的颗粒物能进入到人的呼吸道系统，对人体健康造成危害。而总悬浮颗粒物是很粗糙的颗粒物，其中有很大一部分是PM10到PM100之间的颗粒物的浓度，总悬浮颗粒物跟健康生物学联系比较弱。
潘小川认为，该研究报告得出的结果的生物学依据较差，是比较粗糙的生态学研究。同时，他指出，影响寿命的因素很多，其中最重要的因素是食物和营养因素，报告中这方面的数据是缺失的。
对此，陈玉宇解释称，研究之所以采用总悬浮颗粒物，是因为中国测控PM10、PM2.5是近十年来才开始的，上世纪八九十年代的中国缺乏这些颗粒物的测控数据。同样，营养方面的数据不充分，因此研究也没有采用。
但陈玉宇认为这并不影响研究的科学性，他们的研究在统计学上证明了空气污染对人体健康造成了伤害，至于具体伤害到什么部位，该项研究并不需要论证。
这样的质疑也正说明了一个问题，研究的困难在于获得更好的数据。
“中国很多数据零星地分散在各个不同的机构、管理者和政府部门手中，没有被充分利用起来。”陈玉宇说，该项研究在方法上借鉴了经济学方法，但还需更全面、更高质量的数据，这样研究者才能得出更好的估计和判断。
“研究并非针对供暖政策”
中国长期以来的供暖政策受到了不少环境学者的诟病，但陈玉宇指出，此项研究并非评价供暖政策本身。
“供暖是造成污染水平发生跳跃的一个因素。”陈玉宇表示，研究团队想要体现的是，目前的污染水平已带来相当大的健康成本。
李宏彬亦表示，供暖是造成空气污染的一个原因，其实背后还是燃煤效率问题。研究的目的在于给政府一个依据，加大环境治理的投入。比如，根据这项研究估计的空气污染对寿命的影响，一个简单的计算表明，每年整个社会愿意至少支付GDP的两个百分点的支出，会使得空气颗粒悬浮物的浓度降低100微克/立方米。
陈玉宇表示，治理空气污染不仅是中国面临的问题，也是世界上诸多国家需要解决的问题。
他指出，当前对空气和环境问题的研究应回答三大类问题。一，空气污染或环境污染带来的成本是什么？二，空气污染和环境污染产生的根源是什么？三，在前述两个问题解决的基础上，应通过什么样的办法和政策以最低成本去最有效率地治理好空气？
值得注意的是，这三大类问题在中国的研究还远远不够。陈玉宇希望借此次研究报告带动更多的研究者以更好的研究方法，收集更多更高质量的数据，更准确地评估出空气污染对健康的影响。
陈玉宇说，大气中的PM2.5、PM10和TSP(总悬浮颗粒物)大多来自化石燃料燃烧后排放的气体，政府应设计科学的用煤政策，提高燃煤效率，或通过税收方式抑制污染。
在因果关系分析实证方法中，最优选择为随机实验，但是随机实验的时间成本和经济成本都比较高，而在随机实验不可得的情况下，需要考虑使用其它方法。断点回归（Regression Discontinuity）是仅次于随机实验的，能够有效利用现实约束条件分析变量之间因果关系的实证方法之一。但断点回归仅推断断点处的因果关系，并不能推广到其他样本值，故外部有效性受局限。Abadie等（2009）提出合成控制法（Synthetic Control Methods，SCM）来构造处理组的合理对比对象，克服处理组和对照组之间的差异，还可克服RDD不能外推的缺陷。SCM、RDD、DID（双重差分）、DDD（三重差分）等是目前实证分析的流行工具。SCM与PSM（倾向得分匹配）在构造反事实状态（counter factual state）时各有特色，而RDD在政策分析时也具独特功能。
读者的3个提问：
问题1： 非平行趋势的检验与处理方法
方法一、画时间趋势图
如果在政策干预前有多期数据，则可分别画处理组与控制组的时间趋势图（类似于上图），并直观判断这两组的时间趋势是否平行（比如，考察是否存在Ashenfelter's dip）。如果二者大致平行，则可增强对平行趋势假定的信心。然而，即使在政策干预前两组的时间趋势相同，也无法保证二者在干预后的时间趋势也相同（后者本质上不可观测，因为时间效应已与处理效应混合在一起）。另外，如果只有两期数据，则无法使用此法。
方法二、加入更多的控制变量
从上文的讨论可知，非平行趋势可能由于遗漏变量所导致，故在 zit 中加入更多控制变量，或可缓解内生性。但此法在实践中不易实施。
方法三、假设线性时间趋势
如果假设时间趋势为线性函数，则可加入每位个体的时间趋势项：
在具体回归时，加入个体虚拟变量与时间趋势项 t = 1, 2, ... , T 的交互项即可。然而，线性时间趋势毕竟是较强的假定，不一定能成立。故此法也不完全解决问题，但可作为稳健性检验。
方法四、三重差分法
在一定条件下，可通过引入两个控制组，进行三次差分，称为“三重差分法”（difference-in-differences-in-differences，简记DDD），这样可以更好地控制时间趋势的差异，使得平行趋势假定更易成立。
问题2：关于修正的R平方
问题3：三重差分法及运行
双重差分法的关键假设是实验组与对照组的时间效应一样。这个假设只有通过足够长的时间序列数据才能检验。需要指出的一点是，即使干预发生之前两组时间序列一致，也不能保证干预发生后两组时间序列是一致的。有可能在干预发生的同时在实验组或者对照组中又发生了其他影响产出的事件，则干预发生后两组的时间趋势是不一致的。简单的双重差分估计是有偏的。如图3所示，如果对照组是虚线所示，则双重差分估计是无偏的。但如果对照组是上方的实线，则双重差分估计法是有偏的，偏差部分是在
解决这个问题有两个思路：第一个是寻找更多的对照组，把多个对照组加权构造成一个虚拟的对照组，使得虽然每个对照组都与实验组的时间趋势不一样，但加权后的虚拟对照组的时间趋势与实验组的一样。这个方法被称作综合控制法(Synthetic Control Method)。Abadie & Gardeazabal(2003)(14)用这个方法研究了恐怖冲突对经济发展的影响。
解决这个问题的第二个思路是估算出这个因为时间趋势不同而带来的偏差，然后从双重差分结果中减去这个偏差即可。这被称作三重差分法(Difference-in-differences-in-differences, DDD)。三重差分法的思路是，既然两个地区(分别指实验组和对照组)的时间趋势不一样，那么我们可以分别在两个地区寻找一个没有受到干预影响的人群/行业，通过对这两组的双重差分估算出时间趋势的差异，然后再从原来实验组和对照组的双重差分估算值中减去这个时间趋势差异。Gruber(15)就使用了这种方法。
图3 时间趋势差异造成的估计偏差
三重差分法例证：
三重差分法回归及运行命令
正如上面的二重差分法实际上运用的是OLS做的回归，我们之前说过，倍差法是相当于two-way fixed effect model，里面包括个体效应时间效应，而对于一个这样的panel data，我们可以运用LSDV通过添加个体和时间虚拟变量来回归，或者运用demeaned variables回归来消除个体和时间效应，再加上那些交互效应后就可以像其他fixed effect回归一样。二重差分法一般是在同一个省（地区）区分treatment和control组的，而三重差分法则包括另一个未受到政策冲击的省（地区），来区分treatment和control组的，当然三重差分法要稳健得多。
Empirical Methods in Applied Economics Lecture
Jorn-Steffen Pischke
http://econ.lse.ac.uk/staff/spischke/ec524/evaluation3.pdf
@计量经济学圈
记录一个我们生活在其中的时代社会，一个非常具有潜力的深度与客观兼具的大号，囊括的主题如下：经济、社会、历史、新闻、世界、计量工具。