%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Basic Tobit model for horse adoption data
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
tic; % start stop watch
rand('state',37); % set arbitrary seed for uniform draws
randn('state',37); % set arbitrary seed for normal draws
% Open log file
[fid]=fopen('c:\Klaus\AAEC6564\mlab\logs\mod5_tobit_adoption.txt','w');
if fid==-1;
warning('File could not be opened');
else;
disp('File opened successfully');
end;
% load c:\klaus\AAEC6564\mlab\worksp\adoption.txt;
% data=adoption;
% clear adoption;
%
% save c:\klaus\AAEC6564\mlab\worksp\adoption;
% 'ok'
% break
load c:\klaus\AAEC6564\mlab\worksp\adoption;
%2692 winning bets, each for a different horse
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Variable
%
% 1 id running id
% 2 year 2007 or 2008
% 3 winbet winning bet in dollars (min. bet =$125)
% 4 male 1=male
% 5 female 1=female (base category is "gelding", i.e. neutered
% 6 bay 1=color=bay
% 7 black
% 8 brown
% 9 buckskin
% 10 dun
% 11 gray
% 12 pinto
% 13 redroan
% 14 sorrel baseline color category is a mix, called "other"
% 15 age
% 16 training1= horse has received training prior to auction
y=(data(:,3)-125)./100; %interpret as premium over min. bet, in units of $100
n=length(y);
X=[ones(n,1) data(:,4:end)];
% contents of X
% 1 constant
% 2 male 1=male
% 3 female 1=female (base category is "gelding", i.e. neutered
% 4 bay 1=color=bay
% 5 black
% 6 brown
% 7 buckskin
% 8 dun
% 9 gray
% 10 pinto
% 11 redroan
% 12 sorrel baseline color category is a mix, called "other"
% 13 age age
% 14 training1= horse has received training prior to auction
f0=find(y==0);
f1=find(y>0);
l0=length(f0);
l1=length(f1);
pc=l0/n;
fprintf(fid,'percent censored =\t%6.3f\n',pc);
fprintf(fid,'\n');
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% starting values, priors, and tuners
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% general elements
r1=15000; % burn-ins
r2=5000; % keepers
R=r1+r2;
k=size(X,2);
% generic OLS
bols=inv(X'*X)*X'*y;
res=y-X*bols;
s2=(res'*res)/(n-k);
% elements for beta
mu0=zeros(k,1); %diffuse prior for mean of betas
V0=eye(k)*100; % diffuse prior for varcov of beta
betadraw=bols;
% elements for sig2
v0=1/2;
tau0=1/2; % diffuse prior shape and scale
sig2draw=s2; % use OLS variance as starting draw for GS (Tobias uses "1" in his code)
[betamat,sig2mat]=gs_tobit(X,y,k,n,r1,r2,mu0,V0,betadraw,v0,tau0,...
sig2draw,f0,f1);
'GS done'
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Convergence Checks - Geweke's CD diagnostic
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% put all draws together & run diagnostics
allmat=[betamat;sig2mat];
kdiag=klausdiagnostics_greater0(allmat);
fprintf(fid,'total number of iterations =\t%6.0f \n',R);
fprintf(fid,'burn-in iterations =\t%6.0f \n',r1);
% beta stuff
out=kdiag(1:k,:)';
fprintf(fid,'Output table for betas \n\n');
fprintf(fid,'mean\t\tstd\t\tp(>0)\t\tnse\t\tIEF\t\tm*\t\tcd\n\n');
fprintf(fid,'%6.3f\t%6.3f\t%6.3f\t%6.3f\t%6.3f\t%6.3f\t%6.3f\n',out);
fprintf(fid,'\n');
%sig2 stuff
out=kdiag(k+1:end,:)';
fprintf(fid,'Output table for sig2\n\n');
fprintf(fid,'mean\t\tstd\t\tp(>0)\t\tnse\t\tIEF\t\tm*\t\tcd\n\n');
fprintf(fid,'%6.3f\t%6.3f\t%6.3f\t%6.3f\t%6.3f\t%6.3f\t%6.3f\n',out);
fprintf(fid,'\n');
save c:\klaus\AAEC6564\mlab\worksp\mod5_tobit_adoption X betamat sig2mat;
finish = toc/60;
fprintf(fid,'Time elapsed in minutes \n\n');
fprintf(fid,'%6.3f\n',finish);
st=fclose(fid);
if st==0;
disp('File closed successfully');
else;
warning('Problem with closing file');
end;