/* fig3_1.do

   programs required
       binscrn2.ado  (v 1.1.0 or later)
       bvnellip.ado (v 1.1.0 or later)

   data required: hpns.dta
       hypothetical prenatal screening data

   last update:  20 July 2005

*/
version 8
set more off
cap log close

use ${semt_data}hpns,clear
qui expand nobs

/* get the estimates for ln(rDLR+), lnr(rDLR-) and covariance matrix with -binscrn2- */
quietly binscrn2 y d arm

matrix Betaln = r(ln_rDLR)
matrix Cov = r(Cov_ln_rDLR)

/* generate points for confidence ellipse with -bvnellip- */
bvnellip, plotv(elpos elneg) cov(Cov) beta(Betaln) df(10000) sortv(psort) np(500)
gen elpos2 = exp(elpos)
gen elneg2 = exp(elneg)

sort psort
gen lnbetapos = Betaln[1,1] if _n==1
gen lnbetaneg = Betaln[2,1] if _n==1

#delimit ;
line elneg elpos ||
   scatter lnbetaneg lnbetapos, ms(+) msize(medlarge)
   xscale(r(-.4 2.8) titlegap(3)) xlabel(0 1 2)  xtick(-.4(.2)2.8)
   yscale(r(-2.4 .8) titlegap(3)) ylabel(-2 -1 0) ytick(-2.4(.2).8)
   title(Figure 3.1, position(11) justification(left) margin(0 0 2 0)
          bmargin(0 0 0 0) size(3.8))
   subtitle("95% confidence region for ln(rDLR+),ln(rDLR-)",
              justification(left) position(11) margin(0 0 4 0) size(3.2))
   xtitle("ln(rDLR+)", size(3.6))
   ytitle("ln(rDLR-)", size(3.6))
   plotregion(margin(zero))
   aspectratio(1, place(right))
   xsize(3.0) ysize(3.5) fysize(105)
   legend(off)
   ;