Tsay JASA 1998 Threshold VAR

[TomDoan]

These are "replications" for Tsay (1998), "Testing and Modeling Multivariate Threshold Models", J. of American Statistical Assn, vol. 93, no. 443, 1188-1202. Note that the calculations in the paper had a serious programming error with the log likelihoods being computed incorrectly. (I'm fairly sure that instead of using - T/2 log |Sigma|, they used - T/2 log |T Sigma|). This forced the estimates towards breaks that are roughly equally spaced. The calculations here use the proper likelihoods.

The theory in the paper and the US interest rates example is covered in detail in the 2nd Edition of the Structural Breaks and Switching Models e-course

Zip with programs and data files

Detailed description

[Trankied]

Code: Select all

@gridseries(from=-.30,to=.05,n=80,pts=ngrid) rgrid
*
compute bestaic=%na
do i=1,ngrid-19
   do j=i+20,ngrid
      sweep(group=(thresh<rgrid(i))+(thresh<rgrid(j)),var=homo)
      # g3year g3month
      # constant g3year{1 to p} g3month{1 to p}
      compute thisaic=-2.0*%logl+2.0*%nregsystem
      if .not.%valid(bestaic).or.thisaic<bestaic
         compute bestaic=thisaic,bestlower=rgrid(i),bestupper=rgrid(j)
   end do j
end do i
*

I don't know exactly why set i=1,ngrid-19 and j=i+20
how to determine 19,20 and what's it means?

thanks for helping

[Trankied]

thanks for answer
BTW,the code above

Code: Select all

do i=1,ngrid-19
   do j=i+20,ngrid

I would like to know why set i=1,ngrid-19 and j=i+20,ngrid.
Could you please interpret what they means?

thans for answering

[TomDoan]

If you're doing a grid search for two breaks, the second break point has to come after the first. You also want to make sure that the two have a gap between them.
do i=1,ngrid-19
do j=i+20,ngrid

has I as the subscript for the first break point and J as the second. In order to allow for 20 grid points between them, I can't go beyond ngrid-19, and J has to start at I+20.

[Trankied]

I am appreciate your great help. It's very useful for me.

I have another question:
If I wanna estimate two threshold value for three-regime model, first threshold value get by [-0.3,-0.2],the other by [-0.15,0.05]
How should I code on this example?

[TomDoan]

I am appreciate your great help. It's very useful for me.

I have another question:
If I wanna estimate two threshold value for three-regime model, first threshold value get by [-0.3,-0.2],the other by [-0.15,0.05]
How should I code on this example?

The overall lower and upper bounds are already -.3 and +.05 (note the @GridSeries). If you want to put hard restrictions on the two ranges, you just have to figure out the upper bound for I, and the lower bound for J that correspond to grid values of -.2 and -.15 respectively.

[Trankied]

so if I wanna choose lower break by 15%~40% and upper break by 57%~85%,
If there are 100 data available, could I set i=15,40 j=57,85 and n=100 ?
instead of estimating two break separately, I wanna estimate them at the same time and use AIC to determine it.
How should I code it?

[Trankied]

another question:
if I would like to set threshold is delayed X.
I should

Code: Select all

set thresh =X{1} 

or

Code: Select all

set thresh = X{-1}

[TomDoan]

so if I wanna choose lower break by 15%~40% and upper break by 57%~85%,
If there are 100 data available, could I set i=15,40 j=57,85 and n=100 ?
instead of estimating two break separately, I wanna estimate them at the same time and use AIC to determine it.
How should I code it?

The FROM and TO limits on the @GRIDSERIES are already set to exclude the most extreme values since you need a certain amount of data just to estimate the models. I'm not completely sure what you mean by the percentages, but the Tsay coding is using an equally-spaced grid, rather than a grid based upon empirical values. The principal reason for that is that the number of potential regressions for an empirical grid was too large to be feasible at the time the paper was written.

[TomDoan]

another question:
if I would like to set threshold is delayed X.
I should

Code: Select all

set thresh =X{1} 

or

Code: Select all

set thresh = X{-1}

RATS uses x{n} for lag n, so you want x{1}. x{-1} would be for lead 1.

[yizhuorui]

Dear tom,
I choose a SUR instruction to pick the threshold since my model is a near-VAR.(In fact, i identified error correction term through the Cats and restricted some variable to weekly exogeneity,and added with some dummy variable ,and then every equation become diffrent ,so i have to choose to SUR instruction instead of sweep instruction)

Code: Select all

equation wtii wti
# constant  exchangerate{1 to p} CPI{1 to p} output{1 to p} wti{1 to p} monetary{1 to p} dummywti1 dummywti2 dummywti3 dummywti4 dummywti5
equation exchangeratee exchangerate
# constant  exchangerate{1 to p} CPI{1 to p} output{1 to p} wti{1 to p} monetary{1 to p} dummyexchangerate1 dummyexchangerate2 dummyexchangerate3
equation cpii cpi
# constant  exchangerate{1 to p} CPI{1 to p} output{1 to p} wti{1 to p} monetary{1 to p}  dummycpi1 dummycpi2 $
dummycpi3 dummycpi4 dummycpi5 dummycpi6 dummycpi7 dummycpi8 dummycpi9 dummycpi10 dummycpi11 dummycpi12
equation outputt output
# constant  exchangerate{1 to p} CPI{1 to p} output{1 to p} wti{1 to p} monetary{1 to p} ect dummyoutput1 dummyoutput2 $
dummyoutput3 dummyoutput4 dummyoutput5 dummyoutput6 dummyoutput7 dummyoutput8 dummyoutput9 dummyoutput10 dummyoutput11 $
dummyoutput12 dummyoutput13 dummyoutput14 dummyoutput15 dummyoutput16 dummyoutput17 dummyoutput18
equation monetaryy monetary
# constant  exchangerate{1 to p} CPI{1 to p} output{1 to p} wti{1 to p} monetary{1 to p} dummymonetary1 dummymonetary2
group yzr wtii exchangeratee cpii outputt monetaryy
set thresh = wti{7}
@gridseries(from=-0.083106944,to=0.060975559,n=147,pts=ngrid) rgrid
set aic 1 ngrid = 0.0
*
compute bestaic=%na
set aic 1 ngrid = 0.0
set bic 1 ngrid = 0.0
set cic 1 ngrid = 0.0
*
compute rtest=rgrid(i)
sur(model=yzr,noprint,nosigma,smpl=thresh>rtest) 

compute bic(i)=-2.0*%logl+2.0*71
sur(model=yzr,noprint,nosigma,smpl=thresh<rtest) 

compute cic(i)=-2.0*%logl+2.0*71
compute aic(i) = bic(i)+bic(i)
if i==1.or.aic(i)<bestaic
compute bestaic=aic(i),bestbreak=rtest
end do i
*
*
scatter(footer=$
"AIC Values for Interest Rates Regression vs Break Point")
# rgrid aic
disp "Best Break is" bestbreak "with AIC" bestaic

However ,the output demonstrated “## REG13. Singular Regressions - Check for Collinearity among Rows 1 to 1.”。what should i do to fix this error?(when i didn't add some dummy variable,the result is the same.) I attached with my data

[TomDoan]

Without a full program, I can't really tell, though I assume you want ECT{1}, not current ECT, unless you defined ECT with a lag to start. I'm not really understanding the point of dozens of dummies with a threshold model. The dummies and the threshold model are two rather distinct ways to deal with structural breaks and they don't necessarily work together well. You may very well be running into SUR's that can't be estimated due to the lag of data points in a partition.

[Bach]

Hi Tom,

I am using the Tsay threhold program to determine the threshold for my 3 variable VAR ( using one of the Variable in the model, output) as the threshold. I am however getting an error message. Does the TSAY code allow such configuration ? Do you see any issue with the code below ?

Code: Select all

   *
compute sstart = 2000:3
compute send   = 2014:1
*

set dgdp = gdp - gdp{1}
set dspen = spen - spen{1}
set dtax = tax- tax{1}

set dgdp100 =  100*dgdp
set dspen100 = 100*dspen
set dtax100 =  100*dtax

@dfunit(lags=12) dgdp
@dfunit(lags=12) dspen
@dfunit(lags=12) dtax


spgraph(vfields=3,$
  footer="Figure 3. Time Plots of Growth Series of Tunisian Fiscal Variables")
graph(vlabel="Output")
# dgdp100
graph(vlabel="Government Spending")
# dspen100
graph(vlabel="Tax")
# dtax
spgraph(done)
*
@VARLagSelect(lags=12,crit=aic)
# dspen dgdp dtax
*
compute p =2
compute k =3
*
do d=1,3
   dofor m0 = 50 100
      set thresh = dgdp100{d}
      *
      rls(noprint,order=thresh,condition=m0) dspen100 / rrspen
      # constant dspen100{1 to p} dgdp100{1 to p} dtax100{1 to p}
      rls(noprint,order=thresh,condition=m0) dtax100 / rrtax
      # constant dspen100{1 to p} dgdp100{1 to p} dtax100{1 to p}

      rls(noprint,order=thresh,condition=m0) dgdp100 / rrgdp
      # constant dspen100{1 to p} dgdp100{1 to p} dtax100{1 to p}

      * We need to exclude the conditioning observations, so we generate
      * the series of ranks of the threshold variable over the
      * estimation range.
      *
      order(ranks=rr) thresh %regstart() %regend()
      *
      linreg(noprint,smpl=rr>m0) rrspen / wrspen
      # constant dspen100{1 to p} dgdp100{1 to p} dtax100{1 to p}
      linreg(noprint,smpl=rr>m0) rrtax / wrtax
      # constant dspen100{1 to p} dgdp100{1 to p} dtax100{1 to p}

      linreg(noprint,smpl=rr>m0) rrgdp / wrgdp
      # constant dspen100{1 to p} dgdp100{1 to p} dtax100{1 to p}

      ratio(mcorr=%nreg,degrees=k*%nreg,noprint)
      # rrspen rrtax
      # wrspen wrtax

      # rrgdp rrspen
      # wrgdp wrspen

      # rrgdp rrtax
      # wrgdp wrtax

      disp "D=" d "m0=" m0 @16 "C(d)=" *.## %cdstat @28 "P-value" #.##### %signif
   end dofor m0
end do d
*
* Evaluate the AIC across a grid of threshold settings
*
set thresh = dgdp{1}
@gridseries(from=-4,to=3,n=300,pts=ngrid) rgrid
set aic 1 ngrid = 0.0
*
compute bestaic=%na
*
do i=1,ngrid
   compute rtest=rgrid(i)
   sweep(group=thresh<rtest,var=homo)
   # dspen100 dgdp100 dtax100
   # constant dspen100{1 to p} dgdp100{1 to p} dtax100{1 to p}
   compute aic(i)=-2.0*%logl+2.0*%nregsystem
   if i==1.or.aic(i)<bestaic
      compute bestaic=aic(i),bestbreak=rtest
end do i
*
scatter(footer="AIC Values for Interest Rates Regression vs Break Point")
# rgrid aic
disp "Best Break is" bestbreak "with AIC" bestaic
*

***** /////***** to estimate double break
set thresh = dgdp{1}
*
* This is a slightly different set of grids than used in the paper.
*
@gridseries(from=-4,to=5,n=80,pts=ngrid) rgrid
*
compute bestaic=%na
do i=1,ngrid-19
   do j=i+20,ngrid
      sweep(group=(thresh<rgrid(i))+(thresh<rgrid(j)),var=homo)
       # dspen100 dgdp100 dtax100
      # constant dspen100{1 to p} dgdp100{1 to p} dtax100{1 to p}
      compute thisaic=-2.0*%logl+2.0*%nregsystem
      if .not.%valid(bestaic).or.thisaic<bestaic
         compute bestaic=thisaic,bestlower=rgrid(i),bestupper=rgrid(j)
   end do j
end do i
*
disp "Best double break is at" bestlower "and" bestupper "with AIC" bestaic
       

Alternatively, if I now use the TsayTest procedure to do the same, would the next instruction produce the same accurate result ( or should I remove the dgdp ( which is the threshold) from
# constant dspen100{1 to p} dgdp100{1 to p} dtax100{1 to p} ? as I am getting d=1 is significant under the latter scenario but insignificant if on the ocntrary I include the dgdp in the equation .

Thanks

Code: Select all

 
do d=1,4
  set thresh = dgdp{d}
  @tsaytest(title="Threshold Test with Delay="+d,thresh=thresh) dgdp100
 # constant dspen100{1 to p} dgdp100{1 to p} dtax100{1 to p}
end do d
*
* And with reversed ordering
*
do d=1,4
  set thresh = -dgdp{d}
  @tsaytest(title="Threshold Test Reversed with Delay="+d,thresh=thresh) dgdp
  # constant dspen100{1 to p} dgdp100{1 to p} dtax100{1 to p}
end do d
  

[TomDoan]

Tsay can be extended to three (or more). You just need to do two sets of three series on the RATIO instruction, just adding the corresponding series for the third variable.

[Bach]

Thanks Tom for the quick response. I changed the code as suggested which now reads as follows;

Code: Select all

  
      ratio(mcorr=%nreg,degrees=k*%nreg,noprint)
      # rrspen rrgdp rrtax
      # wrspen wrgdp wrtax
 

However , I am now getting this error message '

Code: Select all

The Error Occurred At Location 1698, Line 29 of loop/block
## REG6. Matrix VCV-1 Is Singular. Check for Collinearity in First 3 series
The Error Occurred At Location 1698, Line 29 of loop/block
## REG6. Matrix VCV-2 Is Singular. Check for Collinearity in First 3 series
The Error Occurred At Location 1698, Line 29 of loop/block
## M1. MCORR = 7 is Illegal. Value should be between 1 and 1
D= 1 m0= 50    C(d)= NA    P-value   NA
## SR10. Missing Values And/Or SMPL Options Leave No Usable Data Points
The Error Occurred At Location 1118, Line 20 of loop/block

Any idea what the cause of the above might be ?
Thanks