王鵬超 韓立彬



〔摘要〕多時點雙重差分法具有準自然試驗特征,可以相對干凈地識別因果效應,廣泛應用于與政策評估相關的研究中,但必須重視其可能存在的估計偏差問題。本文總結了多時點雙重差分法存在的問題和相應的解決措施。通過梳理最新文獻發現,多時點雙重差分法回歸系數識別的是組別—時間處理效應的加權平均,而非受處理個體的平均處理效應。在異質性處理效應下,多時點雙重差分法估計系數有偏,嚴重時估計系數符號會與真實系數符號相反。目前文獻上提出的解決措施可以歸結為一個診斷方法和三類解決方法。其中,診斷方法為Goodman-Bacon的系數分解定理,三類解決方法分別是加總方法、兩步回歸法和堆疊型雙重差分法。
〔關鍵詞〕雙重差分法(DID);多時點雙重差分法;異質性處理效應;組別—時間處理效應
中圖分類號:F064.1 ? ?文獻標識碼:A ? ?文章編號:1008-4096(2023)02-0027-13
一、問題的提出
計量經濟學可信性革命推動了實證經濟學進展,因果推斷成為實證經濟學研究的顯學。2021年,諾貝爾經濟學獎授予Card、Angrist和Imbens三位學者,表彰Card對勞動經濟學領域的實證貢獻,以及Angrist和Imbens對因果推斷方法的貢獻。這充分肯定了因果推斷方法在經濟學中的應用與發展。雙重差分法(Difference-In-Difference,DID)作為應用最廣的因果推斷方法之一,可以相對干凈地識別因果效應,在政策評估中受到國內外學者青睞。本文統計了2005—2020年使用DID方法的中文期刊文章數量。DID已成為國內實證研究者進行學術研究的重要工具。
根據政策實施時點的不同,DID一般可分為單時點DID(Staggered DID)和多時點DID(Multiple DID)。然而學界對多時點DID識別的系數含義與正確性卻較少討論。單時點DID識別的是受處理個體的平均處理效應(Average Treatment Effect on the Treated,ATT),多時點DID識別的是否也是受處理個體的平均處理效應?多時點DID估計系數是否有偏?現有文獻并未過多討論?!禔merican Economic Review》2020年第9期,Chaisemartin和D'Haultfoeuille[1]探討了多時點DID存在的問題,《Journal of Econometrics》在2021年第2期發布了“處理效應”專題,其中3篇文章與多時點DID識別直接相關。表明學術界對這一方法存在問題的高度關注。
最新研究發現,多時點DID估計系數識別的并不是受處理個體的平均處理效應,而是組別—時間處理效應的加權平均。當存在異質性處理效應時,估計系數有偏[1-5]。Goodman-Bacon[2]認為,多時點DID估計系數可分解為多個單時點DID系數的加權平均,權重與每個單時點DID的樣本份額和解釋變量方差相關,且都為正值。然而,部分單時點DID把早接受處理組作為晚接受處理組的對照組,在異質性處理效應下,這部分系數可能為負,從而總體估計系數會存在較大偏差。Chaisemartin和D'Haultfoeuille[1]、Borusyak和Jaravel[3]以及Borusyak等[4]認為,組別—時間處理效應為正,但部分權重為負,導致最終估計結果有偏。雖然事件研究法可以將不同處理時點轉化為處理時點一致的相對時點,但Sun和Abraham[6]證明,事件研究法設定中的每一相對時點系數不僅與該相對時點系數相關,還與回歸方程中其他相對時點系數及被剔除在方程之外的相對時點系數相關。在異質性處理效應下,利用相對時點系數大小檢驗平行趨勢假定是否滿足也會存在問題。
針對多時點DID存在的問題,學者們提出了不同的解決方法,本文將其歸結為一個診斷方法和三類解決方法。其中,診斷方法為Goodman-Bacon[2]的系數分解定理,該方法用于診斷估計系數的偏差程度。第一類解決方法為Sun和Abraham[6]、Callaway和Sant'Anna[7]提出的加總方法,即分別估計每個時期每個組別平均處理效應,再將其加總得到所有受處理個體的平均處理效應;第二類解決方法為Gardner[8]和Borusyak和Jaravel等[3]提出的兩步回歸法;第三類解決方法為堆疊型DID(stacked DID)[9],將每一政策時點前后一段時期內的處理組和干凈的對照組形成一個數據集,之后把所有的數據集堆疊并進行回歸。
二、DID方法的基本原理
為闡明DID的基本原理,考慮包含2個組別和2個時期的2×2 DID情形。組別包括一個處理組(treat_i=1)和一個對照組(treat_i=0),時期包括政策前(post_t=0)和政策后(post_t=1)。政策效果通過式(1)進行估計:
Y_it=β_0 treat_i×post_t+β_1 treat_i+β_2 post_t+ε_it (1)
其中,β_0為政策評估所關注系數,識別的是受處理個體的平均處理效應。這一系數可以用式(2)和式(3)加以分解說明:
β ?_0 ?=[E(Y_it (1)|treat_i=1,post_t=1)-E(Y_it (0)|treat_i=1,post_t=0)]-[E(Y_it (0)|treat_i=0,post_t=1)-E(Y_it (0)|treat_i=0,post_t=0)] (2)
(=[E(Y_it (1)|treat_i=1,post_t=1)-E(Y_it (0)|treat_i=1,post_t=1)]+[E(Y_it (0)|treat_i=1,post_t=1)@ ? -E(Y_it (0)|treat_i=1,post_t=0)]-[E(Y_it (0)|treat_i=0,post_t=1)-E(Y_it (0)|treat_i=0,post_t=0)]) (3)
其中,Y_it (1)和Y_it (0)是潛在結果。式(2)是處理組事后與事前均值的差異減去對照組事后與事前均值的差異,式(3)在式(2)上各加減一項E(Y_it (0)|treat_i=1,post_t=1)。式(3)中第一項表示事后所有受處理個體處理效應的均值,第二項是處理組事后假若未經處理的結果均值減去處理組事前的結果均值,第三項同理第二項。若滿足平行趨勢假定(第二項與第三項相減為0),得到β ?_0={E(Y_it (1)-Y_it (0)|treat_i=1,post_t=1)},β ?_0識別的是受處理個體的平均處理效應。
單個政策時點情形下,為靈活地估計回歸系數,常用個體固定效應λ_i和時間固定效應η_t代替式(1)中的treat_i和post_t。在多個政策時點情形下,若將D_it表示為處理狀態,處理組在政策后受到影響為1,否則為0,由于post_t與個體i和時間t同時相關,多時點DID的回歸方程如式(4)所示:
Y_it=β_0 D_it+λ_i+η_t+ε_it (4)
在多時點DID實際應用中,由于存在多個處理組,無法直接通過對比處理組和對照組結果均值的時間變化,以此檢驗是否滿足平行趨勢假定,因而事件研究法通常被當作替代方法。這一方法不僅可以檢驗事前是否滿足平行趨勢假定,還可以觀察事后政策效果的動態變化。事件研究法的模型設定為式(5):
Y_(i,t)=∑_(l=-k)^(-2)β_l^0 ?D_(i,t)^l+∑_(l=0)^Lβ_l^1 ?D_(i,t)^l+λ_i+η_t+ε_(i,t) (5)
其中,(-k,L)是相對時點l的范圍,D_(i,t)^l是每一相對時點l是否接受處理,接受處理為1,未接受處理為0。實證中常剔除-1期這一相對時點作為基準,每一相對時點系數都表示為相對-1期這一時點系數的大小。假若政策前每一相對時點系數β_l^0都無法拒絕系數為零的假設,則滿足平行趨勢假定,政策后每一相對時點系數β_l^1反映的是政策效果隨時間的變化。式(4)和式(5)的識別策略被廣泛用于漸進性試點和多期試點政策研究中。
三、多時點DID存在的問題
最新研究指出在異質性處理效應下多時點DID估計結果有偏。Baker等[10]利用蒙特卡羅方法,分析了不同處理效應下多時點DID估計系數存在偏差的六種情況,并分別進行了模擬。 圖1和圖2為模擬1—模擬3的結果。由圖1可知,模擬1和模擬2得到的估計系數圍繞真實系數呈正態分布,表明在單時點DID情況下,無論處理效應是否隨時間變化,估計系數都是無偏的。在存在多個處理時點,處理效應不隨時間和處理組組別變化時,模擬3得到的估計系數依然無偏。
圖3和圖4為模擬4—模擬6的結果。由圖3可知,在模擬4—模擬6中,多時點DID估計系數與真實系數存在偏差,偏差不斷增大,并且在模擬6中估計系數符號與真實系數符號相反。
圖3 模擬4—模擬6處理組和控制組結果均值的時間路徑
異質性處理效應既包括處理效應隨不同處理組組別發生變化,也包括同一個處理組在時間維度發生變化。無論處理效應是否隨時間變化,單時點DID都不存在估計系數有偏問題。然而,由于早接受處理組在回歸中被作為晚接受處理組的對照組,導致多時點DID估計系數有偏,尤其在異質性處理效應下,偏差會更大。不僅如此,估計系數識別的也不是受處理個體的平均處理效應,而是組別—時間處理效應的加權平均[1-5]。關于加權處理,Goodman-Bacon[2]對多時點DID系數的加權給出了直觀解釋。 考慮個體為N時期為T的平衡面板數據,假設存在2個政策時點:k和l,假設存在3個組別:早接受處理組k、晚接受處理組l和未受處理組U。
(0,k)時期為PRE(k),(k,l)時期為MID(k,l),(l,T)時期為POST(l)。在不考慮控制變量時,DID回歸方程為式(6):
y_it=β^DD D_it+λ_i+η_t+ε_it (6)
通過分解β ?^DD,證明該系數可表示為4個2×2 DID系數的加權平均[2],如式(7)所示:
β ?^DD=s_ku β ?_ku^(2×2)+s_lu β ?_lu^(2×2)+s_kl^k β ?_kl^(2×2,k)+s_kl^l β ?_kl^(2×2,l) (7)
其中,s_ku=((n_k+n_U )^2 V ?_kU^D)/V ?^D ,s_lu=((n_l+n_U )^2 V ?_lU^D)/V ?^D , s_kl^k=(((n_k+n_l)(1-D ?_l))^2 V ?_kl^(D,k))/V ?^D ,s_kl^l=(((n_k+n_l)D ?_k )^2 V ?_kl^(D,l))/V ?^D 表示每個系數的權重,n_(j,j∈{k,l,U})是組別j樣本量占總體樣本量比重,D ?_(i,i∈{k,l})是(i,T)該時期占總時期T的比重,V ?^D是總體方差,V ?_kU^D,V ?_lU^D,V ?_kl^(D,k),V ?_kl^(D,l)是每個2×2 DID中D_it的方差。每一權重都由每組樣本份額的平方和每組方差與總體方差之比兩部分構成,4個權重之和為1。β ?_ku^(2×2),β ?_lu^(2×2),β ?_kl^(2×2,k),β ?_kl^(2×2,l)分別是早接受處理組k和未受處理組U在(0,T)的系數、晚接受處理組l和未受處理組U在(0,T)的系數、早接受處理組k和晚接受處理組l在(0,l)的系數、晚接受處理組l和早接受處理組k在(k,T)的系數。 其中,β ?_ku^(2×2),β ?_lu^(2×2),β ?_kl^(2×2,k)含義是處理組事后與事前y均值的差異減去對照組事后與事前y均值的差異。β ?_kl^(2×2,l)是將早接受處理組k作為對照組,含義為晚接受處理組l事后和事前y均值的差異減去早接受處理組k作為對照組事后與事前y均值的差異。
由式(2)至(3)可知,β ?_ku^(2×2),β ?_lu^(2×2),β ?_kl^(2×2,k)分別是受處理個體平均處理效應+(處理組的時間趨勢-對照組的時間趨勢),在滿足平行趨勢假定條件下,系數識別的是受處理個體的平均處理效應。但β ?_kl^(2×2,l)不同,β ?_kl^(2×2,l)是受處理個體平均處理效應+(晚接受處理組l的時間趨勢-早接受處理組k作為對照組的時間趨勢)-早接受處理組k處理效應的時間趨勢。若早接受處理組的處理效應隨時間變化較大,β ?_kl^(2×2,l)可能為負值并最終影響總體系數β ?^DD的估計。 因此,當N趨于無窮時,式(7)可以進一步表示為式(8):
plim_(N→∞) β ?^DD=β^DD=VWATT+VWCT-ΔATT (8)
其中,VWATT(Variance-Weighted Average Treatment Effect on the Treated)是方差加權ATT,VWCT(Variance-Weighted Common Trends)是方差加權的平行趨勢,ΔATT是早接受處理組處理效應的時間趨勢。在滿足平行趨勢假定和處理效應恒定(即ΔATT=0)條件下,多時點DID的估計系數為VWATT。當存在異質性處理效應時,多時點DID估計系數就會產生偏差。
更具體地,若將每個2×2 DID表述為是每個組,當每個組組內處理效應恒定,組間各自的處理效應也相同時,如圖1模擬3所示,由于權重之和為1,VWATT就是ATT;當每個組組內處理效應恒定,組間各自的處理效應不同時,如圖3的模擬4所示,估計系數不再是ATT,而是表現為VWATT,權重與每個組的樣本份額和處理變量方差相關;當處理效應隨時間和組別變化時,如圖3的模擬5和模擬6所示,由于式(8)中還需減去ΔATT,因而估計系數會存在較大偏差。更為嚴重時,估計系數符號會與真實系數符號相反,無法準確識別政策效果。
學者們對多時點DID估計系數有偏的解釋存在差別。Goodman-Bacon[2]認為,所有2×2 DID系數的權重為正,但部分系數符號為負,導致多時點DID估計系數存在偏差。Chaisemartin和D'Haultfoeuille[1]、Borusyak和Jaravel[3]以及Borusyak等[4]卻認為,多時點DID系數是組別—時間處理效應的加權平均,組別—時間處理效應為正,但部分權重為負,導致最終估計結果有偏。本文以Chaisemartin和D'Haultfoeuille[1]的研究來解釋這類觀點??紤]城市g-年份t-企業i層面的數據結構。假定每個城市每年至少有一家企業,處理發生在城市層面,城市受處理后所有企業也同樣受到處理。令Y_(i,c,t) (1)和Y_(i,c,t) (0)是企業的潛在結果,D_(c,t)是城市的受處理狀態,D_(c,t)為1表示城市c在年份t受到處理,否則為0,D_(i,c,t)為企業受處理狀態。N_1=∑_(i,c,t)D_(i,c,t) 是所有受處理企業的觀測值數量,所有受處理企業的平均處理效應為:Δ^TR=1/N_1 ?∑_((i, c, t): D_(c, t)=1)[Y_(i,c,t) (1)-Y_(i,c,t) (0)] 。令δ^TR=E(Δ^TR),δ^TR是回歸估計所要識別的ATT。將N_(c,t)是城市—時間組cell(c,t)的觀測值數量,每個城市—時間組的平均處理效應定義為:Δ_(c,t)=1/N_(c,t) ?∑_(i=1)^(N_(c,t))[Y_(i,c,t) (1)-Y_(i,c,t) (0)] 。
δ^TR等價于所有受處理城市—時間組平均處理效應Δ_(c,t)的加權平均,權重為每個城市—時間組樣本與所有受處理企業樣本之比,即式(9):
Δ^TR=E[∑_((c, t): D_(c,t)=1)N_(c,t)/N_1 ?Δ_(c,t) ] (9)
但是在滿足平行趨勢假定的雙向固定效應模型下,實際得到的估計系數為式(10):
β_fe=E[∑_((c,t): D_(c,t)=1)N_(c,t)/N_1 ?w_(c,t) Δ_(c,t) ] ,
w_(c,t)=ε_(c,t)/(∑_((c, t): D_(c,t)=1)N_(c,t)/N_1 ?ε_(c,t) ) (10)
其中,w_(c,t)等式中的ε_(c,t)由D_(c,t)=α+γ_c+λ_t+ε_(c,t)得到。 N_(c,t)/N_1 ?w_(c,t)是每個受處理城市—時間組平均處理效應的權重,由城市—時間組的樣本份額,以及ε_(c,t)與ε_(c,t)均值之比兩部分構成。式(9)和式(10)表明,β_fe與Δ^TR并不相等,估計系數識別的不是ATT,而是組別—時間處理效應的加權平均。部分組別—時間處理效應權重為負時,多時點DID的估計系數有偏。
考慮2個組別和3個時點的DID,第1組在第3時點受到處理,第2組在第2時點受到處理。ε_(c,t)可由式ε_(c,t)=D_(c,t)-D ?_c-D ?_t+D ? ??得到,其中D ?_g是g城市企業受處理狀態變量的均值,D ?_t是t年企業受處理狀態變量的均值,D ? ??是所有企業在所有年份受處理狀態變量的均值。不同組別和不同時間的ε_(c,t)分別是ε_1,3=1/6,ε_2,2=1/3,ε_2,3=-1/6。ε_(c,t)均值為(ε_1,3+ε_2,2+ε_2,3)×1/3=1/9,w_(c,t)分別是w_1,3=3/2, w_2,2=3,w_2,3=-3/2,β_fe是β_fe=1/2 E(Δ_1,3)+E(Δ_2,2)-1/2 E(Δ_2,3),如圖5所示。
由圖5可知,第2組在第3時點處理效應的權重為負,導致多時點DID估計系數有偏。負權重的來源可以通過Goodman-Bacon[2]的分解定理進一步解釋。β_fe可分解為2個2×2 DID系數的加權平均,即β_fe=(DID_1+DID_2)/2。在滿足平行趨勢假定前提下,DID_1和DID_2分別為{(DID_1=[E(Y_2,2)-E(Y_2,1)]-[E(Y_1,2)-E(Y_1,1)]@DID_2=[E(Y_1,3)-E(Y_1,2)]-[E(Y_2,3)-E(Y_2,2)] )},DID_1識別的是ATT,因而DID_1=E(Δ_2,2),DID_2表示為DID_2=E(Δ_1,3)-[E(Δ_2,3)-E(Δ_2,2)],即ATT減去早接受處理組處理效應的時間趨勢,最終可得系數為β_fe=1/2 E(Δ_1,3)+E(Δ_2,2)-1/2 E(Δ_2,3),表明負的權重也是源于早接受處理組被作為晚接受處理組的對照組。在相同時間點上,受處理時間越長的城市,其處理效應的權重也越有可能為負[1]。當負權重非常大時,會導致估計結果產生大的偏差。
Sun和Abraham[6]研究表明,事件研究法可以解決圖3中模擬4和模擬5存在的問題,但無法解決模擬6存在的問題,原因在于每一相對時點上的處理效應不僅與自身時點處理效應相關,還與回歸中其他相對時點及被剔除在等式之外的其他相對時點處理效應相關??紤]個體為N,時期為T+1的樣本。Y_(i,t)為結果變量,Y_(i,e+l)和Y_(i,e+l)^∞分別為個體i在時點t的可觀測結果和未接受處理的反事實結果。D_(i,t)為二值變量,表示個體受處理狀態。將受處理個體i首次受到處理的時點表示為E_i=min{t |D_(i,t)=1},個體未接受處理的時點表示為E_i=∞,隊列(cohort)e為首次接受處理時點E_i相同的所有個體集合。相對時點為l,l=t-e。Sun和Abraham[6]將受處理隊列e在相對時點l的平均處理效應(Cohort-Specific Average Treatment Effect on the Treated,CATT)定義為式(11):
CATT_(e,l)=E[Y_(i,e+l)-Y_(i,e+l)^∞ |E_i=e] (11)
Sun和Abraham[6]基于這一定義得出了主要結論。假設T=3,存在3個隊列。 相對時點l,l∈{-3,-2,-1,0,1,2}?;貧w中常剔除第-1期作為基期,在該例中,由于沒有未受處理隊列,還需額外再剔除1期,最終剔除第-3期和第-1期。 以μ_(-2)表示第-2期系數,在滿足平行趨勢假定條件下,經證明μ_(-2)可分解為式(12)至式(14):
μ_(-2)=∑_(e∈{2,3})ω_(e,-2) ?CATT_(e,-2) (12)
+∑_(e∈{1,2,3})ω_(e,0) ?CATT_(e,0)+∑_(e∈{1,2})ω_(e,1) ?CATT_(e,1)+∑_(e∈{1})ω_(e,2) ?CATT_(e,2) (13)
+∑_(e∈{1,2,3})ω_(e,-1) ?CATT_(e,-1)+∑_(e∈{3})ω_(e,-3) ?CATT_(e,-3) (14)
其中,ω是CATT權重,式(12)是只和第-2期相關CATT的加權平均,所有權重之和為1;式(13)是回歸等式中除第-2期之外其他時期CATT的加權平均,所有權重之和為0;式(14)是被剔除在回歸等式之外其他時期CATT的加權平均,所有權重之和為-1。
若處理前不存在預期性行為,即個體不會在處理前受到影響,當l<0時,對于任意隊列e,CATT_(e,l)=0,式(12)和式(14)都為0。在滿足無預期性行為假設下,存在兩類情形:(1)假設不同隊列處理效應的時間趨勢相同,如圖3的模擬5所示,對于相對時點l的任意隊列e,CATT_(e,l)=CATT_l,由于權重之和為0,所以式(13)為0,最終第-2期系數為0;(2)假設不同隊列處理效應的時間趨勢不同,如圖3的模擬6所示,第-2期系數非0,表明這一相對時點的系數有偏。同理,分解事前其他相對時點和事后相對時點系數也會存在偏誤。在現實案例中,不同隊列更多表現為異質性處理效應。這種情況下利用相對時點系數大小判別是否滿足平行趨勢假定,以及檢驗處理效應的動態變化是存在問題的。
四、多時點DID估計有偏的解決方法
針對在異質性處理效應下多時點DID存在的問題,本文將解決方法總結為:診斷方法,也就是系數分解定理,三類解決方法分別為加總方法、兩步回歸法及堆疊型DID。
(一)系數分解定理
由Goodman-Bacon[2]的DID系數分解可知,當存在2個處理時點加一個對照組時,多時點DID的估計系數表示為4個2×2 DID系數的加權平均。當存在K個不同的政策時點加一個對照組時,多時點DID的估計系數可表示為K×K個系數的加權平均。通過觀察分解的不同系數大小和系數權重,即可診斷多時點DID存在的偏誤多大程度會影響最終估計結果。
關于系數分解定理。Stevenson和Wolfers[11]研究了1969—1985年美國部分州推行的無過錯離婚法案對婦女自殺率的影響,結果表明無過錯離婚法案實行后女性的自殺率有所降低。在不考慮任何控制變量情形下,多時點DID估計系數可以分解為156個2×2 DID系數的加權平均。
圖中所有三角和叉號表示的系數符號為負,且系數無偏,與整體估計系數的符號一致。菱形和圓圈表示的系數會令整體估計有偏,因為這部分樣本把早接受處理組作為晚接受處理組的對照組,并且所有菱形表示的系數符號為正,與整體估計系數的符號相反。菱形和圓圈兩部分所占比重較高,最終會嚴重低估法案實施后對女性自殺率的影響。因此,通過分解系數大小與系數所占權重可以診斷多時點DID的估計偏差。無過錯離婚法案對女性自殺率回歸系數的分解如圖6所示。
圖6 無過錯離婚法案對女性自殺率回歸系數的分解
Goodman-Bacon[2]的系數分解定理并非解決方法,回歸中加入控制變量后得到的圖形也無法區分早接受處理組和晚接受處理組和晚處理組和早處理組的分解,但由于實證文章中通常匯報不加控制變量的回歸結果,Baker等[10]認為此方法具有一般適用性,可以用于診斷只包括雙向固定效應的回歸結果偏差。
(二)加總方法
Sun和Abraham[6]對事件研究法處理多時點DID時存在的問題提出了相應解決方法,即求得每一相對時點下所有隊列的CATT_(e,l),用樣本份額對這一時點下所有CATT_(e,l)加權平均。Sun和Abraham[6]將其提出的估計量稱之為交互加權估計量(Interaction-Weighted Estimator,IW Estimator),分三步求得。第一步,利用雙向固定效應回歸估計CATT_(e,l),回歸等式為式(15):
Y_(i,t)=λ_i+η_t+∑_(e?C)∑_(l≠-1)δ_(e,l) ?(1{E_i=e}?D_(i,t)^l )+ε_(i,t) (15)
其中,C表示對照組,D_(i,t)^l表示是否為相對時點l,1{E_i=e}表示是否為隊列e。與每一相對時點l作對比的是相對時點-1,與每一組隊列e作對比的是對照組C.因此δ_(e,l)識別的是每一相對時點l下每一隊列e的平均處理效應CATT_(e,l)。 第二步,估計每個隊列e在相對時點l的權重Pr{E_i=e|E_i∈[-l,T-l]},該權重含義為:每個隊列e在相對時點l下的樣本占相對時點l里所有隊列樣本的比重。 第三步,每一相對時點的交互加權估計量表示為v ?_l=∑_eδ ?_(e,l) ?(Pr) ?{E_i=e|E_i∈[-l,T-l]}。這一方法可以得到每一相對時點l準確的估計系數。
類似地,Callaway和Sant'Anna[7]提出的解決方法依賴于定義的組別—時點平均處理效應(Group-time Average Treatment Effect),將組別—時點平均處理效應加總,即可得到總的處理效應。組別—時點平均處理效應定義為式(16):
ATT(g,t)=E[Y_t (g)-Y_t (0)|G_g=1] (16)
其中,g是個體首次受到處理的時點,G是首次接受處理時點為g的個體集合。G_g表示是否是首次接受處理時點為g的組別G。Y_t (g)是首次接受處理時點為g的組別在時點t的結果,Y_t (0)是這一組別在t這一時點假若未經處理的潛在結果。Callaway和Sant'Anna[7]考慮了兩類對照組,一類是將從未接受處理個體作為對照組,另一類是將從未接受處理個體以及尚未接受處理個體作為對照組。在滿足無預期性行為假設和平行趨勢假設條件下,若其他變量都不影響個體受處理狀態, 則可以利用兩個不同的對照組分別求得組別—時點平均處理效應,如式(17):
(ATT〖(g,t)〗_unc^nev=E[Y_t-Y_(g-1) |G_g=1]-E[Y_t-Y_(g-1) |C=1]@ATT〖(g,t)〗_unc^ny=E[Y_t-Y_(g-1) |G_g=1]-E[Y_t-Y_(g-1) |D_t=0] ) (17)
其中,nev表示將未接受處理個體作為對照組(C=1),ny表示將從未接受處理個體以及未接受處理個體作為對照組(D_t=0)。unc表示不考慮控制變量,g-1是首次接受處理時點g的前一時點。如表1所示,假設存在4個個體,6個時期,個體1和個體2未接受處理,個體3在t=3接受處理,個體4在t=4接受處理,表中數字表示個體的結果Y。根據上式可得:{ATT〖(3,4)〗_unc^nev=(8-4)-((5-3)+(4-2))/2=2@ATT〖(3,4)〗_unc^ny=(8-4)-((5-3)+(4-2)+(7-5))/3=2)}。
然而,個體受處理狀態通常是非隨機的,在考慮控制變量情形下,三種方法可以獲得ATT(g,t):結果回歸方法(Outcome Regression,OR)、逆概率加權方法(Inverse Probability Weighting,IPW)、雙重穩健方法(Doubly Robust,DR)。三種方法得到的ATT(g,t)在識別上相同,若需要進一步做統計推斷,雙重穩健方法會更加穩?。?]。
依靠求得的ATT(g,t),通過不同形式加總,可以得到總體的平均處理效應、處理效應的動態變化、不同處理組組間處理效應、累積的平均處理效應。 相應表達式如式(18)至式(21):
θ_es (e)=∑_(g∈G)1{g+e≤T}P(G=g|G+e≤T)ATT(g,g+e) (18)
式(18)為每一相對時點e的系數。其中,1{g+e≤T}是指示函數,G是所有個體首次接受處理的時點集合。以相對時點e=-2為例,θ_es (-2)表示為θ_es (-2)=1{1≤6}P{G=3|G-2≤6}ATT(3,3-2)+1{3≤6}P{G=5|G-2≤6}ATT(5,3-2)=P{G=3|G≤8}ATT(3,1)+P{G=5|G≤8}ATT(5,3)。因此,每一相對時點e的系數為所有不同組別G=g在相對時點ATT(g,t)的加權平均,權重為每個組別的樣本份額。
在[e,e^']數據結構平衡的相對時期里,每個相對時點e的系數為式(19):
θ_es^bal (e;e^')=∑_(g∈G)1{g+e^'≤T}P(G=g|G+e^'≤T)ATT(g,g+e) (19)
其中,θ_es^bal (e;e^')表示為在e到e^'范圍內相對時點e的系數。相對時點e在受處理時點附近范圍內,會包含所有受處理隊列。當相對時點e^'距離受處理時點較遠時,可能只包含部分受處理隊列。由于所含隊列數量不同,其系數大小難以解釋,因而可以保留包含所有受處理隊列的相對時期。
每個組別G=g┴?的處理效應為組別內所有受處理時點ATT(g┴?,t)的均值,如式(20)所示:
θ_sel (g┴?)=1/(T-g┴?+1) ∑_(t=g┴?)^TA TT(g┴?,t) (20)
總體系數為每個組別G=g┴?處理效應θ_sel (g┴?)的加權平均,權重為每個組別的樣本份額,如式(21)所示:
θ_sel^O=∑_(g∈G)θ_sel (g┴?)P(G=g|G≤T) (21)
上述兩種方法均是先獲得組別—時點平均處理效應,然后加總得到總的處理效應。兩者的不同之處在于:第一,Sun和Abraham[6]定義的是相對時點上不同隊列的平均處理效應,可以加總得到每一相對時點系數,但Callaway和Sant'Anna[7]定義的是正常時點上不同隊列的平均處理效應,將其加總可得相對時點的平均處理效應、不同組別的平均處理效應、總體的平均處理效應,適用更多情形;第二,Sun和Abraham[6]的解決方法沒有考慮控制變量對結果的影響,反觀Callaway和Sant'Anna[7]提出的獲得ATT(g,t)的三種方法都依賴于控制變量,因而更具有一般性。
(三)兩步回歸法
Gardner[8]開發了兩階段回歸法(Two-Stage Regression Approach),與Borusyak等[4]開發的方法類似,這是因為Borusyak等[4]認為Gardner[8]構造的估計量采取了“插補”形式(Imputation Form),且通過兩步可以獲得無偏的估計系數。本文將這兩種方法統稱為兩步回歸法。
Gardner[8]首先對多時點DID存在的問題給出直觀解釋,然后在此基礎上提出通過兩步法解決多時點DID存在的問題。考慮個體i和時間t層面的數據結構,將受處理時點相同的個體歸為g組,g∈{0,1,...,G},g=0是對照組,若g=2,是處理時點發生在t=2時的個體集合。g組受處理后的時期為p,p∈{0,1,...,P},p=0是受處理之前,若p=2,是組別g=2接受處理后的時期為2。Y_gpit、Y_1gpit和Y_0gpit分別是個體可觀測結果、處理組潛在結果和對照組潛在結果。D_(g,p)表示組別g在時期p是否接受處理。
當存在兩個時期和兩個組別(包括一個處理組和一個對照組)時,2×2 DID回歸等式為Y_gpit=λ_g+γ_p+β_gp D_gp+ε_gpit,系數β_gp識別的是受處理個體的平均處理效應,即β_gp=β_11=E(Y_1gpit-Y_0gpit |D_gp=1),回歸等式用均值形式表示為式(22):
E(Y_gpit |g,p,D_gp)=λ_g+γ_p+β_gp D_gp (22)
當存在多個處理組時,受處理個體的平均處理效應為E(β_gp |D_gp=1)=E(Y_1gpit-Y_0gpit |D_gp=1),即多個處理組平均處理效應的均值。多時點DID回歸等式用均值形式表示,如式(23)所示:
E(Y_gpit |g,p,D_gp)=λ_g+γ_p+E(β_gp |D_gp=1)D_gp+[β_gp-E(β_gp |D_gp=1)] D_gp (23)
其中,E(β_gp |D_gp=1)是需要識別的受處理個體的平均處理效應,β_gp是回歸得到的估計系數。由前文可知,β_gp估計的是多個處理組處理效應的加權平均,β_gp與E(β_gp |D_gp=1)并不相等。當僅有一個處理組時,兩項相減為0。當處理效應為同質時,β_gp等于受處理個體平均處理效應,兩項相減也為0。因此,最后一項可看作擾動項,其會影響多時點DID的最終估計結果。
Gardner[8]提出通過兩階段回歸法解決多時點DID存在的問題。第一階段,保留D_gp=0的樣本,用Y_gpit對λ_g和γ_p進行回歸;第二階段,在全樣本里將Y_gpit-λ ?_g-γ ?_t對D_gp回歸,最終估計系數識別的是受處理個體的平均處理效應。該方法的優點在于易于操作,并且第一階段回歸中也可加入控制變量。此外,該方法還可以拓展到事件研究法的應用中,在第二階段將D_gp變為如式(5)中的相對時點虛擬變量即可。雖然這一方法通過回歸得到無偏的估計系數,但在統計推斷時還需對標準誤進行調整,為此原文作者提供了Stata命令包did2s供實證研究者使用。
相比于Gardner[8]的兩階段回歸法,Borusyak等[4]提出的方法更為直觀。若處理組和對照組事前的差異可以完全由個體固定效應和時間固定效應捕捉,在控制雙向固定效應后,對照組可以看作是處理組的反事實狀態。政策后某一時點t上,處理組個體結果減去對照組個體結果,兩者的差值便是受處理個體i在時間t的處理效應,通過對這一處理效應加權平均,即可得到所有受處理個體的平均處理效應。類似地,若處理組和對照組事前的差異還受其他變量影響,通過控制這些變量也可得到所有受處理個體的平均處理效應。
該方法可以通過兩步得到多時點DID的無偏估計系數,本文以只考慮固定效應情形加以說明。第一步,保留沒有受到處理的樣本,回歸估計處理組個體未接受處理的潛在結果Y_it (0),Y_it (0)=λ_i+η_t+ε_it。第二步,對于處理組樣本,令Y ?_it (0)=λ ?_i+η ?_t,個體i在時間t的處理效應表示為τ ?_it=Y_it-Y ?_it (0)。最終的總體估計系數為所有受處理個體處理效應的加權平均,權重為樣本份額。
(四)堆疊型DID
Cengiz等[9]在研究最低工資對就業的影響時,使用了堆疊型DID。堆疊型DID是指在受處理時點前后(-j,k)范圍內,為受處理時點相同的隊列尋找干凈的對照組并形成一個數據集,數據集中包括這段時期內受處理個體的樣本、從未接受處理個體的樣本和尚未接受處理個體的樣本,之后將所有數據集堆疊并進行回歸?;貧w方程設定為式(24):
Y_mit=∑_(l=-j)^k〖β_l D_mit^l 〗+λ_mi+η_mt+ε_mit (24)
其中,m是數據集,i是個體固定效應,t是時間固定效應,D_mit^l表示是否為相對時點l,λ_mi為數據集—個體聯合固定效應,η_mt是數據集—時間聯合固定效應?;貧w系數β_l的大小可以觀測處理效應的動態變化。雖然這一方法較為直觀,但并沒有相應的理論證明。
總結而言,Goodman-Bacon[2]的分解定理應用于診斷多時點DID系數的偏差。本文比較了不考慮控制變量時的不同解決方法。 在這一數據生成設定下,這幾類方法分別得到的相對時點估計系數相近,與真實系數也較為接近,然而不同解決方法仍存在差異。雖然堆疊型DID在實證研究中得到了應用,但缺乏理論證明,且該方法需要手工合成數據,較為繁瑣。Sun和Abraham[6]提出的方法只在事件研究法下適用,也沒有考慮加入控制變量的情形。Callaway和Sant'Anna[7]、Gardner[8]、Borusyak等[4]提出的方法既可以估計多時點DID系數,也可以應用于事件研究法,回歸中還可以加入控制變量,因而更具有一般性。
五、結論與建議
多時點DID識別的不是受處理個體的平均處理效應,而是組別—時間處理效應的加權平均。特別是,如果存在異質性處理效應,多時點DID的估計系數與真實系數會存在偏差。依據文獻給出的解決思路,本文將解決方法歸納為一個診斷方法和三類解決方法。
隨著DID在經濟學實證研究中的廣泛應用,學者有必要了解多時點DID問題的應對方式,以保證實證中估計系數的一致性?;诖?,本文給出以下四點建議:第一,重視對政策背景和研究設計的討論,清晰地闡述不同組別受政策影響的可能方向;第二,在處理組較少時,可以通過觀察不同處理組和對照組的時間趨勢,以此檢驗異質性處理效應對結果的影響及平行趨勢假定是否滿足;第三,如果是面板數據結構,可利用系數分解定理評估系數偏差大?。坏谒?,從實踐角度出發,在利用多時點DID作為識別策略時,研究者應至少在以上三類解決方法中選擇其一,加強實證結果的可靠性和穩健性。
參考文獻:
[1] ?CHAISEMARTIN C D, D'HAULTFOEUILLE X. Two-way fixed effects estimators with heterogeneous treatment effects [J]. American economic review, 2020, 110(9): 2964-2996.
[2] ?GOODMAN-BACON A. Difference-in-differences with variation in treatment timing [J]. Journal of econometrics, 2021, 225(2): 254-277.
[3] ?BORUSYAK K, JARAVEL X. Revisiting event study designs [R]. SSRN working paper, 2017.
[4] ?BORUSYAK K, JARAVEL X, SPIESS J. Revisiting event study designs: robust and efficient estimation [R].CEPR
discussion papers, 17247.
[5] ?ATHEY S, IMBENS G W. Design-based analysis in difference-in-differences settings with staggered adoption [J]. Journal of econometrics, 2022, 226(1): 62-79.
[6] ?SUN L, ABRAHAM S. Estimating dynamic treatment effects in event studies with heterogeneous treatment effects [J]. Journal of econometrics, 2021, 225(2): 175-199.
[7] ?CALLAWAY B, SANT'ANNA P H C. Difference-in-differences with multiple time periods [J]. Journal of econometrics, 2021, 225(2): 200-230.
[8] ?GARDNER J. Two-stage differences in differences [R]. Working paper, 2021.
[9] ?CENGIZ D, DUBE A, LINDNER A, et al. The effect of minimum wages on low-wage jobs [J]. The quarterly journal of economics, 2021, 134(3): 1405-1454.
[10] ?BAKER A C, LARCKER D F, WANG C C Y. How much should we trust staggered difference-in-differences estimates? [J]. Journal of financial economics, 2022, 144(2): 370-395.
[11] ?STEVENSON B, WOLFERS J. Bargaining in the shadow of the law: divorce laws and family distress [J]. The quarterly journal of economics, 2006, 121(1): 267-288.
[12] ?SANT'ANNA P H C, ZHAO J. Doubly robust difference-in-differences estimators [J]. Journal of econometrics, 2020, 219(1): 101-122.
Potential Problems and Solutions of Staggered Difference-in-Difference Approach
WANG Peng-chao,HAN Li-bin
(Economic and Social Development Research Institute, Dongbei University of Finance and Economics,
Dalian 116025, China)
Summary:As one of the mainstream causal inference methods, difference-in-difference approach has the characteristics of quasi-natural experiment and can relatively exogenously identify the causal effect, so it is favored by scholars at home and abroad. The interpretation of treatment effect estimated by difference-in-difference approach with a single treatment period is well known, but there are only a few studies discussing the interpretation and accuracy of treatment effect estimated ?by staggered difference-in-difference approach. Recently, some latest studies discuss these problems in detail. This paper, ?by means of sorting such literature, summarizes ?the interpretation of treatment effect,potential problems, and corresponding ?solutions of staggered difference-in-difference approach.
The latest literature shows that the coefficient estimated by difference-in-difference approach with a single treatment period is unbiased regardless of the heterogeneity ?effect. ?Staggered difference-in-difference approach identifies ?the weighted average of ?different ?group-time treatment effects. ?The estimated ?coefficient ?is unbiased ?with ?homogeneity ?treatment effect but biased with heterogeneous treatment effect. Because some early-treated groups are taken as the control groups of the late-treatment groups, the estimated coefficients of this part are negative and finally result in the bias of the ?aggregated coefficient. In severe cases, the symbols of both estimated coefficient and real coefficient are opposite. According to the solutions given by the latest literature, the methods to solve the bias of estimated coefficient can be divided into ?'a ?diagnosis method' and 'three kinds of solutions'. The diagnosis method is Goodman-Bacon decomposition theorem. It is used to diagnose the degree of bias by estimating the sizes and weights of different group-time treatment effect. The first is the aggregation method including two ways, all of which are to find comparable control groups for each treatment group and estimate each group-time treatment effect. Then the unbiased estimated coefficient can be obtained by averaging all the group-period effects weighted ?by the sample share. The second ?is ?two-step regression method involving two ways. ?The reason ?for ?uniformly ?terming the two ways ?as ?this ?name lies ?in that they are similar in the solutions and ?the ?unbiased coefficient ?can ?be ?gained ?by ?two ?steps. ?The ?third ?is ?the ?stacked difference-in-difference approach. It aims to find comparable control group for cohorts with the same treatment period and form a data set. This dataset includes the samples of ?treated group, never-treated ?group, ?and not-yet-treated ?group. Then it is to stack all ?the data sets ?and regress ?an ?augmented difference-in-difference specification.
With the wide application of staggered difference-in-difference approach in empirical research of economics, it ?is necessary for empirical researchers to know how to deal with the problems of this method.
Key words: difference-in-difference; staggered difference-in-difference; heterogeneous treatment effect; group-time treatment effect
(責任編輯:李明齊)
收稿日期:2022-11-12
基金項目:國家自然科學基金青年項目“土地資源配置對人力資本空間分布的影響研究:理論、機制與對策”(72003020)
作者簡介:王鵬超(1996—),男,山西晉城人,博士研究生,主要從事區域和城市經濟研究。E-mail:pengchaowang1996@163.com
韓立彬(1988—),男,山東臨沂人,副教授,博士,主要從事區域和城市經濟研究。E-mail:hanlibin@dufe.edu.cn