clear all
set more off, permanently
set matsize 11000, permanently

*set directory
global dir = 

*set ado path
adopath ++


capture: log close
log using "$dir/Extrapolation/log.log", replace


/********************/
/* Data preparation */
/********************/


/* Load data */

use "$dir/Data/analysis_data.dta", clear
desc

/* Generate variables */

foreach x in base_eng base_math {
	gen `x'2 = `x'^2
	gen `x'3 = `x'^3
}

foreach x in eng4 math4 eng7 math8 {
	gen `x'_2 = `x'^2
	gen `x'_3 = `x'^3
}

gen eng4_math4 = eng4*math4
gen eng7_math8 = eng7*math8

gen all = 1
gen male = (female == 0)

egen isee_mq_score = rowmean(math_score quan_score)
egen isee_rv_score = rowmean(verb_score read_score)

/* Sample restrictions */

tab private, miss
keep if private != 1

tab yrs_in_exam_cur, miss
keep if yrs_in_exam_cur == 0 | yrs_in_exam_cur == .

keep if app_grade != .

keep if (app_grade == 7 & schyear >= 1999 & schyear <= 2005) | (app_grade == 9 & schyear >= 2001 & schyear <= 2007)

/* Drop redundant variables */

drop bls_comp bla_comp obr_comp bls_rank bla_rank obr_rank gpa verb quan read math 
drop female_a black_a white_a asian_a hisp_a other_a gpa_score isee_score verb_score quan_score read_score math_score 
drop obr_cutoff bla_cutoff bls_cutoff obr_gencomp bla_gencomp bls_gencomp obr_stdcomp bla_stdcomp bls_stdcomp 
drop obr_aacutoff_comp obr_aacutoff_gencomp bla_aacutoff_comp bla_aacutoff_gencomp bls_aacutoff_comp bls_aacutoff_gencomp
drop *_level *raw pcorr_* *barrons_* *attend_* *persist_* *graduate_* *apnbr* *apscr* *apstem*
drop inexam* inbps* in9sch* graderep bpsyrs innext* prop_* distance_* peer2_* peergap_* apprank_* peer_q*
drop *_sat *_satv *_satm *_satw *_psat *_psatr *_psatm *_psatw
drop *_eng4 *_math4 *_eng7 *_math7 *_eng8 *_math8 *_scie8 *_sciehs *_biohs *_physhs *_chemhs *_techhs
drop risk_set_prefs

/* Prepare an estimate file */

est clear
gen temp = 1
quietly: reg temp
esttab using "$dir/Extrapolation/Tables/estimates.csv", nogaps compress nonotes nonumbers replace
drop temp


/**********************/
/* RD plots/estimates */
/**********************/


/* RD plots */
		
local xmin = -20
local xstep = 5
local xmax = 20

local ymin = 0
local ystep = 0.5
local ymax = 2

foreach grade in 7 9 {
	foreach outcome in offer to peer_math peer_eng math10 eng10 {
	
		local xtitle = "Standardized Running Variable"

		if ("`outcome'" == "offer") {
			local title = "Offer Rate"
			local ytitle = "Offer Rate"
			
			local ymin = 0
			local ystep = 0.25
			local ymax = 1
		}
		else if ("`outcome'" == "to") {
			local title = "Enrollment Rate"
			local ytitle = "Enrollment Rate"

			local ymin = 0
			local ystep = 0.25
			local ymax = 1
		}
		else if ("`outcome'" == "peer_eng") {
			local title = "Peer English"
			local ytitle = "Peer Mean"

			local ymin = -0.5
			local ystep = 0.5
			local ymax = 2.5
		}
		else if ("`outcome'" == "peer_math") {
			local title = "Peer Math"
			local ytitle = "Peer Mean"
			
			local ymin = -0.5
			local ystep = 0.5
			local ymax = 2.5
		}
		else if ("`outcome'" == "eng10") {
			local title = "10th Grade ELA"
			local ytitle = "Standardized Score"
			
			local ymin = 0
			local ystep = 0.5
			local ymax = 2
		}
		else if ("`outcome'" == "math10") {
			local title = "10th Grade Math"
			local ytitle = "Standardized Score"
			
			local ymin = 0
			local ystep = 0.5
			local ymax = 2
		}
		
		foreach sch in obr bla bls {
			preserve
			if ("`outcome'" == "math10" | "`outcome'" == "eng10") {
				keep if sample_`grade'_`outcome' == 1
			}
			keep if experiment_`sch' == 1 & disc_`sch' == 1 & app_grade == `grade'
			count 
			
			if ("`outcome'" == "offer" | "`outcome'" == "to") {
				local y = "`outcome'_`sch'"
			}
			else {
				local y = "`outcome'"
			}
			
			local runvar_controls = "runvar_`sch' intrunvar_`sch'"

			if ("`sch'" == "obr") {
				local subtitle = "O'Bryant"
			}
			else if ("`sch'" == "bla") {
				local subtitle = "Latin Academy"
			}
			else if ("`sch'" == "bls") {
				local subtitle = "Latin School"
			}
			
			capture: rdob_mod2 `y' runvar_`sch', c(0) dm(1)
			if (_rc == 0) {
				local h_opt = r(h_opt)
			}
			else {
				local h_opt = 10
			}
			gen lambda = (1 - abs(runvar_`sch'/`h_opt'))*(abs(runvar_`sch'/`h_opt') <= 1)				
			lpoly `y' runvar_`sch' if above_`sch' == 0, kernel(triangle) bw(`h_opt') degree(1) gen(pred0) at(runvar_`sch') nograph
			lpoly `y' runvar_`sch' if above_`sch' == 1, kernel(triangle) bw(`h_opt') degree(1) gen(pred1) at(runvar_`sch') nograph
			gen pred = pred0 if above_`sch' == 0
			replace pred = pred1 if above_`sch' == 1
			
			egen temp_x = cut(runvar_`sch'), at(-100(1)100)
			replace temp_x = temp_x  + (1/2)
			bysort temp_x: egen temp_y = mean(`y')
			bysort temp_x: gen temp_i = _n
			
			keep if runvar_`sch' >= `xmin' & runvar_`sch' <= `xmax'
			
			graph twoway (scatter temp_y temp_x if temp_i == 1) ///
				(line pred runvar_`sch' if above_`sch' == 0, lcolor(black) lpattern(solid) lwidth(medthick) sort) ///
				(line pred runvar_`sch' if above_`sch' == 1, lcolor(black) lpattern(solid) lwidth(medthick) sort), ///
				legend(off) xtitle("`xtitle'") ytitle("`ytitle'") ///
				subtitle(`subtitle') xline(0) xlabel(`xmin'(`xstep')`xmax') ylabel(`ymin'(`ystep')`ymax') graphregion(color(white) fcolor(white)) plotregion(color(white) fcolor(white)) plotregion(style(none)) saving(temp_`outcome'_`sch', replace)

			restore
		}	
	}
	
	graph combine temp_offer_obr.gph temp_offer_bls.gph temp_to_obr.gph temp_to_bls.gph, ///
		cols(2) imargin(medium) iscale(0.5) ycommon xcommon t1title("Offer Rate", pos(12) ring(0) size(small)) t2title("Enrollment Rate", pos(0) ring(0) size(small)) ///
		graphregion(color(white) fcolor(white)) plotregion(color(white) fcolor(white))  plotregion(style(none)) 
	graph export "$dir/Extrapolation/Plots/nonparametric_rd_plot_offer_enrollment_`grade'.eps", replace
	
	graph combine temp_peer_math_obr.gph temp_peer_math_bls.gph temp_peer_eng_obr.gph temp_peer_eng_bls.gph, ///
		cols(2) imargin(medium) iscale(0.5) ycommon xcommon t1title(Peer Math, pos(12) ring(0) size(small)) t2title(Peer English, pos(0) ring(0) size(small)) ///
		graphregion(color(white) fcolor(white)) plotregion(color(white) fcolor(white))  plotregion(style(none)) 
	graph export "$dir/Extrapolation/Plots/nonparametric_rd_plot_peer_mean_`grade'.eps", replace
	
	graph combine temp_math10_obr.gph temp_math10_bls.gph temp_eng10_obr.gph temp_eng10_bls.gph, ///
		cols(2) imargin(medium) iscale(0.5) ycommon xcommon t1title(10th Grade Math, pos(12) ring(0) size(small)) t2title(10th Grade ELA, pos(0) ring(10) size(small)) ///
		graphregion(color(white) fcolor(white)) plotregion(color(white) fcolor(white))  plotregion(style(none)) 
	graph export "$dir/Extrapolation/Plots/nonparametric_rd_plot_`grade'.eps", replace
}

/* RD estimates */

foreach grade in 7 9 {
	foreach outcome in math10 eng10 {	
		est clear
		foreach sch in obr bla bls {
		
			local fixed_effects = "i.risk_set i.app_prefs i.yr_`outcome'"
			local se = "robust"
			
			preserve
			quietly: keep if sample_`grade'_`outcome' == 1 & experiment_`sch' == 1 & disc_`sch' == 1
			
			egen clusterid = group(yr_`outcome' sch_`outcome')
			
			quietly: xi: reg `outcome' above_`sch' runvar_`sch' intrunvar_`sch' runvar2_`sch' intrunvar2_`sch' runvar3_`sch' intrunvar3_`sch' `fixed_effects', `se'
			eststo rd_`grade'_`outcome'_`sch'_par	

			rdob_mod2 `outcome' runvar_`sch', c(0) dm(1)
			local h_opt = r(h_opt)
			gen lambda = (1 - abs(runvar_`sch'/`h_opt'))*(abs(runvar_`sch'/`h_opt') <= 1)
			quietly: xi: reg `outcome' above_`sch' runvar_`sch' intrunvar_`sch' `fixed_effects' [aw = lambda], `se'
			eststo rd_`grade'_`outcome'_`sch'_nonpar
			
			restore
		}
		local options = `"b(3) se(3) mtitle nogaps compress nonotes nonumbers star(* 0.1 ** 0.05 *** 0.01)"'
		esttab rd_*_par rd_*_nonpar using "$dir/Extrapolation/Tables/estimates.csv", `options' keep(above) rename(above_obr above above_bla above above_bls above) append
	}
}


/******************************************/
/* Parametric/Dong & Lewbel extrapolation */
/******************************************/

		
local xmin = -20
local xstep = 5
local xmax = 20

local ymin = 0
local ystep = 0.5
local ymax = 2

local ylabel = "0(0.5)2"

local xtitle = "Standardized Running Variable"
local ytitle = "Standardized Score"

foreach grade in 7 9 {
	foreach outcome in math10 eng10 {
			
		foreach sch in obr bls {
		
			local fixed_effects = "i.risk_set i.app_prefs i.yr_`outcome'"
			local se = "robust"
		
			forval dl = 0/1 {
				forval i = 1/3 {
					preserve
					keep if experiment_`sch' == 1 & sample_`grade'_`outcome' == 1 & disc_`sch' == 1

					if (`i' == 2 & "`sch'" == "obr") {
						local title = "O'Bryant"
					}
					else if (`i' == 2 & "`sch'" == "bls") {
						local title = "Latin School"
					}
					else {
						local title = `"."", color("scheme background")"'
					}
					
					if (`i' == 1) {
						orthpoly runvar_`sch' if `outcome' != . & above_`sch' == 0 , degree(3) gen(orthog1 orthog2 orthog3)
						reg `outcome' orthog1 if above_`sch' == 0, `se'
						eststo orthg0_`grade'_`outcome'_`sch'_1
						drop orthog1 orthog2 orthog3
						orthpoly runvar_`sch' if `outcome' != . & above_`sch' == 1, degree(3) gen(orthog1 orthog2 orthog3)
						reg `outcome' orthog1 if above_`sch' == 1, `se'
						eststo orthg1_`grade'_`outcome'_`sch'_1
						drop orthog1 orthog2 orthog3
						
						local runvar_controls = "runvar_`sch' intrunvar_`sch'"
						local subtitle = "Linear"
					}
					else if (`i' == 2) {
						orthpoly runvar_`sch' if `outcome' != . & above_`sch' == 0 , degree(3) gen(orthog1 orthog2 orthog3)
						reg `outcome' orthog1 orthog2 if above_`sch' == 0, `se'
						eststo orthg0_`grade'_`outcome'_`sch'_2
						drop orthog1 orthog2 orthog3
						orthpoly runvar_`sch' if `outcome' != . & above_`sch' == 1, degree(3) gen(orthog1 orthog2 orthog3)
						reg `outcome' orthog1 orthog2 if above_`sch' == 1, `se'
						eststo orthg1_`grade'_`outcome'_`sch'_2
						drop orthog1 orthog2 orthog3
						
						local runvar_controls = "runvar_`sch' runvar2_`sch' intrunvar_`sch' intrunvar2_`sch'"
						local subtitle = "Quadratic"
					}
					else if (`i' == 3) {
						orthpoly runvar_`sch' if `outcome' != . & above_`sch' == 0 , degree(3) gen(orthog1 orthog2 orthog3)
						reg `outcome' orthog1 orthog2 orthog3 if above_`sch' == 0, `se'
						eststo orthg0_`grade'_`outcome'_`sch'_3
						drop orthog1 orthog2 orthog3
						orthpoly runvar_`sch' if `outcome' != . & above_`sch' == 1, degree(3) gen(orthog1 orthog2 orthog3)
						reg `outcome' orthog1 orthog2 orthog3 if above_`sch' == 1, `se'
						eststo orthg1_`grade'_`outcome'_`sch'_3
						drop orthog1 orthog2 orthog3
						
						local runvar_controls = "runvar_`sch' runvar2_`sch' runvar3_`sch' intrunvar_`sch' intrunvar2_`sch' intrunvar3_`sch'"
						local subtitle = "Cubic"
					}
					
					if ("`sch'" == "obr") {
						local pred0_color = "none"
						local pred1_color = "red"
					}
					else if ("`sch'" == "bls") {
						local pred0_color = "red"
						local pred1_color = "none"
					}
					
					if (`dl' == 0) {
					
						quietly: reg `outcome' above_`sch' `runvar_controls'
						predict pred	
						quietly: reg `outcome' above_`sch' `runvar_controls' if above_`sch' == 0
						predict pred0
						quietly: reg `outcome' above_`sch' `runvar_controls' if above_`sch' == 1
						predict pred1
						
						foreach s in 1 5 10 15 {
							if ("`sch'" == "obr") { 
								local s1 = -`s'
								local s2 = (-`s')^2
								local s3 = (-`s')^3
							}
							else {
								local s1 = `s'
								local s2 = `s'^2
								local s3 = `s'^3
							}
							quietly: reg `outcome' above_`sch' `runvar_controls' `fixed_effects', `se'
							local n = e(N)
							if (`i' == 1) {
								lincom above_`sch' + `s1'*intrunvar_`sch'
							}
							else if (`i' == 2) {
								lincom above_`sch' + `s1'*intrunvar_`sch' + `s2'*intrunvar2_`sch'
							}
							else if (`i' == 3) {
								lincom above_`sch' + `s1'*intrunvar_`sch' + `s2'*intrunvar2_`sch' + `s3'*intrunvar3_`sch'
							}									
							mat b = (`r(estimate)')
							mat V = (`r(se)'^2)
							mat colnames b = est
							mat rownames b = est
							mat colnames V = est
							mat rownames V = est
							eret post b V, obs(`n')
							eststo x_`grade'_`outcome'_`sch'_par`i'_`s'	
						}
					}
					
					if (`dl' == 1) {
						rdob_mod2 `outcome' runvar_`sch' if disc_`sch' == 1, c(0) dm(1)
						local h_opt = r(h_opt)
						gen lambda = (1 - abs(runvar_`sch'/`h_opt'))*(abs(runvar_`sch'/`h_opt') <= 1)				
						lpoly `outcome' runvar_`sch' if above_`sch' == 0, kernel(triangle) bw(`h_opt') degree(`i') gen(pred0) at(runvar_`sch') nograph
						lpoly `outcome' runvar_`sch' if above_`sch' == 1, kernel(triangle) bw(`h_opt') degree(`i') gen(pred1) at(runvar_`sch') nograph
						gen pred = pred0 if above_`sch' == 0
						replace pred = pred1 if above_`sch' == 1
						drop pred0 pred1
						quietly: reg `outcome' above_`sch' `runvar_controls' [aw = lambda] if above_`sch' == 0
						predict pred0
						quietly: reg `outcome' above_`sch' `runvar_controls' [aw = lambda] if above_`sch' == 1
						predict pred1
						
						foreach s in 1 5 10 15 {
							capture: drop lambda_dl
							if ("`sch'" == "obr") {
								local s1 = -`s'
								local s2 = (-`s')^2
								local s3 = (-`s')^3
								gen lambda_dl = (1 - abs((runvar_`sch' + `s'*(1 - above_`sch'))/`h_opt'))*(abs((runvar_`sch' + `s'*(1 - above_`sch'))/`h_opt') <= 1)
							}
							else {
								local s1 = `s'
								local s2 = `s'^2
								local s3 = `s'^3
								gen lambda_dl = (1 - abs((runvar_`sch' - `s'*above_`sch')/`h_opt'))*(abs((runvar_`sch' - `s'*above_`sch')/`h_opt') <= 1)
							}
							quietly: reg `outcome' above_`sch' `runvar_controls' `fixed_effects' [aw = lambda_dl], `se'
							local n = e(N)
							if (`i' == 1) {
								lincom above_`sch' + `s1'*intrunvar_`sch'
							}
							else if (`i' == 2) {
								lincom above_`sch' + `s1'*intrunvar_`sch' + `s2'*intrunvar2_`sch'
							}
							else if (`i' == 3) {
								lincom above_`sch' + `s1'*intrunvar_`sch' + `s2'*intrunvar2_`sch' + `s3'*intrunvar3_`sch'
							}									
							mat b = (`r(estimate)')
							mat V = (`r(se)'^2)
							mat colnames b = est
							mat rownames b = est
							mat colnames V = est
							mat rownames V = est
							eret post b V, obs(`n')
							eststo x_`grade'_`outcome'_`sch'_dl`i'_`s'	
						}
					}

					egen temp_x = cut(runvar_`sch'), at(-100(1)100)
					replace temp_x = temp_x  + (1/2)
					bysort temp_x: egen temp_y = mean(`outcome')
					bysort temp_x: gen temp_i = _n
					
					keep if runvar_`sch' >= `xmin' & runvar_`sch' <= `xmax'
					replace pred0 = . if pred0 < `ymin' | pred0 > `ymax'
					replace pred1 = . if pred1 < `ymin' | pred1 > `ymax'
					
					graph twoway (scatter temp_y temp_x if temp_i == 1) ///
						(lpoly pred runvar_`sch' if above_`sch' == 0, bw(10) degree(3) lcolor(black) lpattern(solid) lwidth(medthick)) ///
						(lpoly pred runvar_`sch' if above_`sch' == 1, bw(10) degree(3) lcolor(black) lpattern(solid) lwidth(medthick)) ///
						(lpoly pred1 temp_x if above_`sch' == 0, bw(10) degree(3) lcolor(red) lpattern(shortdash) lwidth(medthick) lcolor(`pred1_color')) ///
						(lpoly pred0 temp_x if above_`sch' == 1, bw(10) degree(3) lcolor(red) lpattern(shortdash) lwidth(medthick) lcolor(`pred0_color')), ///
						legend(order(2 4) size(small) label(2 "Fitted") label(4 "Extrapolated")) xtitle("`xtitle'") ytitle("`ytitle'") ///
						title(`title') subtitle(`subtitle') xline(0) xlabel(`xmin'(`xstep')`xmax') ylabel(`ymin'(`ystep')`ymax') ///
						graphregion(color(white) fcolor(white)) plotregion(color(white) fcolor(white)) plotregion(style(none)) saving(temp_`sch'_dl`dl'_`i', replace)

					restore
				}
			}
		}
		grc1leg temp_obr_dl0_1.gph temp_obr_dl0_2.gph temp_obr_dl0_3.gph temp_bls_dl0_1.gph temp_bls_dl0_2.gph temp_bls_dl0_3.gph, cols(3) imargin(tiny) iscale(0.5) ycommon xcommon ///
			graphregion(color(white) fcolor(white)) plotregion(color(white) fcolor(white)) plotregion(style(none))
		graph export "$dir/Extrapolation/Plots/parametric_extrapolation_`grade'_`outcome'.eps", replace
		grc1leg temp_obr_dl1_1.gph temp_obr_dl1_2.gph temp_obr_dl1_3.gph temp_bls_dl1_1.gph temp_bls_dl1_2.gph temp_bls_dl1_3.gph, cols(3) imargin(tiny) iscale(0.5) ycommon xcommon ///
			graphregion(color(white) fcolor(white)) plotregion(color(white) fcolor(white)) plotregion(style(none))
		graph export "$dir/Extrapolation/Plots/dong_lewbel_extrapolation_`grade'_`outcome'.eps", replace

		local options = `"b(3) se(3) mtitle nogaps compress nonotes nonumbers star(* 0.1 ** 0.05 *** 0.01)"'
		esttab orthg*_`outcome'_*_1  using "$dir/Extrapolation/Tables/estimates.csv", `options' append
		esttab orthg*_`outcome'_*_2  using "$dir/Extrapolation/Tables/estimates.csv", `options' append
		esttab orthg*_`outcome'_*_3  using "$dir/Extrapolation/Tables/estimates.csv", `options' append
		esttab x_`grade'_`outcome'_*_par1* x_`grade'_`outcome'_*_dl1* using "$dir/Extrapolation/Tables/estimates.csv", `options' append
		esttab x_`grade'_`outcome'_*_par2* x_`grade'_`outcome'_*_dl2* using "$dir/Extrapolation/Tables/estimates.csv", `options' append
		esttab x_`grade'_`outcome'_*_par3* x_`grade'_`outcome'_*_dl3* using "$dir/Extrapolation/Tables/estimates.csv", `options' append
	}
}


/*********************/
/* CIA extrapolation */
/*********************/


/* Check common support on covariates */

quietly: tab app_prefs, gen(app_prefs_)

foreach grade in 7 9 {
	foreach sch in obr bls {
		preserve
		quietly: keep if app_grade == `grade' & experiment_`sch' == 1 & disc_`sch' == 1
		quietly: keep if abs(runvar_`sch') <= 20
		
		mat results = J(1, 4, .)
		foreach x of varlist eng4 math4 eng7 math8 female black white asian hisp free lep special app_prefs_* {
			quietly: su `x' if above_`sch' == 0
			mat row = (r(min), r(max))
			quietly: su `x' if above_`sch' == 1
			mat row = (row, r(min), r(max))
			mat rownames row = "`x'"
			mat results = (results\ row)
		}
		mat results = results[2..., 1...]
		mat colnames results = "min_0" "max_0" "min_1" "max_1"
		disp "`grade' - `sch'"
		mat list results
		
		restore
	}
}

drop app_prefs_*

/* Check common support on pscore */

foreach grade in 7 9 {

	local controls1 = "female black white asian hisp"
	local controls2 = "free lep special"
	if (`grade' == 7) {
		local controls3 = "eng4 math4"
	}
	else {
		local controls3 = "eng4 math4 eng7 math8"	
	}

	foreach outcome in math10 eng10 {
		foreach window in 20 15 10 {
			foreach sch in obr bls {
			
				local ytitle = "Frequency"
				local xtitle = "Propensity Score"
			
				if ("`sch'" == "obr") {
					local fixed_effects = "i.risk_set i.app_prefs"
				}
				else if ("`sch'" == "bls") {
					local fixed_effects = "i.risk_set"
				}
			
				if (`window' == 20 & "`sch'" == "obr") {
					local title = "O'Bryant"
				}
				else if (`window' == 20 & "`sch'" == "bls") {
					local title = "Latin School"
				}
				else {
					local title = `"."", color("scheme background")"'
				}
				
				if ("`sch'" == "obr") {
					local ylabel = `"-150 "150" -100 "100" -50 "50" 0 "0" 50 "50" 100 "100" 150 "150""'
					local ymax = 150
				}
				else if ("`sch'" == "bls") {
					local ylabel = `"-75 "75" -50 "50" -25 "25" 0 "0" 25 "25" 50 "50" 75 "75""'
					local ymax = 75
				}
							
				preserve
				quietly: keep if sample_`grade'_`outcome' == 1 & experiment_`sch' == 1 & disc_`sch' == 1
				quietly: keep if abs(runvar_`sch') <= `window'
				egen temp = rowmiss(`outcome' `controls1' `controls2' `controls3') 			
				quietly: keep if temp == 0
				drop temp
				
				xi: logit above_`sch' `controls1' `controls2' `controls3' `fixed_effects'
				local n = e(N)
				predict pred
				egen temp = cut(pred), at(0(0.1)1)
				egen temp_x = group(temp)
				replace temp_x = temp_x - 0.5
				drop temp
				
				bysort temp_x: egen temp = count(temp_x) if above_`sch' == 0
				bysort temp_x: egen temp_y0 = mean(temp)
				drop temp				
				quietly: count if temp_x != . & above_`sch' == 0
				replace temp_y0 = -temp_y0
				quietly: su temp_y0
				if (`r(min)' < -`ymax') {
					local text = -`r(min)'
				}
				else {
					local text = ""
				}
				replace temp_y0 = -`ymax' if temp_y0 < -`ymax' & temp_y0 != .
				
				bysort temp_x: egen temp = count(temp_x) if above_`sch' == 1
				bysort temp_x: egen temp_y1 = mean(temp)
				drop temp
				count if temp_x != . & above_`sch' == 1
				replace temp_y1 = temp_y1
				
				bysort temp_x: gen temp_i = _n

				graph twoway (bar temp_y1 temp_x if temp_i == 1) ///
					(bar temp_y0 temp_x if temp_i == 1) ///
					, xtitle("`xtitle'") ytitle("`ytitle'") legend(label(1 "Treated") label(2 "Untreated") size(small)) ///
					xlabel(0 "0" 2 "0.2" 4 "0.4" 6 "0.6" 8 "0.8" 10 "1") ylabel(`ylabel') text(-71 0.5 "`text'", size(small)) ///
					title(`title') /*subtitle("Window = +/- `window' (N = `n')")*/ ///
					graphregion(color(white) fcolor(white)) plotregion(color(white) fcolor(white)) plotregion(style(none)) saving(temp_`sch'_`window', replace)
				
				restore
			}
		}
*		grc1leg temp_obr_20.gph temp_obr_15.gph temp_obr_10.gph temp_bls_20.gph temp_bls_15.gph temp_bls_10.gph, cols(3) imargin(tiny) iscale(0.5) graphregion(color(white) fcolor(white)) plotregion(color(white) fcolor(white)) plotregion(style(none))
		grc1leg temp_obr_20.gph temp_bls_20.gph, cols(2) imargin(tiny) iscale(0.5) graphregion(color(white) fcolor(white)) plotregion(color(white) fcolor(white)) plotregion(style(none))
		graph export "$dir/Extrapolation/Plots/common_support_plot_`grade'_`outcome'.eps", replace
	}
}

/* CIA test */

foreach grade in 7 9 {

	local controls1 = "female black white asian hisp"
	local controls2 = "free lep special"
	if (`grade' == 7) {
		local controls3 = "eng4 math4"		
	}
	else {
		local controls3 = "eng4 math4 eng7 math8"	
	}
	
	foreach outcome in math10 eng10 {
	
		if ("`outcome'" == "eng10") {
			local title = "10th Grade ELA"
		}
		else if ("`outcome'" == "math10") {
			local title = "10th Grade Math"
		}
	
		local ytitle = "Standardized Score (Residual)"
		local xtitle = "Standardized Running Variable"
				
		foreach window in 20 15 10 {
			est clear
			foreach sch in obr bls {
			
				if ("`sch'" == "obr") {
					local subtitle = "O'Bryant"
				}
				else if ("`sch'" == "bls") {
					local subtitle = "Latin School"
				}
			
				local fixed_effects = "i.risk_set i.app_prefs i.yr_`outcome'"
				local se = "robust"
				
				preserve
				quietly: keep if sample_`grade'_`outcome' == 1 & experiment_`sch' == 1 & disc_`sch' == 1
				quietly: keep if abs(runvar_`sch') <= `window'
				egen temp = rowmiss(`outcome' `controls1' `controls2' `controls3') 			
				quietly: keep if temp == 0
				
				egen clusterid = group(yr_`outcome' sch_`outcome')
				
				forval k = 0/1 {
					quietly: xi: reg `outcome' runvar_`sch' `fixed_effects' if above_`sch' == `k', `se'
					eststo citest1_`grade'_`outcome'_`sch'_`k'_`window'
					quietly: eststo: xi: reg `outcome' runvar_`sch' `controls1' `controls2' `controls3' `fixed_effects' if above_`sch' == `k', `se'
					eststo citest2_`grade'_`outcome'_`sch'_`k'_`window'
			
					quietly: eststo: xi: reg `outcome' `controls1' `controls2' `controls3' `fixed_effects' if above_`sch' == `k', `se'
					predict res_y`k', res							
				}
				
				gen res_y = (1 - above_`sch')*res_y0 + above_`sch'*res_y1
				egen temp_x = cut(runvar_`sch'), at(-100(1)100)
				replace temp_x = temp_x  + (1/2)
				bysort temp_x: egen temp_y = mean(res_y)
				bysort temp_x: gen temp_i = _n
				
				rdob_mod2 res_y runvar_`sch', c(0) dm(1)
				local h_opt = r(h_opt)
				
				graph twoway (scatter temp_y temp_x if temp_i == 1) ///
					(lpoly res_y runvar_`sch' if above_`sch' == 0, kernel(triangle) degree(1) bw(`h_opt') lcolor(black) lpattern(solid) lwidth(thick)) ///
					(lpoly res_y runvar_`sch' if above_`sch' == 1, kernel(triangle) degree(1) bw(`h_opt') lcolor(black) lpattern(solid) lwidth(thick)) ///
					, ylab(-1(0.5)1) xline(0) yline(0, lcolor(black) lpattern(shortdash) lwidth(medthin)) subtitle(`subtitle') xtitle("`xtitle'") ytitle("`ytitle'") legend(off) ///
					graphregion(color(white) fcolor(white)) plotregion(color(white) fcolor(white)) plotregion(style(none)) saving(temp_`outcome'_`sch'_`window', replace)
					
				restore
			}
			local options = `"b(3) se(3) mtitle nogaps compress nonotes nonumbers star(* 0.1 ** 0.05 *** 0.01)"'
			esttab citest1_* citest2_* using "$dir/Extrapolation/Tables/estimates.csv", `options' keep(runvar) rename(runvar_obr runvar runvar_bls runvar) append		
		}
	}
	graph combine temp_math10_obr_20.gph temp_math10_bls_20.gph temp_eng10_obr_20.gph temp_eng10_bls_20.gph, ///
		cols(2) imargin(medium) iscale(0.5) ycommon xcommon t1title(10th Grade Math, pos(12) ring(0) size(small)) t2title(10th Grade ELA, pos(0) ring(10) size(small)) ///
		graphregion(color(white) fcolor(white)) plotregion(color(white) fcolor(white)) plotregion(style(none))
	graph export "$dir/Extrapolation/Plots/citest_plot_`grade'_20.eps", replace
}

/* CIA estimates - sharp */

capture: program drop ciest
program ciest, rclass
	version 11.0
	syntax varlist(fv), outcome(string) sch(string) method(string) trim(real)
	
	noisily {
	
		rdob_mod2 `outcome' runvar_`sch', c(0) dm(1)
		local h_opt = r(h_opt)
		gen lambda = (1 - abs(runvar_`sch'/`h_opt'))*(abs(runvar_`sch'/`h_opt') <= 1)		
		quietly: xi: reg `outcome' above_`sch' runvar_`sch' intrunvar_`sch' i.risk_set i.app_prefs i.yr_`outcome' [aw = lambda]
		local delta_rd = _b[above_`sch']
		drop lambda
				
		gen temp = (runvar_`sch' > -5 & runvar_`sch' <= 5)
		su temp if `outcome' != .
		local q = r(mean)
		capture noisily: xi: logit temp `varlist' if `outcome' != ., iterate(100)
		if (_rc == 0 & `e(converged)' == 1) {
			predict pred_q if `outcome' != .
		}
		drop temp
		else {
			gen pred_q = .		
		}
			
		if ("`method'" == "kline") {
			reg `outcome' `varlist' if above_`sch' == 0
			local n0 = e(N)
			predict pred0, xb
			reg `outcome' `varlist' if above_`sch' == 1
			local n1 = e(N)
			predict pred1, xb
			
			gen effect = pred1 - pred0
			su effect if above_`sch' == 0
			local effect0 = r(mean)
			su effect if above_`sch' == 1
			local effect1 = r(mean)
			
			gen w = pred_q/`q'
			su w
			replace w = w/r(mean)
			gen temp = w*effect
			su temp
			local delta_cia = r(mean)
			drop w effect temp
			
			local diff_cia_rd = `delta_cia' - `delta_rd'			
		}
		else if ("`method'" == "hir") {
			capture noisily: xi: logit above_`sch' `varlist' if `outcome' != ., iterate(100)
			if (_rc == 0 & `e(converged)' == 1) {
				predict pred_p if `outcome' != .
				replace pred_p = . if (pred_p < `trim' | pred_p > 1 - `trim')
				
				su pred_p if above_`sch' == 0
				local n0 = r(N)
				su pred_p if above_`sch' == 1
				local n1 = r(N)
				su above_`sch'
				local p = r(mean)
				
				gen w0 = (1 - above_`sch')/(1 - pred_p)
				gen w1 = above_`sch'/pred_p
				
				gen w00 = w0*((1 - pred_p)/(1 - `p'))
				gen w10 = w1*((1 - pred_p)/(1 - `p'))
				
				gen w01 = w0*(pred_p/`p')
				gen w11 = w1*(pred_p/`p')
				
				gen w0q = w0*(pred_q/`q')
				gen w1q = w1*(pred_q/`q')
				
				foreach w in w00 w10 w01 w11 w0q w1q {
					su `w' if `outcome' != .
					replace `w' = `w'/r(mean)
				}
		
				gen temp = (w10 - w00)*`outcome'
				su temp
				local effect0 = r(mean)
				drop temp
				
				gen temp = (w11 - w01)*`outcome'
				su temp
				local effect1 = r(mean)
				drop temp
				
				gen temp = (w1q - w0q)*`outcome'
				su temp
				local delta_cia = r(mean)
				drop temp
				
				local diff_cia_rd = `delta_cia' - `delta_rd'
				
				drop pred_p w0 w1 w00 w10 w01 w11 w0q w1q
			}
			else {
				local n0 = .
				local n1 = .
				local effect0 = .
				local effect1 = .
				local diff_cia_rd = .
			}
		}
		drop pred_q
	}
	return scalar effect0 = `effect0'
	return scalar effect1 = `effect1'
	return scalar n0 = `n0'
	return scalar n1 = `n1'
	return scalar diff_cia_rd = `diff_cia_rd'
end

foreach grade in 7 9 {

	local controls1 = "female black white asian hisp"
	local controls2 = "free lep special"
	if (`grade' == 7) {
		local controls3 = "eng4 math4"
	}
	else {
		local controls3 = "eng4 math4 eng7 math8"	
	}
	
	est clear
	foreach outcome in math10 eng10 {
		foreach window in 20 {
			foreach method in kline hir {
				foreach sch in obr bls {
					noisily {
						if ("`method'" == "kline" | "`sch'" == "obr") {
							local fixed_effects = "i.risk_set i.app_prefs i.yr_`outcome'"
						}
						else {
							local fixed_effects = "i.risk_set"
						}
						
						preserve
						keep if sample_`grade'_`outcome' == 1 & experiment_`sch' == 1 & disc_`sch' == 1
						keep if abs(runvar_`sch') <= `window'
						egen temp = rowmiss(`outcome' `controls1' `controls2' `controls3') 			
						keep if temp == 0
						drop temp
						
						if ("`method'" == "hir5") {
							local tmethod = "hir"
							local trim = 0.05
						}	
						else if ("`method'" == "hir10") {
							local tmethod = "hir"
							local trim = 0.10
						}
						else {
							local tmethod = "`method'"
							local trim = 0
						}
						
						ciest `controls1' `controls2' `controls3' `fixed_effects', outcome(`outcome') sch(`sch') method(`tmethod') trim(`trim')
						local est0 = `r(effect0)'
						local est1 = `r(effect1)'
						local n0 = `r(n0)'
						local n1 = `r(n1)'
						local diff_cia_rd_est = `r(diff_cia_rd)'

						if ("`method'" == "kline") {
						
							local ytitle = "Standardized Score"
							local xtitle = "Standardized Running Variable"
						
							lpoly pred0 runvar_`sch', degree(1) kernel(triangle) generate(fit0) at(runvar_`sch') nograph
							lpoly pred1 runvar_`sch', degree(1) kernel(triangle) generate(fit1) at(runvar_`sch') nograph
							
							egen temp_x = cut(runvar_`sch'), at(-100(1)100)
							replace temp_x = temp_x  + (1/2)
							bysort temp_x: egen temp_y = mean(`outcome')
							bysort temp_x: egen temp_fit0 = mean(fit0)
							bysort temp_x: egen temp_fit1 = mean(fit1)
							bysort temp_x: egen temp_pred0 = mean(pred0)
							bysort temp_x: egen temp_pred1 = mean(pred1)
							bysort temp_x: gen temp_i = _n
							
							if ("`sch'" == "obr") {
								local pred0_color = "none"
								local pred1_color = "red"
								local subtitle = "O'Bryant"
							}
							else if ("`sch'" == "bls") {
								local pred0_color = "red"
								local pred1_color = "none"
								local subtitle = "Latin School"
							}
							
							if ("`outcome'" == "math10") {
								local title = "10th Grade Math"
							}
							else if ("`outcome'" == "eng10") {
								local title = "10th Grade ELA"
							}
						
							graph twoway (scatter temp_y temp_x if abs(temp_x) <= `window' & temp_i == 1) ///
								(scatter temp_pred1 temp_x if abs(temp_x) <= `window' & runvar_`sch' < 0 & temp_i == 1, msymbol(x) mcolor(`pred1_color')) ///
								(scatter temp_pred0 temp_x if abs(temp_x) <= `window' & runvar_`sch' >= 0 & temp_i == 1, msymbol(x) mcolor(`pred0_color')) ///
								(lpoly fit0 runvar_`sch' if abs(runvar_`sch') <= `window' & runvar_`sch' < 0, degree(1) kernel(triangle) lcolor(black) lpattern(solid) lwidth(medthick)) ///
								(lpoly fit1 runvar_`sch' if abs(runvar_`sch') <= `window' & runvar_`sch' >= 0, degree(1) kernel(triangle) lcolor(black) lpattern(solid) lwidth(medthick)) ///
								(lpoly fit1 runvar_`sch' if abs(runvar_`sch') <= `window' & runvar_`sch' < 0, degree(1) kernel(triangle) lpattern(shortdash) lwidth(medthick) lcolor(`pred1_color')) ///
								(lpoly fit0 runvar_`sch' if abs(runvar_`sch') <= `window' & runvar_`sch' >= 0, degree(1) kernel(triangle) lpattern(shortdash) lwidth(medthick) lcolor(`pred0_color')) ///
								, legend(order(4 6) size(small) label(4 "Fitted") label(6 "Extrapolated")) xtitle("`xtitle'") ytitle("`ytitle'") ///
								subtitle(`subtitle') xline(0) xlabel(-20(5)20) graphregion(color(white) fcolor(white)) plotregion(color(white) fcolor(white)) plotregion(style(none)) saving(temp_`outcome'_`method'_`window'_`sch', replace)
							capture: drop fit0 fit1
							
							drop pred* temp*
						}		
						
						save "$dir/Extrapolation/bstemp.dta", replace
						set seed 123
						mat bs_est = J(1, 3, .)
						
						forval b = 1/500 {				
							bsample			
							
							quietly: ciest `controls1' `controls2' `controls3' `fixed_effects', outcome(`outcome') sch(`sch') method(`tmethod') trim(`trim')
							mat bs_est = (bs_est\ (`r(effect0)', `r(effect1)', `r(diff_cia_rd)'))	
							
							use "$dir/Extrapolation/bstemp.dta", clear
						}
						mat bs_est = bs_est[2..., 1...]
						mat list bs_est
						
						svmat bs_est
						su bs_est1
						local se0 = `r(sd)'
						su bs_est2
						local se1 = `r(sd)'
						su bs_est3
						local diff_cia_rd_se = `r(sd)'
						disp `diff_cia_rd_est'
						disp `diff_cia_rd_se'
						local tstat = `diff_cia_rd_est'/`diff_cia_rd_se'
						if (`tstat' == . ){
							local tstat = 0
						}
						
						if ("`sch'" == "obr") {
							mat b = (`est0')
							mat V = (`se0'^2)
						}
						else {
							mat b = (`est1')
							mat V = (`se1'^2)
						}
						mat colnames b = est
						mat rownames b = est
						mat colnames V = est
						mat rownames V = est
						eret post b V
						estadd scalar n0  = `n0'
						estadd scalar n1  = `n1'
						estadd scalar tstat  = `tstat'
						eststo ciest_`grade'_`outcome'_`sch'_`window'_`method'
											
						restore
					}
				}
			}
		}
	}
	
	grc1leg temp_math10_kline_20_obr.gph temp_math10_kline_20_bls.gph temp_eng10_kline_20_obr.gph temp_eng10_kline_20_bls.gph, ///
		cols(2) imargin(medium) iscale(0.5) ycommon xcommon t1title(10th Grade Math, pos(12) ring(0) size(small)) t2title(10th Grade ELA, pos(0) ring(10) size(small)) ///
		graphregion(color(white) fcolor(white)) plotregion(color(white) fcolor(white)) plotregion(style(none))
	graph export "$dir/Extrapolation/Plots/cia_sharp_plot_`grade'_kline_20.eps", replace
				
	local options = `"b(3) se(3) mtitle nogaps compress nonotes nonumbers star(* 0.1 ** 0.05 *** 0.01)"'

	esttab *_20_kline *_20_hir using "$dir/Extrapolation/Tables/estimates.csv", `options' keep(est) scalars(n0 n1 tstat) sfmt(0 0 3) noobs append
*	esttab *_15_kline *_15_hir using "$dir/Extrapolation/Tables/estimates.csv", `options' keep(est) scalars(n0 n1 tstat) sfmt(0 0 3) noobs append
*	esttab *_10_kline *_10_hir using "$dir/Extrapolation/Tables/estimates.csv", `options' keep(est) scalars(n0 n1 tstat) sfmt(0 0 3) noobs append

*	esttab *_20_hir5 *_20_hir10 using "$dir/Extrapolation/estimates.csv", `options' keep(est) scalars(n0 n1) sfmt(0) noobs append
*	esttab *_15_hir5 *_15_hir10 using "$dir/Extrapolation/estimates.csv", `options' keep(est) scalars(n0 n1) sfmt(0) noobs append
*	esttab *_10_hir5 *_10_hir10 using "$dir/Extrapolation/estimates.csv", `options' keep(est) scalars(n0 n1) sfmt(0) noobs append
}


/* ***HERE***


/* CIA estimates - fuzzy */

capture: program drop ciest_fuzzy
program ciest_fuzzy, rclass
	version 11.0
	syntax varlist(fv), outcome(string) treatment(string) sch(string) method(string) trim(real)
	
	noisily {
		if ("`method'" == "kline") {
			reg `outcome' `varlist' if above_`sch' == 0
			local n0 = e(N)
			predict pred_rf0 if `outcome' != ., xb
			reg `outcome' `varlist' if above_`sch' == 1
			local n1 = e(N)
			predict pred_rf1 if `outcome' != ., xb
			
			reg `treatment' `varlist' if above_`sch' == 0
			predict pred_fs0 if `outcome' != ., xb
			reg `treatment' `varlist' if above_`sch' == 1
			predict pred_fs1 if `outcome' != ., xb
			
			gen effect_rf = pred_rf1 - pred_rf0
			gen effect_fs = pred_fs1 - pred_fs0

			su effect_rf if above_`sch' == 0
			local rf0 = r(mean)
			su effect_rf if above_`sch' == 1
			local rf1 = r(mean)
			
			su effect_fs if above_`sch' == 0
			local fs0 = r(mean)
			su effect_fs if above_`sch' == 1
			local fs1 = r(mean)
			
			local late0 = `rf0'/`fs0'
			local late1 = `rf1'/`fs1'		
			
			drop pred_* effect_*		
		
		}
		else if ("`method'" == "hir") {
			capture noisily: xi: logit above_`sch' `varlist' if `outcome' != ., iterate(100)
			if (_rc == 0 & `e(converged)' == 1) {
				predict pred_p if `outcome' != .
				replace pred_p = . if (pred_p < `trim' | pred_p > 1 - `trim')
				
				su pred_p if above_`sch' == 0
				local n0 = r(N)
				su pred_p if above_`sch' == 1
				local n1 = r(N)
				su above_`sch'
				local p = r(mean)
				
				gen w0 = (1 - above_`sch')/(1 - pred_p)
				gen w1 = above_`sch'/pred_p
				
				gen w00 = w0*((1 - pred_p)/(1 - `p'))
				gen w10 = w1*((1 - pred_p)/(1 - `p'))
				
				gen w01 = w0*(pred_p/`p')
				gen w11 = w1*(pred_p/`p')
				
				foreach w in w00 w10 w01 w11 {
					su `w' if `outcome' != .
					replace `w' = `w'/r(mean)
				}
		
				gen temp = (w10 - w00)*`outcome'
				su temp
				local rf0 = r(mean)
				drop temp
				
				gen temp = (w11 - w01)*`outcome'
				su temp
				local rf1 = r(mean)
				drop temp
				
				gen temp = (w10 - w00)*`treatment'
				su temp
				local fs0 = r(mean)
				drop temp
				
				gen temp = (w11 - w01)*`treatment'
				su temp
				local fs1 = r(mean)
				drop temp
				
				local late0 = `rf0'/`fs0'
				local late1 = `rf1'/`fs1'				
				
				drop pred_p w0 w1 w00 w10 w01 w11
			}
			else {
				local n0 = .
				local n1 = .
				local rf0 = .
				local rf1 = .
				local fs0 = .
				local fs1 = .
				local late0 = .
				local late1 = .
			}
		}
		else if ("`method'" == "abadie") {
			capture: xi: logit above_`sch' `varlist' if `outcome' != ., iterate(100)
			if (_rc == 0 & `e(converged)' == 1) {
				predict pred_p if `outcome' != .
				replace pred_p = . if (pred_p < `trim' | pred_p > 1 - `trim')
				
				su pred_p if above_`sch' == 0
				local n0 = r(N)
				su pred_p if above_`sch' == 1
				local n1 = r(N)
				su above_`sch'
				local p = r(mean)

				gen w = (above_`sch' - pred_p)/((1 - pred_p)*pred_p)
				gen w0 = w*((1 - pred_p)/(1 - `p'))
				gen w1 = w*(pred_p/`p')
				
				gen kappa = 1 - to_`sch'*(1 - above_`sch')/(1 - pred_p) - (1 - to_`sch')*above_`sch'/pred_p				
				gen kappa0 = kappa*((1 - pred_p)/(1 - `p'))
				gen kappa1 = kappa*(pred_p/`p')
				
				gen w0_y = w0*`outcome'
				gen w1_y = w1*`outcome'
				
				su w0_y
				local rf0 = r(mean)
				su w1_y
				local rf1 = r(mean)
				
				su kappa0
				local fs0 = r(mean)
				su kappa1
				local fs1 = r(mean)
				
				local late0 = `rf0'/`fs0'
				local late1 = `rf1'/`fs1'
				
				drop pred_p w w0 w1 kappa kappa0 kappa1 w0_y w1_y
			}
			else {
				local n0 = .
				local n1 = .
				local rf0 = .
				local rf1 = .
				local fs0 = .
				local fs1 = .
				local late0 = .
				local late1 = .
			}
		}
	}
	return scalar late0 = `late0'
	return scalar late1 = `late1'
	return scalar rf0 = `rf0'
	return scalar rf1 = `rf1'
	return scalar fs0 = `fs0'
	return scalar fs1 = `fs1'
	return scalar n0 = `n0'
	return scalar n1 = `n1'
end

foreach grade in /*7*/ 9 {

	local controls1 = "female black white asian hisp"
	local controls2 = "free lep special"
	if (`grade' == 7) {
		local controls3 = "eng4 math4"
	}
	else {
		local controls3 = "eng4 math4 eng7 math8"	
	}
	
	
	foreach method in kline1 hir1 kline2 hir2 abadie {
		foreach window in 20 {
			est clear
			foreach outcome in math10 eng10 {
				foreach sch in obr bls {
						
					disp "`outcome' - `window' - `method' - `sch'"
				
					noisily {
					
						if ("`method'" == "kline1" | "`method'" == "hir1") {
							local treatment = "to_`sch'"
						}
						else {
							local treatment = "`sch'_yrs_`outcome'"
						}

						if ("`method'" == "kline1" | "`method'" == "kline2" | "`sch'" == "obr") {
							local fixed_effects = "i.risk_set i.app_prefs i.yr_`outcome'"
						}
						else {
							local fixed_effects = "i.risk_set"
						}
						
						if ("`method'" == "abadie5") {
							local tmethod = "abadie"
							local trim = 0.05
						}	
						else if ("`method'" == "abadie10") {
							local tmethod = "abadie"
							local trim = 0.10
						}
						else if ("`method'" == "kline1" | "`method'" == "kline2") {
							local tmethod = "kline"
							local trim = 0
						}
						else if ("`method'" == "hir1" | "`method'" == "hir2") {
							local tmethod = "hir"
							local trim = 0
						}
						else {
							local tmethod = "`method'"
							local trim = 0
						}
						
						preserve
						keep if sample_`grade'_`outcome' == 1 & experiment_`sch' == 1 & disc_`sch' == 1
						keep if abs(runvar_`sch') <= `window'
						egen temp = rowmiss(`outcome' `controls1' `controls2' `controls3') 			
						keep if temp == 0
						drop temp
						
						ciest_fuzzy `controls1' `controls2' `controls3' `fixed_effects', outcome(`outcome') treatment(`treatment') sch(`sch') method(`tmethod') trim(`trim')
						local late0 = `r(late0)'
						local late1 = `r(late1)'
						local fs0 = `r(fs0)'
						local fs1 = `r(fs1)'
						local rf0 = `r(rf0)'
						local rf1 = `r(rf1)'
						local n0 = `r(n0)'
						local n1 = `r(n1)'
						
						save "$dir/Extrapolation/bstemp.dta", replace
						set seed 123
						mat bootstrap_est = J(1, 6, .)
						forval b = 1/500 {	
							bsample			
							quietly: ciest_fuzzy `controls1' `controls2' `controls3' `fixed_effects', outcome(`outcome') treatment(`treatment') sch(`sch') method(`tmethod') trim(`trim')
							mat bootstrap_est = (bootstrap_est\ (`r(late0)', `r(late1)', `r(fs0)', `r(fs1)', `r(rf0)', `r(rf1)'))	
							use "$dir/Extrapolation/bstemp.dta", clear
						}
						svmat bootstrap_est
						su bootstrap_est1
						local late0_se = `r(sd)'
						su bootstrap_est2
						local late1_se = `r(sd)'
						su bootstrap_est3
						local fs0_se = `r(sd)'
						su bootstrap_est4
						local fs1_se = `r(sd)'
						su bootstrap_est5
						local rf0_se = `r(sd)'
						su bootstrap_est6
						local rf1_se = `r(sd)'
						
						if ("`sch'" == "obr") {
							mat b = (`late0')
							mat V = (`late0_se'^2)
						}
						else {
							mat b = (`late1')
							mat V = (`late1_se'^2)
						}
						mat colnames b = est
						mat rownames b = est
						mat colnames V = est
						mat rownames V = est
						eret post b V
						estadd scalar n0  = `n0'
						estadd scalar n1  = `n1'
						eststo gr`grade'`outcome'`sch'`window'`method'late
						
						if ("`sch'" == "obr") {
							mat b = (`fs0')
							mat V = (`fs0_se'^2)
						}
						else {
							mat b = (`fs1')
							mat V = (`fs1_se'^2)
						}
						mat colnames b = est
						mat rownames b = est
						mat colnames V = est
						mat rownames V = est
						eret post b V
						estadd scalar n0  = `n0'
						estadd scalar n1  = `n1'
						eststo gr`grade'`outcome'`sch'`window'`method'fs
						
						if ("`sch'" == "obr") {
							mat b = (`rf0')
							mat V = (`rf0_se'^2)
						}
						else {
							mat b = (`rf1')
							mat V = (`rf1_se'^2)
						}
						mat colnames b = est
						mat rownames b = est
						mat colnames V = est
						mat rownames V = est
						eret post b V
						estadd scalar n0  = `n0'
						estadd scalar n1  = `n1'
						eststo gr`grade'`outcome'`sch'`window'`method'rf
						
						restore
					}
				}
			}
			local options = `"b(3) se(3) mtitle nogaps compress nonotes nonumbers star(* 0.1 ** 0.05 *** 0.01)"'
			esttab *`window'`method'fs *`window'`method'rf *`window'`method'late using "$dir/Extrapolation/Tables/estimates.csv", `options' keep(est) scalars(n0 n1) sfmt(0) noobs append
		}
	}
}

/*

/* Conditional Effect Ignorability estimates - sharp */

capture: program drop ceiest
program ceiest, rclass
	version 11.0
	syntax varlist(fv), outcome(string) sch(string) /*fe(string)*/
	
	noisily {
	
		foreach x of varlist `varlist' {
			disp "`x'"
			gen temp_above_`sch'_`x' = above_`sch'*`x'
			gen temp_runvar_`sch'_`x' = runvar_`sch'*`x'
			gen temp_intrunvar_`sch'_`x' = intrunvar_`sch'*`x'
		}
		

		reg `outcome' `varlist' if above_`sch' == 0
		local n0 = e(N)
		reg `outcome' `varlist' if above_`sch' == 1
		local n1 = e(N)
	
*		reg `outcome' above_`sch' runvar_`sch' intrunvar_`sch' `varlist' /*`fe'*/ temp_above_*
		reg `outcome' above_`sch' runvar_`sch' intrunvar_`sch' `varlist' `fe' temp_above_* temp_runvar_* 
*		reg `outcome' above_`sch' runvar_`sch' intrunvar_`sch' `varlist' `fe' temp_above_* temp_runvar_* temp_intrunvar_*

*		rdob_mod2 `outcome' runvar_`sch' if disc_`sch' == 1, c(0) dm(1)
*		local h_opt = r(h_opt)
*		gen lambda = (1 - abs(runvar_`sch'/`h_opt'))*(abs(runvar_`sch'/`h_opt') <= 1)
*		reg `outcome' above_`sch' runvar_`sch' intrunvar_`sch' `varlist' `fe' temp_above_* [aw = lambda]
*		reg `outcome' above_`sch' runvar_`sch' intrunvar_`sch' `varlist' `fe' temp_above_* temp_runvar_* [aw = lambda]
*		reg `outcome' above_`sch' runvar_`sch' intrunvar_`sch' `varlist' `fe' temp_above_* temp_runvar_* temp_intrunvar_* [aw = lambda]
*		drop lambda

		gen effect = _b[above_`sch']
		foreach x in `varlist' {
			replace effect = effect + _b[temp_above_`sch'_`x']*`x'
		}
	
		su effect if above_`sch' == 0
		local effect0 = r(mean)
		su effect if above_`sch' == 1
		local effect1 = r(mean)

		drop temp_* effect

	}
	return scalar effect0 = `effect0'
	return scalar effect1 = `effect1'
	return scalar n0 = `n0'
	return scalar n1 = `n1'
end

capture: tab risk_set, gen(risk_set_)
capture: tab app_prefs, gen(app_prefs_)
capture: tab yr_math10, gen(yr_math10_)
capture: tab yr_eng10, gen(yr_eng10_)

foreach grade in 7 9 {

	local controls1 = "female black white asian hisp"
	local controls2 = "free lep special"
	if (`grade' == 7) {
		local controls3 = "eng4 math4"
	}
	else {
		local controls3 = "eng4 math4 eng7 math8"	
	}
	
	est clear
	foreach outcome in math10 eng10 {
		foreach window in 20 /*15 10*/ {
			foreach sch in obr bls {
				noisily {
					local fixed_effects = "risk_set_* app_prefs_* yr_`outcome'_*"
					
					preserve
					keep if sample_`grade'_`outcome' == 1 & experiment_`sch' == 1 & disc_`sch' == 1
					keep if abs(runvar_`sch') <= `window'
					egen temp = rowmiss(`outcome' `controls1' `controls2' `controls3') 			
					keep if temp == 0
					drop temp
					
					ceiest `controls1' `controls2' `controls3' `fixed_effects', outcome(`outcome') sch(`sch') /*fe(`fixed_effects')*/
					local est0 = `r(effect0)'
					local est1 = `r(effect1)'
					local n0 = `r(n0)'
					local n1 = `r(n1)'
					
					ciest `controls1' `controls2' `controls3' `fixed_effects', outcome(`outcome') sch(`sch') method(kline) trim(0)
					local diff_cia_cei_est0 = `est0' - `r(effect0)'
					local diff_cia_cei_est1 = `est1' - `r(effect1)'
					drop pred*
					
					save "$dir/Extrapolation/bstemp.dta", replace
					set seed 123
					mat bs_est = J(1, 4, .)
					forval b = 1/500 {	
						bsample			
						
						ceiest `controls1' `controls2' `controls3' `fixed_effects', outcome(`outcome') sch(`sch') /*fe(`fixed_effects')*/
						local bs_est0 = `r(effect0)'
						local bs_est1 = `r(effect1)'
						
						ciest `controls1' `controls2' `controls3' `fixed_effects', outcome(`outcome') sch(`sch') method(kline) trim(0)
						local bs_diff_cia_cei_est0 = `est0' - `r(effect0)'
						local bs_diff_cia_cei_est1 = `est1' - `r(effect1)'
						drop pred*
						
						mat bs_est = (bs_est\ (`bs_est0', `bs_est1', `bs_diff_cia_cei_est0', `bs_diff_cia_cei_est1'))
						
						use "$dir/Extrapolation/bstemp.dta", clear
					}
					mat bs_est = bs_est[2..., 1...]
					svmat bs_est
					
					su bs_est1
					local se0 = `r(sd)'
					su bs_est2
					local se1 = `r(sd)'

					su bs_est3
					local diff_cia_cei_se0 = `r(sd)'
					su bs_est4
					local diff_cia_cei_se1 = `r(sd)'				
					
					if ("`sch'" == "obr") {
						mat b = (`est0')
						mat V = (`se0'^2)
						local tstat = `diff_cia_cei_est0'/`diff_cia_cei_se0'
					}
					else {
						mat b = (`est1')
						mat V = (`se1'^2)					
						local tstat = `diff_cia_cei_est1'/`diff_cia_cei_se1'
					}
					mat colnames b = est
					mat rownames b = est
					mat colnames V = est
					mat rownames V = est
					eret post b V
					estadd scalar n0  = `n0'
					estadd scalar n1  = `n1'
					estadd scalar tstat = `tstat'
					eststo ceiest_`grade'_`outcome'_`sch'_`window'
										
					restore
				}
			}
		}
	}
				
	local options = `"b(3) se(3) mtitle nogaps compress nonotes nonumbers star(* 0.1 ** 0.05 *** 0.01)"'

	esttab *_20 using "$dir/Extrapolation/Tables/estimates.csv", `options' keep(est) scalars(n0 n1 tstat) sfmt(0 0 3) noobs append
*	esttab *_15 using "$dir/Extrapolation/Tables/estimates.csv", `options' keep(est) scalars(n0 n1 tstat) sfmt(0 0 3) noobs append
*	esttab *_10 using "$dir/Extrapolation/Tables/estimates.csv", `options' keep(est) scalars(n0 n1 tstat) sfmt(0 0 3) noobs append
}

drop risk_set_* app_prefs_* yr_math10_* yr_eng10_*

*/


*/ ***HERE***


capture: log close