*! rd 2.9: package to implement regression discontinuity design *! author Austin Nichols * rd 2.9 30sep2016 fix undocumented bingraph option; use ineq implied by strineq and only graph 1 obs per bin * rd 2.8 9may2012 add check on OB calc failure e.g. r(2000) and divide Z by 10 until no failure * rd 2.7 14dec2011 add -mbnum- option to allow numlist of actual numeric values for bandwidths rather than multiples of default bandwidth * rd 2.6 12sep2011 fixed ddens option * rd 2.5 7sep2011 double reported bandwidth of rec kernel to match IK, add -covcoefdiff- option to not constrain covar() vars to have same coefs on either side of Z0 * rd 2.4 11jun2011 fixed strineq option, added rec kernel, added weights to new estim proc, and added cluster() and covar() options * rd 2.3 7jun2011 added e(w) for compatibility with rd_obs output (and help file) * rd 2.2 24may2011 added option strineq to assume X_T==1 iff Z>0 * rd 2.1 23may2011 switched back to X==1 iff Z>=0 * rd 2.0 20mar2011 shifted to weighted linear regression, -suest-, and IV instead of -lpoly- by default * rd moved to rd_obs (mar 2011) to reflect shift in estimation strategy * OB calc now (nov 2010) based on rdob.ado from http://www.economics.harvard.edu/faculty/imbens/software_imbens * rd 1.0.6 15oct2010 fixed bug where n+2 exceeds N in z grid * rd 1.0.5 25aug2010 fixed bug in noSCatter option * rd 1.0.4 17dec2009 added Imbens and Kalyanaraman approach to bw calc as default; http://www.nber.org/papers/w14726 * rd 1.0.3 15nov2007 fixed bug in graph command, made ddens optionally off * rd 1.0.2 5nov2007 first public release prog rd, eclass loc version "2.9 30sep2016" if replay() eret di else { loc vmult 1 loc lpoly "locpoly" loc k "tri" loc bw "w" if (_caller()<10) version 8.2 else { version 10 loc lpoly "lpoly" loc k "k(tri)" loc bw "bw" } syntax varlist(min=2 max=3) [if] [in] [aw fw pw/] [, s(string) Mult(real 1) GENerate(str) at(str) FORCElinear CLuster(varlist) /* */ n(integer 0) GRaph BINgraph binvar(varlist max=1) noSCatter Width(real 0) BWidth(real 0) Kernel(str) COVar(varlist numeric) /* */ defbw Degree(integer 1) mbw(numlist >0 integer) z0(real 0) lwald NOIsily SUPPresse format(string) STRineq COVCoefdiff /* */ Placebo DDens Tdisc x(varlist numeric) SCOpt(str) LINEopt(str) bdep bdopt(string) rcapopt(string) OXline level(real 95) MBNum(numlist >=0) *] marksample touse if !inrange(`level'*100,1,9999)|(mod(`level'*100,1)!=0) loc level 95 if `"`format'"'=="" loc format %4.3g conf format `format' * fix "above condition" and "below condition" depending whether strict ineq assumed or not if "`strineq'"=="" { loc abc ">=0" loc bec "<0" } else { loc abc ">0" loc bec "<=0" } tempname beta if "`bingraph'"!="" loc graph "graph" if "`graph'"!="" { cap `lpoly' * any Stata that does not recognize -lpoly- for graphing must download -locpoly- if _rc==199 { cap net from http://www.stata-journal.com/software/sj6-4 net inst st0053_3 } } if "`generate'`at'"!="" { di as error "Do not specify " as res "generate " as err "or " as res "at " as err "options." di as err "Instead use " as res "s(stubname) " as err "syntax to generate new variables." error 198 } if `bwidth'!=0 & `width'!=0 { di as error "Cannot specify both " as res "width " as err "and " as res "bwidth " as err "options." error 198 } if "`kernel'"=="" loc kernel="triangle (default)" if strpos("`kernel'","tri")==0&strpos("`kernel'","rec")==0 { di as res "rd" as err " allows only " as res "tri" as err "angle or " as res "rec" as err "tangle kernels; use " as res "rd_obs " as err "for old options." error 198 } if strpos("`kernel'","rec")>0 loc k: subinstr loc k "tri" "rec" if "`lpoly'"=="locpoly" & "`exp'"!="" { preserve tempvar wtvar dropme qui { gsort -`touse' g double `wtvar'=`exp'-`exp'[_n-1] if _n>1 & `touse' su `wtvar', meanonly replace `wtvar'=`exp'/r(min) loc N=_N expand `wtvar' g byte `dropme'=_n>`N' } } else if "`exp'"!="" loc wt="[aw=`exp']" if "`exp'"=="" loc exp 1 tokenize `varlist' loc was2="`2'" if "`3'"=="" { loc nolwald=1 loc 3="`2'" tempvar 2 qui g byte `2'=(`3'>=`z0') loc was2="" di as res "Two variables specified; treatment is " di as res "assumed to jump from zero to one at Z=`z0'. " di as res _n " Assignment variable Z is `3'" di as res " Treatment variable X_T unspecified" di as res " Outcome variable y is `1'" _n } else { di as res "Three variables specified; jump in treatment " di as res "at Z=`z0' will be estimated. Local Wald Estimate" di as res "is the ratio of jump in outcome to jump in treatment." di as res _n " Assignment variable Z is `3'" di as res " Treatment variable X_T is `2'" di as res " Outcome variable y is `1'" _n } tempvar z obs close loc was3="`3'" tempvar 3 qui g double `3'=(`was3'-`z0') qui compress `3' la var `3' "`: var lab `was3''" loc w=0 if `bwidth'!=0 loc w=`bwidth' else if `width'!=0 loc w=`width' loc divide 1 if `w'==0 & "`mbnum'"=="" { loc obfail 1 loc divide 1 while `obfail'==1 { qui { loc obfail 0 tempvar Y1 Y2 temp1 temp2 temp3 temp4 qui su `3' if `touse' `wt', d loc Sx = r(sd) loc N = r(N) loc h1 = 1.84*`Sx'*(`N'^(-1/5)) qui su `1' if -`h1'<=`3' & `3'`bec' & `touse' `wt', meanonly loc Yh1n = r(mean) loc Nh1n = r(N) qui su `1' if `3'`abc' & `3'<=`h1' & `touse' `wt', meanonly loc Yh1p = r(mean) loc Nh1p = r(N) loc fxc = (`Nh1p'+`Nh1n')/(2*`N'*`h1') g double `Y1' =(`1'-`Yh1n')^2 su `Y1' if -`h1'<=`3' & `3'`bec' & `touse' `wt', meanonly loc Y1sum =r(sum) g double `Y2' =(`1'-`Yh1p')^2 su `Y2' if `3'`abc' & `3'<=`h1' & `touse' `wt', meanonly loc Y2sum=r(sum) loc sigma2c=(`Y1sum'+`Y2sum')/(`Nh1p'+`Nh1n') su `3' if `3'`abc' & `touse' `wt', d scalar medXp = r(p50) loc Np=r(N) su `3' if `3'`bec' & `touse' `wt', d scalar medXn = r(p50) loc Nn=r(N) g byte `temp1'=(`3'`abc') if !mi(`3') g double `temp2'=`3' g double `temp3'=`3'^2 g double `temp4'=`3'^3 regress `1' `temp1' `temp2' `temp3' `temp4' if `3' >= scalar(medXn) & `3' <= scalar(medXp) & `touse' `wt' loc m3c=6*_b[`temp4'] loc h2p=3.56*((`sigma2c'/(`fxc'*max((`m3c')^2, 0.01)))^(1/7))*(`Np'^(-1/7)) loc h2n=3.56*((`sigma2c'/(`fxc'*max((`m3c')^2, 0.01)))^(1/7))*(`Nn'^(-1/7)) regress `1' `temp2' `temp3' if 0<=`3' & `3'<=`h2p' & `touse' `wt' loc m2pc=2*_b[`temp3'] loc N2p=e(N) regress `1' `temp2' `temp3' if -`h2n'<=`3' & `3'`bec' & `touse' `wt' loc m2nc=2*_b[`temp3'] loc N2n=e(N) loc rp=(720*`sigma2c')/(`N2p'*((`h2p')^4)) loc rn=(720*`sigma2c')/(`N2n'*((`h2n')^4)) if strpos("`k'","rec")>0 loc CK = 5.4/2 else loc CK = 3.4375 loc w= `CK'*(((2*`sigma2c')/(`fxc'*(((`m2pc'-`m2nc')^2)+(`rp'+`rn'))))^(1/5))*`N'^(-1/5) loc w100=`mult'*`w' loc repw = cond(strpos("`k'","rec")>0, `w100'*2, `w100') } if `obfail'==1 { qui replace `3'=`3'/10 loc divide=`divide'*10 } } } su `3' if `touse', meanonly loc minz=r(min) loc maxz=r(max) if `minz'>0 | `maxz'<0 { di as err _n "Assignment variable Z should have cutoff at `z0'" di as err "But range of assignment variable Z does not include `z0':" if `z0'!=0 su `was3' if `touse' else su `3' if `touse' error 198 } if "`mbw'"=="" & "`mbnum'"=="" loc mbw "100 50 200" else if "`mbnum'"=="" { loc j100 "100" loc mbw: list uniq mbw loc mbw: list mbw - j100 loc mbw "100 `mbw'" } if "`mbw'"=="" & "`mbnum'"!="" { gettoken w mbnum: mbnum loc mbw 100 foreach v in `mbnum' { loc mbw "`mbw' `=round(`v'/`w'*100)'" } loc w100=`mult'*`w' loc w=`mult'*`w' } if "`cluster'"!="" { tempvar clustv egen `clustv'=group(`cluster') loc clv cluster(`clustv') } if "`w100'"=="" loc w100=`mult'*`w' if "`graph'"!="" { di as res "Command used for graph: " as txt "`lpoly'" as res "; Kernel used: " as txt "`kernel'" if `n'==0 { * calculate steps for 50 pts on the right loc step=(`maxz')/50 if floor(-(`minz')/(`maxz'/50))>3 loc N=50+floor(-(`minz')/(`maxz'/50)) else loc N=54 cap set obs `N' qui replace `touse'=0 if mi(`touse') loc n=50 } else { * calculate steps for `n' pts loc step=(`maxz'-`minz')/(`n'-1) cap set obs `n' loc N=`n' * reset `n' to n pts on the right loc n=floor(`maxz'/`step') } if `w'==0 & "`defbw'"=="" { sort `3' g `close'=sum(`touse'*(`3'>0)) qui g long `obs'=_n if `touse' su `obs' if `close'>0, meanonly loc ub=`3'[`=r(min)+30'] loc lb=`3'[`=r(min)-30'] loc w=max(`ub',abs(`lb')) qui count if (`3'`abc' & `3'<`step'*2) & `touse' if r(N)>40 loc w=`step'*2 qui count if (`3'`bec' & `3'>-`step'*2) & `touse' if r(N)>40 loc w=`step'*2 } if "`w100'"=="" loc w100=`mult'*`w' qui { g `z'=(_n-1)*`step' in 1/`=`n'+1' replace `z'=-(_n-`=`n'+1')*`step' in `=`n'+2'/`N' } la var `z' "Assignment variable relative to cutoff" if "`binvar'"=="" loc binvar="`z'" loc opt "at(`z') nogr `k' deg(`degree') `options'" if "`defbw'"!="" { `lpoly' `1' `3' `wt' if `touse', gen(`f0') `opt' loc w100=`mult'*r(width) } foreach i of loc mbw { loc lw=`i'/100*`w100' if "`i'"=="100" loc i tempvar i`i'f0 i`i'f1 i`i'g0 i`i'g1 loc b="`bw'(`lw')" loc opt "at(`z') nogr `k' `b' deg(`degree') `options'" qui { if "`ddens'"!="" { tempvar df0 df1 kdensity `3' if `3'`abc' & `touse' `wt', at(`z') nogr `k' `b' gen(`df1') replace `df1'=. if `z'<0 kdensity `3' if `3'`bec' & `touse' `wt', at(`z') nogr `k' `b' gen(`df0') replace `df0'=. if `z'>0 su `touse' if `3'`abc' `wt' loc sumw1=r(sum_w) su `touse' if `3'`bec' `wt' loc sumw0=r(sum_w) replace `df1'=`df1'*`sumw1'/(`sumw1'+`sumw0') replace `df0'=`df0'*`sumw0'/(`sumw1'+`sumw0') loc dd`i'=`df1'[1]-`df0'[1] `noisily' di as txt "Bandwidth: " as res `lw' _c `noisily' di as txt "; Jump Estimate for Density: " _c `noisily' di as res (`df1'[1]-`df0'[1]) if "`graph'"!="" { tempvar kd kdensity `3' if `touse' `wt', at(`z') gen(`kd') nogr line `kd' `df1' `df0' `z', sort leg(off) name(ddens`i', replace) ti("Jump in Density of Assignment Variable" "Bandwidth=`lw'") } } foreach lhs in `x' `1' { `lpoly' `lhs' `3' `wt' if `3'`bec' & `touse', gen(`i`i'f0') `opt' replace `i`i'f0'=. if `z'>0 `lpoly' `lhs' `3' `wt' if `3'`abc' & `touse', gen(`i`i'f1') `opt' replace `i`i'f1'=. if `z'<0 loc numerat=(`i`i'f1'[1]-`i`i'f0'[1])*`divide' `noisily' di as txt "Bandwidth: " as res `lw' _c `noisily' di as txt "; Jump Estimate for `lhs': " _c `noisily' di as res (`i`i'f1'[1]-`i`i'f0'[1]) if "`lhs'"!="`1'" { cap loc rown: rownames `beta' mat `beta'=nullmat(`beta') \ `numerat' mat rownames `beta'= `rown' `lhs'`i' } if "`graph'"!="" { if "`scatter'"=="" { if "`bingraph'"!="" { tempvar only1perbin bys `binvar' (`3'): g byte `only1perbin'=_n==ceil(_N/2) tempvar y g `y'=. levelsof `binvar' if `touse', loc(bins) loc binT: word count `bins' forv bint=1/`binT' { loc lastbin loc nextbin cap loc lastbin: word `=`bint'-1' of `bins' if "`lastbin'"=="" loc lastbin=c(mindouble) cap loc nextbin: word `=`bint'+1' of `bins' if "`nextbin'"=="" loc nextbin=. su `lhs' `wt' if `3'`abc'`lastbin' & `3'`bec'`nextbin' &`touse' replace `y'=r(mean) if `3'`abc'`lastbin' & `3'`bec'`nextbin' &`touse'&`only1perbin' } } else loc y="`lhs'" loc sp "sc `y' `3' if `touse', mc(gs14) `scopt' ||" } tw `sp' line `i`i'f0' `i`i'f1' `z', lw(thick thick) lp(l -) leg(off) ti("`:var la `lhs''" "Bandwidth `lw'") name(`lhs'`i', replace) `lineopt' } if "`s'"!="" { ren `i`i'f0' `s'`i'`lhs'0 ren `i`i'f1' `s'`i'`lhs'1 } else drop `i`i'f0' `i`i'f1' } if "`lwald'"=="" & "`nolwald'"=="" { `lpoly' `2' `3' `wt' if `3'`bec' & `touse', gen(`i`i'g0') `opt' replace `i`i'g0'=. if `z'>0 `lpoly' `2' `3' `wt' if `3'`abc' & `touse', gen(`i`i'g1') `opt' replace `i`i'g1'=. if `z'<0 if "`graph'"!="" { if "`scatter'"=="" { if "`bingraph'"!="" { tempvar y g `y'=. levelsof `binvar' if `touse', loc(bins) loc binT: word count `bins' forv bint=1/`binT' { loc lastbin loc nextbin cap loc lastbin: word `=`bint'-1' of `bins' if "`lastbin'"=="" loc lastbin=c(mindouble) cap loc nextbin: word `=`bint'+1' of `bins' if "`nextbin'"=="" loc nextbin=. su `2' `wt' if `3'`abc'`lastbin' & `3'`bec'`nextbin' &`touse' replace `y'=r(mean) if `3'`abc'`lastbin' & `3'`bec'`nextbin' &`touse'&`only1perbin' } } else loc y="`2'" loc sp "sc `y' `3' if `touse', mc(gs14) `scopt' ||" } tw `sp' line `i`i'g0' `i`i'g1' `z', lw(thick thick) lp(l -) leg(off) ti("`:var la `2''" "Bandwidth `lw'") name(`2'`i', replace) `lineopt' } loc denomin=(`i`i'g1'[1]-`i`i'g0'[1])*`divide' if "`s'"!="" { ren `i`i'g0' `s'`i'`2'0 ren `i`i'g1' `s'`i'`2'1 } } else { loc denomin=1 if "`s'"!="" & "`was2'"!="" { g byte `i`i'g0'=0 if `z'<=0 g byte `i`i'g1'=1 if `z'`abc' ren `i`i'g0' `s'`i'`was2'0 ren `i`i'g1' `s'`i'`was2'1 } } if "`lwald'"=="" & "`nolwald'"=="" { * cap loc rown: rownames `beta' * mat `beta'=nullmat(`beta') \ `numerat' * mat rownames `beta'= `rown' numer`i' * cap loc rown: rownames `beta' * mat `beta'=nullmat(`beta') \ `denomin' * mat rownames `beta'= `rown' denom`i' } cap loc rown: rownames `beta' * mat `beta'=nullmat(`beta') \ `=`numerat'/`denomin'' * mat rownames `beta'= `rown' lwald`i' * cap loc rown: rownames `beta' if "`dd`i''"!="" { mat `beta'=nullmat(`beta') \ `dd`i'' mat rownames `beta'= `rown' ddens`i' } } loc w`i'=`lw' if strpos("`k'","rec")>0 loc w`i'=`lw'*2 di as txt "Bandwidth: " as res `lw' _c di as txt "; loc Wald Estimate: " _c di as res `numerat'/`denomin' } if "`s'"!="" ren `z' `s'`3' cap compress `s'* cap drop if `dropme' if "`suppresse'"!="" { matrix `beta'=`beta'' qui count if `touse' loc N=r(N) tempvar ee qui g byte `ee'=`touse' eret post `beta', esample(`ee') ereturn scalar N = `N' } else { * cap mat drop `beta' loc rown } } if "`suppresse'"=="" { tempname var diag tempvar above zbelow zabove g byte `above'=(`3'`abc'0) if `touse' g double `zabove'=(`above')*(`3') g double `zbelow'=(1-`above')*(`3') if "`covcoefdiff'"!="" & "`covar'"!="" { foreach v of loc covar { tempvar a`v' g double `a`v''=(`above')*(`v') g double `b`v''=(1-`above')*(`v') loc acovar `acovar' `a`v'' loc bcovar `bcovar' `b`v'' } } if "`covcoefdiff'"=="" & "`covar'"!="" { loc acovar `covar' loc bcovar } foreach i of loc mbw { loc w`i'=`i'/100*`w100' if "`i'"=="100" loc i noi di as res "Estimating for bandwidth `w`i''" tempvar kwt if strpos("`k'","rec")==0 g double `kwt'=max(0,`w`i''-abs(`3'))*(`exp') else g double `kwt'=(`exp')*((-`w`i''<=(`3'))&((`3')<`w`i'')) * double reported bandwidth for rect kernel to match IK parametrization if strpos("`k'","rec")>0 loc w`i'=`w`i''*2 foreach lhs in `x' `1' { qui if !("`lwald'"=="" & "`nolwald'"=="") { reg `lhs' `above' `zabove' `zbelow' `acovar' `bcovar' [pw=`kwt'], `clv' loc numerat=_b[`above'] loc var1=_se[`above']^2 mat `var'=nullmat(`var') \ `var1' cap loc rown: rownames `beta' mat `beta'=nullmat(`beta') \ `numerat' if "`lhs'"=="`1'" mat rownames `beta'= `rown' lwald`i' else mat rownames `beta'= `rown' `lhs'`i' cap loc rown: rownames `beta' } qui if "`lwald'"=="" & "`nolwald'"=="" { tempname v suv sub rb rv dtmp r1 r2 mat `rb'=J(1,3,.) mat `rv'=J(3,3,0) reg `lhs' `above' `zabove' `zbelow' `acovar' `bcovar' [aw=`kwt'] est sto `r1' reg `2' `above' `zabove' `zbelow' `acovar' `bcovar' [aw=`kwt'] su `2' [aw=`kwt'], mean loc bad1 0 if !inrange(_b[`above']+_b[_cons],r(min),r(max))|!inrange(_b[_cons],r(min),r(max)) { noi di as err "A predicted value of treatment at cutoff lies outside feasible range;" loc bad1=1 if "`forcelinear'"!="" { noi di as err "interpret local Wald estimates with extreme caution." } else { noi di as err "switching to local mean smoothing for treatment discontinuity." reg `2' `above' `acovar' `bcovar' [aw=`kwt'] } } est sto `r2' suest `r1' `r2', `clv' if _se[`above']==0 { di as err "using -ivreg- to estimate due to failure of -suest-" if `bad1'==0 { reg `lhs' `above' `zabove' `zbelow' `acovar' `bcovar' [pw=`kwt'], `clv' mat `rb'[1,1]=_b[`above'] mat `rv'[1,1]=_se[`above']^2 reg `2' `above' `zabove' `zbelow' `acovar' `bcovar' [pw=`kwt'], `clv' mat `rb'[1,2]=_b[`above'] mat `rv'[2,2]=_se[`above']^2 ivreg `lhs' (`2'=`above') `zabove' `zbelow' `acovar' `bcovar' [aw=`kwt'], `clv' mat `rb'[1,3]=_b[`2'] mat `rv'[3,3]=_se[`2']^2 } else { reg `lhs' `above' `3' `acovar' `bcovar' [pw=`kwt'] mat `rb'[1,1]=_b[`above'] mat `rv'[1,1]=_se[`above']^2 reg `2' `above' `3' `acovar' `bcovar' [pw=`kwt'] mat `rb'[1,2]=_b[`above'] mat `rv'[2,2]=_se[`above']^2 ivreg `lhs' (`2'=`above') `3' `acovar' `bcovar' [pw=`kwt'], `clv' mat `rb'[1,3]=_b[`2'] mat `rv'[3,3]=_se[`2']^2 } } else { mat `sub'=e(b) mat `rb'[1,1]=`sub'[1,"`r1'_mean:`above'"],`sub'[1,"`r2'_mean:`above'"] mat `suv'=e(V) mat `rv'[1,1]=`suv'["`r1'_mean:`above'","`r1'_mean:`above'"],`suv'["`r1'_mean:`above'","`r2'_mean:`above'"]\ `suv'["`r1'_mean:`above'","`r2'_mean:`above'"],`suv'["`r2'_mean:`above'","`r2'_mean:`above'"] nlcom [`r1'_mean]`above'/ [`r2'_mean]`above' mat `rb'[1,3]=r(b) mat `rv'[3,3]=r(V) } cap conf mat `diag' if _rc==0 { loc posit=rowsof(`diag')+1 loc siz=rowsof(`diag')+rowsof(`rv') mat `dtmp'=J(`siz',`siz',0) mat `dtmp'[1,1]=`diag' mat `dtmp'[`posit',`posit']=`rv' mat `diag'=`dtmp' } else { mat `diag'=`rv' } cap loc rown: rownames `beta' mat `beta'=nullmat(`beta') \ `rb'' if "`lhs'"=="`1'" mat rownames `beta'= `rown' numer`i' denom`i' lwald`i' else mat rownames `beta'= `rown' `lhs'numer`i' `lhs'denom`i' `lhs'`i' cap loc rown: rownames `beta' } } } if !("`lwald'"=="" & "`nolwald'"=="") mat `diag'=diag(`var') matrix `beta'=`beta'' mat rownames `diag'= `rown' mat colnames `diag'= `rown' mat `diag'=`vmult'*`diag' qui count if `touse' loc N=r(N) eret post `beta' `diag', esample(`touse') ereturn scalar N = `N' } ereturn loc depvar "`1'" ereturn loc cmd "rd" ereturn loc version "`version'" ereturn loc rdversion "`version'" foreach i of loc mbw { if "`i'"=="100" ereturn scalar w`i'=`w100' else ereturn scalar w`i'=`w`i'' } cap ereturn scalar w=`w100' eret di cap drop if `dropme' if "`bdep'"!="" { tempvar bwid est lbci ubci loc estn 1 qui { g `bwid'=. g `est'=. g `lbci'=. g `ubci'=. foreach i of local mbw { if "`oxline'"!="" & "`i'"=="100" loc ox "xli(`w100')" if "`i'"=="100" loc i replace `bwid'=`w`i'' in `estn' replace `est'=_b[lwald`i'] in `estn' replace `lbci'=_b[lwald`i']-invnormal(1/2-`level'/200)*_se[lwald`i'] in `estn' replace `ubci'=_b[lwald`i']+invnormal(1/2-`level'/200)*_se[lwald`i'] in `estn' loc xl `"`xl' `w`i'' "`:di `format' `w`i'''" "' loc estn=`estn'+1 } tw rcap `lbci' `ubci' `bwid', `rcapopt'||sc `est' `bwid', xla(`xl') leg(lab(1 "CI") lab(2 "Est")) xti(Bandwidth) yti(Estimated effect) `ox' `bdopt' } } } end