% Program to solve Hands-on 1.2
% written by Marcel Bluhm in January 2012
% prepare machine
clc;
clear all;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Load data %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
[num, txt, raw]= xlsread('DataPWT',3);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Estimate Equation (1.10) using POLS %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% - form dependent variable
Y=num(:,[1;2;4]);
Y(Y(:,1)==1970,:)=[]; % eliminate first period grgdp
% - form independent variable
X=[num(:,[1;2]),ones(size(num,1),1),num(:,3)];% constant and independent variable
X(X(:,1)==2005,:)=[]; % eliminate last period rgdpl
% - estimate Equation (1.2)
theta=(X(:,3:4)'*X(:,3:4))^(-1)*X(:,3:4)'*Y(:,3);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Obtain the pooled OLS residuals %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
v=[Y(:,[1;2]),Y(:,3)-X(:,3:4)*theta]; % error terms
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Compute consistent estimator of sig_v %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
sig_v=1/(size(v,1)-size(X,2))*(v(:,3)'*v(:,3));
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Compute consistent estimator of sig_c %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
vv=0; %initialize variable
for i=1:num(end,2)
% obtain current cross-section
vi=v(v(:,2)==i,3);
for t=1:size(vi,1)-1
for s=t+1:size(vi,1)
vv=vv+vi(t)*vi(s);
end
end
end
N=num(end,2); % number of cross-sections
T=sum(numel(v(v(:,2)==i)),1); % number of observations per cross-section
sig_c=1/(N*T*(T-1)/2-size(X,2))*vv;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Form Omega and obtain parameter vector of interest %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
sig_u=sig_v-sig_c;
omega=sig_u*eye(T)+sig_c*ones(T,1)*ones(T,1)';
p1=zeros(2,2); % initialize first part of Equation (1.11)
p2=zeros(2,1); % initialize second part of Equation (1.11)
for i=1:N
Xi=X(X(:,2)==i,3:4);
Yi=Y(Y(:,2)==i,3);
p1=p1+Xi'*omega*Xi;
p2=p2+Xi'*omega*Yi;
end
beta_re=p1^(-1)*p2;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Compute error terms and avar(beta_re) %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
v_=[Y(:,1:2),Y(:,3)-X(:,3:4)*beta_re]; % year id v_
p1=zeros(2,2); % initialize first part of Equation (1.12)
p2=zeros(2,2); % initialize second part of Equation (1.12)
for i=1:N
Xi=X(X(:,2)==i,3:4);
vi_=v_(v_(:,2)==i,3);
p1=p1+Xi'*omega^(-1)*Xi;
p2=p2+Xi'*omega^(-1)*(vi_*vi_')*omega^(-1)*Xi;
end
avar_beta_re=p1^(-1)*p2*p1^(-1);
t_stat=beta_re./(diag(avar_beta_re).^(.5));
result=[beta_re,t_stat];
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Test whether the model actually contains an unobserved effect %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
nm=0; % initialize numerator of Equation (13)
den=0; % initialize denominator of Equation (14)
for i=1:N
% obtain current cross-section
vi=v(v(:,2)==i,3);
temp=0;
for t=1:T-1
for s=t+1:T
nm=nm+vi(t)*vi(s);
temp=temp+vi(t)*vi(s);
%den=den+(vi(t)*vi(s))^2;
end
end
den=den+temp^2;
end
TS=nm/den^(0.5);