*
* @GregoryHansen(options) depvar start end
* # other "x" variables
*
* Implements the Gregory-Hansen cointegration test ("Residual-based
* tests for cointegration in models with regime shifts" Journal of
* Econometrics (1996)). The cointegrating regression is allowed to have
* a level break, a full structural break, or a trend with a level break.
* The break point t0 is unknown and determined by finding the minimum
* value for the ADF statistic. The number of lags of the change in the
* residual used in computing the ADF statistic can either be input, or
* can be selected automatically using AIC, BIC or general-to-specific
* pruning by t-tests. All of the automatic lag selections are applied to
* the base model, not the individual break points.
*
* Options:
*   MODEL=[CONSTANT]/TREND/FULLBREAK
*     CONSTANT includes an intercept and a level shift dummy.
*     TREND includes an intercept, trend and a level shift dummy.
*     FULLBREAK includes no trend, but dummies for the intercept and the slopes
*
*   LAGMETHOD=[INPUT]/AIC/BIC/TTEST
*     Selects the method for deciding the number of additional lags. If INPUT, the
*     number of lags given by the MAXLAGS option is used. If AIC, the AIC-minimizing
*     value between 0 and MAXLAGS is used; if BIC, it's the BIC-minimizing value, and
*     if t-test, the number of lags for which the last included lag has a marginal
*     significance level less than the cutoff given by the SIGNIF option.
*
*   MAXLAGS=number of additional lags (LAGMETHOD=INPUT) or the maximum number of lags
*     to consider (other LAGMETHOD's) [number of observations**.25]
*
*   SIGNIF=cutoff significance level for LAGMETHOD=TTEST [.10]
*
*   SMPL=SMPL series with entries to skip
*
*   GRAPH/[NOGRAPH]
*     GRAPH requests a graph of the unit-root statistics with the different break points
*
*   PI=fraction of data range to skip at either end when examining possible break points [.15]
*
* Variables Defined:
*  %cdstat is the unit root test statistic
*  %%breakpoint is the entry at which %cdstat is achieved.
*
* Revision Schedule:
*  09/2003 Written by Tom Doan, Estima
*  10/2004 Correction made to independent variables in regression
*
procedure GregoryHansen depvar start end
type series depvar
type integer start end
*
option choice  model     1 constant trend fullbreak
option integer maxlags
option choice  lagmethod 1 input aic bic ttest
option real    signif    .10
option series  smpl
option switch  graph
option real    pi     .15
*
local integer maxlag bestlag startl endl lag znobs tstart tend m i time
local series trend breakint breakstats e de
local real icmult ic icmin crit01 crit05
local vect[integer] others xlist
local vect[series] xdummies
local rect c01 c05
*
* These are the .01 and .05 critical values. The rows are the number of
* x variables, and the columns are the model
*
compute c01=||-5.13,-5.45,-5.47|-5.44,-5.80,-5.97|-5.77,-6.05,-6.51|-6.05,-6.36,-6.92||
compute c05=||-4.61,-4.99,-4.95|-4.92,-5.29,-5.50|-5.28,-5.57,-6.00|-5.56,-5.83,-6.41||
*
* Figure out the range available for the series and one lag
*
enter(varying) xlist
compute m=%rows(xlist)
if model==3
   dim xdummies(m)

inquire(reglist) startl<<start endl<<end
# depvar{0 1} xlist
set trend startl endl = t
if model==2 {
   linreg(noprint) depvar / e
   # constant trend xlist
}
else {
   linreg(noprint) depvar / e
   # constant xlist
}
set de = e-e{1}
*
* Determine default value of maxlag
*
compute znobs=endl-startl+1
if .not.%defined(maxlags)
   compute maxlag=fix(znobs**.25)
else
   compute maxlag=maxlags

if lagmethod==1
   compute bestlag=maxlag
else {
   cmom(smpl=smpl) startl+maxlag endl
   # constant e{1} de{0 to maxlag}
   compute icmult=%if(lagmethod==2,2.0/%nobs,log(%nobs)/%nobs)
   linreg(cmom,noprint) de
   # constant e{1}
   compute ic=log(%rss/%nobs)+icmult*%nreg
   compute icmin=ic,bestlag=0
   do lag=maxlag,1,-1
      linreg(cmom) de
      # constant e{1} de{1 to lag}
      if lagmethod==4 {
         if %ttest(%tstats(2+lag),%ndf)<=signif {
            compute bestlag=lag
            break
         }
      }
      else {
         compute ic=log(%rss/%nobs)+icmult*%nreg
         if ic<icmin
            compute icmin=ic,bestlag=lag
      }
   end do lag
}
compute znobs=endl-startl+1-bestlag
compute tstart=startl+bestlag+fix(pi*znobs),tend=endl-fix(pi*znobs)
set breakstats tstart tend = 0.0
compute others=xlist
if model==2 {
   enter(varying) others
   # others trend
}
else
if model==3 {
   enter(varying) others
   # others xdummies
}
*
do time=tstart,tend
   set breakint   startl endl = t>=time
   if model==3 {
      do i=1,m
         set xdummies(i) startl endl = xlist(i){0}*(t>=time)
      end do i
   }
   linreg(noprint,smpl=smpl) depvar startl endl e
   # constant breakint others
   set de = e-e{1}
   if bestlag>0 {
      linreg(noprint,smpl=smpl) de startl+bestlag endl
      # e{1} de{1 to bestlag}
   }
   else {
      linreg(noprint,smpl=smpl) de startl+bestlag endl
      # e{1}
   }
   compute breakstats(time) = %tstats(1)
end do time
*
if graph {
   graph(header="Gregory-Hansen Cointegration Tests")
   # breakstats tstart tend
}
ext(noprint) breakstats startl+bestlag endl
compute %cdstat=%minimum,%%breakpoint=%minent
disp
disp "Gregory-Hansen Cointegration Tests"
disp "Variables"
disp %l(depvar)
do i=1,m
   disp %l(xlist(i))
end do i
compute crit01=c01(m,model),crit05=c05(m,model)

if model==1
   disp "Break in Intercept. No Trend"
else
if model==2
   disp "Break in Intercept. Trend Included"
else
   disp "Full Structural Break"

disp(hold) "Including" bestlag "Lags of Difference Selected by "
if lagmethod==1
   disp "User"
else
if lagmethod==2
   disp "AIC"
else
if lagmethod==3
   disp "BIC"
else
if lagmethod==4
   disp "Truncation at Longest Lag with P-value<" signif

disp "Minimum T-Statistic" %cdstat "at" %datelabel(%%breakpoint)
disp "Critical Values are 1%" *.## crit01 "and 5%" crit05
end

