*
* @HTUnit( options ) series start end
*
* Implements tests from Harris and Tzavalis(1999), "Inference for unit
* roots in dynamic panels where the time dimension is fixed". Journal
* of Econometrics, vol 91, pp 201–226.
*
* This has a null of unit root, with asymptotics assuming large N, fixed
* T.
*
* Requires version 7.30 or later of RATS due to the TITLE option on
* REPORT(ACTION=DEFINE). That just affects the title of the window if
* you reload it with the Window--Restore Report operation. If you remove
* that, it will work with 7.1 or later.
*
* Parameters:
*
*   series     = series to analyze
*   start end  = range to use in performing tests [by default, maximum range of series]
*
* Options:
*   DET=NONE/[CONSTANT]/TREND
*      Individual-specific deterministic components to remove
*
*   SMPL=series with zeros in entries to be skipped [not used]
*
*   TITLE=title for output ["Harris-Tzavalis Test: Series ..."]
*   [PRINT]/NOPRINT show test output
*
* Variables Defined:
*   %CDSTAT - the Z statistic
*   %SIGNIF - the average DF statistic
*   %NGROUP - number of individuals
*
* Revision Schedule:
*  07/2010 Written by Tom Doan, Estima
*
procedure HTUnit series start end
type series 	series
type integer 	start end
*
option choice  det      2 none constant trend
option switch  print  	1
option string  title
option series  smpl
*
local vect[int] 	detvars
local integer		startl endl starti endi
local series      itrend
local vect[ser]	resids
local real        mustar sigstar
local report      htreport
local series      icount
local real        pt nt wrongsize
local string      ltitle
local string      descript
*
if .not.%defined(series) {
   disp "###Syntax: @HTUNIT( options ) >>>>>SERIES<<<<<"
   return
}
inquire(series=series) startl<<start endl<<end
compute starti=%indiv(startl),endi=%indiv(endl)
panel(icount=1,smpl=smpl,compress) series startl endl icount
sstats(max) starti endi icount>>pt
sstats      starti endi icount==pt>>nt
sstats(max) starti endi %if(icount<>pt.and.icount<>0,t,0)>>wrongsize
if wrongsize>0 {
   disp "@HTUNIT requires balanced data set."
   disp "Individual "+fix(wrongsize)+" has "+fix(icount(fix(wrongsize)))+" data points"
   return
}
*
* Create the individual trend series in case we need it.
*
set itrend = %period(t)
*
* Create the set of regressors and calculate the variance and
* bias adjustment.
*
if det==1 {
   dim detvars(0)
   compute mustar=0.0
   compute sigstar=2.0/(pt*(pt-1))
   compute descript="None"
}
else
if det==2 {
   compute detvars=||constant||
   compute mustar=3.0/(pt+1)
   compute sigstar=3*(17.0*pt^2-20*pt+17)/(5.0*(pt-1)*(pt+1)^3)
   compute descript="Constant"
}
else {
   compute detvars=||constant,itrend||
   compute mustar=7.5/(pt+2)
   compute sigstar=15.0*(193*pt^2-728*pt+1147)/(112*(pt+2)^3*(pt-2))
   compute descript="Constant + Trend"
}
*
* Sweep out the individual specific components
*
sweep(group=%indiv(t),smpl=smpl,series=resids) startl endl
# series
# detvars
*
* Copy the residuals out of the VECT[SERIES] into %resids.
*
set %resids = resids(1)
*
* Do the regression of %resids on its lag.
*
linreg(noprint) %resids
# %resids{1}
*
* Compute the HT statistic and its significance level (one tailed).
*
compute %cdstat=((%beta(1)-1.0)+mustar)*sqrt(nt/sigstar)
compute %signif=%cdf(%cdstat)
compute %ngroup=fix(nt)
*
if %defined(title)
   compute ltitle=title
else
   compute ltitle="Harris-Tzavalis Test: Series "+%l(series)

report(use=htreport,action=define,title=ltitle)
report(use=htreport,atrow=1,atcol=1,span) ltitle
report(use=htreport,atrow=2,atcol=1,span) "Test has fixed T, large N asymptotics"
report(use=htreport,atrow=3,atcol=1,span) "Null is rho(i)=1. Alternative is rho(i)==rho<>1"
report(use=htreport,atrow=4,atcol=1,span) "Individual Specific Components: "+descript
report(use=htreport,atrow=6,atcol=1) "N"      fix(nt)
report(use=htreport,atrow=7,atcol=1) "T"      fix(pt)
report(use=htreport,atrow=8,atcol=1) "Rho" 	 %beta(1)
report(use=htreport,atrow=9,atcol=1) "Z" 	    %cdstat
report(use=htreport,atrow=10,atcol=1) "Signif" %signif
report(use=htreport,action=format,width=8,atcol=2)
if print
   report(use=htreport,action=show)
end

