/****************************************************************************
* Donohue and Wolfers (2006) 
* Figure 4
* Distribution of Regression-Estimated Effects Across States
****************************************************************************/


clear
set more off
set matsize 800
set mem 250m
cd "..."

use "StatePanel.dta", clear

* Death penalty Regression Code
xi: reg pc_mur legal i.st i.year [w=popul] if sample_ds==1, cluster(st)
xi: reg pc_mur legal rpc_inc ur ipolice nonwhite age15to19 age20to24 i.st i.decade [w=popul] if sample_ds==1, cluster(st)
xi: reg pc_mur legal rpc_inc ur ipolice nonwhite age15to19 age20to24 i.st i.year [w=popul] if sample_ds==1, cluster(st)
gen ex10=0
for num 1/10: replace ex10=ex10+lX.execs
replace ex10=1 if ex10>0
gen active=legal*ex10
gen passive=legal*(ex10==0)
xi: reg pc_mur active passive rpc_inc ur ipolice nonwhite age15to19 age20to24 i.st i.year [w=popul] if sample_ds==1, cluster(st)
gen adopt=1 if legal==1 & l.legal==0 & l.legal~=. & year>=1960
gen abolish=1 if legal==0 & l.legal==1 & l.legal~=. & year>=1960
sort st year
gen n=_n
egen ad=rank(n) if adopt==1
egen ab=rank(n) if abolish==1
sort stid year
for num 1/46: gen abX=1 if ab==X \ replace abX=l.abX if l.abX~=. \ replace abX=0 if abX==.
for num 1/41: gen adX=1 if ad==X \ replace adX=l.adX if l.adX~=. \ replace adX=0 if adX==.

xi: reg pc_mur ab1-ab46 ad1-ad41 rpc_inc ur ipolice nonwhite age15to19 age20to24 i.st i.year [w=popul] if sample_ds==1, cluster(st)
lincom (ab1+ab2+ab3+ab4+ab5+ab6+ab7+ab8+ab9+ab10+ab11+ab12+ab13+ab14+ab15+ab16+ab17+ab18+ab19+ab20+ab21+ab22+ab23+ab24+ab25+ab26+ab27+ab28+ab29+ab30+ab31+ab32+ab33+ab34+ab35+ab36+ab37+ab38+ab39+ab40+ab41+ab42+ab43+ab44+ab45+ab46)/46
lincom (ad1+ad2+ad3+ad4+ad5+ad6+ad7+ad8+ad9+ad10+ad11+ad12+ad13+ad14+ad15+ad16+ad17+ad18+ad19+ad20+ad21+ad22+ad23+ad24+ad25+ad26+ad27+ad28+ad29+ad30+ad31+ad32+ad33+ad34+ad35+ad36+ad37+ad38+ad39+ad40+ad41)/41

gen se=.
gen beta=.
for num 1/46: replace beta=_b[abX] if ab==X & _se[abX]>0 \ replace se=_se[abX] if ab==X
for num 1/41: replace beta=_b[adX] if ad==X & _se[adX]>0 \ replace se=_se[adX] if ad==X
gen ciu=beta+invttail(e(df_r), .025)*se
gen cil=beta-invttail(e(df_r), .025)*se
gen str20 ci="("+string(cil,"%9.1f")+", "+string(ciu,"%9.1f")+")"
tomode statename, by(st) replace
sort adopt statename 
list statename year beta ci adopt if abolish==1 | adopt==1, clean
gen wt=1/(se*se)
for any abolish adopt: reg beta if X==1 \ reg beta if X==1 [w=wt] \ reg beta if X==1 [w=popul] \ summ beta if beta>0 & X==1 \ summ beta if beta<0 & X==1


* Actual code to create Figure 4 begins here
kdensity beta if adopt==1 , xline(0, lcolor(black)) title(Death Penalty Reinstatement) xlabel(-6(2)6) xtitle("Estimated Effect on Homicide Rate" "Annual murders per 100,000 people") name(adopt, replace)
kdensity beta if abolish==1 , xline(0, lcolor(black)) title(Death Penalty Abolition) xlabel(-6(2)6) xtitle("Estimated Effect on Homicide Rate" "Annual murders per 100,000 people") name(abolish, replace)
graph combine adopt abolish, ycommon xcommon title("Distribution of Regression-Estimated Effects Across States") note("Kernel density estimates using Epanechnikov kernel") xsize(10) ysize(7.5)