estm.m

Rice University的Robin C. Sickles professor开发的专门用于paneldata model test and estimation 的program

% Main Program for PDE.m
% Written by Wonho Song, September 2003
% Updated, February 2004
% E-mail: [e-mail], [e-mail]


function []=estm(y,x,n,g_para,pss_para,kss_para,EST);  

[nt,p]=size(x);  t=nt/n;

 trunc1 = g_para(1);
 trunc2 = g_para(2);
   tmtr = g_para(3);
out_fig = g_para(4:6);
out_tab = g_para(7:9);
firmeff = g_para(10);

     NB = pss_para(1);
 gr_pss = pss_para(2:4);
pss_out = pss_para(5);
     p1 = pss_para(6);

     LS = kss_para(1);
 gr_kss = kss_para(2:4);



tt=(1:t)';
tm=kron(ones(n,1),tt);    % linear time trend

if tmtr==1; x_tr=[tm x]; p1=p1+1; [nt,p]=size(x_tr); 
else; x_tr=x;
end;


if p1>=p; 
disp(' ')
disp('p1 must be less than the number of variables p.')
disp(' ')
break;
end;

if sum(EST)==0;
disp(' ')
disp('No estimator is chosen. Please select at least one estimator.')
disp(' ')
break;
end;

if trunc1>0&trunc1<1/n;
disp(' ')
disp('trunc1 is too small. Please choose larger value. ')
disp(' ')
break;
end;

if trunc2>0&trunc2<1/n;
disp(' ')
disp('trunc2 is too small. Please choose larger value. ')
disp(' ')
break;
end;


disp(' ')
disp('Please be patient!')
disp('It may take long time to compute estimates if the data size is large.')
disp(' ')


%************************************************************
%    Finding bandwidth for PSS estimators using bootstrap   *
%************************************************************

[s1,s2w,s2g,s3]=bootpss(y,x,n,NB,gr_pss,pss_out,p1,EST);
 

Label=[]; Label1=[]; Label2=[]; COEF=[]; SE=[]; EFFECTS=[]; EFFICIENCY=[]; 
AVGEFF=[]; metacurve=[];  allcurves=[]; lgd=[]; lgdtv=[]; R2all=[]; INDEFF=[];

disp(' ')
disp('Computing... ')
disp(' ')


%************************************************
%       TRADITIONAL PANEL DATA ESTIMATION       *
%************************************************

if EST(1)==1;

Label=[Label '   FIX    ' ];
Label1=[Label1 '   FIX    ' ];
Label2=[Label2 '   FIX    ' ];

[b_with,se_with,sigf,ef,FIX,EFIX,R2w]=panelf(y,x_tr,n);

INDEFF=[INDEFF FIX];

if trunc1>0;
[temp2,J]=sort(EFIX);
mmm=J(ceil([1:trunc1*n n-trunc1*n+1:n]));
EFIX(mmm)=[];
[n2,temp]=size(EFIX);
else; n2=n;
end;


FIX=kron(FIX,ones(t,1));
EFIX=kron(EFIX,ones(t,1));

COEF=[COEF b_with];
SE=[SE se_with];

if trunc1==0&trunc2==0;
EFFECTS=[EFFECTS FIX];
EFFICIENCY=[EFFICIENCY EFIX];
end;

AVGEFF=[AVGEFF mean(EFIX)];

VFIX=reshape(EFIX,t,n2);      VFIX2=mean(VFIX')';
if firmeff==1; save VFIX; end;


R2all=[R2all R2w];

allcurves=[allcurves VFIX2];
lgd=[lgd; 'FIX  '];

disp('FIX done')

elseif EST(1)==2;

Label=[Label '   FIX   ' '    RND    '];
Label1=[Label1 '   FIX   ' '    RND    '];
Label2=[Label2 '   FIX   ' '    RND    '];

[b_with,se_with,sigf,ef,b_gls,se_gls,sigr,er,hw,FIX,RND,EFIX,ERND,R2w,R2r]=panel(y,x_tr,n);

INDEFF=[INDEFF FIX RND];

if trunc1>0;
[temp2,J]=sort(EFIX);
mmm=J(ceil([1:trunc1*n n-trunc1*n+1:n]));
EFIX(mmm)=[];
[nw2,temp]=size(EFIX);

[temp2,J]=sort(ERND);
mmm=J(ceil([1:trunc1*n n-trunc1*n+1:n]));
ERND(mmm)=[];
[ng2,temp]=size(ERND);
else; nw2=n; ng2=n;
end;


FIX=kron(FIX,ones(t,1));
RND=kron(RND,ones(t,1));

EFIX=kron(EFIX,ones(t,1));
ERND=kron(ERND,ones(t,1));

COEF=[COEF b_with b_gls];
SE=[SE se_with se_gls];

if trunc1==0&trunc2==0;
EFFECTS=[EFFECTS FIX RND];
EFFICIENCY=[EFFICIENCY EFIX ERND];
end;

AVGEFF=[AVGEFF mean(EFIX) mean(ERND)];

VFIX=reshape(EFIX,t,nw2);      VFIX2=mean(VFIX')';
VRND=reshape(ERND,t,ng2);      VRND2=mean(VRND')';

if firmeff==1; save VFIX; save VRND; end;

R2all=[R2all R2w R2r];

allcurves=[allcurves VFIX2 VRND2];
lgd=[lgd; 'FIX  ';'RND  '];

disp('FIX and RND done')

end;



%************************************************
%              Hausman-Taylor                   *
%************************************************

if EST(2) ==1;

Label=[Label '   HT     '];
Label1=[Label1 '   HT     '];
Label2=[Label2 '   HT     '];

[b_ht,se_ht,HT,EHT,R2h]=ht(y,x_tr,p1,n);

INDEFF=[INDEFF HT];

if trunc1>0;
[temp2,J]=sort(EHT);
mmm=J(ceil([1:trunc1*n n-trunc1*n+1:n]));
EHT(mmm)=[];
[n2,temp]=size(EHT);
else; n2=n;
end;


HT=kron(HT,ones(t,1));
EHT=kron(EHT,ones(t,1));

COEF=[COEF b_ht];
SE=[SE se_ht];

if trunc1==0&trunc2==0;
EFFECTS=[EFFECTS HT];
EFFICIENCY=[EFFICIENCY EHT];
end;

AVGEFF=[AVGEFF mean(EHT)];

VHT=reshape(EHT,t,n2);        VHT2=mean(VHT')';

if firmeff==1; save VHT; end;

R2all=[R2all R2h];

allcurves=[allcurves  VHT2];
lgd=[lgd; 'HT   '];

disp('HT done')

end;


%************************************************
%                  PSS1                         *
%************************************************

if EST(3)==1;

Label=[Label '  PSS1    '];
Label1=[Label1 '  PSS1    '];
Label2=[Label2 '  PSS1    '];

[b_pss1,se_pss1,PSS1,EPSS1,R2p1]=pss1(y,x_tr,p1,n,s1);

INDEFF=[INDEFF PSS1];

if trunc1>0;
[temp2,J]=sort(EPSS1);
mmm=J(ceil([1:trunc1*n n-trunc1*n+1:n]));
EPSS1(mmm)=[];
[n2,temp]=size(EPSS1);
else; n2=n;
end;


PSS1=kron(PSS1,ones(t,1));
EPSS1=kron(EPSS1,ones(t,1));

COEF=[COEF b_pss1];
SE=[SE se_pss1];

if trunc1==0&trunc2==0;
EFFECTS=[EFFECTS PSS1];
EFFICIENCY=[EFFICIENCY EPSS1];
end;

AVGEFF=[AVGEFF mean(EPSS1)];

VPSS1=reshape(EPSS1,t,n2);    VPSS12=mean(VPSS1')';

if firmeff==1; save VPSS1; end;

R2all=[R2all R2p1];

allcurves=[allcurves  VPSS12];
lgd=[lgd; 'PSS1 '];

disp('PSS1 done')

end;


%************************************************
%                   PSS2                        *
%************************************************

if EST(4)==1;

Label=[Label '  PSS2W  ' '   PSS2G  '];
Label1=[Label1 '  PSS2W  ' '   PSS2G  '];
Label2=[Label2 '  PSS2W  ' '   PSS2G  '];

[b_pss2w,se_pss2w,temp1,temp2,PSS2W,temp4,EPSS2W,temp5,R2p2w,temp6]=pss2(y,x_tr,n,s2w);

[temp3,temp,b_pss2g,se_pss2g,temp4,PSS2G,temp5,EPSS2G,temp6,R2p2g]=pss2(y,x_tr,n,s2g);

INDEFF=[INDEFF PSS2W PSS2G];

if trunc1>0;
[temp2,J]=sort(EPSS2W);
mmm=J(ceil([1:trunc1*n n-trunc1*n+1:n]));
EPSS2W(mmm)=[];
[nw2,temp]=size(EPSS2W);

[temp2,J]=sort(EPSS2G);
mmm=J(ceil([1:trunc1*n n-trunc1*n+1:n]));
EPSS2G(mmm)=[];
[ng2,temp]=size(EPSS2G);
else; nw2=n; ng2=n;
end;


PSS2W=kron(PSS2W,ones(t,1));
EPSS2W=kron(EPSS2W,ones(t,1));

PSS2G=kron(PSS2G,ones(t,1));
EPSS2G=kron(EPSS2G,ones(t,1));

COEF=[COEF b_pss2w b_pss2g];
SE=[SE se_pss2w se_pss2g];

if trunc1==0&trunc2==0;
EFFECTS=[EFFECTS PSS2W PSS2G];
EFFICIENCY=[EFFICIENCY EPSS2W EPSS2G];
end;

AVGEFF=[AVGEFF mean(EPSS2W) mean(EPSS2G)];

VPSS2W=reshape(EPSS2W,t,nw2);  VPSS2W2=mean(VPSS2W')';
VPSS2G=reshape(EPSS2G,t,ng2);  VPSS2G2=mean(VPSS2G')';

if firmeff==1; save VPSS2W; save VPSS2G; end;

R2all=[R2all R2p2w R2p2g];

allcurves=[allcurves  VPSS2W2 VPSS2G2];
lgd=[lgd; 'PSS2W';'PSS2G'];

disp('PSS2 done')

end;


%************************************************
%                     PSS3                      *
%************************************************

if EST(5) ==1;

Label=[Label '   PSS3   '];
Label1=[Label1 '   PSS3   '];
Label2=[Label2 '   PSS3   '];

[b_pss3,b_pss3l,se_pss3,se_pss3l,PSS3,EPSS3,betai,gammai,se_beta,se_gamma,R2p3]=pss3(y,x_tr,n,s3);

INDEFF=[INDEFF PSS3];

if trunc1>0;
[temp2,J]=sort(EPSS3);
mmm=J(ceil([1:trunc1*n n-trunc1*n+1:n]));
EPSS3(mmm)=[];
[n2,temp]=size(EPSS3);
else; n2=n;
end;


PSS3=kron(PSS3,ones(t,1));
EPSS3=kron(EPSS3,ones(t,1));

COEF=[COEF b_pss3];
SE=[SE se_pss3];

if trunc1==0&trunc2==0;
EFFECTS=[EFFECTS PSS3];
EFFICIENCY=[EFFICIENCY EPSS3];
end;

AVGEFF=[AVGEFF mean(EPSS3)];

VPSS3=reshape(EPSS3,t,n2);    VPSS32=mean(VPSS3')';

if firmeff==1; save VPSS3; end;


R2all=[R2all R2p3];

allcurves=[allcurves  VPSS32];
lgd=[lgd; 'PSS3 '];

disp('PSS3 done')

end;


%***********************************************
%                  CSS                         *
%***********************************************

if EST(6) ==1;

Label=[Label '   CSSW   '];
Label1=[Label1 '   CSSW   '];
Label2=[Label2 '   CSSW   '];

[b_cssw,se_cssw,CSSW,R2cw]=cssw(y,x,n);

INDEFF=[INDEFF mean(CSSW)'];

if trunc1>0;

mmm=[];
for i=1:t
[temp2,J]=sort(CSSW(i,:));
mmm=[mmm J(ceil([1:trunc1*n n-trunc1*n+1:n]))];
end;
I=unique(mmm);  
CSSW(:,I)=[];
[t2,n2]=size(CSSW);
nt2=t2*n2;

else; nt2=nt; n2=n;
end;


ECSSW=[]; 
for i=1:t
ECSSW=[ECSSW;exp(CSSW(i,:)-max(CSSW(i,:)))];
end;

CSSW=reshape(CSSW,nt2,1);
ECSSW=reshape(ECSSW,nt2,1);

if tmtr==1; b_cssw=[0;b_cssw]; se_cssw=[0;se_cssw]; end;

COEF=[COEF b_cssw ];
SE=[SE se_cssw ];

if trunc1==0&trunc2==0;
EFFECTS=[EFFECTS CSSW];
EFFICIENCY=[EFFICIENCY ECSSW];
end;

AVGEFF=[AVGEFF mean(ECSSW)];

VCSSW=reshape(ECSSW,t,n2);  VCSSW2=mean(VCSSW')';

if firmeff==1; save VCSSW; end;

R2all=[R2all R2cw];


metacurve=[metacurve VCSSW2];
allcurves=[allcurves VCSSW2];
lgd=[lgd; 'CSSW '];
lgdtv=[lgdtv; 'CSSW '];

disp('CSSW done')