/*

MHEGY: Seasonal Unit Root Tests for Monthly Data

 @MHEGY(options) series start end

   MHEGY performs the HEGY test for seasonal unit roots
   for a series of monthly data. It uses the procedure
   suggested by J.J. Beaulieu and J.A. Miron (J. of
   Econometrics 55, 1993, pp. 305-28). See their article
   for details on the test.


 Parameters

   series     The series, which is tested for seasonal
              cointegration.

   start end  Range to use for the test. If you have not
              set a SMPL, this defaults to the largest
              range of series. Note that these
              parameters define the range of data to use.
              The range for the regression might be
              smaller due to lags. See below.


 Options

   DET=NONE/[CONSTANT]/TREND
     This option defines the deterministic components
     (beneath a possible det. seasonal) of the test
     regression. Note that with DET=TREND a linear trend
     and a constant are included.

   DSEAS/[NODSEAS]
     This option determines, whether a deterministic
     seasonal is included in the test regression beneath
     the tested unit root. I strongly suggest to include
     a deterministic seasonal, because experience has
     shown, that most economic time series have
     deterministic (and no unit root) seasonality.

   CRIT=AIC/[BIC]
     With this option you select the criterion, which
     determines the choice of the lag length for the
     test regression.

   MAXLAG=max. # of lags [12]
     The maximum number of lags for the test regression.


 Notes

   Before you SOURCE this procedure, you have to SOURCE
   the LAGSELECT procedure (written by Norm Morin,
   available on www.estima.com). I use this wonderful
   procedure to choose the lag length for the test
   regression.

   If you want to perform the original HEGY test for a
   quarterly time series, use the procedure HEGY.SRC by
   Suliman Al-Turki. It's also available on
   www.estima.com.


 Example

   calendar 1960 1 12
   allocate 12 1997:12
   open data m.prn
   data(format=prn,org=obs,missing=0) 1960:1 1997:4 tu0800
   log tu0800 / m
   source(noecho) lagselec.src
   source(noecho) mhegy.src
   @mhegy(det=trend) m


 Version
   1.0 1999/06/07
   2.0 March 27, 2003 -- Fixed filters to run from 1-12, not 0-11. By Estima

 Author
   Ulrich Leuchtmann
   INVESCO Asset Management, Frankfurt, Germany
   email ulrich_leuchtmann@fra.invesco.com

*/
procedure mhegy series start end
type series series
type integer start end
*
option choice  det     2 none constant trend
option switch  dseas   1
option choice  crit    2 aic bic
option integer maxlag 12
*
local integer startl endl nobs span lags ii i
local real fac F0304 F0506 F0708 F0910 F1112
local vector tpi
local series y1 y2 y3 y4 y5 y6 y7 y8 y9 y10 y11 y12 y13
local series d1 d2 d3 d4 d5 d6 d7 d8 d9 d10 d11 d12 trend
dimension tpi(12)


*

inquire(series=series) startl>>start endl>>end
comp nobs = endl - startl + 1
display ' '
display 'HEGY-Test on Seasonal Integration (Monthly Version)'
display ' for series ' %label(series) 'from' %datelabel(startl) 'to' %datelabel(endl)
inquire(seasonal) span
if (span.NE.12) {
 display 'Use procedure MHEGY for monthly data only!'
 display ' '
 return
}
display ' Number of observations:' nobs '  Regression based on' (nobs-12) 'observations.'
display ' '
*
set trend = t
if dseas {
 seasonal  d1 startl endl 12  1
 seasonal  d2 startl endl 12  2
 seasonal  d3 startl endl 12  3
 seasonal  d4 startl endl 12  4
 seasonal  d5 startl endl 12  5
 seasonal  d6 startl endl 12  6
 seasonal  d7 startl endl 12  7
 seasonal  d8 startl endl 12  8
 seasonal  d9 startl endl 12  9
 seasonal d10 startl endl 12 10
 seasonal d11 startl endl 12 11
 seasonal d12 startl endl 12 12
}
*
filter series (startl+12) endl y1
#  1 to 12
# 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0 1.0
filter series (startl+12) endl y2
#  1 to 12
# -1.0  1.0 -1.0  1.0 -1.0  1.0 -1.0  1.0 -1.0  1.0 -1.0  1.0
filter series (startl+12) endl y4
#  1 to 12
#  1.0  0.0 -1.0  0.0  1.0  0.0 -1.0  0.0 -1.0  0.0 -1.0  0.0
set y3 (startl+12) endl = y4{1}
comp fac = sin(2.0*%pi/3.0)
filter series (startl+12) endl y6
#  1 to 12
#  fac -fac  0.0  fac -fac  0.0  fac -fac  0.0  fac -fac  0.0
set y5 (startl+12) endl = -(sqrt(3.0)/3.0) * (y6 + 2.0 * y6{1})
filter series (startl+12) endl y8
#  1 to 12
#  fac  fac  0.0 -fac -fac  0.0  fac  fac  0.0 -fac -fac  0.0
set y7 (startl+12) endl = -(sqrt(3.0)/3.0) * (y8 - 2.0 * y8{1})
filter series (startl+12) endl y10
#  1 to 12
#  0.5 -fac  1.0 -fac  0.5  0.0 -0.5  fac -1.0  fac -0.5  0.0
set y9 (startl+12) endl = -(sqrt(3.0) * y10 + 2.0 * y10{1})
filter series (startl+12) endl y12
#  1 to 12
#  0.5  fac  1.0  fac  0.5  0.0 -0.5 -fac -1.0 -fac -0.5  0.0
set y11 (startl+12) endl = -(sqrt(3.0) * y12 - 2.0 * y12{1})
filter series (startl+12) endl y13
#  0 12
#  1.0 -1.0
*
if (det.eq.1) {
   if dseas {
      display ' Deterministic components: deterministic seasonals.'
      @lagselect(crit=crit,det=none,exog,$
                 maxlag=maxlag,noprint) y13 (startl+12) endl
      # d1 d2 d3 d4 d5 d6 d7 d8 d9 d10 d11
      comp lags = %imax(lagaic,lagbic)
      if lags {
         linreg(noprint) y13 (startl+12) endl
         # y1 y2 y3 y4 y5 y6 y7 y8 y9 y10 y11 y12 $
           d1 d2 d3 d4 d5 d6 d7 d8 d9 d10 d11 y13{1 to lags}
       }
      else {
         linreg(noprint) y13 (startl+12) endl
         # y1 y2 y3 y4 y5 y6 y7 y8 y9 y10 y11 y12 $
           d1 d2 d3 d4 d5 d6 d7 d8 d9 d10 d11
       }
    }
   else {
      display ' Deterministic components: none.'
      @lagselect(crit=crit,det=none,noexog,$
                 maxlag=maxlag,noprint) y13 (startl+12) endl
      comp lags = %imax(lagaic,lagbic)
      if lags {
         linreg(noprint) y13 (startl+12) endl
         # y1 y2 y3 y4 y5 y6 y7 y8 y9 y10 y11 y12 $
           y13{1 to lags}
       }
      else {
         linreg(noprint) y13 (startl+12) endl
         # y1 y2 y3 y4 y5 y6 y7 y8 y9 y10 y11 y12
       }
    }
 }
else if (det.eq.2) {
   if dseas {
      display ' Deterministic components: constant, deterministic seasonals.'
      @lagselect(crit=crit,det=constant,exog,$
                 maxlag=maxlag,noprint) y13 (startl+12) endl
      # d1 d2 d3 d4 d5 d6 d7 d8 d9 d10 d11
      comp lags = %imax(lagaic,lagbic)
      if lags {
         linreg(noprint) y13 (startl+12) endl
         # y1 y2 y3 y4 y5 y6 y7 y8 y9 y10 y11 y12 $
           constant d1 d2 d3 d4 d5 d6 d7 d8 d9 d10 d11 y13{1 to lags}
       }
      else {
         linreg(noprint) y13 (startl+12) endl
         # y1 y2 y3 y4 y5 y6 y7 y8 y9 y10 y11 y12 $
           constant d1 d2 d3 d4 d5 d6 d7 d8 d9 d10 d11
       }
    }
   else {
      display ' Deterministic components: constant.'
      @lagselect(crit=crit,det=constant,noexog,$
                 maxlag=maxlag,noprint) y13 (startl+12) endl
      comp lags = %imax(lagaic,lagbic)
      if lags {
         linreg(noprint) y13 (startl+12) endl
         # y1 y2 y3 y4 y5 y6 y7 y8 y9 y10 y11 y12 $
           constant y13{1 to lags}
       }
      else {
         linreg(noprint) y13 (startl+12) endl
         # y1 y2 y3 y4 y5 y6 y7 y8 y9 y10 y11 y12 $
           constant
       }
    }
 }
else {
   if dseas {
      display ' Deterministic components: constant, trend, deterministic seasonals.'
      @lagselect(crit=crit,det=trend,exog,$
                 maxlag=maxlag,noprint) y13 (startl+12) endl
      # d1 d2 d3 d4 d5 d6 d7 d8 d9 d10 d11
      comp lags = %imax(lagaic,lagbic)
      if lags {
         linreg(noprint) y13 (startl+12) endl
         # y1 y2 y3 y4 y5 y6 y7 y8 y9 y10 y11 y12 $
           constant trend d1 d2 d3 d4 d5 d6 d7 d8 d9 d10 d11 y13{1 to lags}
       }
      else {
         linreg(noprint) y13 (startl+12) endl
         # y1 y2 y3 y4 y5 y6 y7 y8 y9 y10 y11 y12 $
           constant trend d1 d2 d3 d4 d5 d6 d7 d8 d9 d10 d11
       }
    }
   else {
      display ' Deterministic components: constant, trend.'
      @lagselect(crit=crit,det=trend,noexog,$
                 maxlag=maxlag,noprint) y13 (startl+12) endl
      comp lags = %imax(lagaic,lagbic)
      if lags {
         linreg(noprint) y13 (startl+12) endl
         # y1 y2 y3 y4 y5 y6 y7 y8 y9 y10 y11 y12 $
           constant trend y13{1 to lags}
       }
      else {
         linreg(noprint) y13 (startl+12) endl
         # y1 y2 y3 y4 y5 y6 y7 y8 y9 y10 y11 y12 $
           constant trend
       }
    }
 }
*
do ii=1,12
 comp tpi(ii) = %beta(ii)/sqrt(%seesq*%xx(ii,ii))
end do
exclude(noprint)
# y3 y4
compute F0304 = %cdstat
exclude(noprint)
# y5 y6
compute F0506 = %cdstat
exclude(noprint)
# y7 y8
compute F0708 = %cdstat
exclude(noprint)
# y9 y10
compute F0910 = %cdstat
exclude(noprint)
# y11 y12
compute F1112 = %cdstat
*
if (crit.eq.1)
{
  display ' Lag selection using the AIC criterion.'
  display ' Minimum AIC at lag' lagaic
}
else
{
  display ' Lag selection using the BIC criterion.'
  display ' Minimum BIC at lag' lagbic
}
display ' '
display '              Cycles        t-Tests               '
display ' Frequency   per Year   pi_odd   pi_even   F-Test '
display '**************************************************'
display @2 '  0  ' @13 '   0   '  @24 ###.## tpi(1)  @33 '   -  '       @43 '   -  '
display @2 ' pi  ' @13 '   6   '  @24 '   -  '       @33 ###.## tpi(2)  @43 '   -  '
display @2 ' pi/2' @13 '3 and 9'  @24 ###.## tpi(3)  @33 ###.## tpi(4)  @43 ###.## F0304
display @2 '2pi/3' @13 '8 and 4'  @24 ###.## tpi(5)  @33 ###.## tpi(6)  @43 ###.## F0506
display @2 ' pi/3' @13 '2 and 10' @24 ###.## tpi(7)  @33 ###.## tpi(8)  @43 ###.## F0708
display @2 '5pi/6' @13 '7 and 5'  @24 ###.## tpi(9)  @33 ###.## tpi(10) @43 ###.## F0910
display @2 ' pi/6' @13 '1 and 11' @24 ###.## tpi(11) @33 ###.## tpi(12) @43 ###.## F1112
display ' '
*
end mhegy










