*! tvgc v0.01 JOtero 01may2020 *! modified from NEWtvgc0.ado cfb 18nov2020 *! v1.1 cfb 11dec2020 min #boot, validate prefix, graph options *! v1.2 cfb 23dec2020 added preserve/restore to enable if qualifier *! v1.3 cfb 24dec2020 added merge to retain prefix variables if selected *! v1.4 cfb 17jan2021 added restab *! v1.5 cfb 01mar2021 adjust handling of notitle *! v1.5 cfb 24sep2022 remove mata clear *! v1.5 cfb 18jan2024 requires v14 for runiformint() capture program drop tvgc // DISABLE mata: mata clear program tvgc, rclass version 14 syntax varlist(min=2 numeric ts) [if] [in] [, /// TREND /// PREfix(string) /// P(integer -1) /// D(integer -1) /// WINdow(integer -1) /// MATRIX /// ROBUST /// BOOT(integer -1) /// SEED(integer -1) /// SIZEcontrol(integer 12) /// GRAPH EPS PDF /// noTITLE noPRINT RESTAB] marksample touse qui tsset loc tsf `r(tsfmt)' preserve qui keep if `touse' loc tvar `r(timevar)' markout `touse' `tvar' quietly tsreport if `touse' if r(N_gaps) { display in red "sample may not contain gaps" exit } if `p'==-1 { local p = 2 } if `d'==-1 { local d = 1 } if `d' > 0 { loc la "LA-" } loc vce 0 if "`robust'" != "" { loc vce 1 } tempvar en trd quietly gen `en' = _n quietly gen `trd' = sum(`touse') local lastobs = `trd'[_N] su `tvar' if `touse', mean loc t1 = strofreal(`r(min)',"`tsf'") loc t2 = strofreal(`r(max)',"`tsf'") if `window'>`lastobs' { display in red "initial window exceeds available observations" exit } if `window'<=0 { local wwid = floor(0.2*`lastobs') } else if `window'>0 { local wwid = `window' } // BS logic // if ( `boot'>0) { local bootrepl = cond(`boot'==-1, 199, max(20,`boot')) // } if `seed'!=-1 { local seednum = `seed' set seed `seednum' } local depvar : word 1 of `varlist' local xvars : list varlist - depvar local numvars : word count `varlist' local numxvars : word count `xvars' local case = cond("`trend'" == "", 1, 2) // test variable creation if selected loc mmrg = cond("`prefix'" == "", 0, 1) if "`prefix'" != "" { forvalues i = 2/`numvars' { local excli: word `i' of `varlist' confirm new var `prefix'forward_`excli' confirm new var `prefix'rolling_`excli' confirm new var `prefix'recursive_`excli' } } su `en' if `trd'>0 & !mi(`trd') & `touse', mean loc first = `r(min)' loc last = `r(max)' - `wwid' loc full = `r(max)' * generate the variables in the lag-augmented part if (`d'>0) { local pp1 = `p' + 1 local ppd = `p' + `d' } local ppd1 = `p' + `d' + 1 local numvars : word count `varlist' forvalues i = 2/`numvars' { mata: mat_`i' = J(`last'+1,`last'+1,.) } loc rhsv local lagvarlistx local lagvarlistz forvalues i = 1/`numvars' { tempvar tmpvarx`i' local var`i' : word `i' of `varlist' qui gen `tmpvarx`i'' = `var`i'' if (`i' > 1) { loc rhsv "`rhsv' `var`i''" } } forvalues ii = 1/`p' { forvalues i = 1/`numvars' { tempvar lagtmpvarx`i'_`ii' qui gen `lagtmpvarx`i'_`ii'' = l`ii'.`tmpvarx`i'' local lagvarlistx `lagvarlistx' `lagtmpvarx`i'_`ii'' } } if (`d'>0) { forvalues ii = `pp1'/`ppd' { forvalues i = 1/`numvars' { tempvar lagtmpvarz`i'_`ii' qui gen `lagtmpvarz`i'_`ii'' = l`ii'.`tmpvarx`i'' local lagvarlistz `lagvarlistz' `lagtmpvarz`i'_`ii'' } } } loc en = `full'-`first'+1 forv t = `first'/`full' { loc wo = `t' + `wwid' - 1 forv tt = `wo'/`full' { mata: tvgcprep("`varlist'","`lagvarlistx'","`lagvarlistz'","`trd'",`case',`t',`tt',`ppd1',`numvars',`p',`d',`vce') forvalues i = 2/`numvars' { local var`i' : word `i' of `varlist' mata: mat_`i'[`t',(`tt'-`wwid'+1)] = Wstat[`i'] } } } local fr = `full'-`last'-1 if "`trend'" == "" { display as res _n "Time-varying `la'VAR Granger causality test, `t1' - `t2'" } else if "`trend'" == "trend" { display as res _n "Time-varying `la'VAR Granger causality test including trend, `t1' - `t2'" } // set trace on mata: gcres = J(`=`numvars'-1',3, .) loc wn =`sizecontrol' * 6 tvgc_boot `varlist', `trend' p(`p') d(`d') boot(`bootrepl') window(`wn') sizecontrol(`sizecontrol') mata: tvgc_cv() forvalues i = 2/`numvars' { local excli: word `i' of `varlist' loc i1 = `i'-1 mata: tvgcres(`i1', mat_`i') if "`prefix'" != "" { qui gen `prefix'forward_`excli' = . qui gen `prefix'rolling_`excli' = . qui gen `prefix'recursive_`excli' = . loc vn `prefix'forward_`excli' `prefix'rolling_`excli' `prefix'recursive_`excli' mata: st_view(v`i'=.,((`fr'+1)::`full'),tokens(st_local("vn"))) mata: v`i'[.,.]= mat_`i'[1,.]', diagonal(mat_`i'), colmax(mat_`i')' // set trace on // graph only available with prefix if "`graph'" != "" { loc t1g "Forward expanding Wald test" loc t2g "Rolling Wald test" loc t3g "Recursive expanding Wald test" lab var `tvar' " " loc j 0 tempvar l90 l95 qui g `l90'=. qui g `l95'=. // di in r "title | `title' " foreach v of local vn { loc j = `j'+1 loc rv : word `i' of `varlist' loc cv90 = cv90[`i1',`j'] loc cv95 = cv95[`i1',`j'] qui replace `l90' = `cv90' if !mi(`v') qui replace `l95' = `cv95' if !mi(`v') // loc t1txt "" // loc t2txt "" if "`title'" == "notitle" { tsline `v' `l90' `l95' if !mi(`v'), ylab(,angle(0) labs(small)) /// scheme(s2mono) graphregion(color(white)) yti(" ") /// xlab(#6, labs(small)) name(`v', replace) legend(off) } if "`title'" != "notitle" { tsline `v' `l90' `l95' if !mi(`v'), ylab(,angle(0) labs(small)) /// scheme(s2mono) graphregion(color(white)) yti(" ") /// xlab(#6, labs(small)) /// title("`t`j'g' for `depvar' G-caused by `rv', `t1' - `t2'", size(medium)) /// subtitle("with 90th (--) and 95th (-) percentiles of bootstrapped test statistics", size(small)) /// name(`v', replace) legend(off) // loc t1txt "`t`j'' for `depvar' G-caused by `rv', `t1' - `t2'" // loc t2txt "with 90th (--) and 95th (-) percentiles of bootstrapped test statistics" } if "`eps'" != "" { qui gr export `v'.eps, replace } if "`pdf'" != "" { qui gr export `v'.pdf, replace } } } // end of graph logic } // end of prefix logic if "`matrix'" != "" { mata: st_matrix("m_`excli'",mat_`i') return matrix m_`excli' = m_`excli' } } // end loop over RHS vars loc cn Max_Wald_forward Max_Wald_rolling Max_Wald_recursive mata: st_matrix("gcres",gcres) mat rownames gcres = `rhsv' mat colnames gcres = `cn' matlist gcres, for(%18.3f) tit("TVGC `robust' test statistics for H0: `depvar' is GC") mat gcres0 = gcres return matrix gcres = gcres mata: st_local("T", strofreal(rows(mat_2))) foreach c in 90 95 99 { mat rownames cv`c' = `rhsv' mat colnames cv`c' = `cn' matlist cv`c', for(%18.3f) tit("`c'th percentile of test statistics [`bootrepl' replications]") mat gccv`c'0 = cv`c' // matlist gccv`c'0 return matrix gccv`c' = cv`c' } if "`restab'" == "restab" { // output of 95, 99 loc nrow = (`numvars'-1) * 3 mat restab = J( `nrow',3, .) loc k 1 forv i=2/`numvars' { forv j=1/3 { mat restab[`k',`j'] = gcres0[`=`i'-1',`j'] mat restab[`=`k'+1',`j'] = gccv950[`=`i'-1',`j'] mat restab[`=`k'+2',`j'] = gccv990[`=`i'-1',`j'] } loc k = `k' + 3 } // matlist restab, format(%9.3f) // from outtable.tex tempname hh s1 s2 s3 file open `hh' using restab.tex, write replace file write `hh' "\begin{table}[htbp]" _n local hl "\hline" loc fmt %9.3f local l "l" local align "c" loc nc 4 forv i=2/`nc' { local l "`l'`align'" } file write `hh' "\begin{tabular}{`l'}" _n "\toprule" _n file write `hh' "Variable & Max Wald FE & Max Wald RO & Max Wald RE \\ [2pt] \midrule" _n loc k 1 forv i=1/`=`numvars'-1' { loc vn: word `i' of `rhsv' file write `hh' "\$`vn'\$ & " `fmt' (restab[`k',1]) " & " `fmt' (restab[`k',2]) " & " `fmt' (restab[`k',3]) " \\" _n sca `s1' = "("+strofreal(restab[`=`k'+1',1],"`fmt'")+")" sca `s2' = "("+strofreal(restab[`=`k'+1',2],"`fmt'")+")" sca `s3' = "("+strofreal(restab[`=`k'+1',3],"`fmt'")+")" file write `hh' " & " (`s1') " & " (`s2') " & " (`s3') " \\" _n sca `s1' = "["+strofreal(restab[`=`k'+2',1],"`fmt'")+"]" sca `s2' = "["+strofreal(restab[`=`k'+2',2],"`fmt'")+"]" sca `s3' = "["+strofreal(restab[`=`k'+2',3],"`fmt'")+"]" file write `hh' " & " (`s1') " & " (`s2') " & " (`s3') " \\ [3pt] " _n loc k=`k'+3 } file write `hh' "\bottomrule " _n "\end{tabular}" _n file write `hh' "\end{table}" _n file close `hh' } // end restab return local cmd "tvgc" return scalar T = `T' return scalar p = `p' return scalar d = `d' return scalar bootrepl = `bootrepl' return scalar window = `wwid' return scalar sizecontrol = `sizecontrol' return local depvar `depvar' return local rhsvars `rhsv' return local tstart `t1' return local tend `t2' tempfile mm qui save "`mm'" restore // merge in created variables if selected if `mmrg' == 1 { qui merge 1:1 `tvar' using "`mm'" drop _merge } end mata: void tvgcres(real rhs, real matrix gctest) { external real matrix gcres gcres[rhs,1] = max(gctest[1, ]) gcres[rhs,2] = max(diagonal(gctest)) gcres[rhs,3] = max(gctest) } void tvgc_cv() { // each triple of columns in bsmat refer to test stats for one RHS variable external real matrix bsmat rcv90 = mm_quantile(bsmat, 1, 0.90) rcv95 = mm_quantile(bsmat, 1, 0.95) rcv99 = mm_quantile(bsmat, 1, 0.99) st_matrix("cv90",colshape(rcv90,3)) st_matrix("cv95",colshape(rcv95,3)) st_matrix("cv99",colshape(rcv99,3)) } void tvgcprep(string scalar varlist, string scalar lagvarlistx, string scalar lagvarlistz, string scalar trend, real kase, real t, real tt, real ppd1, real numvars, real p, real d, real vce) { external real matrix Y, X, Z, Wstat external real scalar n st_view(Y=., ., tokens(varlist)) Y = Y[t..tt,.] Y = Y[|ppd1,1 \ rows(Y),cols(Y)|] st_view(X=., ., tokens(lagvarlistx)) X = X[t..tt,.] X = X[|ppd1,1 \ rows(X),cols(X)|] if (d>0) { st_view(Z=., ., tokens(lagvarlistz)) Z = Z[t..tt,.] Z = Z[|ppd1,1 \ rows(Z),cols(Z)|] } st_view(trd=., ., tokens(trend)) trd=trd[t..tt,.] trd = trd[|ppd1,1 \ rows(trd),cols(trd)|] n = rows(Y) iota=J(n,1,1) if (kase==1) { tau=iota } else { tau= trd, iota } if (d==0) { X = X, tau } else { X = X, Z, tau } Wstat=J(numvars,1,.) y = Y[.,1] XpX = quadcross(X, X) XpXi = invsym(XpX) bj = XpXi*quadcross(X, y) ej = y - X*bj ejsq = ej:^2 sigmasq = quadsum(ejsq)/n if (vce==0) { omega = sigmasq :* XpXi } else { omega = XpXi * quadcross(X, ejsq, X) * XpXi } _makesymmetric(omega) for(j=2;j<=numvars;j++){ R=J(p,rows(bj),0) R[1,j]=1 if (p>1) { for(k=2;k<=p;k++) { R[k,(j+(k-1)*numvars)]=1 } } Wstat[j] = ((bj'*R')*invsym(R*omega*R')*(R*bj)) } } end