*! genstest v1.0.6 28feb2020 *! authors Zachary Flynn *! 1.0.6: updated to work with GMM in Stata 16 cap program drop genstest program genstest, rclass syntax [anything] [if] [in] [pweight/] [, null(string) trim(real 0.15)/* */ INSTruments(string) WMATrix(string) cluster(varname)/* */ winitial(string) ci(string) /* */ nuisS varS center DERIVative(string)*/* */ sb init(numlist) gmm(string) test(string)/* */ sb_lags(string)/* */ stab small igmm twostep] if ("`nuisS'"=="" & "`varS'"=="" & "`sb'"=="sb") { di in red "Note: using the nuisS and/or varS options will decrease computation time for the single-break tests." } if ("`e(cmd)'"=="gmm" & "`anything'"=="") { if ("`test'"!="") { local nnull = wordcount("`test'") local sexp = "`e(sexp_1)'" tokenize "`test'" local i 1 while (`i'<=wordcount("`test'")) { while ( regexm("`sexp'","{``i''[:_cons]*[=]*[0-9\.]*}") ) { local sexp = regexr("`sexp'","{``i''[:_cons]*[=]*[0-9\.]*}","<``i''>") } local i = `i' + 1 } } else if ("`anything'" != "") { local sexp = "`anything'" local nnull = wordcount("`null'") } else { local sexp = "`e(sexp_1)'" local nnull = wordcount("`null'") } while (regexm("`sexp'", ":_cons")) { local sexp = regexr("`sexp'", ":_cons", "") } local first = "e(cmdline)" local comma = "" } else { ParseEquation "zero" "`anything'" local nnull = wordcount("`r(names)'") local first = "anything" if ("`test'"!="" & `nnull'==0) { local nnull = wordcount("`test'") local expres = "`anything'" local first = "expres" tokenize "`test'" local i 1 while (`i' <= wordcount("`test'")) { while ( regexm("`expres'","{``i''[=]*[0-9\.]*}") ) { local expres = regexr("`expres'","{``i''[=]*[0-9\.]*}","<``i''>") } local i = `i' + 1 } } local comma = "," } if ("`ci'"!="") { ParseCi `ci' local alpha = `r(alpha)' local single = "`r(single)'" local allpv = "`r(allpv)'" local rnge = "`r(range)'" local points = "`r(points)'" local tograph = "`r(autograph)'" local i = 1 tempname Sci aveSci expSci supSci qllSci tempname avestabSci expstabSci supstabSci qllstabSci mata `Sci' = J(1, `nnull'+1, .) if ("`sb'" == "sb") { mata `aveSci' = J(1, `nnull'+1, .) mata `expSci' = J(1, `nnull'+1, .) mata `supSci' = J(1, `nnull'+1, .) } mata `qllSci' = J(1, `nnull'+1, .) if ("`stab'" == "stab") { if ("`sb'" == "sb") { mata `avestabSci' = J(1, `nnull'+1, .) mata `expstabSci' = J(1, `nnull'+1, .) mata `supstabSci' = J(1, `nnull'+1, .) } mata `qllstabSci' = J(1, `nnull'+1, .) } if (`nnull' == 1) { tokenize `rnge' local start = `1' local stop = `2' tokenize `points' local stepsize = (`stop'-`start')/`1' local cur = `start' local dispci = 1 local crits = "" local critsstab = "" while (`cur' <= `stop') { genstest_internal ``first'' `if' `in' `pweight' `comma' null(`cur') trim(`trim')/* */ inst(`instruments') wmatrix(`wmatrix') cluster(`cluster') winitial(`winitial')/* */ `nuisS' `varS' `center' deriv(`derivative')/* */ `sb' init(`init') test(`test') gmm(`sexp')/* */ sb_lags(`sb_lags') `allpv' critsstab(`critsstab') crits(`crits') `options' `stab' `small'/* */ alpha(`alpha') `igmm' if ("`allpv'"=="allpv") { mata `Sci' = `Sci' \ (`r(par)', `cur') if ("`sb'"=="sb") { mata `aveSci' = `aveSci' \ (`r(pavear)', `cur') mata `expSci' = `expSci' \ (`r(pexpar)', `cur') mata `supSci' = `supSci' \ (`r(psupar)', `cur') } mata `qllSci' = `qllSci' \ (`r(pqllar)', `cur') if ("`stab'"=="stab") { mata `qllstabSci' = `qllstabSci' \ (`r(pqllstabar)', `cur') if ("`sb'"=="sb") { mata `avestabSci' = `avestabSci' \ (`r(pavestabar)', `cur') mata `supstabSci' = `supstabSci' \ (`r(psupstabar)', `cur') mata `expstabSci' = `expstabSci' \ (`r(pexpstabar)', `cur') } } } else { local crits = "`r(crits)'" local critsstab = "`r(critsstab)'" mata `Sci' = `Sci' \ (`r(ar)', `cur') if ("`sb'"=="sb") { mata `aveSci' = `aveSci' \ (`r(avear)', `cur') mata `expSci' = `expSci' \ (`r(expar)', `cur') mata `supSci' = `supSci' \ (`r(supar)', `cur') } mata `qllSci' = `qllSci' \ (`r(qllar)', `cur') if ("`stab'"=="stab") { mata `qllstabSci' = `qllstabSci' \ (`r(qllstabar)', `cur') if ("`sb'"=="sb") { mata `avestabSci' = `avestabSci' \ (`r(avestabar)', `cur') mata `supstabSci' = `supstabSci' \ (`r(supstabar)', `cur') mata `expstabSci' = `expstabSci' \ (`r(expstabar)', `cur') } } } local cur = `cur' + `stepsize' } } else if (`nnull'==2) { local dispci = 0 tokenize `rnge' local start1 = `1' local stop1 = `2' local start2 = `3' local stop2 = `4' tokenize `points' local step1 = (`stop1'-`start1')/`1' + ((`stop1'-`start1') <= 0) local step2 = (`stop2'-`start2')/`2' + ((`stop2'-`start2') <= 0) local cur1 = `start1' local cur2 = `start2' local crits = "" while (`cur1' <= `stop1') { while (`cur2' <= `stop2') { genstest_internal ``first'' `if' `in' `pweight' `comma' null(`cur1' `cur2') trim(`trim')/* */ inst(`instruments') wmatrix(`wmatrix') cluster(`cluster') winitial(`winitial')/* */ `nuisS' `varS' `center' deriv(`derivative')/* */ `sb' init(`init') test(`test') gmm(`sexp') alpha(`alpha')/* */ sb_lags(`sb_lags') `allpv' crits(`crits') critsstab (`critsstab') `options' `stab' `small' `igmm' if ("`allpv'"=="allpv") { mata `Sci' = `Sci' \ (`r(par)', `cur1', `cur2') if ("`sb'"=="sb") { mata `aveSci' = `aveSci' \ (`r(pavear)', `cur1', `cur2') mata `expSci' = `expSci' \ (`r(pexpar)', `cur1', `cur2') mata `supSci' = `supSci' \ (`r(psupar)', `cur1', `cur2') } mata `qllSci' = `qllSci' \ (`r(pqllar)', `cur1', `cur2') if ("`stab'"=="stab") { if ("`sb'"=="sb") { mata `avestabSci' = `avestabSci' \ (`r(pavestabar)', `cur1', `cur2') mata `expstabSci' = `expstabSci' \ (`r(pexpstabar)', `cur1', `cur2') mata `supstabSci' = `supstabSci' \ (`r(psupstabar)', `cur1', `cur2') } mata `qllstabSci' = `qllstabSci' \ (`r(pqllstabar)', `cur1', `cur2') } } else { local crits = "`r(crits)'" local critsstab = "`r(critsstab)'" mata `Sci' = `Sci' \ (`r(ar)', `cur1', `cur2') if ("`sb'"=="sb") { mata `aveSci' = `aveSci' \ (`r(avear)', `cur1', `cur2') mata `expSci' = `expSci' \ (`r(expar)', `cur1', `cur2') mata `supSci' = `supSci' \ (`r(supar)', `cur1', `cur2') } mata `qllSci' = `qllSci' \ (`r(qllar)', `cur1', `cur2') if ("`stab'"=="stab") { if ("`sb'"=="sb") { mata `avestabSci' = `avestabSci' \ (`r(avestabar)', `cur1', `cur2') mata `expstabSci' = `expstabSci' \ (`r(expstabar)', `cur1', `cur2') mata `supstabSci' = `supstabSci' \ (`r(supstabar)', `cur1', `cur2') } mata `qllstabSci' = `qllstabSci' \ (`r(qllstabar)', `cur1', `cur2') } } local cur2 = `cur2' + `step2' } local cur1 = `cur1' + `step1' local cur2 = `start2' } } local paramnames = "`r(names)'" if ("`allpv'"=="") { tokenize `crits' local Scrit = invchi2tail(r(kz)-r(nest), `alpha') mata `Sci' = select (`Sci', `Sci'[.,1]:<=`Scrit') if ("`sb'"=="sb") { mata `aveSci' = select (`aveSci', `aveSci'[.,1]:<=`1') mata `expSci' = select (`expSci', `expSci'[.,1]:<=`2') mata `supSci' = select (`supSci', `supSci'[.,1]:<=`3') mata `qllSci' = select (`qllSci', `qllSci'[.,1]:<=`4') } else { mata `qllSci' = select (`qllSci', `qllSci'[.,1]:<=`1') } if ("`stab'"=="stab") { tokenize `critsstab' if ("`sb'"=="sb") { mata `avestabSci' = select (`avestabSci', `avestabSci'[.,1]:<=`1') mata `expstabSci' = select (`expstabSci', `expstabSci'[.,1]:<=`2') mata `supstabSci' = select (`supstabSci', `supstabSci'[.,1]:<=`3') mata `qllstabSci' = select (`qllstabSci', `qllstabSci'[.,1]:<=`4') } else { mata `qllstabSci' = select (`qllstabSci', `qllstabSci'[.,1]:<=`1') } } } mata st_matrix("`Sci'", `Sci') if ("`sb'"=="sb") { mata st_matrix ("`aveSci'", `aveSci') mata st_matrix ("`expSci'", `expSci') mata st_matrix ("`supSci'", `supSci') } mata st_matrix("`qllSci'", `qllSci') if ("`stab'"=="stab") { if ("`sb'"=="sb") { mata st_matrix ("`avestabSci'", `avestabSci') mata st_matrix ("`expstabSci'", `expstabSci') mata st_matrix ("`supstabSci'", `supstabSci') } mata st_matrix("`qllstabSci'", `qllstabSci') } if (`nnull'==2 & "`tograph'"=="autograph") { cap confirm matrix `qllSci' if (!_rc) { dograph "`qllSci'" "qLL-S" "qllS" "`paramnames'" `alpha' "`allpv'" } cap confirm matrix `Sci' if (!_rc) { dograph "`Sci'" "S" "S" "`paramnames'" `alpha' "`allpv'" } if ("`sb'" == "sb") { cap confirm matrix `aveSci' if (!_rc) { dograph "`aveSci'" "ave-S" "aveS" "`paramnames'" `alpha' "`allpv'" } cap confirm matrix `expSci' if (!_rc) { dograph "`expSci'" "exp-S" "expS" "`paramnames'" `alpha' "`allpv'" } cap confirm matrix `supSci' if (!_rc) { dograph "`supSci'" "sup-S" "supS" "`paramnames'" `alpha' "`allpv'" } } if ("`stab'"=="stab") { cap confirm matrix `qllstabSci' if (!_rc) { dograph "`qllstabSci'" "qLL-stab-S" "qllstabS" "`paramnames'" `alpha' "`allpv'" } if ("`sb'" == "sb") { cap confirm matrix `avestabSci' if (!_rc) { dograph "`avestabSci'" "ave-stab-S" "avestabS" "`paramnames'" `alpha' "`allpv'" } cap confirm matrix `expstabSci' if (!_rc) { dograph "`expstabSci'" "exp-stab-S" "expstabS" "`paramnames'" `alpha' "`allpv'" } cap confirm matrix `supstabSci' if (!_rc) { dograph "`supstabSci'" "sup-stab-S" "supstabS" "`paramnames'" `alpha' "`allpv'" } } } } genstestgraph `Sci', alpha(`alpha') `allpv' local outSci = "`r(output)'" cap return matrix Sci = `Sci' if ("`sb'" == "sb") { genstestgraph `aveSci', alpha(`alpha') `allpv' local outaveSci = "`r(output)'" genstestgraph `expSci', alpha(`alpha') `allpv' local outexpSci = "`r(output)'" genstestgraph `supSci', alpha(`alpha') `allpv' local outsupSci = "`r(output)'" cap return matrix aveSci = `aveSci' cap return matrix expSci = `expSci' cap return matrix supSci = `supSci' } genstestgraph `qllSci', alpha(`alpha') `allpv' local outqllSci = "`r(output)'" cap return matrix qllSci = `qllSci' if ("`stab'"=="stab") { genstestgraph `qllstabSci', alpha(`alpha') `allpv' cap return matrix qllstabSci = `qllstabSci' local outqllstabSci = "`r(output)'" if ("`sb'"=="sb") { genstestgraph `avestabSci', alpha(`alpha') `allpv' cap return matrix avestabSci = `avestabSci' local outavestabSci = "`r(output)'" genstestgraph `expstabSci', alpha(`alpha') `allpv' cap return matrix expstabSci = `expstabSci' local outexpstabSci = "`r(output)'" genstestgraph `supstabSci', alpha(`alpha') `allpv' cap return matrix supstabSci = `supstabSci' local outsupstabSci = "`r(output)'" } } } else { local dispci = 0 } genstest_internal ``first'' `if' `in' `pweight' `comma' null(`null') trim(`trim')/* */ inst(`instruments') wmatrix(`wmatrix') cluster(`cluster') winitial(`winitial')/* */ `nuisS' `varS' `center' deriv(`derivative')/* */ `sb' init(`init') test(`test') gmm(`sexp')/* */ sb_lags(`sb_lags') alpha(`alpha') allpv `options' `stab' `small' `igmm' local nest = r(nest) local kz = r(kz) local pS = `r(par)' if ("`sb'"=="sb") { local paveS = `r(pavear)' local pexpS = `r(pexpar)' local psupS = `r(psupar)' } local pqllS = `r(pqllar)' local st1 = `r(ar)' if ("`sb'"=="sb") { local st2 = `r(avear)' local st3 = `r(expar)' local st4 = `r(supar)' } local st5 = `r(qllar)' local names = r(names) local N = r(N) return scalar nest = `nest' return scalar kz = `kz' return scalar S = `st1' if ("`sb'"=="sb") { return scalar aveS = `st2' return scalar expS = `st3' return scalar supS = `st4' } return scalar qllS = `st5' return scalar N = `N' if ("`stab'" == "stab") { local pqllstabS = `r(pqllstabar)' local qllstabS = `r(qllstabar)' return scalar pqllstabS = `pqllstabS' return scalar qllstabS = `qllstabS' if ("`sb'" == "sb") { local pavestabS = `r(pavestabar)' local avestabS = `r(avestabar)' return scalar pavestabS = `pavestabS' return scalar avestabS = `avestabS' local pexpstabS = `r(pexpstabar)' local expstabS = `r(expstabar)' return scalar pexpstabS = `pexpstabS' return scalar expstabS = `expstabS' local psupstabS = `r(psupstabar)' local supstabS = `r(supstabar)' return scalar psupstabS = `psupstabS' return scalar supstabS = `supstabS' } } local null = "`r(null)'" local linelength = 15 local ll2 = 11 local ll3 = 8 local ll4 = 20 if (`dispci'==1) { local mkdiv = `"di as txt "{hline `linelength'}{c +}{hline `ll2'}{c +}{hline `ll3'}{c +}{hline `ll4'}""' local mktop = `"di as txt "{hline `linelength'}{c TT}{hline `ll2'}{c TT}{hline `ll3'}{c TT}{hline `ll4'}""' local mktitle = `"di as txt "{lalign 15:Test}{c |}{center 11:Statistic}{c |}{center 7:P-value }{c |} CI (alpha=`alpha')""' local mkbot = `"di as txt "{hline `linelength'}{c BT}{hline `ll2'}{c BT}{hline `ll3'}{c BT}{hline `ll4'}""' } else { local mktop = `"di as txt "{hline `linelength'}{c TT}{hline `ll2'}{c TT}{hline `ll3'}""' local mktitle = `"di as txt "{lalign 15:Test}{c |}{center 11:Statistic}{c |}{center 7:P-value }""' local mkdiv = `"di as txt "{hline `linelength'}{c +}{hline `ll2'}{c +}{hline `ll3'}""' local mkbot = `"di as txt "{hline `linelength'}{c BT}{hline `ll2'}{c BT}{hline `ll3'}""' } `mktop' `mktitle' `mkdiv' di as txt "{lalign 15:S}{c |}" _continue di "{txt} " %-9.6f `st1' " {c |}" _continue return scalar pS = `pS' di "{txt} " %-4.3f `pS' _continue if (`dispci'==1) { di " {c |} `outSci'" } else { di } di as txt "{lalign 15:qLL-S}{c |}" _continue di " " %-9.6f `st5' " {c |}" _continue return scalar pqllS = `pqllS' di " " %-4.3f `pqllS' _continue if (`dispci'==1) { di " {c |} `outqllSci'" } else { di } if ("`sb'"=="sb") { di as txt "{lalign 15:ave-S}{c |}" _continue di "{txt} " %-9.6f `st2' " {c |}" _continue return scalar paveS = `paveS' di "{txt} " %-4.3f `paveS' _continue if (`dispci'==1) { di " {c |} `outaveSci'" } else { di } di as txt "{lalign 15:exp-S}{c |}" _continue di "{txt} " %-9.6f `st3' " {c |}" _continue return scalar pexpS = `pexpS' di "{txt} " %-4.3f `pexpS' _continue if (`dispci'==1) { di " {c |} `outexpSci'" } else { di } di as txt "{lalign 15:sup-S}{c |}" _continue di "{txt} " %-9.6f `st4' " {c |}" _continue return scalar psupS = `psupS' di "{txt} " %-4.3f `psupS' _continue if (`dispci'==1) { di " {c |} `outsupSci'" } else { di } } if ("`stab'" == "stab") { di as txt "{lalign 15:qLL-stab-S}{c |}" _continue di " " %-9.6f `qllstabS' " {c |}" _continue di " " %-4.3f `pqllstabS' _continue if (`dispci'==1) { di " {c |} `outqllstabSci'" } else { di } if ("`sb'"=="sb") { di as txt "{lalign 15:ave-stab-S}{c |}" _continue di " " %-9.6f `avestabS' " {c |}" _continue di " " %-4.3f `pavestabS' _continue if (`dispci'==1) { di " {c |} `outavestabSci'" } else { di } di as txt "{lalign 15:exp-stab-S}{c |}" _continue di " " %-9.6f `expstabS' " {c |}" _continue di " " %-4.3f `pexpstabS' _continue if (`dispci'==1) { di " {c |} `outexpstabSci'" } else { di } di as txt "{lalign 15:sup-stab-S}{c |}" _continue di " " %-9.6f `supstabS' " {c |}" _continue di " " %-4.3f `psupstabS' _continue if (`dispci'==1) { di " {c |} `outsupstabSci'" } else { di } } } `mkbot' if ("`names'"!="") { di "{res}Tested null hypothesis vector: <`names'> = <" _continue tokenize "`null'" local i 1 while (`i' <= wordcount("`names'")) { if ("`null'" == "zero") { di " 0.000" _continue } else { di " " %-7.3f ``i'' _continue } local i = `i' + 1 } } di ">" di "{res}Number of Instruments - Included Instruments: `kz' - `nest'" di "{res}Number of Observations: `N'" end cap program drop dograph program dograph args ciname title name params alpha allpv tempname var1 var2 tmp if ("`allpv'"=="allpv") { mata `tmp' = cishell (select (`ciname', `ciname'[.,1] :>= `alpha')) } else { mata `tmp' = cishell(`ciname') } mata st_matrix ("`var1'", `tmp'[.,1]) mata st_matrix ("`var2'", `tmp'[.,2]) qui svmat `var1' qui svmat `var2' qui sum `var2' local min2 = r(min) - r(sd) local max2 = r(max) + r(sd) qui sum `var1' local var1N = r(N) local min1 = r(min) - r(sd) local max1 = r(max) + r(sd) tokenize `params' twoway (line `var2' `var1' if mod(_n, 2)==0, lcolor(black)) (line `var2' `var1' if mod(_n,2)==1, lcolor(black)) /* */ (line `var2' `var1' if _n <= 2, lcolor(black)) (line `var2' `var1' if _n >= `var1N'-1, lcolor(black)), /* */ name("`name'", replace) xscale(range(`min1' `max1')) yscale(range(`min2' `max2')) title("`title'") /* */ xtitle("`1'") ytitle("`2'") legend(off) end cap program drop ParseCi program ParseCi, rclass syntax anything [, alpha(real 0.05) points(numlist) single allpv AUtograph ] return local alpha = `alpha' return local autograph = "`autograph'" local i 1 local rpoints = "" if (wordcount("`points'")>2) { noi di "Confidence Interval only supports up to two parameters at this time." } else if (wordcount("`points'")==2) { tokenize "`points'" local rpoints = "`1' `2'" local i = `i' + 1 } else { local rpoints = "`points'" } return local points = "`rpoints'" if (wordcount("`anything'")>4) { noi di "Confidence Interval only supports up to two parameters at this time." } else if (wordcount("`anything'")==4) { tokenize "`anything'" local rrange = "`1' `2' `3' `4'" } else if (wordcount("`anything'")==2) { tokenize "`anything'" local rrange = "`1' `2'" } else { noi di "Confidence Interval Range incorrectly specified." } return local allpv = "`allpv'" return local single = "`single'" return local range = "`rrange'" end cap program drop genstestgraph program genstestgraph, rclass syntax name(name=ci) [, alpha(real 0.05)* allpv] tempname res if ("`allpv'"=="allpv") { mata `res' = genstest_display (`ci', `alpha', 1) } else { mata `res' = genstest_display (`ci', `alpha', 0) } mata st_local ("tmp", `res') return local output = "`tmp'" end capture program drop genstest_internal program genstest_internal, rclass syntax [anything] [if] [in] [pweight/] [, null(string) trim(real 0.15)/* */ INSTruments(string) WMATrix(string) cluster(varname)/* */ alpha(real 0.05) winitial(string) ci(string) /* */ nuisS varS center DERIVative(string)*/* */ sb init(numlist) /* */ gmm(string) test(string) /* */ sb_lags(string) allpv crits(numlist) critsstab(numlist) stab small igmm twostep] if ("`sb'"=="sb") { local sb = "" } else { local sb = "sb" } local yes = 1 local no = 0 local tau = `trim' local cbar = 10 if ("`null'"=="last") { local sexp = "`gmm'" ParseEquation "zero" "`sexp'" local names = "`r(names)'" tempname b nummat mat `b' = e(b) tokenize "`names'" local nn = wordcount("`names'") local i = 1 local null = "" while (`i' <= `nn') { mat `nummat' = `b'[1,"``i'':_cons"] local num = `nummat'[1,1] local null = "`null' `num'" local i = `i' + 1 } } marksample touse tempvar wvar if ("`weight'"!="") { qui gen `wvar' = `exp' if ("`weight'" == "fweight") { local wtexp = `"[`weight'=`exp']"' } else if ("`weight'" == "pweight") { local wtexp = `"[aweight=`exp']"' } else { local wtexp = "" } qui sum `wvar' if `touse' `wtexp', meanonly if ("`weight'"=="fweight") { local wf = 1 } else { local wf = r(N)/r(sum_w) } } else { qui gen `wvar' = 1 local wf = 1 } qui replace `wvar' = `wf'*`wvar' if ("`null'"=="") { local null = "zero" } tempname spec mata `spec' = genstestSpec() mata `spec'.center = `no' if ("`center'"=="center") { mata `spec'.center = `yes' } mata `spec'.sb_center = `spec'.center /* Set Estimation Type */ if ("`igmm'"=="igmm") { mata `spec'.twostep = `no' mata `spec'.sb_twostep = `no' } else { mata `spec'.twostep = `yes' mata `spec'.sb_twostep = `yes' } /* Set initial weighting matrix */ if ("`winitial'" == "identity") { mata `spec'.identity =`yes' } else { mata `spec'.identity = `no' } mata `spec'.sb_identity = `spec'.identity /* Parse Expression */ if ("`gmm'" == "") { local expression = "`anything'" gl GENSTEST_GMM_POST = 0 } else { local expression = "`gmm'" gl GENSTEST_GMM_POST = 1 } ParseEquation "`null'" "`expression'" local equation = "`r(equation)'" local nest = `r(nest)' local names = "`r(names)'" local est_names = "`r(est_names)'" local orignull = "`null'" tempname eq local i 1 if ("`init'"!="") { tokenize "`init'" } while (`i' <= `nest') { if ("`init'"!="") { scalar GENSTEST_EST`i' = ``i'' } else { scalar GENSTEST_EST`i' = 0 } local i = `i' + 1 } qui gen `eq' = `equation' if `touse' markout `touse' `eq' mata `spec'.nest = `nest' mata `spec'.hasderiv = `no' if ("`derivative'"!="") { ParseDeriv `nest' "`derivative'" "`options'" "`r(est_names)'" "`null'" "`names'" mata `spec'.hasderiv = `yes' local i 1 while (`i' <= `nest') { local deriv`i' = "`r(deriv`i')'" local i = `i' + 1 } } else { mata `spec'.hasderiv = `no' } if ("`init'"!="" & `nest'==wordcount("`init'")) { tokenize "`init'" local i 2 mata `spec'.init = `1' while (`i'<=`nest') { mata `spec'.init = `spec'.init, ``i'' local i = `i' + 1 } } else { mata `spec'.init = J(1,`nest',0) } ********************Parse wmatrix options*********************************** tokenize "`wmatrix'" ParseWMatrix `1', `2' lags(`3') sb_lags(`sb_lags') local wmtype = `r(wmtype)' mata `spec'.wmtype = `wmtype' local hc 1 local hac 2 local clst 3 if (`wmtype'==`hc') { mata `spec'.adjf = `r(adjf)' } else if (`wmtype'==`hac') { mata `spec'.l_opt = `r(l_opt)' mata `spec'.l_auto = `r(l_auto)' mata `spec'.l_nlags = `r(l_nlags)' mata `spec'.sbl_opt = `r(sbl_opt)' mata `spec'.sbl_auto = `r(sbl_auto)' mata `spec'.kernel = `r(k_type)' } else if (`wmtype' == `clst') { local cluster = "`r(cluster)'" } if ("`small'"=="small") { mata `spec'.small = 1 } else { mata `spec'.small = 0 } cap confirm number `tau' if (_rc) { noi di in red "Tau must be a number." exit } else { local k 1 while (`k'<=4) { local j = `k'*.05 if (`tau'==`j') { continue, break } local k = `k' + 1 } if (`k'==5) { noi di in red "Tau must be 0.05, 0.10, 0.15, or 0.20." exit } } cap confirm existence "`instruments'" if (_rc) { cap confirm existence `e(inst_1)' if (!_rc) { tempvar cons qui gen `cons' = 1 local inst = "`e(inst_1)'" local inst = regexr("`instruments'", "_cons", "`cons'") } else { cap confirm variable `e(insts)' if (!_rc) { local inst = "`e(insts)'" } else { noi di in red "There must be at least one instrument." exit } } } local nest = `nest' tempvar cons ParseInst `instruments' local instruments = "`r(inst)'" local constant = "`r(constant)'" if ("`constant'" == "constant") { qui gen `cons' = 1 local instruments = "`r(inst)' `cons'" } local kz : word count `instruments' tokenize "`instruments'" local instruments = "" local i 1 while (`i'<=`kz') { tempvar z`i' qui gen `z`i'' = `wvar'*``i'' if `touse' local instruments = "`instruments' `z`i''" local i = `i' + 1 markout `touse' `z`i'' } tempname z mata `z' = st_data(.,"`instruments'",0) ********************Set Up Clusters (if applicable)************************ cap confirm variable `cluster' tempname clst cl tempname clustvar if (!_rc) { qui gen `clustvar' = `cluster' if `touse' mata `clst' = st_data(.,"`clustvar'", 0) mata `cl' = gencluster(`clst') mata `spec'.cluster = `cl' mata `spec'.sb_clst = 0 // We do full-sample AR first mata `spec'.wmtype = 3 } else { if ("`cluster'"!="") { noi di "Cluster variable does not exist." exit } else { mata `cl' = 0 } } **************S Computation************ tempname gam m_touse phi ar y tempname large mata `spec'.fs_wmatrix = 0 mata `spec'.instr = `z' if (`nest' > 0) { mata `m_touse' = J(rows(`z'),1,1) mata `gam' = `spec'.matgmm(`m_touse') mata `y' = `spec'.ufunc(`gam') mata `phi' = invsym(`spec'.wm1) } else { mata `y' = `spec'.ufunc(0) mata `gam' = 0 mata `phi' = `spec'.PHI(`y', `spec'.instr) mata `spec'.wm1 = invsym (`phi') } tempname xx dd mata `xx' = cross (`z',`z') mata `dd' = `spec'.coeff(`y',`spec'.instr) mata `ar' = `dd''*(1/rows(`z'))*`xx'*`spec'.wm1*`xx'*`dd' if ("`nuisS'"=="nuisS") { mata `spec'.gam = `gam' mata `spec'.nuisS = `yes' } else { mata `spec'.gam = 0 mata `spec'.nuisS = `no' } if ("`varS'" == "varS") { mata `spec'.fs_wmatrix = `spec'.wm1 mata `spec'.varS = `yes' } else { mata `spec'.fs_wmatrix = 0 mata `spec'.varS = `no' } mata `spec'.fs_gam = `gam' **************Single-Break Tests******************* if ("`sb'"=="") { tempname genar mata `genar' = `spec'.wmtype < 3 ?`spec'.arsinglebreak(`tau') : `spec'.clarsinglebreak(`z',`tau') } ************QLL Tests************************* tempname qll qllar pv1 pv pvstab1 pvstab mata `qll' = `spec'.qllstab(`y',`z',`phi') mata `qllar' = (`cbar'/(1+`cbar'))*(`ar') + `qll' if ("`allpv'"=="allpv") { if ("`sb'"=="") { mata `pv1' = pvalue (`genar' \ `qllar', `tau', `nest', `kz', 1)' } else { mata `pv1' = pvalue (`qllar', `tau', `nest', `kz', 0)' } vectonum `pv1' local pv = "`r(num)'" } tempname arvalue mata st_numscalar("`arvalue'",`ar') if ("`stab'" == "stab" & "`allpv'"=="allpv") { if ("`sb'"=="") { mata `pvstab1' = pvalue ((`genar' \ `qllar') - J( 4, 1, `ar'), `tau', `kz', `kz', 1)' } else { mata `pvstab1' = pvalue (`qllar' - J(1, 1, `ar'), `tau', `kz', `kz', 0)' } vectonum `pvstab1' local pvstab = "`r(num)'" } *************Return Values of Interest************************* tempname avear supar expar qllarstat if ("`sb'"=="") { mata st_numscalar("`avear'",`genar'[1,1]) mata st_numscalar("`expar'",`genar'[2,1]) mata st_numscalar("`supar'", `genar'[3,1]) } mata st_numscalar("`qllarstat'",`qllar'[1,1]) if ("`sb'" == "") { return local avear = `avear' return local supar = `supar' return local expar = `expar' } return local qllar = `qllarstat' if ("`stab'"=="stab") { return local qllstabar = `qllarstat' - `arvalue' if ("`sb'" == "") { return local avestabar = `avear' - `arvalue' return local supstabar = `supar' - `arvalue' return local expstabar = `expar' - `arvalue' } } return scalar kz = `kz' return scalar nest = `nest' mata st_numscalar ("r(N)", rows(`spec'.instr)) return local null = "`null'" return scalar ar = `arvalue' if ("`allpv'"=="allpv") { local arpval = chi2tail(`kz'-`nest',`arvalue') return local par = `arpval' } return local names = "`names'" if ("`allpv'"=="allpv") { tokenize `pv' if ("`sb'"=="") { return local pavear = `1' return local pexpar = `2' return local psupar = `3' } if ("`sb'"=="") { local qllpv = `4' } else { local qllpv = `1' } return local pqllar = `qllpv' if ("`stab'"=="stab") { tokenize `pvstab' if ("`sb'"=="") { return local pavestabar = `1' return local pexpstabar = `2' return local psupstabar = `3' local qllpv = `4' } else { local qllpv = `1' } return local pqllstabar = `qllpv' } } tempname N mata st_numscalar("`N'",rows(`z')) return scalar N = `N' tempname cv cvstab if ("`crits'"=="" & "`allpv'"=="") { if ("`sb'"=="") { mata `cv' = chkstat(`kz', `nest', `tau', `alpha', 1)' if ("`stab'"=="stab") { mata `cvstab' = chkstat(`kz', `kz', `tau', `alpha', 1)' } } else { mata `cv' = chkstat(`kz', `nest', `tau', `alpha', 0)' if ("`stab'"=="stab") { mata `cvstab' = chkstat(`kz', `kz', `tau', `alpha', 0)' } } vectonum `cv' return local crits = "`r(num)'" if ("`stab'"=="stab") { vectonum `cvstab' return local critsstab = "`r(num)'" } } else { return local crits = "`crits'" return local critsstab = "`critsstab'" } local i 1 while (`i'<=`nest') { scalar drop GENSTEST_EST`i' local i = `i' + 1 } mata mata drop GENSTEST_* scalar drop GENSTEST_REALLY_TEMP end cap program drop vectonum program vectonum, rclass ****Convert Mata vector to Stata numlist***** args vec local result = "" tempname N z mata st_numscalar("`N'",cols(`vec')) local i 1 while (`i'<=`N') { mata st_numscalar("`z'", `vec'[1,`i']) local result = "`result' " + string(`z') local i = `i' + 1 } return local num = "`result'" end cap program drop ParseEquation program ParseEquation, rclass **** Parse equations into useable form **** args null expression if ("`null'"!="zero") { tokenize "`null'" local length = wordcount("`null'") } local equation = "`expression'" local names = "" local i 1 while (1) { local pos1 = strpos("`equation'","<") if (`pos1' == 0) { continue, break } local substring = substr("`equation'",`pos1',.) local pos2 = strpos("`substring'",">") local name`i' = substr("`substring'",2,`pos2'-2) while(regexm("`equation'","<`name`i''>") == 1) { if ("`null'"!="zero") { local equation = regexr("`equation'","<`name`i''>","``i''") } else { local equation = regexr("`equation'","<`name`i''>","0") } } local names = "`names' `name`i''" local i = `i' + 1 } local i 1 local est_names = "" while (1) { local pos1 = strpos("`equation'","{") local substring = substr("`equation'",`pos1',.) local pos2 = strpos("`substring'","}") if (`pos1' == 0) { continue, break } local est_name`i' = substr("`substring'",2,`pos2'-2) while (regexm("`equation'","{`est_name`i''}") == 1) { local equation = regexr("`equation'","{`est_name`i''}", "GENSTEST_EST`i'") } local est_names = "`est_names' `est_name`i''" local i = `i' + 1 } return local nest = `i' - 1 return local est_names = "`est_names'" return local names = "`names'" return local equation = "`equation'" end cap program drop ParseDerivHelper program ParseDerivHelper, rclass syntax [anything], DERIVative(string)* est_names(string) [null_names(string)] tokenize "`est_names'" local i 1 local match 0 while (`i' <= wordcount("`est_names'")) { if (regexm("`derivative'", "/``i''")) { local match 1 continue, break } local i = `i' + 1 } if ($GENSTEST_GMM_POST == 1) { local i 1 tokenize "`null_names'" while (`i' <= wordcount("`null_names'")) { while (regexm("`derivative'","{``i''}")) { local derivative = regexr("`derivative'","{``i''}","<``i''>") } local i = `i' + 1 } } tokenize `est_names' local i = 1 while (`i' <= wordcount("`est_names'")) { while (regexm("`derivative'", "{``i''}")) { local derivative = regexr("`derivative'", "{``i''}", "GENSTEST_EST`i'") } local i = `i' + 1 } if (`match'!=0) { return local deriv = "`derivative'" } else { return local deriv = "" } return local oderivs = "`options'" end cap program drop ParseDeriv program ParseDeriv, rclass args nest strwderiv1 derivs names null null_names tempname nderiv pos_st pos_end mata `nderiv' = "`derivs'" ParseDerivHelper, derivative(`strwderiv1') `derivs' est_names("`names'") null_names("`null_names'") tempname rwderiv1 if ("`r(deriv)'"!= "") { mata `rwderiv1' = "`r(deriv)'" local i = 2 } else { local i = 1 } while (`i'<=`nest') { cap ParseDerivHelper, `r(oderivs)' est_names("`names'") null_names("`null_names'") if (!_rc) { if ("`r(deriv)'"!="") { tempname rwderiv`i' mata `rwderiv`i'' = "`r(deriv)'" } } if ("`r(deriv)'"!="") { local i = `i' + 1 } } local i 1 local j 1 tokenize "`names'" while (`i' <= `nest') { tempname deriv`i' local cur_name = "``i''" while (`j' <= `nest') { tempname iscorr eqpos mata `iscorr' = strpos(`rwderiv`j'',"/`cur_name'") mata st_local("found",strofreal(`iscorr')) if (`found'!=0) { mata `eqpos' = strpos(`rwderiv`j'', "=") mata `deriv`i'' = substr(`rwderiv`j'',`eqpos'+1,.) continue, break } local j = `j' + 1 } local j 1 local i = `i' + 1 } local i 1 local j 1 tokenize "`names'" while (`i' <= `nest') { mata st_local("stderiv`i'",`deriv`i'') ParseEquation "`null'" "`stderiv`i''" return local deriv`i' = "`r(equation)'" local i = `i' + 1 } end cap program drop ParseInst program ParseInst, rclass syntax anything [, noCONStant] if ("`constant'" == "") { return local constant = "constant" } else { return local constant = "nc" } unab inst : `anything', name (inst()) return local inst = "`inst'" end cap program drop ParseWMatrix program ParseWMatrix, rclass syntax [anything(name=wmatrix)] [, NWest BArtlett ANDerson GAllant PARzen /* */ quadraticspectral qs lags(string) sb_lags(string) Andrews*] local unadj 0 local hc 1 local hac 2 local clst 3 local yes 1 local no 0 if ("`quadraticspectral'"=="") { local quadraticspectral = "`qs'" } if ("`wmatrix'"=="" | "`wmatrix'"=="robust") { local wmtype = `hc' local adjf = 1 } else if ("`wmatrix'"=="cluster") { local cluster = "`options'" local wmtype = `clst' } else if ("`wmatrix'"=="unadjusted") { local wmtype = `unadj' } else if ("`wmatrix'"=="hac") { local wmtype = `hac' } else if (strmatch("`wmatrix'","hc?")) { local wmtype = `hc' local adjf = substr("`wmatrix'",3,3) } if (`wmtype'==`hac') { if ("`nwest'"!="" | "`bartlett'"!="") { local k_type = 0 } else if ("`anderson'"!="" | "`quadraticspectral'"!="" | "`andrews'"!="") { local k_type = 1 } else if ("`gallant'"!="" | "`parzen'"!="") { local k_type = 2 } if ("`lags'"=="" | "`lags'"=="auto" | "`lags'"=="automatic") { local l_auto = `yes' local l_opt = `no' local l_nlags = `no' } else if ("`lags'"=="opt" | "`lags'"=="optimal") { local l_opt = `yes' local l_auto = `no' local l_nlags = `no' } else { local l_nlags = `lags' local l_opt = `no' local l_auto = `no' } if ("`sb_lags'" == "") { if (`l_nlags' == `yes') { local sbl_auto = `yes' local sbl_opt = `no' local sbl_nlags = `no' } else { local sbl_nlags = `no' local sbl_opt = `l_opt' local sbl_auto = `l_auto' } } else { if ("`sb_lags'"=="" | "`sb_lags'"=="auto" | "`sb_lags'"=="automatic") { local sbl_auto = `yes' local sbl_opt = `no' } else if ("`sb_lags'"=="opt" | "`sb_lags'"=="optimal") { local sbl_opt = `yes' local sbl_auto = `no' } else { local sbl_opt = `no' local sbl_auto = `no' } } } return local wmtype = `wmtype' if (`wmtype' == `hac') { return local l_opt = `l_opt' return local l_auto = `l_auto' return local l_nlags = `l_nlags' return local sbl_opt = `sbl_opt' return local sbl_auto = `sbl_auto' return local k_type = `k_type' } else if (`wmtype' == `hc') { return local adjf = `adjf' } else if (`wmtype' == `clst') { return local cluster = "`cluster'" } end mata: string scalar genstest_display (real matrix ci, real scalar alpha, real scalar allpv) { if (cols (ci) == 3) { return ("") } if (rows (ci) == 0) { return ("Rejected Grid") } if (allpv) { acc = select (ci, ci[.,1] :>= alpha) if (rows (acc) <= 1) { return ("Rejected Grid") } } else { acc = ci } real scalar c0, c1 c0 = round (min (acc[.,2]), 0.001) c1 = round (max (acc[.,2]), 0.001) r = "[" + strofreal(c0) + ", " + strofreal(c1) + "]" return (r) } real matrix cishell (real matrix ci) { uni = uniqrows (ci[.,2]) out = J(1, 2, .) for (i=1; i<=rows(uni); i++) { cur = select (ci, ci[.,2] :== uni[i, 1]) min1 = min (cur[.,3]) max1 = max (cur[.,3]) out = out \ (uni[i], min1) out = out \ (uni[i], max1) } out = out[2..rows(out),.] return(out) } real matrix pvalue (real matrix value, real scalar tau, real scalar kx, real scalar kz, real scalar sb) { fh = _fopen (st_strscalar ("c(sysdir_plus)") + "crit_values.matrix", "r") if (fh < 0) { fh = fopen("crit_values.matrix", "r") } rn = kz*(kz+1)*0.5 + kz - kx + 1 stat_cur = sb ? 4 : 1 pval = J(1,1,.) cv_cur = fgetmatrix (fh) row_cur = cv_cur[rn, 2..cols(cv_cur)] vec_cur = select (row_cur, row_cur:>value[stat_cur,1]) pval = pval \ cols(vec_cur)/1000 stat_cur = 1 if (sb==1) { k = ceil(tau / 0.05 - 1) for (j=1; j<=(3*k); j++) { cv_cur = fgetmatrix (fh) } cv_cur = fgetmatrix (fh) row_cur = cv_cur[rn, 2..cols(cv_cur)] vec_cur = select (row_cur, row_cur:>value[stat_cur,1]) stat_cur = stat_cur + 1 pval = pval \ cols(vec_cur)/1000 cv_cur = fgetmatrix (fh) row_cur = cv_cur[rn, 2..cols(cv_cur)] vec_cur = select (row_cur, row_cur:>value[stat_cur,1]) stat_cur = stat_cur + 1 pval = pval \ cols(vec_cur)/1000 cv_cur = fgetmatrix (fh) row_cur = cv_cur[rn, 2..cols(cv_cur)] vec_cur = select (row_cur, row_cur:>value[stat_cur,1]) stat_cur = stat_cur + 1 pval = pval \ cols(vec_cur)/1000 } pval = pval[2..rows(pval),1] if (sb) { pval = pval[2..rows(pval),1] \ pval[1,1] } fclose (fh) return (pval) } real matrix chkstat (real scalar kz, real scalar kx, real scalar tau, real scalar a, real scalar sb) { clm = 1000-a*1000+1 rn = kz*(kz+1)*0.5 + kz - kx + 1 fh = _fopen (st_strscalar ("c(sysdir_plus)") + "crit_values.matrix", "r") if (fh < 0) { fh = fopen ("crit_values.matrix", "r") } stat_cur = 1 crit = J(1,1,.) cv_cur = fgetmatrix (fh) stat_cur = stat_cur + 1 crit = crit \ cv_cur[rn, clm] if (sb==1) { k = ceil(tau / 0.05 - 1) for (j=1; j<=(3*k); j++) { cv_cur = fgetmatrix (fh) } cv_cur = fgetmatrix (fh) stat_cur = stat_cur + 1 crit = crit \ cv_cur[rn, clm] cv_cur = fgetmatrix (fh) stat_cur = stat_cur + 1 crit = crit \ cv_cur[rn, clm] cv_cur = fgetmatrix (fh) stat_cur = stat_cur + 1 crit = crit \ cv_cur[rn, clm] } crit = crit[2..rows(crit),1] if (sb) { crit = crit[2..rows(crit),1] \ crit[1,1] } fclose (fh) return (crit) } class genstestSpec { real scalar wmtype, adjf, nest real scalar small real scalar l_opt, l_auto, l_nlags real scalar sbl_opt, sbl_auto, sbl_nlags real scalar inwmatrix, cur_t, varS, nuisS real scalar center, sb_center, hasderiv real scalar identity, twostep real scalar sb_twostep, sb_identity real matrix instr, fs_wmatrix, fs_vcov, wm1, wm2 real matrix init, gam, fs_gam pointer matrix cluster pointer matrix c1, c2 // cluster breaks real scalar sb_clst scalar kernel real matrix cur_v // For optimal HAC real matrix PHI() real matrix arsinglebreak() real matrix clarsinglebreak() real scalar qllstab() real matrix coeff() real matrix homoskedastic() real matrix HC() real matrix HAC() pointer matrix VARfilter() real matrix clustervar() real scalar optlag() real scalar sigmaj() real scalar nwest() real scalar andrews() real scalar gallant() void objf() void objfss() real rowvector twostep() real rowvector igmm() real matrix matgmm() real matrix ufunc() real matrix nuisfunc() } real matrix genstestSpec::PHI(real matrix u, real matrix z) { if (fs_wmatrix!=0) { phi = fs_wmatrix } else { if (wmtype==0) { phi = homoskedastic(u,z) } if (wmtype==1) { phi = HC(u,z) } if (wmtype==2) { phi = HAC(u,z) } if (wmtype==3) { phi = clustervar(z,u) } } return(phi) } real matrix genstestSpec::arsinglebreak(real scalar tau) { real scalar i, T, kz, start, finish T = rows(instr) kz = cols(instr) start = floor(T*tau) finish = floor(T*(1-tau)) ars = J(1,1,.) l_auto = sbl_auto l_opt = sbl_opt for(i=start;i<=finish;i++) { touse = J(i,1,1) \ J(T-i,1,0) if ((nest > 0) & (gam==0)) { nuis = matgmm(touse) } else if (nest > 0) { nuis = gam } else { nuis = 0 } ufull = ufunc(nuis) u1 = select(ufull,touse) z1 = select(instr,touse) x1 = cross (z1,z1) d1 = coeff(u1,z1) ar1 = d1'*(1/rows(z1))*x1*wm1*x1*d1 u2 = select(ufull,!touse) z2 = select(instr,!touse) x2 = cross (z2,z2) d2 = coeff (u2,z2) ar2 = d2'*(1/rows(z2))*x2*wm2*x2*d2 temp = ar1 + ar2 ars = ars \ J(1,1,temp) } ars = ars[2..rows(ars),1] arave = sum(ars)/rows(ars) arsup = max(ars) arexp = 2*ln(sum(exp(0.5*ars))/rows(ars)) genar = arave \ arexp \ arsup return(genar) } real scalar genstestSpec::qllstab(real matrix Y, real matrix Z, real matrix phi) { T = rows(Z) zk = cols(Z) oneszk = J(1,zk,1) U = Y*oneszk O = L = . eigensystem(invsym(phi),O,L) sqrtphi = O*sqrt(diag(L))*O' sqrtphi = Re(sqrtphi) V = (U:*Z)*sqrtphi D = I(T) D = D[2..T,1..T] D = D', J(T,1,0) D = I(T) - D r = 1 - (10/T) Rcol = r :^ (0::(T-1)) ones = J(1,1,1), J(1,(T-1),0) R = Toeplitz(Rcol,ones) W = R*(D*V) rtilde = r*Rcol B = invsym(rtilde'*rtilde)*rtilde'*W ehatmat = W - rtilde*B ehatmat = ehatmat :^ 2 tssre = sum(ehatmat) vt = J(T,1,1) C = invsym(vt'*vt)*vt'*V vhatmat = V - vt*C vhatmat = vhatmat :^ 2 tssrv = sum(vhatmat) qllarstab = Re(tssrv-r*tssre) return(qllarstab) } real matrix genstestSpec::coeff(real matrix y, real matrix x) { real matrix delta delta = invsym(cross(x,x))*cross(x,y) return(delta) } real matrix genstestSpec::homoskedastic(real matrix u, real matrix x) { N = rows(x) kx = cols(x) varres = (1/(N-kx))*(u'u) Phi = varres*(cross(x,x)*(1/N)) return(Phi) } real matrix genstestSpec::HC(real matrix u, real matrix x) { N = rows(x) kx = cols(x) if (!(adjf==1)) { hat = x*invsym(cross(x,x))*x' h = (diagonal(hat))' } e = (u :^ 2)' if (adjf==1) { e = (N/(N-kx))*e } if ((1= 0.97) { Stilde[i,j] = 0.97 } else if (S[i,j] <= -0.97) { Stilde[i,j] = -0.97 } else { Stilde[i,j] = S[i,j] } } } AA = U*diag(Stilde)*Vt imA = I(k)-AA mMhat = Y-X*AA' results = &mMhat \ &imA return(results) } real scalar genstestSpec::optlag (real matrix h, real scalar mrate) { real vector sigma n = rows(h) k = cols(h) w = J(k,1,1) f = h :* w' m = floor(20*(n/100)^mrate) sigma = J(1,1,.) sgm = 0 for(j=0;j<=m;j++) { for(i=j+1;i<=n;i++) { sgm = sgm + f[i,.]*(f[(i-j),.]') } sigma = sigma, (sgm/n) } sigma = sigma[2..length(sigma)] s0 = sigma[1] + 2*sum(sigma[2..length(sigma)]) s1 = 2*sum((range(1,m,1)') :* sigma[2..length(sigma)]) s2 = 2*sum((range(1,m,1) :^ 2)' :* sigma[2..length(sigma)]) qrate = mrate == (2/9) ? 1/(2 + 1) : 1/(4+1) if (mrate == (2/9)) { c = floor(1.1447 * (((s1/s0)^2)^qrate)) } else if (mrate == (2/25)) { c = floor(1.3221 * (((s2/s0)^2)^qrate)) } else if (mrate == (4/25)) { c = floor(2.6614*((s2/s0)^2)^qrate) } lags = floor(c*(n^(qrate))) if (mrate > (2/25)) { lags = round(lags) < m ? round(lags) : m } else { lags = lags < m ? lags : m } return(lags) } real scalar genstestSpec::sigmaj(real scalar j, real matrix h, real matrix w, real scalar m, real scalar T) { accum = 0 for (t=j+2;t<=T;t++) { accum = accum + (w'*h[t,.]')*(w'*h[t-j,.]') } return (accum/T) } real scalar genstestSpec::nwest(real scalar l, real scalar maxlag) { if (l_auto) { maxlag = floor(4*(cur_t/100)^(2/9)) l_nlags = maxlag } else if (l_opt && maxlag==0) { maxlag = optlag(cur_v,2/9) l_nlags = maxlag } else if (!l_opt) { maxlag = l_nlags } if (l > maxlag) { l = maxlag + 1 } wt = 1 - (l/(maxlag+1)) return(wt) } real scalar genstestSpec::andrews(real scalar l, real scalar maxlag) { if (l_auto) { maxlag = 4 * (cur_t/100)^(2/25) l_nlags = maxlag } else if (l_opt && maxlag==0) { maxlag = optlag(cur_v,2/25) l_nlags = maxlag } else if (!l_opt) { maxlag = l_nlags } z = l/(maxlag+1) th = 6.0*pi()*z/5.0 if (z == 0) { wt = 1 } else { wt = 3.0*(sin(th)/th - cos(th))/(th^2) } return(wt) } real scalar genstestSpec::gallant(real scalar l, real scalar maxlag) { if (l_auto) { maxlag = 4 * (cur_t/100)^(4/25) l_nlags = maxlag } else if (l_opt && maxlag==0) { maxlag = optlag(cur_v,4/25) l_nlags = maxlag } else if (!l_opt) { maxlag = l_nlags } z = l/(maxlag+1) if (z <= .5) { wt = 1 - 6*z^2 + 6*z^3 } else if ((z >= 0.5) & (z <= 1)) { wt = 2*(1-z)^3 } else { wt = 0 } return(wt) } void genstestSpec::objf(todo,c,crit,g,H) { z = instr u = ufunc(c) m = z'*u - (1/rows(z))*z'u w = wm1 crit = 0.5*(m'*w*m) if (todo>=1) { g = u'*z*w*z'*nuisfunc(c) } } void genstestSpec::objfss(todo,c,crit,g,H) { external GENSTEST_TOUSE z = instr z1 = select(z, GENSTEST_TOUSE) z2 = select(z, !GENSTEST_TOUSE) u = ufunc(c) u1 = select(u, GENSTEST_TOUSE) u2 = select(u, !GENSTEST_TOUSE) m1 = z1'*u1 m2 = z2'*u2 crit = 0.5*(m1'*wm1*m1 + m2'*wm2*m2) if (todo==1) { v = nuisfunc(c) v1 = select(v,GENSTEST_TOUSE) v2 = select(v,!GENSTEST_TOUSE) g = u1'*z1*wm1*z1'*v1 + u2'*z2*wm2*z2'*v2 } } void optim_wrapper(todo,c,class genstestSpec M,crit,g,H) { M.objf(todo,c,crit,g,H) } void optimss_wrapper(todo,c,class genstestSpec M,crit,g,H) { M.objfss(todo,c,crit,g,H) } real rowvector genstestSpec::twostep() { external GENSTEST_TOUSE z = instr z1 = select(z, GENSTEST_TOUSE) z2 = select(z, !GENSTEST_TOUSE) if (z1 == z) { wm1 = identity ? I(cols(z)) : invsym(cross(z1,z1)) S = optimize_init() optimize_init_evaluator(S, &optim_wrapper()) optimize_init_argument(S,1,this) optimize_init_conv_maxiter(S,25) optimize_init_tracelevel(S,"none") optimize_init_which(S, "min") optimize_init_params(S, init) if (hasderiv) { optimize_init_evaluatortype(S,"d1") } else { optimize_init_evaluatortype(S,"d0") } p = optimize(S) } else { p = fs_gam } u = ufunc(p) u1 = select(u,GENSTEST_TOUSE) if (sb_clst == 1) { origcluster = cluster cluster = c1 wm1 = invsym(PHI(u, z)) } else { wm1 = invsym(PHI(u1,z1)) } if (z1!=z) { u2 = select(u,!GENSTEST_TOUSE) if (sb_clst == 1) { cluster = c2 wm2 = invsym(PHI(u, z)) } else { wm2 = invsym(PHI(u2,z2)) } } else { wm2 = 0 } S = optimize_init() if (z1 == z) { optimize_init_evaluator(S, &optim_wrapper()) } else { optimize_init_evaluator(S, &optimss_wrapper()) } optimize_init_argument(S,1,this) optimize_init_conv_maxiter(S,25) optimize_init_tracelevel(S, "none") optimize_init_which(S, "min") optimize_init_params(S, init) if (hasderiv) { optimize_init_evaluatortype(S,"d1") } else { optimize_init_evaluatortype(S,"d0") } p = optimize(S) if (sb_clst == 1) { cluster = origcluster } return(p) } real rowvector genstestSpec::igmm() { external GENSTEST_TOUSE z = instr z1 = select(z, GENSTEST_TOUSE) z2 = select(z, !GENSTEST_TOUSE) if (z1 == z) { wm1 = identity ? I(cols(z)) : invsym(cross(z1,z1)) S = optimize_init() optimize_init_evaluator(S, &optim_wrapper()) optimize_init_argument(S,1,this) optimize_init_conv_maxiter(S,25) optimize_init_tracelevel(S,"none") optimize_init_which(S, "min") optimize_init_params(S, init) if (hasderiv) { optimize_init_evaluatortype(S,"d1") } else { optimize_init_evaluatortype(S,"d0") } p = optimize(S) } else { p = fs_gam } p1 = p :+ 1 while (norm (p1-p) > 0.00001) { u = ufunc(p) u1 = select(u,GENSTEST_TOUSE) if (sb_clst == 1) { origcluster = cluster cluster = c1 wm1 = invsym(PHI(u, z)) } else { wm1 = invsym(PHI(u1,z1)) } if (z1!=z) { u2 = select(u,!GENSTEST_TOUSE) if (sb_clst == 1) { cluster = c2 wm2 = invsym(PHI(u, z)) } else { wm2 = invsym(PHI(u2,z2)) } } else { wm2 = 0 } S = optimize_init() if (z1 == z) { optimize_init_evaluator(S, &optim_wrapper()) } else { optimize_init_evaluator(S, &optimss_wrapper()) } optimize_init_argument(S,1,this) optimize_init_conv_maxiter(S,25) optimize_init_tracelevel(S, "none") optimize_init_which(S, "min") optimize_init_params(S, init) if (hasderiv) { optimize_init_evaluatortype(S,"d1") } else { optimize_init_evaluatortype(S,"d0") } p1 = p p = optimize(S) if (sb_clst == 1) { cluster = origcluster } } return(p) } real matrix genstestSpec::matgmm(real matrix touse) { external GENSTEST_TOUSE GENSTEST_TOUSE = touse if (fs_wmatrix!=0) { wm1 = fs_wmatrix wm2 = fs_wmatrix S = optimize_init() if (touse == J(rows(touse),1,1)) { optimize_init_evaluator(S, &optim_wrapper()) } else { optimize_init_evaluator(S, &optimss_wrapper()) } optimize_init_argument(S,1,this) optimize_init_conv_maxiter(S,25) optimize_init_tracelevel(S, "none") optimize_init_which(S, "min") optimize_init_params(S, init) if (hasderiv) { optimize_init_evaluatortype(S,"d1") } else { optimize_init_evaluatortype(S,"d0") } p = optimize(S) } else if (twostep == 1) { p = twostep() } else { p = igmm() } return(p) } real matrix genstestSpec::ufunc(real rowvector b) { p = cols(b) for(i=1;i<=p;i++) { st_numscalar("GENSTEST_REALLY_TEMP",i) stata("local i = GENSTEST_REALLY_TEMP") call = st_macroexpand("GENSTEST_EST`" + "i" + "'") st_numscalar(call,b[i]) } (void) st_addvar("double", name=st_tempname()) call = st_macroexpand("qui replace " + name + " = `" + "wvar" + "'*(`" + "equation"/* */ + "') " + "if " + "`" + "touse" + "'") stata(call) u = st_data(.,name,0) st_dropvar(name) return(u) } void reParseExpression (real matrix null) { for (i=1;i<=cols(null);i++) { st_local ("new_null_temp", null[i,1]) stata ("local result = `" + "result" + "' " + "`" + "new_null_temp" + "'") } } real matrix genstestSpec::nuisfunc(real rowvector b) { real matrix g p = cols(b) for(i=1;i<=p;i++) { st_numscalar("GENSTEST_REALLY_TEMP",i) stata("local i = GENSTEST_REALLY_TEMP") call = st_macroexpand("GENSTEST_EST`" + "i" + "'") st_numscalar(call,b[i]) } for(i=1;i<=p;i++) { st_numscalar("GENSTEST_REALLY_TEMP",i) stata("local i = GENSTEST_REALLY_TEMP") (void) st_addvar("double", name=st_tempname()) call = st_macroexpand("qui replace " + name + " = `" + "wvar'*(`" + "deriv`" + "i" + "'')" + " if " + "`" + "touse" + "'") stata(call) temp = st_data(.,name,0) if (i==1) { g = temp } else { g = g, temp } st_dropvar(name) } return(g) } pointer matrix gencluster(real matrix cluster) { unique = uniqrows(cluster) Nc = rows(unique) result = J(1,1,&.) for(i=1;i<=Nc;i++) { p = getindexbyv(unique[i,1],cluster) result = result \ p } result = result[2..rows(result),1] return(result) } pointer getindexbyv(real scalar value, real matrix cluster) { /*Creates list of indicies where cluster variable has certain value*/ N = rows(cluster) result = J(1,1,.) for(i=1;i<=N;i++) { if(cluster[i,1] == value) { result = result \ i } } result = result[2..rows(result),1] p = &result return(p) } real matrix clustering(real scalar dimension, cluster) { /*Creates a touse vector for a number of clusters*/ N = rows(cluster) touse = J(dimension,1,0) for(i=1;i<=N;i++) { temp = *cluster[i,1] for(j=1;j<=rows(temp);j++) { k = temp[j,1] touse[k,1] = 1 } } return(touse) } real matrix genstestSpec::clarsinglebreak(real matrix Z, real scalar tau) { origcluster = cluster Nc = rows(cluster) begin = floor(tau*Nc) if (begin == 0) { begin = 1 } finish = floor((1-tau)*Nc) ars = J(1,1,.) sb_clst = 1 for(i=begin;i<=finish;i++) { c1 = origcluster[1..i,1] c2 = origcluster[i+1..rows(origcluster), 1] touse = clustering(rows(Z),c1) if ((nest > 0) & (gam==0)) { nuis = matgmm(touse) } else if (nest > 0) { nuis = gam } else { nuis = 0 } ufull = ufunc(nuis) u1 = select(ufull,touse) z1 = select(instr,touse) zz1 = cross(z1,z1) cluster = c1 d1 = coeff(u1,z1) v1 = zz1*invsym(clustervar(Z, ufull))*zz1 u2 = select(ufull,!touse) z2 = select(instr,!touse) cluster = c2 d2 = coeff(u2,z2) zz2 = cross (z2,z2) v2 = zz2*invsym(clustervar(Z,ufull))*zz2 temp = (1/rows(z1))*d1'*v1*d1 + (1/rows(z2))*d2'*v2*d2 ars = ars \ J(1,1,temp) } ars = ars[2..rows(ars),1] avear = sum(ars)/rows(ars) supar = max(ars) expar = 2*ln(sum(exp(0.5*ars))/rows(ars)) genar = avear \ expar \ supar sb_clst = 0 cluster = origcluster return(genar) } real matrix genstestSpec::clustervar(real matrix z, real matrix u) { Nc = rows(cluster) T = rows(z) vcov = J(cols(z),cols(z),0) touse = J(T,1,0) resid = u adjf = 1 for(i=1;i<=Nc;i++) { temp = *(cluster[i,1]) for(j=1;j<=rows(temp);j++) { tn = temp[j,1] touse[tn,1] = 1 } ut = select(u,touse) zt = select(z,touse) psi = colsum((ut :* zt)) vcov = vcov + cross(psi,psi) touse = J(T,1,0) } return(vcov*(1/T)) } end