*! plausexog: Estimating bounds with a plausibly exogenous exclusion restriction *! Version 3.1.2 July 01, 2020 @ 15:03:27 *! Author: Damian Clarke (application of code and ideas of Conley et al., 2012) *! Much of the heart of this code comes from the Conley et al implementation *! Contact: damian.clarke@usach.cl /* version highlights: 0.1.0: UCI and LTZ working. No graphing 0.2.1: Graphing added 1.0.0: Completed beta testing. Output complete. SSC v 1.0.0 1.0.1: Minor change to graph label options. 1.1.0: Weighting incorporated 1.2.0: Bug fix: very long names passed through syntax 2.0.0: Now Allowing for all arbitrary distributions with simulation algorithm 2.0.1: Graph issue when direct effect is negative for UCI (github issue #2) 2.1.0: Updating for Stata version 11.0 (requires mvtnorm) 3.0.0: Allows use in IC and Small Stata with simulations. Simplification of syntax. 3.1.0: Bug fix on confidence interval calculation (norm vs t) Chapman/Batini 3.1.2: Dropping observations with missing X or Z for LTZ (comment TH/LT) */ cap program drop plausexog program plausexog, eclass version 11.0 #delimit ; syntax anything(name=0 id="variable list") [if] [in] [fweight pweight aweight iweight] [, grid(real 2) gmin(numlist) gmax(numlist) level(real 0.95) omega(numlist min=1) mu(numlist min=1) GRAph(varlist) GRAPHOMega(numlist min=2 max=22) graphmu(numlist min=2 max=22) graphdelta(numlist) VCE(string) DISTribution(string) seed(numlist min=1 max=1) iterations(integer 5000) * ] ; #delimit cr preserve ******************************************************************************** *** (1) Unpack arguments, check for valid syntax, general error capture ******************************************************************************** local 0: subinstr local 0 "=" " = ", count(local equal) local 0: subinstr local 0 "(" " ( ", count(local lparen) local 0: subinstr local 0 ")" " ) ", count(local rparen) tokenize `0' local method `1' macro shift if "`method'"!="uci"&"`method'"!="ltz"&"`method'"!="upwci" { dis as error "Method of estimation must be specified." dis "Re-specify using uci, ltz or upcwi (see help file for more detail)" exit 200 } if `equal'!=1|`lparen'!=1|`rparen'!=1 { dis as error "Specification of varlist is incorrect." dis as error "Ensure that syntax is: method yvar [exog] (endog=iv), [opts]" exit 200 } if `level'<=0|`level'>=1 { dis as error "Confidence level was requested as `level'" dis as error "The confidence level must be between 0 and 1. default is level(0.95)" exit 200 } local yvar `1' macro shift local varlist1 while regexm(`"`1'"', "\(")==0 { local varlist1 `varlist1' `1' macro shift } local varlist2 while regexm(`"`1'"', "=")==0 { local var=subinstr(`"`1'"', "(", "", 1) local varlist2 `varlist2' `var' macro shift } local varlist_iv while regexm(`"`1'"', "\)")==0 { local var=subinstr(`"`1'"', "=", "", 1) local varlist_iv `varlist_iv' `var' macro shift } foreach list in varlist1 varlist2 varlist_iv { fvexpand ``list'' local `list' `r(varlist)' } local allout `varlist1' `varlist2' constant local allexog `varlist1' `varlist_iv' constant local count1 : word count `varlist1' local count2 : word count `varlist2' local count_iv : word count `varlist_iv' local count_all : word count `allout' local count_exog : word count `allexog' local countmin : word count `gmin' local countmax : word count `gmax' if `count2'>`count_iv' { dis as error "Specify at least as many instruments as endogenous variables" exit 200 } if "`method'"=="uci" { if `countmin'!=`count_iv'|`countmax'!=`count_iv' { dis as error "You must define as many gamma values as instrumental variables" dis "If instruments are believed to be valid, specify gamma=0 for gmin and gmax" exit 200 } foreach item in min max { local count=1 foreach num of numlist `g`item'' { local g`count'`item'=`num' local ++count } } } if "`method'"=="ltz" & length("`distribution'")==0 { if length("`omega'")==0|length("`mu'")==0 { dis as error "For ltz, omega and mu values must be provided" exit 200 } local count_omega: word count `omega' local count_mu : word count `mu' if `count_omega'!=`count_iv' { dis as error "Variance assumptions for each plausibly exogenous IV should be provided" dis as error "`count_iv' IVs are included, so ensure that omega has `count_iv' elements" exit 200 } if `count_mu'!=`count_iv' { dis as error "Mean assumptions for each plausibly exogenous IV should be provided" dis as error "`count_iv' IVs are included, so ensure that mu has `count_iv' elements" exit 200 } } local derr1 "Simulation-based estimates require the user-written ado mvtnorm." local derr2 "To use this method, please first install mvtnorm from the SSC" local derr3 "(ssc install mvtnorm)." if length("`distribution'")!=0 { cap which rmvnormal if _rc!= 0 { dis as error "`derr1' `derr2' `derr3'" error 200 } if "`method'"!="ltz" { dis as error "The distribution option can only be specified with ltz" error 200 } if length("`omega'")!=0|length("`mu'")!=0 { dis as error "Omega and mu should not be used with the distribution option" error 200 } local distribution: subinstr local distribution "," " , ", all local distcnt : list sizeof distribution local jj=1 foreach j of numlist 1(1)`distcnt' { local dist`jj': word `j' of `distribution' local dist`jj': subinstr local dist`jj' "," "" if "`dist`jj''"!="" local ++jj } local expD = 2+2*`count_iv' local expS = 2+1*`count_iv' local cntD = 2*`count_iv' local cntS = 1*`count_iv' local derr1 "If specifying a distribution with" local derrb "and `count_iv' instrument(s)" local derr2 "parameter(s) must be specified" local accept "normal, uniform, chi2, poisson, t, gamma, special" if "`dist1'"=="normal" { if `jj'!=`expD' { dis as error "`derr1' normal `derrb', `cntD' `derr2' (mean and std dev)." exit 200 } foreach ivn of numlist 1(1)`count_iv' { local ivl = `ivn'*2 local ivh = `ivn'*2+1 local gammaCall`ivn' rnormal(`dist`ivl'', `dist`ivh'') } } else if "`dist1'"=="uniform" { if `jj'!=`expD' { dis as error "`derr1' uniform `derrb', `cntD' `derr2' (min and max)." exit 200 } foreach ivn of numlist 1(1)`count_iv' { local ivl = `ivn'*2 local ivh = `ivn'*2+1 local gammaCall`ivn' `dist`ivl''+(`dist`ivh''-`dist`ivl'')*runiform() } } else if "`dist1'"=="chi2" { if `jj'!=`expS' { dis as error "`derr1' chi2 `derrb', `cntS' `derr2' (degrees of freedom)." exit 200 } foreach ivn of numlist 1(1)`count_iv' { local ivl = `ivn'+1 if `dist`ivl''<1 { dis as error "At least 1 degree of freedom must be specified for chi2" exit 200 } local gammaCall`ivn' rchi2(`dist`ivl'') } } else if "`dist1'"=="poisson" { if `jj'!=`expS' { dis as error "`derr1' poisson `derrb', `cntS' `derr2' (distribution mean)." exit 200 } foreach ivn of numlist 1(1)`count_iv' { local ivl = `ivn'+1 if `dist`ivl''<1 { dis as error "At least a mean of 1 must be specified for poisson" exit 200 } local gammaCall`ivn' rpoisson(`dist`ivl'') } } else if "`dist1'"=="t" { if `jj'!=`expS' { dis as error "`derr1' t `derrb', `cntS' `derr2' (degrees of freedom)." exit 200 } foreach ivn of numlist 1(1)`count_iv' { local ivl = `ivn'+1 if `dist`ivl''<1 { dis as error "At least 1 degree of freedom must be specified for t" exit 200 } local gammaCall`ivn' rt(`dist`ivl'') } } else if "`dist1'"=="gamma" { if `jj'!=`expD' { dis as error "`derr1' gamma `derrb', `cntD' `derr2' (shape and scale)." exit 200 } foreach ivn of numlist 1(1)`count_iv' { local ivl = `ivn'*2 local ivh = `ivn'*2+1 if `dist`ivl''<=0|`dist`ivh''<=0 { dis as error "The shape and scale parameter for gamma must be > 0" exit 200 } local gammaCall`ivn' rgamma(`dist`ivl'',`dist`ivh'') } } else if "`dist1'"=="special" { local sperr1 "To define your own distribution you must specify" local sperr2 "valid variable with the empirical distribution for each IV" if `jj'!=`expS' { dis as error "`sperr1' one `sperr2'" exit 200 } foreach ivn of numlist 1(1)`count_iv' { local ivl = `ivn'+1 cap sum `dist`ivl'' if _rc!=0 { dis as error "`sperr1' a `sperr2'" exit 200 } } } else { dis as error "The distribution option can only specify: `accept'" error 200 } } dis "Estimating Conely et al.'s `method' method" dis "Exogenous variables: `varlist1'" dis "Endogenous variables: `varlist2'" dis "Instruments: `varlist_iv'" ******************************************************************************** *** (2) Estimate model under assumption of gamma=0 ******************************************************************************** local cEx : word count `varlist1' local cEn : word count `varlist2' local cIV : word count `varlist_iv' *dis "Trying to run intial reg with `cEx' exog vars `cEn' endog and `cIV' insts" qui ivregress 2sls `yvar' `varlist1' (`varlist2'=`varlist_iv') `if' `in' /* */ [`weight' `exp'], vce(`vce') qui estimates store __iv if length("`graph'")!=0 { local CI = -invnormal((1 - `level')/2) local lcomp = _b[`graph'] - `CI'*_se[`graph'] local ucomp = _b[`graph'] + `CI'*_se[`graph'] } ******************************************************************************** *** (3) Union of Confidence Intervals approach (uci) *** Here we are creating a grid and testing for each possible gamma combo: *** ie - {g1min,g2min,g3min}, {g1max,g2min,g3min}, ..., {g1max,g2max,g3max} *** This conserves much of the original (Conley et al.) code, which does *** this in quite a nice way. ******************************************************************************** if "`method'"=="uci" { local points=1 foreach gnum of numlist 1(1)`points' { local cRatio = `gnum'/`points' foreach item in min max { local count=1 foreach num of numlist `g`item'' { local g`count'`item'l=`num'*`cRatio' local ++count } } ************************************************************************** *** (3a) Iterate over iter, which is each possible combination of gammas ************************************************************************** local iter=1 while `iter' <= (`grid'^`count_iv') { local R=`iter'-1 local w=`count_iv' **Create weighting factor to grid gamma. If grid==2, gamma={max,min} while `w'>0 { local a`w' = floor(`R'/(`grid'^(`w'-1))) local R = `R'-(`grid'^(`w'-1))*`a`w'' local gamma`w' = `g`w'minl' + ((`g`w'maxl'-`g`w'minl')/(`grid'-1))*`a`w'' local --w } tempvar Y_G qui gen `Y_G'=`yvar' local count=1 foreach Z of local varlist_iv { qui replace `Y_G'=`Y_G'-`Z'*`gamma`count'' local ++count } *********************************************************************** *** (3b) Estimate model based on assumed gammas, memoize conf intervals *********************************************************************** qui ivregress 2sls `Y_G' `varlist1' (`varlist2'=`varlist_iv') `if' /* */ `in' [`weight' `exp'], vce(`vce') *********************************************************************** *** (3c) Check if variable is not dropped (ie dummies) and results *********************************************************************** mat b2SLS = e(b) mat cov2SLS = e(V) local vars_final local counter=0 foreach item in `e(exogr)' `e(instd)' _cons { if _b[`item']!=0|_se[`item']!=0 { local vars_final `vars_final' `item' local ++counter } } ereturn scalar numvars = `counter' mat b2SLSf = J(1,`counter',.) mat se2SLSf = J(`counter',`counter',.) tokenize `vars_final' foreach num of numlist 1(1)`counter' { mat b2SLSf[1,`num']=_b[``num''] mat se2SLSf[`num',`num']=_se[``num''] } mat CI = -invnormal((1 - `level')/2) mat ltemp = vec(b2SLSf) - CI*vec(vecdiag(se2SLSf)) mat utemp = vec(b2SLSf) + CI*vec(vecdiag(se2SLSf)) *********************************************************************** *** (3d) Check if CI from model is lowest/highest in union (so far) *********************************************************************** foreach regressor of numlist 1(1)`counter' { if `iter'==1 { local l`regressor'=. local u`regressor'=. } local l`regressor' = min(`l`regressor'',ltemp[`regressor',1]) local u`regressor' = max(`u`regressor'',utemp[`regressor',1]) } local ++iter } if `gnum'==`points' { dis in yellow _newline dis "Conley et al (2012)'s UCI results" _col(55) "Number of obs = " e(N) dis in yellow in smcl "{hline 78}" dis "Variable" _col(13) "Lower Bound" _col(29) "Upper Bound" dis in yellow in smcl "{hline 78}" tokenize `vars_final' foreach regressor of numlist 1(1)`counter' { if `l`regressor''>`u`regressor''{ dis in green "``regressor''" _col(13) `u`regressor'' _col(29) `l`regressor'' } else { dis in green "``regressor''" _col(13) `l`regressor'' _col(29) `u`regressor'' } foreach var of local varlist2 { if `"`var'"'==`"``regressor''"' { if `l`regressor''>`u`regressor''{ ereturn scalar lb_`var'=`u`regressor'' ereturn scalar ub_`var'=`l`regressor'' } else { ereturn scalar lb_`var'=`l`regressor'' ereturn scalar ub_`var'=`u`regressor'' } } } } dis in yellow in smcl "{hline 78}" } } if length("`graph'")!=0 { if `count_iv'>1 { dis as error "Graphing with UCI only supported with 1 plausibly exogenous variable" exit 200 } local points=7 matrix __graphmat = J(`points',4,.) if (`g1max'<=0&`g1min'<=0)|(`g1max'>=0&`g1min'>=0) { local range = `g1max'-`g1min' local jj=0 while `jj'<=`points'-1 { local cRatio = (`jj'/(`points'-1))*`range' local gammaC = `cRatio'+`g1min' local ++jj tempvar Y_G qui gen `Y_G'=`yvar' qui replace `Y_G'=`Y_G'-`varlist_iv'*`gammaC' qui ivregress 2sls `Y_G' `varlist1' (`varlist2'=`varlist_iv') `if' /* */ `in' [`weight' `exp'], vce(`vce') local CI = -invnormal((1 - `level')/2) local ltemp = _b[`graph'] - `CI'*_se[`graph'] local utemp = _b[`graph'] + `CI'*_se[`graph'] if `jj'==1 { local l`regressor'=`lcomp' local u`regressor'=`ucomp' } local l`regressor' = min(`l`regressor'',`ltemp') local u`regressor' = max(`u`regressor'',`utemp') matrix __graphmat[`jj',1]=`gammaC' matrix __graphmat[`jj',3]=`l`regressor'' matrix __graphmat[`jj',4]=`u`regressor'' } } else { local range = `g1max'-0 local jj=3 while `jj'<=`points'-1 { local cRatio = (`jj'-3)/(`points'-4)*`range' local gammaC = `cRatio'+0 local ++jj tempvar Y_G qui gen `Y_G'=`yvar' qui replace `Y_G'=`Y_G'-`varlist_iv'*`gammaC' qui ivregress 2sls `Y_G' `varlist1' (`varlist2'=`varlist_iv') `if' /* */ `in' [`weight' `exp'], vce(`vce') local CI = -invnormal((1 - `level')/2) local ltemp = _b[`graph'] - `CI'*_se[`graph'] local utemp = _b[`graph'] + `CI'*_se[`graph'] if `jj'==4 { local l`regressor'=`lcomp' local u`regressor'=`ucomp' } local l`regressor' = min(`l`regressor'',`ltemp') local u`regressor' = max(`u`regressor'',`utemp') matrix __graphmat[`jj',1]=`gammaC' matrix __graphmat[`jj',3]=`l`regressor'' matrix __graphmat[`jj',4]=`u`regressor'' } local range = `g1min'-0 local jj=0 while `jj'<=`points'-5 { local cRatio = (`jj'+1)/(`points'-4)*`range' local gammaC = 0+`cRatio' local ++jj tempvar Y_G qui gen `Y_G'=`yvar' qui replace `Y_G'=`Y_G'-`varlist_iv'*`gammaC' qui ivregress 2sls `Y_G' `varlist1' (`varlist2'=`varlist_iv') `if' /* */ `in' [`weight' `exp'], vce(`vce') local CI = -invnormal((1 - `level')/2) local ltemp = _b[`graph'] - `CI'*_se[`graph'] local utemp = _b[`graph'] + `CI'*_se[`graph'] if `jj'==1 { local l`regressor'=`lcomp' local u`regressor'=`ucomp' } local l`regressor' = min(`l`regressor'',`ltemp') local u`regressor' = max(`u`regressor'',`utemp') matrix __graphmat[`jj',1]=`gammaC' matrix __graphmat[`jj',3]=`l`regressor'' matrix __graphmat[`jj',4]=`u`regressor'' } } } } ******************************************************************************** *** (5) Local to Zero approach (ltz) ******************************************************************************** if "`method'"=="ltz" { tempvar const qui gen `const'=1 if length("`distribution'")!=0&length("`graph'")!=0 { dis as error "Graphing is not enabled with the distribution option" dis as error "For graphing and LTZ, omega and mu must be used" exit 200 } ***************************************************************************** *** (5a) Remove any colinear elements to ensure that matrices are invertible *** For the case of the Z vector, this requires running the first stage regs ***************************************************************************** if `count1'!=0 { unab testvars: `varlist1' local usevars1 foreach var of local testvars { cap dis _b[`var'] if _rc!=0 continue if _b[`var']!=0 local usevars1 `usevars1' `var' } } unab testvars: `varlist2' local usevars2 foreach var of local testvars { cap dis _b[`var'] if _rc!=0 continue if _b[`var']!=0 local usevars2 `usevars2' `var' } foreach var of varlist `yvar' `varlist_iv' `usevars2' `usevars1' { qui drop if `var'==. } ***************************************************************************** *** (5b) Form moment matrices: ZX and ZZ ***************************************************************************** mat ZX = J(1,`count_all',.) mat ZZ = J(1,`count_exog',.) unab allvars: `varlist_iv' `usevars1' `const' tokenize `allvars' mat vecaccum a = `1' `usevars2' `usevars1' `if' `in' mat ZX = a while length("`2'")!= 0 { mat vecaccum a = `2' `usevars2' `usevars1' `if' `in' mat ZX = ZX\a macro shift } tokenize `allvars' mat vecaccum a = `1' `varlist_iv' `usevars1' `if' `in' mat ZZ = a while length("`2'")!= 0 { mat vecaccum a = `2' `varlist_iv' `usevars1' `if' `in' mat ZZ = ZZ\a macro shift } scalar s1=rowsof(ZZ) scalar s2=rowsof(ZX) if length("`distribution'")==0 { tempname omega_in mu_in matrix `omega_in'= J(s1,s1,0) matrix `mu_in'= J(s1,1,0) local vv=1 foreach val of numlist `mu' { matrix `mu_in'[`vv',1] =`val' local ++vv } local vv=1 foreach val of numlist `omega' { matrix `omega_in'[`vv',`vv'] =`val' local ++vv } } ***************************************************************************** *** (5c) Form estimates if non-normal distribution *** Coded based on the algorithm described on page 265 of REStat **************************************************************************** if length("`distribution'")!=0 { if c(flavor)=="Small" { dis "The distribution option should be used with care with Small Stata" dis "In this case bounds are calculated based on a maximum of 100 simulations" dis "It may be preferable to use a Guassian (exact) prior via mu() and omega()" local iterations=100 } else if c(flavor)=="IC"&c(SE)==0&c(MP)==0 { set matsize 800 local iterations=800 dis "The distribution option should be used with care with Stata IC" dis "In this case bounds are calculated based on a maximum of 800 simulations" dis "It may be preferable to use a Guassian (exact) prior via mu() and omega()" } else { set matsize 10000 } if length("`seed'")!= 0 { set seed `seed' } if "`dist1'"=="special" { local distvars foreach ivn of numlist 1(1)`count_iv' { local dv = `ivn'+1 local distvars `distvars' `dist`dv'' } mkmat `distvars', matrix(specialgamma) nomissing matrix mnsp = rowsof(specialgamma) local nsp = mnsp[1,1] matrix gammaDonors = J(`iterations',`count_iv',.) foreach ivn of numlist 1(1)`count_iv' { local gd=1 while `gd'<`iterations' { matrix gammaDonors[`gd',`ivn']=specialgamma[ceil(runiform()*`nsp'),`ivn'] local ++gd } } } qui estimates restore __iv qui estat vce matrix varcovar = r(V) mata: st_matrix("betas", select(st_matrix("e(b)"), st_matrix("e(b)") :!=0)) matrix mnvars = colsof(betas) local nvars = mnvars[1,1] local betas local zeros foreach num of numlist 1(1)`nvars' { local betas `betas' `=betas[1,`num']' local zeros `zeros' 0 } **matrix gamma = J(`nvars',1,0) matrix gamma = J(`count_exog',1,0) matrix A = inv(ZX'*inv(ZZ)*ZX)*ZX' dis "Simulating `iterations' iterations. This may take a moment." foreach num of numlist 1(1)`nvars' { matrix betasSim = J(`iterations',`nvars',.) } local iter = 1 while `iter' <= `iterations' { foreach ivn of numlist 1(1)`count_iv' { if "`dist1'"=="special" local gammaCall`ivn' = gammaDonors[`iter',`ivn'] matrix gamma[`ivn',1]=`gammaCall`ivn'' } matrix F = A*gamma qui rmvnormal, mean(`zeros') sigma(varcovar) matrix gsims = r(rmvnormal) foreach num of numlist 1(1)`nvars' { local input = gsims[1,`num'] matrix betasSim[`iter',`num'] = `input'+F[`num',1] } local ++iter } foreach num of numlist 1(1)`nvars' { mata : st_matrix("betasSim", sort(st_matrix("betasSim"), `num')) local lbci = 1-(1-`level')/2 local ubci = (1-`level')/2 local l`num' = betas[1,`num']-betasSim[round(`iterations'*`lbci'),`num'] local u`num' = betas[1,`num']-betasSim[round(`iterations'*`ubci'),`num'] *local lowerbound = betasSim[round(`iterations'*0.025),`num'] *local upperbound = betasSim[round(`iterations'*0.975),`num'] *dis "Bound for variable `num' is [`lowerbound',`upperbound']" } local Cint "Conley et al (2012)'s LTZ results" dis in yellow _newline dis "`Cint' (Non-Gaussian)" _col(55) "Number of obs = " e(N) dis in yellow in smcl "{hline 78}" dis "Variable" _col(13) "Lower Bound" _col(29) "Upper Bound" dis in yellow in smcl "{hline 78}" local vars_final local counter=0 foreach item in `e(instd)' `e(exogr)' _cons { if _b[`item']!=0|_se[`item']!=0 { local vars_final `vars_final' `item' local ++counter } } tokenize `vars_final' foreach regressor of numlist 1(1)`nvars' { dis in green "``regressor''" _col(13) `l`regressor'' _col(29) `u`regressor'' foreach var of local varlist2 { if `"`var'"'==`"``regressor''"' { ereturn scalar lb_`var'=`l`regressor'' ereturn scalar ub_`var'=`u`regressor'' } } } dis in yellow in smcl "{hline 78}" exit } ***************************************************************************** *** (5di) Form augmented var-covar and coefficient matrices for graphing ***************************************************************************** if length("`graph'")!=0 { local Nomegas : word count `graphomega' local Nmus : word count `graphmu' local Gerr1 "graphmu and graphomega" if `Nomegas'==0|`Nmus'==0 { dis as error "When specifing graph and ltz, the `Gerr1' options are required" exit 272 } if `Nomegas'!=`Nmus' { dis as error "`Gerr1' must take the same number of elements" exit 272 } matrix __graphmatLTZ = J(`Nomegas',4,.) local countdelta : word count `graphdelta' if `countdelta'!=0 { local j=1 foreach num of numlist `graphdelta' { matrix __graphmatLTZ[`j',1]=`num' local ++j } } tempname omegaC muC tokenize `graphmu' local j=1 foreach item in `graphomega' { matrix `omegaC'=J(s1,s1,0) matrix `omegaC'[1,1]=`item' matrix `muC'=J(s1,1,0) matrix `muC'[1,1]=``j'' qui estimates restore __iv mat Vc = e(V) + inv(ZX'*inv(ZZ)*ZX)*ZX' * `omegaC' * ZX*inv(ZX'*inv(ZZ)*ZX) mat bc = e(b) - (inv(ZX'*inv(ZZ)*ZX)*ZX' * `muC')' ereturn post bc Vc matrix CI = -invnormal((1-`level')/2) if `countdelta'==0 { scalar delta=`omegaC'[1,1] matrix __graphmatLTZ[`j',1]=delta } matrix __graphmatLTZ[`j',2]=_b[`graph'] matrix __graphmatLTZ[`j',3]=_b[`graph']-_se[`graph']*CI matrix __graphmatLTZ[`j',4]=_b[`graph']+_se[`graph']*CI local ++j } } ***************************************************************************** *** (5dii) Form augmented var-covar and coefficient matrices (see appendix) ***************************************************************************** qui estimates restore __iv mat V1 = e(V) + inv(ZX'*inv(ZZ)*ZX)*ZX' * `omega_in' * ZX*inv(ZX'*inv(ZZ)*ZX) mat b1 = e(b) - (inv(ZX'*inv(ZZ)*ZX)*ZX' * `mu_in')' ***************************************************************************** *** (5e) Determine lower and upper bounds ***************************************************************************** mat CI = -invnormal((1-`level')/2) mat lb = b1 - vecdiag(cholesky(diag(vecdiag(V1))))*CI mat ub = b1 + vecdiag(cholesky(diag(vecdiag(V1))))*CI dis _newline dis "Conley et al. (2012)'s LTZ results" _col(55) "Number of obs = " e(N) local lev=`level'*100 set level `lev' ereturn post b1 V1 ereturn display local CI = -invnormal((1 - `level')/2) foreach var of local varlist2 { ereturn scalar lb_`var'=_b[`var']-`CI'*_se[`var'] ereturn scalar ub_`var'=_b[`var']+`CI'*_se[`var'] } set level 95 } ******************************************************************************** *** (6) Visualise as per Conely et al., (2012) Figures 1-2 ******************************************************************************** if length("`graph'")!=0 { if "`method'"=="uci" { svmat __graphmat if length("`options'")==0 { local options ytitle("Estimated Beta for `graph'") scheme(s1color) /* */ xtitle("{&delta}") title("Union of Confidence Interval Approach") /* */ note("Methodology described in Conley et al. (2012)") /* */ legend(order(1 "Upper Bound (UCI)" 2 "Lower Bound (UCI)")) } sort __graphmat1 twoway line __graphmat3 __graphmat1, lpattern(dash) lcolor(black) || /// line __graphmat4 __graphmat1, lpattern(dash) lcolor(black) /// `options' } if "`method'"=="ltz" { svmat __graphmatLTZ if length("`options'")==0 { local options ytitle("{&beta}") xtitle("{&delta}") /* */ title("Local to Zero Approach") /* */ note("Methodology described in Conley et al. (2012)") /* */ legend(order(1 "Point Estimate (LTZ)" 2 "CI (LTZ)")) } twoway line __graphmatLTZ2 __graphmatLTZ1, || /// line __graphmatLTZ3 __graphmatLTZ1, lpattern(dash) lcolor(black) || /// line __graphmatLTZ4 __graphmatLTZ1, lpattern(dash) lcolor(black) /// `options' } cap drop __graphmat* } ******************************************************************************** *** (7) Clean up ******************************************************************************** estimates drop __iv restore end