capture program drop switch_oprobit version 10 program define switch_oprobit, eclass sortpreserve *!version 1.0.0 cagregory *!version 1.1.0 update to fix bugs with standard errors of treatment fx if replay() { if ("`e(cmd)'" != "switch_oprobit") error 301 Replay `0' } else Estimate `0' end program define Estimate, eclass syntax anything(equalok) [if] [in] [pweight iweight fweight], treat(varlist) /// [Robust CLuster(varname) level(integer 95)] *parse gettoken depvar2 indvar2: anything gettoken depvar1 indvar1: treat *local indvar2 "`indvar2' `depvar1'" *options mlopts mlopts, `options' //ml options if "`cluster'" != "" { local clopt cluster(`cluster') } if "`weight'" != "" { tempvar wvar quietly gen double `wvar' `exp' local wgt "[`weight'=`wvar']" local awgt "[aw=`wvar']" } if "`level'"!="" { tempname cilevel global cilevel = `level' } *mark sample marksample touse markout `touse' `depvar1' markout `touse' `depvar2' markout `touse' `indvar1' markout `touse' `indvar2' markout `touse' `cluster' *get startvals qui levelsof(`depvar2') if `touse' local t: word count `r(levels)' global ncut = `t'-1 *treatment eq: probit tempname b_trt startmat ystar di in gr _newline "Estimating Treatment equation" probit `depvar1' `indvar1' `wgt' if `touse', nocoef qui predict `ystar', xb mat `b_trt' = e(b) *outcome eq: ordered probit, treated di in gr _newline "Estimating Ordered Outcome Equation for Treatment Group" oprobit `depvar2' `indvar2' `wgt' if `depvar1' & `touse' , nocoef tempname b b_out_t cuts_t zstar_t qui predict `zstar_t', xb mat `b' = e(b) local j = colsof(`b') mat `b_out_t' =`b'[1,1..`j'-$ncut] *cutpoints mat `cuts_t'= `b'[1,(`j'-$ncut)+1..`j'] *outcome eq: ordered probit, untreated di in gr _newline "Estimating Ordered Outcome Equation for Untreated Group" oprobit `depvar2' `indvar2' `wgt' if `depvar1'==0 & `touse' , nocoef tempname b_u b_out_u cuts_u zstar_u qui predict `zstar_u', xb mat `b_u' = e(b) local j = colsof(`b_u') mat `b_out_u' =`b_u'[1,1..`j'-$ncut] *cutpoints mat `cuts_u'= `b_u'[1,(`j'-$ncut)+1..`j'] *rho qui corr `zstar_t' `ystar' local rho1_init = min(r(rho),.5) qui corr `zstar_u' `ystar' local rho0_init = min(r(rho),.5) *starting values mat `startmat' = `b_trt', `b_out_u', `b_out_t', `rho0_init', `rho1_init', `cuts_u', `cuts_t' *mat list `startmat' local cutpts1 local cutpts0 forvalues i = 1/$ncut { local cutpts1 "`cutpts1' /cut_1`i'" local cutpts0 "`cutpts0' /cut_0`i'" } qui levelsof(`depvar2') global nchoices: word count `r(levels)' *displaying ancilliary parameters local athrho diparm(atanh_rho0) local athrho "`athrho' diparm(atanh_rho1) diparm(__sep__)" local rho diparm(atanh_rho0, tanh label("rho0")) local rho "`rho' diparm(atanh_rho1, tanh label("rho1")) diparm(__sep__)" local di_cuts1 local di_cuts0 forv i = 1/$ncut { local di_cuts1 "`di_cuts1' diparm(cut_1`i')" local di_cuts0 "`di_cuts0' diparm(cut_0`i')" } local di_cuts "`di_cuts0' `di_cuts1' diparm(__sep__)" *full model di in gr _newline "Estimating Full Model" ml model lf0 switchoprobit_work ("`depvar1'": `depvar1' = `indvar1') /// ("`depvar2'_0": `depvar2' = `indvar2', noconstant) /// ("`depvar2'_1": `depvar2' = `indvar2', noconstant) /// /atanh_rho0 /atanh_rho1 `cutpts0' `cutpts1' if `touse' /// `wgt' /// , /// title("Ordered Probit Switching Regression") /// init(`startmat', copy) /// search(off) /// technique(bhhh) /// nonrtol /// maximize /// `clopt' /// `mlopts' /// `robust' /// `di_cuts' /// `athrho' /// `rho' *ml init `startmat', copy *ml check *ml maximize, nonrtol difficult ereturn scalar cuts = $ncut ereturn scalar k_aux = 2*($ncut+2) ereturn local treatx "`indvar1'" ereturn local rhsout "`indvar2'" ereturn local cmd "switch_oprobit" ereturn local depvar "`depvar1' `depvar2'_0 `depvar2'_1 " ereturn local predict "switch_oprobit_p" *ereturn display switchoprobit_replay end //end parsing and calling command program define switchoprobit_replay ml display, level($cilevel) neq(3) noomitted tempname lrtest pval qui test _b[atanh_rho1:_cons]=_b[atanh_rho0:_cons] local `lrtest' = round(r(chi2),.01) local `pval' = round(r(p),.0001) di in gr `"{p} LR test of distinct regimes = ``lrtest''. Probability of identical regimes for treated and untreated states ``pval''.{p_end}"' end //end parsing and calling command