/*The model is D=I(v+x'B+e>0), E(ez)=0, v=s'A+eta, 
  eta homoskedastic normal, eta independent of e,s,x. 
  s is union of regressors x and instruments z*/ 

*.simpbin1 dv, bin() exog() very()
*set trace on

cap pro drop simpbin1
pro def simpbin1
version 8.0

local varlist "required existing min(2) max(2)"
local if "optional"
local options BIN(string) EXOG(string) VERY(string)"
parse "`*'"
                           
preserve

tokenize "`varlist'"                  
di in g "Regression of `2' on `very'"
qui reg `2' `very'

tempvar eta sigma
predict `eta' if e(sample), resid
gen `sigma'=sqrt(e(rss)/(_N))

tempvar f
gen `f'=normden(`eta', 0, `sigma')     /*assume v model error is normal,estimate density*/
zrnd `f'                               /* round near zero density away from zero*/
   
gen `1'star=(`1'-(`2'>=0))/`f'                          

di in g "                "
di in g "Based on conditional normal density of `2', f(`2'|.)"
di in g "                "
di in g "`1'star is constructed as (`1'-(`2'>=0))/ f(`2'|.)"

ivreg `1'star (`bin'=`exog'),robust
drop `1'star

di in g "Note the Std. Err. and Conf. Interval are wrong! "
di  "                              "
di in g "Bootstrap command can be used to get correct Std Err. and Conf. Interval."
di "                                      "
di in y `"Try it, :)) -  bs "simpbin1 varlist" _b, size() level(), reps() saving() replace"'

end