***Stata code for Stanford workshop
***Afternoon session on calculating bounds

*estimating marginal distributions (slide 22)
scalar p_10 = .35
scalar p_11 = .65
scalar p_00 = .49
scalar p_01 = .51

*estimating stratum membership distribution (slide 23)
scalar pi_not = p_10
scalar pi_aot = p_01
scalar pi_iot = 1 - p_10 - p_01

*reporting estimated stratum membership distribution
display pi_not
display pi_iot
display pi_aot

*estimting marginal outcomes (slide 29)
scalar theta_10 = .38
scalar theta_11 = .89
scalar theta_00 = .34
scalar theta_01 = .86

*estimating phi_1 (slide 26)
scalar phi_1 = pi_iot / (pi_iot + pi_aot)

*estimating phi_0 
scalar phi_0 = pi_not / (pi_not + pi_iot)

******************************************
*bounds on margin of theta_11 (slide 31)
******************************************

*setting possible values for theta_1_aot & calculating associated value of theta_1_iot
matrix theta_1_aot = J(101, 1, .)
matrix theta_1_iot = J(101, 1, .)


forvalues n = 1(1)101 {

	matrix theta_1_aot[`n',1] = (`n' - 1) * .01
	matrix theta_1_iot[`n',1] = (1/phi_1) * theta_11 - ((1 - phi_1) / phi_1) * theta_1_aot[`n',1]
	}
	
*check the matrices
matrix list theta_1_aot
matrix list theta_1_iot
	

*translating matrices to variables for analysis
svmat theta_1_aot, names(theta_1_aot)
svmat theta_1_iot, names(theta_1_iot)

*scatterplot of theta_1_iot and theta_1_aot
scatter theta_1_iot theta_1_aot


*scatterplot of theta_1_iot and theta_1_aot
scatter theta_1_iot theta_1_aot || scatter theta_1_iot theta_1_aot if inrange(theta_1_iot, 0, 1) & inrange(theta_1_aot, 0, 1), mcolor(red) legend(off)


*boundary values for theta_1_iot & theta_1_aot 
sum theta_1_iot theta_1_aot if inrange(theta_1_iot, 0, 1) & inrange(theta_1_aot, 0, 1)


***********************************************
*bounds on margin of theta_00 (slide 32 - 33)
***********************************************

***try it first yourself, but code at the bottom of the do file!

*****************************************************************************
*assuming monotonicity of outcome potential values in IOT stratum (slide 34)
*****************************************************************************

*scatterplot of theta_0_not and theta_0_iot (values in red ruled out)
scatter theta_0_not theta_0_iot || scatter theta_0_not theta_0_iot if inrange(theta_0_iot, 0, .49), mcolor(red)

*boundary values for theta_0_iot & theta_0_not
sum theta_0_not theta_0_iot if inrange(theta_0_iot, 0, .49)



*****************************************************************************
*exploring joint distributions of treatment effects (slide 35)
*****************************************************************************
matrix pce_not = J(52, 50, .)
matrix pce_iot = J(52, 50, .)
matrix pce_aot = J(52, 50, .)

matrix iot_1 = J(52, 1, .)
matrix iot_0 = J(50, 1, .)

*setting possible values for theta_iot_1
forvalues n = 1(1)52 {

	matrix iot_1[`n',1] = .49 + (`n' - 1)*.01
	
	}
*setting possible values for 
forvalues n = 1(1)50 {

	matrix iot_0[`n',1] = (`n' - 1)*.01 
	
	}

*estimating possible values of pce_iot
forvalues i = 1(1)52 {

	forvalues j = 1(1)50 {
	
	*pce_iot
	matrix pce_iot[`i',`j'] = iot_1[`i',1] - iot_0[`j',1]
	
	*pce_not
	matrix pce_not[`i',`j'] = -.096 + .4 * iot_0[`j',1]
	 
	*pct_aot
	matrix pce_aot[`i',`j'] = .274 - .275 * iot_1[`i',1]
	
	}
	}

*setting matsize
set matsize 11000

*make the matrix a vector
matrix pce_iot_v = vec(pce_iot)
matrix pce_not_v = vec(pce_not)
matrix pce_aot_v = vec(pce_aot)

*translate vectors to data columns for analysis
svmat pce_aot_v, names(pce_aot)
svmat pce_iot_v, names(pce_iot)
svmat pce_not_v, names(pce_not)

*calculating descriptive statistics
sum pce*

*calculating descriptive statistics under assumption that pce_iot > pce_not
sum pce* if pce_iot > pce_not



























































***********************************************
*bounds on margin of theta_00 (slide 32 - 33)
***********************************************

*setting possible values for theta_0_iot & calculating associated value of theta_0_not
matrix theta_0_iot = J(101, 1, .)
matrix theta_0_not = J(101, 1, .)

forvalues n = 1(1)101 {

	matrix theta_0_iot[`n',1] = (`n' - 1) * .01
	matrix theta_0_not[`n',1] = (1 / phi_0) * theta_00 - ((1 - phi_0) / phi_0) * theta_0_iot[`n',1]
	}

svmat theta_0_iot, names(theta_0_iot)
svmat theta_0_not, names(theta_0_not)


*scatterplot of theta_0_not and theta_0_iot
scatter theta_0_not theta_0_iot

*boundary values for theta_0_iot & theta_0_not
sum theta_0_not theta_0_iot