************************************************************************
*
*       RATS PROCEDURE FOR MULTIVARIATE GROUP MEAN PANEL FMOLS TESTS
*                     IN HETEROGENEOUS PANELS
*
*  By Peter Pedroni, Williams College: ppedroni@williams.edu
*
*  Includes revisions made by Prof. Pedroni to his original program
*  that allow for common time effects in unbalanced panels, automatic
*  selection of bandwidth for kernels, and automatic adjustment of display
*  for number of regressors
*
*  Revisions:
*   06/04 Converted to procedure by Tom Doan, Estima
*
*   09/05 Changed: compute term1 = tr(G21+V21)
*              To: compute term1 = tr(G12+V21)
*          to match Pedroni's original program
*
* Syntax:
*  @panelfm( options ) start end
*  # list of variables (list dependent variable first)
*
*
* Parameters:
*  start,end  = Range of entries to use. Defaults to all available data.
*               If you set these, do so as follows:
*                start = first observation of the first individual in the desired range
*                end   = last observation of the last individual in the desired range
*
* Options:
*  block = time dimension [by default is taken from CALENDAR(PANELOBS)]
*          (see note below for unbalanced panels)
*  tdum/[notdum]   Use TDUM to subtract out common time effects
*  bvec  = VECTOR  Allows you to supply a vector of hypothesized slope
*                  coefficients [default is zeros]
*
* For unbalanced panels, you should enter into your data set either blanks or the value
* NA for missing values, so that the end result is that each member has uniform length T,
* and the data is stacked into a column of length N*T.
*
*
*                            REFERENCES:                               *
*                                                                      *
*   For Theory and Monte Carlos:                                       *
*                                                                      *
*   Pedroni, Peter (2000) "Fully Modified OLS for Heterogeneous        *
*   Cointegrated Panels," Advances in Econometrics, Vol. 15, 93-130,   *
*   NONSTATIONARY PANELS, PANEL COINTEGRATION AND DYNAMIC PANELS,      *
*   JAI Press.                                                         *
*                                                                      *
*   Working paper version: Pedroni, Peter (1996), Indiana University   *
*   Working Papers in Economics, No. 96-020.                           *
*                                                                      *
*   For general discussion and application:                            *
*                                                                      *
*   Pedroni, Peter (2001) "Purchasing Power Parity Tests in            *
*   Cointegrated Panels," Review of Economics and Statistics, 83,      *
*   727-731.                                                           *
*                                                                      *
*                                                                      *
************************************************************************
*
procedure panelfm  startp endp
type integer startp endp
*
option integer block
option switch  tdum   0
option vector  bveco
*
local integer tperiods start end
*
local vect W
local vec[real] temp tden iFMOLS iFMtstat itemp1 pfg pfgt O21 G21 G12 V21 tvecreal
local rect[real] gammastar gXX ytemp O22 G22 V22
local vec[series] xvec dxvec indfm indfmt dvec rvec g21ser g12ser
local rect[rect] G
local rect[series] g22ser
local vec[series] ytilde
local series u ehat g11ser icount
local vect[int] datavec rhs
local integer nsecs
local integer mlag
local integer i j k m
local integer maxlag q r s
local symm v
local rect g11vec g21vec g12vec g22vec
local rect term1 term2 itemp2
local series ystar tvec
local real o11
local vec bvec
local vec[label] labelvec
local integer unbal
*
* Get list of variables (LHS variable first)
*

enter(varying) datavec
compute m=%rows(datavec)-1
*
* Create list to be used for RHS (that is, datavec w/o the
* first entry, which will be used for the dependent variable)
*
dim rhs(m)
ewise rhs(i)=datavec(i+1)
dim labelvec(m+1)
ewise labelvec(i)=%l(datavec(i))

*
*  Get the block size
*

if .not.%defined(block)
   compute tperiods=2//1-1//1
else
   compute tperiods=block

*
* Get the data range
*

inquire(reglist) start<<startp end<<endp
# datavec

*
* Normalize start to begin at the first entry in an individual's
* record, and end to end at the last entry.
*

compute start=(start-1)/tperiods*tperiods+1
compute end  =((end-1)/tperiods+1)*tperiods
compute nsecs=(end-start+1)/tperiods

*
* Run an OLS mainly to figure out what data points are available.
*

linreg(noprint) datavec(1) start end u
# rhs

* NOTE: If using Version 6, you can uncomment this line
* and the CALENDAR statement below if you want the procedure
* to reset back to your original CALENDAR setting:
*
* compute savecal = %calendar()
*
calendar(panelobs=tperiods)
panel(icount=1,compress) u / icount
dim tvecreal(nsecs)
ewise tvecreal(i) = icount(i)
compute unbal = fix(%minvalue(tvecreal)<%maxvalue(tvecreal))
*
* calendar(recall=savecal)
*


*
* Dimension everything
*

dim temp(m) tden(m) iFMOLS(m) iFMtstat(m) itemp1(m) $
              pfg(m) pfgt(m) O21(m) G21(m) G12(m) V21(m)
dim gammastar(1,m) gXX(m,m) ytemp(1,1) $
              O22(m,m) G22(m,m) V22(m,m)
dim xvec(m) dxvec(m) indfm(m) indfmt(m) dvec(m+1) rvec(m+1) $
               g21ser(m) g12ser(m)
dim G(2,2)
dim g22ser(m,m)

dim ytilde(m+1)

dim bvec(m)

if .not.%defined(bveco)
   compute bvec=%const(0.)
else
   ewise bvec(i)=bveco(i)


                        *** MAIN SOURCE CODE ***

display ' '
display '      Currently computing panel statistics. Please wait. '
display ' '

*
*  Copy either the raw data, or (if tdum) the deviations from time period means
*  means into ytilde.
*

do K=1,m+1
   if tdum==1
      panel(block=tperiods,entry=1.0,time=-1.0) datavec(K) start end ytilde(k)
   else
      set ytilde(k) start end = datavec(k){0}
end do k

ewise pfg(i) = 0.0
ewise pfgt(i) = 0.0

panel(compress,icount=1.0) u start end tvec

do J=1,Nsecs
   compute mlag = fix(%round((4.0*(Tvec(J)/100.0)**(2.0/9.0)),0))
   do K=1,M+1
      move ytilde(k) (j-1)*tperiods+1 j*tperiods dvec(k) 1
      diff dvec(k) 2 tperiods rvec(k)
      diff(center) rvec(k) 2 tperiods rvec(k)
   end do k
   do K=1,M
      set xvec(K) 1 Tperiods = dvec(K+1)
      set dxvec(K) 2 Tperiods = rvec(K+1)
   end do K
   linreg(noprint) dvec(1) 2 Tperiods ehat
   # dvec(2) to dvec(m+1) constant

   set rvec(1) =  ehat

   compute maxlag = mlag
   VCV(matrix=V,noprint) 2 Tperiods
   # rvec(1) to rvec(m+1)
   mcov(lwindow=newey,lags=maxlag,noprint) 2 Tperiods
   # rvec(1) to rvec(m+1)

   compute O11 = (1.0/(Tperiods-1.0))*%CMOM(1,1)
   do q=1,m
      compute O21(q)= (1.0/(Tperiods-1.0))*%CMOM(q+1,1)
      compute V21(q)= V(q+1,1)
      do r=1,m
         compute O22(q,r) = (1.0/(Tperiods-1.0))*%CMOM(q+1,r+1)
         compute V22(q,r) = V(q+1,r+1)
      end do r
   end do q

   compute maxlag = mlag
   dimension W(maxlag)
   do I=1,maxlag
      compute W(I) = 1.0 - I/(maxlag+1.0)
   end do I

   cross(covariances,noprint) rvec(1) rvec(1) 2 Tperiods 1 maxlag g11ser
   do q=1,m
      cross(covariances,noprint) rvec(q+1) rvec(1) 2 Tperiods 1 maxlag g21ser(q)
      cross(covariances,noprint) rvec(1) rvec(q+1) 2 Tperiods 1 maxlag g12ser(q)
      do r=1,m
         cross(covariances,noprint) rvec(q+1) rvec(r+1) 2 Tperiods 1 maxlag g22ser(q,r)
      end do r
   end do q

   make g11vec 1 maxlag maxlag
   # g11ser
   compute G(1,1) = %scalar(tr(W)*g11vec)
   do q=1,m
      make g21vec 1 maxlag maxlag
      # g21ser(q)
      compute G21(q) = %scalar(tr(W)*g21vec)
      make g12vec 1 maxlag maxlag
      # g12ser(q)
      compute G12(q) = %scalar(tr(W)*g12vec)
      do r=1,m
         make g22vec 1 maxlag maxlag
         # g22ser(q,r)
         compute G22(q,r) = %scalar(tr(W)*g22vec)
      end do r
   end do q
   do s=2,Tperiods
      ewise temp(I) = dxvec(I)(s)
      compute ytemp = tr(o21)*inv(o22)*temp
      set ystar s s = dvec(1) - ytemp(1,1)
   end do s

   compute term1 = tr(G12+V21)
   compute term2 = tr(o21)*inv(o22)*(G22+V22)
   compute gammastar = term1 - term2

   linreg(noprint) ystar 2 Tperiods
   # xvec(1) to xvec(m) constant

   ewise gXX(I,J) = %XX(I,J)
   ewise itemp1(I) = %beta(I)
   compute itemp2 = (Tperiods-1.0)*gammastar*gXX
   compute iFMOLS = itemp1 - tr(itemp2)
   ewise tden(I) = sqrt(%XX(I,I)*o11)
   ewise iFMtstat(I) = (iFMOLS(I)-bvec(I))/tden(I)

   do k=1,m
     set indfm(k) J J = iFMOLS(k)
     set indfmt(k) J J = iFMtstat(k)
   end do k

   compute pfg = pfg + iFMOLS
   compute pfgt = pfgt + iFMtstat

   dofor k=ehat rvec(1) ystar
      set k 1 Tperiods = %NA
   end dofor k

end do J

compute pfg = (1.0/Nsecs)*pfg
compute pfgt = (1.0/sqrt(Nsecs))*pfgt


***********************************************************************

    *** DISPLAY INDIVIDUAL AND PANEL RESULTS ***

display ' '
display 'Dependent Variable :' labelvec(1)

display ' '
if tdum == 0 ; { ; display 'common time dummies NOT included' ; }
if tdum == 1 ; { ; display 'common time dummies included' ; }
display ' '
display 'Null vector for t-stats is: ' ; write bvec
display ' '
display '                        INDIVIDUAL FMOLS RESULTS'
display '                        (t-stats in parentheses)'
display ' '

display ' '
display @4 'Member ' @+8 'Variable' @+17 'Coefficient' @+13 't-statistic'
display @2 '***********' @+5 '*************' @+12 '**************' @+10 '**************'
do i=1,Nsecs
   do k=1,m
      display @5 'No.'+i @22 labelvec(k+1) @47 ##.## indfm(k)(i) @71 '(' ##.## indfmt(k)(i) ')'
   end do k
   display @2 '-----------' @+5 '-------------' @+12 '--------------' @+10 '--------------'
end do i
display @2 '***********' @+5 '*************' @+12 '**************' @+10 '**************'
display ' '
display ' '
display ' '

display '                    PANEL GROUP FMOLS RESULTS    '
display ' '
display @20 'Coefficient' @42 't-statistic'
do k=1,m
 display @10 labelvec(k+1) @21 ##.## pfg(k) @44 '(' ##.## pfgt(k) ')'
end do k
display ' '
display ' '
if unbal == 0
{
   display '        Nsecs =' Nsecs ', Tperiods =' Tperiods $
                         ', no. regressors =' M
}
if unbal == 1
{
   display '        Nsecs =' Nsecs ', Tperiods = (unbalanced)' $
                         ', no. regressors =' M
}
display '         ******************************************** '
end












