set mem 400
log using c:\temp\PROGRAM3.log, replace
quietly{
**************************************************************
*                                                            *
*      THIS COMPUTER PROGRAM IS FOR THE PANEL DATA           *
*      SPECIFICATION WHERE THE DISTURBANCES ARE BOTH         *
*      SPATIALLY AND TIME-WISE AUTOCORRELATED                *
*                                                            *
*        Y=ALPHA+X*BETA+U, U=(IT#(RHO*WEIGHTINGMATRIX))*U+E, *
*        E=(eT#MU)+V                                         *
*                                                            *
**************************************************************
**************************************************************
*EXPLANATION OF THE VARIABLES USED IN THE PROGRAM            *
*Y.......VECTOR OF DEPENDENT VARIABLE                        *
*X.......VECTOR OF EXOGENOUS EXPLANATORY VARIABLE            *
*U.......VECTOR OF DISTURBANCES                              *
*MMAT....SPATIAL WEIGHTING MATRIX(REMAINS SAME FOR ALL TIME  *
*E.......VECTOR OF INNOVATIONS                               *
*RHO.....SPATIAL AUTOREGRESSIVE PARAMETER                    *
*AHOLS...OLS ESTIMATE OF ALPHA                               *
*AHGLS...GLS ESTIMATE OF ALPHA                               *
*BHOLS...OLS ESTIMATE OF BETA                                *
*BHGLS...GLS ESTIMATE OF BETA                                *
*RHOHGM..GENERAL MOMENTS ESTIMATOR OF RHOH                   *
*NS......NUMBER OF CROSS SECTIONAL UNITS                     *
*T.......NUMBER OF TIME PERIODS                              *
*NT......TOTAL SAMPLE SIZE                                   *
**************************************************************

**************************************************************
*                      READ THE DATA                         *
**************************************************************

infile m1-m100 using c:\temp\mmat.dat
save c:\temp\mmat.dta, replace

*use c:\temp\mmat
scalar NS=100
scalar T=4
scalar NT=NS*T
set matsize 100     

mkmat m1-m100, mat(MMAT)
mat MMM=trace(MMAT'*MMAT)/_N
mat IMAT=I(NS)
mat TMAT=I(T)
mat JMAT=J(T,T,1/T)

set matsize 400
mat MMAT1=TMAT#MMAT
mat Q0=(TMAT-JMAT)#IMAT
mat Q1=JMAT#IMAT

mat MW=MMAT'*MMAT
mat MWT=trace(MW)/NS
mat MWTT1=MW*(MMAT+MMAT')
mat MWTT=trace(MWTT1)/NS
mat MWW=MMAT*MMAT
mat MWWW=MW*MW
mat MWWWW=MW+MWW
mat MQ0=trace(MWWW)/NS
mat MQ1=trace(MWWWW)/NS
drop m1-m100
**************************************************************
*         IN THE FIRST STEP ESTIMATE THE MODEL BY OLS        *
*              AND GET THE ESTIMATED DISTURBANCES            *
**************************************************************
infile y x1 x2 using c:\temp\var3.dat
merge using c:\temp\mmat
reg y x2 x1, noconstant 
predict r, resid
gen r2=r^2
egen ss=sum(r2)
scalar SHOLS=ss/(NS-1)
scalar AHOLS=_b[x1]
scalar BHOLS=_b[x2]
}

**************************************************************
*       STEP 1:  OLS ESTIMATORS OF REGRESSION PARAMETERS     *              *
**************************************************************
scalar list AHOLS BHOLS SHOLS

quietly{
**************************************************************
*   IN THE SECOND STEP USE THE ESTIMATES OF DISTURBANCES     * 
*          FOR THE GENERAL MOMENTS(GM) ESTIMATOR             *
**************************************************************

****************************************************************
* SET FLAG L FOR CHOICE OF PROCEDURE TO ESTIMATE RHO AND SIGMA *
*                                                              *
* L = 1..  SELECT INITIAL GM ESTIMATORS                        *
*   = 2..  SELECT PARTIALLY WEIGHTED GM ESTIMATORS             *
*   = 3..  SELECT FULLY WEIGHTED GM ESTIMATORS                       *
****************************************************************
scalar L=1
}

gmproc3 r


**************************************************************
* STEP 2, 1st. Part:                                         *
* INITIAL GM ESTIMATES FOR RHO AND THE VARIANCE COMPONENTS   *
**************************************************************
scalar list RHOH SIGV SIG1

**************************************************************
* STEP 2, 2nd. Part:                                         *
* WEIGHTED NLS ESTIMATORS FOR RHO AND THE VARIANCE COMPONENTS*
*   (FOR L=1 this part is ommited and the output should be   *      
*    equal to the one obtained in the 1st. part of step 2)   *
**************************************************************
scalar list RHOGM SIGVV SIG11


quietly{ 
**************************************************************
*    IN THE THIRD STEP USE THE GM ESTIMATES FOR SPATIAL      *
*  AUTOREGRESSIVE PARAMETER, CORRECT FOR SPATIAL CORRELATION *
*  IN DISTURBANCES AND  ESTIMATE THE CORRECTED MODEL BY OLS  *
************************************************************** 

merge using c:\temp\var3
} 
glsproc3 y x1 x2 RhoGM

**************************************************************
*              STEP 3:  FEASIBLE GLS ESTIMATORS              *
**************************************************************
scalar list AHFGLS BHFGLS SHFGLS

log close
clear