qui gen `touse' = e(sample)		

di "Table 4 loop.  MC Reps = `mcreps', BS Reps = `bsreps', Number of states sampled = `numstates'"

do_the_monte_carlo , statelist("`statelist'") sourcedata("`mysourcedata'") savedata("`mcout_`ns''")
		mcreps(`num_mc_reps') bsreps(`num_bs_reps') numstates(`ns')

qui save `savedata' , replace

syntax  , [mcreps(integer 3) bsreps(integer 4) numstates(integer 6) statelist(string) sourcedata(string) savedata(string)]

qui bootstrap _b _se , reps(`bsreps') nodots cluster(pseudo_cluster) idcluster(pseudo_cluster_2) saving(`pairs_bs_out' , double)
	: reg lnwage policy age age2 yrseduc  , cluster(pseudo_cluster_2)

/* calls program "table4_np_boot", which is defined in t34_programs.do.  That program re-does the "deviations from means" and then 
re-runs the regression, to get the betas and (clustered) standard errors. */
qui bootstrap b=r(mybeta) se=r(myse) , reps(`bsreps') nodots cluster(pseudo_cluster) idcluster(pseudo_cluster_2) 
	seed(10101) saving(`pairs_bs_out' , double)
	: table4_np_boot , lhs(`lhs') rhs(`rhs') keyrhs(`keyrhs') clusterid(pseudo_cluster_2)

	ereturn scalar df_m =              `df_m'

*	ereturn scalar rmse =              `rmse'
	ereturn scalar mss =               `mss'

	ereturn scalar r2_a =              `r2_a'

		ereturn scalar N_clus_`clusvar' =  `Gclusvar`c''

return scalar mybeta = `matabeta'

return scalar mybeta = `matabeta'

graph twoway (scatter p_rad_res p_webb_res if numstates == `gg') 
		(line p_webb_res p_webb_res if numstates == `gg') ,
		title(G=`gg') name(g_sc_`gg' , replace)

local se_ar1_clu_bs = r(sd)

local se_clu_bs = r(sd)

gen t_pairs_bs = ((b - `b_fe') / se)

gen t_pairs_bs = ((b - `b') / se)

/* calls program "table4_np_boot", which is defined in t34_programs.do.  That program re-does the "deviations from means" and then 
re-runs the regression, to get the betas and (clustered) standard errors. */
qui bootstrap b=r(mybeta) se=r(myse) , reps(`bsreps') nodots cluster(pseudo_cluster) idcluster(pseudo_cluster_2) 
	seed(10101) saving(`pairs_bs_out' , double)
	: table4_np_boot , lhs(`lhs') rhs(`rhs') keyrhs(`keyrhs') clusterid(pseudo_cluster_2)

di "se rob,"				_column(20) "-"  "," 					_column(35) "-"

di "se CR2,"				_column(20) "-"  "," 					_column(35) "-"

di "se pairs BS,"			_column(20) %8.4f `se_ar1_clu_bs' "," 	_column(35) %6.4f `p_bsse_Tg1'

di "se rob,"				_column(20) %8.4f `se_fe_rob' "," 		_column(35) %6.4f `p_rob_Tnk'

di "se CR2,"				_column(20) %8.4f `se_fe_CR2' "," 		_column(35) %6.4f `p_CR2_Tg1'

di "se pairs BS,"			_column(20) %8.4f `se_fe_clu_bs' "," 	_column(35) %6.4f `p_bsse_Tg1'

di "se rob,"		%8.4f `se_ols_robust'

di "se boot,"		%8.4f `se_ols_bootstrap'

di "se rob,"		"     -   "

di "se boot,"		%8.4f `se_re_bootstrap'

local se_ar1_clu_bs = r(sd)

local se_fe_clu = _se[`regrhs']

local se_fe_CR2 = `se_CR2'

local se_fe_rob = _se[`regrhs']

local se_ols_cluster = _se[policy] 		

local se_ols_robust = _se[policy] 			

local se_re_cluster = _se[policy]

local se_wild_rademacher_restricted = _se[policy]

local se_wild_rademacher_restricted = _se[`key_rhs']

local se_wild_rademacher_unrestricted = _se[`key_rhs']

local se_wild_webb_restricted = _se[`key_rhs']

local se_wild_webb_restricted = _se[`key_rhs']

		    di "          See Cameron, Gelbach & Miller, "

/* Random Effects model with cluster bootstrap for standard errors */
xtreg lnwage policy age age2 yrseduc , re vce(bootstrap , reps(399) seed(10101) )

/* calls program "table4_np_boot", which is defined in t34_programs.do.  That program re-does the "deviations from means" and then 
re-runs the regression, to get the betas and (clustered) standard errors. */
qui bootstrap b=r(mybeta) se=r(myse) , reps(`bsreps') nodots cluster(pseudo_cluster) idcluster(pseudo_cluster_2) 
	seed(10101) saving(`pairs_bs_out' , double)
	: table4_np_boot , lhs(`lhs') rhs(`rhs') keyrhs(`keyrhs') clusterid(pseudo_cluster_2)

/* now make a small subsample.  15% of 20% is 3% of original */
set seed 10102

set seed 10103

set seed `myseed'

set seed 10101

set seed 10104

		    di " -> NOTE: Raw estimated variance matrix was non-positive semi-definite."

do set_directory_macros

do set_directory_macros

do set_directory_macros

local shortcommand = subinstr("`command'","`lhs'","y_wild_rademacher_restricted",.)

local shortcommand = subinstr("`command'","`lhs'","y_wild_webb_restricted",.)

local shortcommand = subinstr("`command'","`lhs'_deviation","y_wild_rademacher_restricted",.)

local shortcommand = subinstr("`command'","`reglhs'","y_wild_rad_restricted_deviation",.)

		di "    If some variable in your options list contains the string 'robust', you will"

generate sort_order = uniform()

generate sort_order = uniform()

generate sort_order = uniform()

sort sort_order

sort sort_order

program define cgmreg, eclass byable(onecall) sortpreserve

/* assign treatment dummy to 1/2 of states (round down) */
gen sortvar = uniform()

sort sortvar

merge 1:m statefip using "`sourcedata'" , assert(match using) keep(match) nogenerate

/* generate state policy variable */
use statefip using "`sourcedata'" , replace

qui use "`sourcedata'" , replace

do_the_monte_carlo , sourcedata("`mysourcedata'") savedata("`mcout_`ns''")
		mcreps(`num_mc_reps') bsreps(`num_bs_reps') numstates(`ns')

/* get the dataset ready for estimation */
qui t3_setup , keylhs("lnwage_sy") sourcedata("${madedat}/CPS_panel") outdata(`main_data')

/* the t3_ols code will partial out state and year fixed effects */
t3_ols , bsreps(`bsreps') sourcedata(`main_data') lhs(lnwage_sy) rhs(policy) keyrhs(policy) yearonly

do_the_monte_carlo , statelist("`statelist'") sourcedata("`mysourcedata'") savedata("`mcout_`ns''")
		mcreps(`num_mc_reps') bsreps(`num_bs_reps') numstates(`ns')

syntax  , [mcreps(integer 3) bsreps(integer 4) numstates(integer 6) sourcedata(string) savedata(string)]