请台湾同胞把程序贴出来，大家讨论

new;
library optmum;

format /m1 /rd 17,7;

n=10460;

load data[n,17]=D:\final.txt;

y=data[1:n,6];

x1=data[1:n,8];x2=data[1:n,9];x3=data[1:n,10];x4=data[1:n,11];x5=data[1:n,12];

z1=data[1:n,14];z2=data[1:n,17];z3=data[1:n,13];z4=data[1:n,7];z5=data[1:n,16];

at=data[1:n,4]; da=data[1:n,2]; disp=data[1:n,15];

T=rows(y);


@=================================================================@


  START=1;  @1993.....@

c1= .1324215   ;               @  Const1  @

c2= -.0000593 ;            @  pe  @
c3=.0123564 ;            @  booktomarket  @
c4=-.1677171 ;            @  cp  @
c5=.1027254 ;            @  leverageratio  @
c6= .0038986 ;            @  dividendratio  @
c7=  0.4344141 ;            @  std1  @

PRMTR_IN=c1~c2~c3~c4~c5~c6~c7;

PRMTR_IN=PRMTR_IN';

@ Maximum Likelihood Estimation @
@==================================================@

{xout,fout,gout,cout}=optmum(&lik_fcn,PRMTR_in);

PRM_FNL=TRANS(xout);             @ Estimated  coefficients, constrained@

output on;
"==FINAL OUTPUT========================================================";

"initial values of prmtr is";
                  trans(prmtr_in)';
"==============================================================";
"Calculating Hessian..... Please be patient!!!!";
                  hout0=hessp(&lik_fcn,xout);
                  hout=inv(hout0);

grdn_fnl=gradfd(&TRANS,xout);
Hsn_fnl=grdn_fnl*hout*grdn_fnl';
SD_fnl =sqrt(diag(Hsn_fnl)); @Standard errors of the estimated coefficients@



"===============================================================";

print "likelihood value is  "; -1*fout;
print "Estimated parameters are:"; xout;
print "Standard errors of parameters are:"; sd_fnl;
end;


@===========The Forcasting Performance ================@

        Cont=PRM_FNL[1,1];
      
         b1=PRM_FNL[2,1];
         b2=PRM_FNL[3,1];
         b3=PRM_FNL[4,1];
         b4=PRM_FNL[5,1];         
         b5=PRM_FNL[6,1];
         Var2=PRM_FNL[7,1]^2;

i=1;Ce=0;SSE=0;SAE=0;
do until i>10460;
ce=y[i,1]-cont-b1*x1[i,1]-b2*x2[i,1]-b3*x3[i,1]-b4*x4[i,1]-b5*x5[i,1];
SSE=SSE+ce^2;
SAE=SAE+abs(ce);
i=i+1;
endo;
M=meanc(y);

RWSSE=sumc((y[1:T,1]-M)^2);
RWSAE=sumc(abs(y[1:T,1]-M));

end;
@========================================================================@
proc lik_fcn(PRMTR1);

local prmtr,ppr,qpr,pr_vl,pr_val,likv,lik_fcn,pr_tr,prob_t,prob_dd,A,lik,
      EN,i,tmp,pr_stt0,pr_stl0,PAI,e,pr_trf,prob_,a0,a1,b0,b1,g1,PRMTR_in,
      g2,pro_,VAR,b2,pr_trf0,u0,u1,VAR_pre,MEAN_pre,b3,b4,b5,gradfd,
      Resd_pre,Resd_filt,VAR_fil,cont,h1,h2,c,d,c1,m,u,grdn_fnl,Hsn_fnl,
      cont_a,cont_u,h1_a,h1_u,h2_a,h2_u,a_a,a_u,xout,hout0,hout,SD_fnl,
      b_a,b_u,g2_a,g2_u,g1_a,g1_u,cor1,cor2,Cont1,Cont2,maxlik,
      COV,var2,ce,sse,sae,e1;

PRMTR=trans(PRMTR1);

LOCATE 15,1; PRMTR';

/* DEFINE PARAMETERS */

        Cont=PRMTR[1,1];

         b1=PRMTR[2,1];
         b2=PRMTR[3,1];
         b3=PRMTR[4,1];
         b4=PRMTR[5,1];         
         b5=PRMTR[6,1];
         Var2=PRMTR[7,1]^2;


  LIKV=0.0;
  i=1;e=0;
  do until i>10460;
     
   VAR=Var2;                                       
   
   PAI=1./SQRT(2.*PI.*VAR);
   e=y[i,1]-cont-b1*x1[i,1]-b2*x2[i,1]-b3*x3[i,1]-b4*x4[i,1]-b5*x5[i,1];  

   e=e./sqrt(var);

   PR_VL=PAI.*EXP(-0.5*e.*e);
                           
   PR_VAL=SUMC(PR_VL);        
   LIK=-1*LN(PR_VAL);

   
LIKV = LIKV+LIK;
i= i+1;
endo;

               LOCATE 2,35;"LIKV=";;LIKV;
retp(LIKV);
endp;
@++++++++++++++++++++++++++++++++++++++++++++++++++++++++@

@++++++++++++++++++++++++++++++++++++++++++++++++++++++++@



@++++++++++++++++++++++++++++++++++++++++++++++++++++++++@




proc TRANS(c0); @ constraining values of reg. coeff.@
         local c1,m,u,d1,d2,d3,d4,d5,d6;

         c1=c0;

  
retp(c1);
endp;