proc (3)=est_gls(x,y,theta,lambda,t,n,names);
/*
est_gls.src is written by Kristian Jönsson (November 12, 2002), Department of Economics, Lund University.
Contact info:  kristian.jonsson@nek.lu.se

See separate documentation and 'Econometric Analysis of Panel Data' 2nd edition, Baltagi (2001).

Note:
This code can be used freely as long as proper reference is given.
No performance guarantee is made.
Bugreports are welcome.
*/
local P_mu,P_lambda;
local y_re,z_re,z,beta,se,t_val,pv,res,df;
local k;
local i;
local mask,fmat,th_method,z_out;
local Z_mu,P,Q;
local betam;
local res_back,y_back;

th_method= "Wallace & Hussain (1969)" | "Amemiya (1971)" | "Swamy    and Arora (1972)" | "Nerlove (1971)";

P_mu=eye(n).*.ones(t,t)/t;     /* Time average for different individuals */
P_lambda=ones(n,n)/n.*.eye(t);  /* Individual average for different time-periods */ 

Z_mu=eye(n).*.ones(T,1);
/* Define P and Q matrices */
/* Choose one way or another to demean the data */
P=Z_mu*inv(Z_mu'*Z_mu)*Z_mu';
Q=eye(n*t)-P;

z=ones(rows(x),1)~x;

betam=zeros(4,cols(z));
res_back=zeros(rows(x),4);
y_back=zeros(rows(x),4);

mask = ones(1,5); mask[1]=0; /*zero for character output and 1 for numeric output */
fmat=zeros(5,3); fmat[1,.]="*.*s"~8~8;
for i (2,5,1);
fmat[i,.]="*.*lf"~12~3; 
endfor;


for i (1,4,1);

y_re=y-theta[i,1]*P_mu*y-theta[i,2]*P_lambda*y+theta[i,3]*P_mu*P_lambda*y;
z_re=z-theta[i,1]*P_mu*z-theta[i,2]*P_lambda*z+theta[i,3]*P_mu*P_lambda*z;

y_back[.,i]=y_re;

beta=y_re/z_re;
betam[i,.]=beta';
k=cols(z_re);
res=y_re-z_re*beta;
res_back[.,i]=res;
se=sqrt(diag(res'*res/(n*t-k)*inv(z_re'*z_re)));
t_val=beta./se;
df=rows(y)-rows(beta);
pv=2*cdftc(abs(t_val),df);
print "**************************************************************************";
format/sa/rds 30,-1;
print $ "Weight method: " th_method[i];
format /m1 /rd 8,4; print "Weight 1 (theta 1): " theta[i,1];
format /m1 /rd 8,4; print "Weight 2 (theta 2): " theta[i,2];
format /m1 /rd 8,4; print "Weight 3 (theta 3): " theta[i,3];
print "---Variance Components---";
format /m1 /rd 8,4; print "Lambda_1 : " lambda[i,1];
format /m1 /rd 8,4; print "Lambda_2 : " lambda[i,2];
format /m1 /rd 8,4; print "Lambda_3 : " lambda[i,3];
format /m1 /rd 8,4; print "Lambda_4 : " lambda[i,4];
print "**************************************************************************";
print "Variable" "        Coef  " "   Std err " "    t-value   " "  p-value     ";
print "********" "        ****  " "   ******* " "    *******   " "  *******     ";
z_out=printfm(names~beta~se~t_val~pv,mask,fmat);
print "**************************************************************************";
print;print;

endfor;
retp(y_back,betam',res_back);
endp;
