%---------------------------------------------;
% ECG 790: Fall 2006, program for Assessed    ;
% Problem set #1 Question 2                   ;
%---------------------------------------------;


%-------------------------------------------;
% Load data (into matrix subset_JIVE)       ;
%-------------------------------------------;

load 'c:\arh\courses\ecg790\2006\computer\Angrist-krueger\assprb106\subset_JIVE.txt';

%-------------------------------------------;
% Select data used in estimation            ;
%-------------------------------------------;
year=subset_JIVE(:,1)-29;
state=subset_JIVE(:,2);
quar=subset_JIVE(:,3);
ed=subset_JIVE(:,4);
logw=subset_JIVE(:,5);

%------------------------------------------;
% Sample size is n                         ;
%------------------------------------------;

[n,ncyear]=size(year);

%------------------------------------------;
% Construct necessary dummy variables      ;
%------------------------------------------;

%------------------------------------------;
% first specify dimensions of matrices to  ;
% contain the dummies                      ;
%------------------------------------------;

% year of birth: 10 years (1930-1939) - drop one to avoid linear dependence;
y=zeros(n,9);

% state of birth: 4 (AK, KS, KY, TN) - drop one to avoid linear dependence;
s=zeros(n,3);

% interaction term quarter of birth x year of birth dummies

% number of quarter is four but drop one to avoid linear dependence;

iqy=3*9;
qy=zeros(n,iqy);

% interaction term quarter of birth x state of birth dummies

iqs=3*3;
qs=zeros(n,iqs);

%----------------------------;
% create dummies             ;
%----------------------------;

for i=1:n
   if year(i)~=10 
   j=year(i);
   y(i,j)=1;
   end
   if (quar(i)~=1) & (year(i)~=10)
   k=quar(i)+(j-1)*3-1;
   qy(i,k)=1;
   end;
   if state(i)~=1
   j=state(i)-1;
   s(i,j)=1;
   end
   if (quar(i)~=1) & (state(i)~=1) 
   k=quar(i)+(j-1)*3-1;
   qs(i,k)=1;
   end;
 end;

%------------------------------------------;
% Collect variables into data matrices     ;
%                                          ;
% X= explanatory variables                 ;
% Z= instrumental varaibles                ;
%------------------------------------------;

X=[ones(n,1),ed, y, s];
Z=[ones(n,1), y, s, qy, qs];

PZ  = Z*inv(Z'*Z)*Z';

%----------------------------;
% (i) 2SLS estimate          ;
%----------------------------;

thetaIV  = (X'*PZ*X)\X'*PZ*logw;

fid = fopen('AKoutput.dat', 'w');
fprintf(fid,'(i) 2SLS estimator of parameter on education: %6.5f\n\n', thetaIV(2,1));
fclose(fid);


%---------------------------;
% (ii) 95% CI               ;
%---------------------------;

%----------------------------------------------;
% Calculate the covariance matrix of           ;
% 2SLS using (2.29) in Hall (2005)             ;
%----------------------------------------------;

e = logw-X*thetaIV;
sigmaIV     = e'*e/n;
varIV       = sigmaIV*inv(X'*PZ*X);

fid = fopen('AKoutput.dat', 'a');
fprintf(fid, '(ii) 95 CI for returns to education parameter: (%6.5f, %6.5f)\n\n',...
    thetaIV(2,1)-1.96*sqrt(varIV(2,2)),thetaIV(2,1)+1.96*sqrt(varIV(2,2)));
fclose(fid);

%------------------------------------------;
% (iii)Overidentifying restrictions test   ;
%------------------------------------------;

J=e'*PZ*e/sigmaIV;
[nrz,q]=size(Z);
[nrx,p]=size(X);
df=q-p;
pval=1-chi2cdf(J,df);

fid = fopen('AKoutput.dat', 'a');
fprintf(fid,'(iii) The overidentifying restrictions test is: %6.5f with p-value %1.4f', J,pval );

fclose(fid);