***************************************************************
// Filename: test_sun_shapiro_cutbins.do
// Author: eddie yu
// Date: 2023-10-16
// Task: Sun and Shapiro Test
//
***************************************************************

clear all
snapshot erase _all
est clear

cd "$output/Work2024"

use "$data/Replication Data/figure6_data_A.dta", clear


* Create Decile bins of PR Exposure
local idx = 1
gen bin_pr_link = `idx' if pr_link_i == 0 
foreach j in 0.5 1 1.5 2 2.5 3 3.5 4 {
	di `j'
	local idx = `idx' + 1
	replace bin_pr_link = `idx' if pr_link_i > `j' - 0.5 & pr_link_i <= `j'
	summ pr_link_i if bin_pr_link == `idx'
}
replace bin_pr_link = `idx' + 1 if pr_link_i > 4


tab bin_pr_link, summ(pr_link_i) // 40% of sample has PR exposure = 0


// Regressions
keep if year == 1990 | year == 1995 | year == 2012
xi i.year|pr_link_i , noomit
drop _IyeaXp*1995

xi i.year*i.bin_pr_link, noomit
drop _IyeaXbin_*_1
drop _IyeaXbin_1995_*

local loopvalues = "2/10"

forval i = `loopvalues' {
	label var _IyeaXbin_2012_`i' "Bin `i'"
}


local loopvalues = "2/10"

forval i = `loopvalues' {
	label var _IyeaXbin_2012_`i' "Bin `i'"
}


forval i = `loopvalues' {
	eststo bin_`i':  reghdfe emp_growth _IyeaX* [aw=wgt] if bin_pr_link == 1 | bin_pr_link == `i', a(i.pid i.year#i.industr ) cl(fips_state indust ) nocons
	lincom _IyeaXbin_2012_`i'
		estadd local beta_pr = string(r(estimate), "%8.3fc")
		estadd local se_pr   = "(" + string(r(se), "%8.3fc") + ")"
		estadd local p_pr    = "[" + string(r(p), "%8.3fc") + "]"
		estadd local ci_pr   = "[" + string(r(lb), "%8.3fc") + "," + string(r(ub), "%8.3fc") + "]"
		estadd local lb_pr   = string(r(lb), "%8.3fc")
		estadd local ub_pr   = string(r(ub), "%8.3fc")
	summ pr_link_i if bin_pr_link == `i', d
		estadd local avg_pr_link = string(r(mean), "%8.3fc")
		estadd local beta_pr_scaled = string( _b[_IyeaXbin_2012_`i'] / r(mean), "%8.3fc")
		estadd local iqr_adj_pr_link = string(r(mean)/$iqr, "%8.3fc")
		estadd local beta_pr_scaled_iqr = string( _b[_IyeaXbin_2012_`i'] / (r(mean)/$iqr ), "%8.3fc")
}


esttab bin_* using sun_shapiro_pr_cutbins_longdiff.tex, drop(*) stats(beta_pr se_pr p_pr ci_pr avg_pr_link iqr_adj_pr_link beta_pr_scaled  beta_pr_scaled_iqr N, /// 
	labels("Long Diff Estimate" "s.e." "p-val" "95\% CI" ///
			  "Mean PR Exposure" "Mean PR Exposure (IQR adjusted)" "Estimate / Mean PR Exposure" "Estimate / IQR adjusted PR Exposure")) /// 
	mlabels(, nodepvars) ///
	nostar nonotes nonum nogap replace label nocons collabels(none) varwidth(40) 
	

	
	
********************************************************************************
* JUNK CODE BELOW
********************************************************************************

* create sample of regression to run DiD w/ TWFE
forval v = 2/10 {
	gen reg_sample`v' = 0
	replace reg_sample`v' = 1 if bin_pr_link == `v' | bin_pr_link == 1
	tab bin_pr_link if reg_sample`v' == 1 // check if we have totally unaffected & each group
}

* iqr
sum pr_link_i if base_emp > 0, d 
scalar iqr = r(p75)-r(p25)


// Graphs
xi i.year|pr_link_i , noomit
drop _IyeaXp*1995

xi bin_pr_link, noomit

forval i = 2/10 {
	preserve
// 	local i = 2
	keep if reg_sample`i'==1 
	tab bin_pr_link
	di `avg_pr_link'
	eststo emp_`i':  reghdfe emp_growth _IyeaX* [aw=wgt], a(i.pid i.year#i.industr) cl(fips_state indust ) nocons
	lincom _IyeaXpr__2012
		estadd local beta_pr = string(r(estimate), "%8.3fc")
		estadd local se_pr   = "(" + string(r(se), "%8.3fc") + ")"
		estadd local p_pr    = "[" + string(r(p), "%8.3fc") + "]"
		estadd local ci_pr   = "[" + string(r(lb), "%8.3fc") + "," + string(r(ub), "%8.3fc") + "]"
		estadd local lb_pr   = string(r(lb), "%8.3fc")
		estadd local ub_pr   = string(r(ub), "%8.3fc")
	summ pr_link_i if pr_link_i > 0, d
		estadd local avg_pr_link = string(r(mean), "%8.3fc")
		estadd local beta_pr_scaled = string( _b[_IyeaXpr__2012] / r(mean), "%8.3fc")
		estadd scalar iqr_adj_pr_link = string(r(mean)/iqr, "%8.3fc")
		estadd local beta_pr_scaled_iqr = string( _b[_IyeaXpr__2012] / iqr_adj_pr_link, "%8.3fc")
// 	est save ster/emp_reg_`i', replace
	restore
}



esttab emp_* using sun_shapiro_pr_cutbins_longdiff.tex, drop(*) stats(beta_pr se_pr p_pr ci_pr avg_pr_link beta_pr_scaled iqr_adj_pr_link beta_pr_scaled_iqr N, /// 
	labels("Long Diff Estimate" "s.e." "p-val" "95\% CI" ///
			  "Mean PR Exposure" "Estimate / Mean PR Exposure" "Mean PR Exposure (IQR adjusted)" "Estimate / IQR adjusted PR Exposure")) /// 
	mtitles("Bin2" "Bin3" "Bin4" "Bin5" "Bin6" "Bin7" "Bin8" "Bin9" "Bin10") ///
	nostar nonotes nonum nogap replace label nocons collabels(none) varwidth(40)