proc (1)=HAC(x,y,n,t,names);
/*
hac.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.
*/

/* Estimates Arellano's HAC VC-matrix */
/* See notes p 51 */
/* Est.Asy.Var(b)=inv(x'x)*sum_i(x_i'e_i*e_i'*x_i)*inv(x'x) */
/* x_i is a kx1 vector containing the i:th observation of the k variables */
/* A degrees of freedom correction is performed to obtain results identical to E-views (Note however vc is derived asy) */

local res;
local S0,temp_sum,var_w;
local nobs,k;
local i;
local z,fmat,mask;
local name;
local Q, x_tr;


Q=eye(n).*.eye(t)-eye(n).*.(ones(t,t)/t)-(ones(n,n)/n).*.eye(t)+(ones(n,n)/n).*.(ones(t,t)/t);
x_tr=Q*x;

res=Q*y-Q*x*((Q*y)/(Q*x));


nobs=rows(x);
k=cols(x);


temp_sum=zeros(k,k);
for i (1,n,1);
temp_sum=temp_sum+x_tr[(i-1)*t+1:(i-1)*t+t,.]'*(res[(i-1)*t+1:(i-1)*t+t,.]*res[(i-1)*t+1:(i-1)*t+t,.]')*x_tr[(i-1)*t+1:(i-1)*t+t,.];
endfor;
S0=temp_sum;

var_w=inv(x_tr'*x_tr)*S0*inv(x_tr'*x_tr);

name=names[2:k+1];
mask = ones(1,2); mask[1]=0; /*zero for character output and 1 for numeric output */
fmat=zeros(2,3); fmat[1,.]="*.*s"~8~8; fmat[2,.]="*.*lf"~12~3;

print;print;print;


print "************ Arellano's HAC SEs in 2-way FE model ************";
z=printfm(name~sqrt(diag(var_w)),mask,fmat);
print "*******************************************************************";

retp(sqrt(diag(var_w)));
endp;