% Program to solve Hands-on 1.5
% written by Marcel Bluhm in January 2012
% prepare machine
clc;
clear all;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Simulate Data %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Note: Since it is a pooled regression, no time or cross-section ids
% necessary
N=100; % number of cross-sections
T=50; % number of time periods
theta_ini=[2;5];
X=[ones(N*T,1),sqrt(5)*randn(N*T,1)];
e=randn(N*T,1);
Y=X*theta_ini+e;
Y(Y<0)=0;
Y(Y>0)=1;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Estimate coefficient vector via ML %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% - Generate starting guess via OLS
theta_start=(X'*X)^(-1)*X'*Y;
%minus_likelihood=pooledprobit(theta_start,Y,X);
%minus_likelihood=pooledprobit(theta,Y,X);
[theta,fval,exitflag,output,s,H]=fminunc('pooledprobit',theta_start,[],Y,X);
avar_theta=inv(H);
t_stat=theta./(diag(avar_theta).^(.5));
result=[theta,t_stat];
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Run MC Simulation with 5000 replications %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
storemat=[];
for i=1:1000
i
e=randn(N*T,1);
Y=X*theta_ini+e;
Y(Y<0)=0;
Y(Y>0)=1;
theta_start=(X'*X)^(-1)*X'*Y;
theta=fminunc('pooledprobit',theta_start,[],Y,X);
storemat=[storemat,theta];
end
result_MC=mean(storemat,2);