经管之家App
让优质教育人人可得
立即打开
[此贴子已经被作者于2008-2-15 0:39:53编辑过]
就是这个问题！
We want to minimize the sse of the regression referenced in the program using a grid of
points on the unknown parameters beta10 to beta40. Again, if you find a program that
yield a result to this nonlinear regression, you may use it.
The optimal scenario is the one where beta10 is close to 0.67.
You can modify the program, improve on it, as long as you do not change the problem at
hand. It is a simplified version of a program I am running on my computer. I just cannot
seem to make it converge to anything.
try the initial grid suggested in the program, then adapt as you get results. Be careful
not to narrow it down too fast. This function has lots of minima, you may miss the
solution. beta10 should definitely be less than 1.
// the problem is: minimize the sum of squared errors of the following regression
// x[.,1]=beta10/(beta20+beta30*x[.,4]+beta40*x[.,5]-x[.,3])+epsilon
// we want an optimal beta10 close to 0.67
new; outwidth 80;
format /mb1 /ros 16,8;
//output file
output file=d:\CGPV05\new.out reset;
screen off;output off;
//input file
n=300; // # of observations
load x[n,7]=d:\CGPV08\gauss\x.txt;
// sum of total squares
sstfirst=x[.,1]'*x[.,1];
KM=30;KM1=15;
sse=zeros(KM^3*KM1,1);r2=zeros(KM^3*KM1,1);
sst=x[.,1]'*x[.,1];
screen off;
vbarhat=zeros(n,1);
ssel=zeros(n,1);
// grid on the parameters
/*results for bbar*/
theta=0.67969667;/* std (0.010407291)*/
cst=10.464800; /* (0.63402387 )*/
vbarz=5.275000; /* (0.19825535 )*/
vbarz2=-0.003109; /* (0.011462653)*/
beta10mi=0.4;beta10ma=0.98;
beta20mi=2.0;beta20ma=2.7;
beta30mi=2.0;beta30ma=2.4;
beta40mi=-6e-3;beta40ma=1e-3;
betj0=zeros(4,KM^3*KM1);beta0=zeros(4,KM^3*KM1);
betj1=zeros(4,KM);
k1=1;num=1;pum=0;
do until k1>KM;
betj1[1,k1]=beta10mi+(k1-1)*(beta10ma-beta10mi)/KM;
k2=1;
do until k2>KM;
betj1[2,k2]=beta20mi+(k2-1)*(beta20ma-beta20mi)/KM;
k3=1;
do until k3>KM1;
/*betj1[3,k3]=beta30mi+(k3-1)*(beta30ma-beta30mi)/KM;*/
k4=1;
do until k4>KM;
betj0[1,num]=betj1[1,k1];
betj0[2,num]=/*beta20mi*/betj1[2,k2];
betj0[3,num]=/*beta30mi//betj1[3,k3]*/beta30mi+(k3-1)*(beta30ma-beta30mi)/KM1;
betj0[4,num]=beta40mi+(k4-1)*(beta40ma-beta40mi)/KM;
// end of the grid
// sum of squared errors is sse, r-square is r2
beta0[.,num]=betj0[.,num];
sse[num]=ssec(beta0[.,num]);
r2[num]=1-sse[num]/sst;
screen off;
num=num+1;
fip:k4=k4+1;
endo;
k3=k3+1;
endo;
k2=k2+1;
endo;
k1=k1+1;
endo;
/* Final Output */
output on;
output file=d:\CGPV08\gauss\sqr_020808.out;
i=1;do until i>KM^3*KM1;
if(beta0[.,i]==zeros(4,1));
sse=1e256;
endif;
i=i+1;
endo;
d2=minindc(sse); bsol2=beta0[.,d2]; // optimal solution is bsol2
/*optimal solution is */
screen on;
print"";print "final output with square v";
print ""; print "KM" KM;
print "";print "KM1 # interval for the linear term" KM1;
print "";
print "SST" sst;
print "SSE" beta0[.,1:num-1]'~sse[1:num-1]~r2[1:num-1];
vbarhat=zeros(n,1);
vbarhat=beta0[2,d2]+beta0[3,d2].*x[.,4]+beta0[4,d2].*x[.,5];
ssel=zeros(n,1);
ssel=(x[.,1]-mbet(beta0[.,d2],x))^2;
reg=zeros(n,1);
reg=mbet(beta0[.,d2],x);
print "*************";
print "data for Q." x[.,1]~x[.,3]~vbarhat~ssel~reg;
print "^^^^^^^^^^^^^";
print "fitted values, z, bbar,Yhat" x[.,4]~x[.,3]~reg;
/* number of violations of vbar>bbar*/
vio=0;
i=1;do until i>n;
if(vbarhat<x[i,3]);vio=vio+1;endif;
i=i+1;
endo;
print "number of violations, number of obs" vio~n;
print""; print "d2" d2;
print "";print "optimal beta" bsol2;
print "";
print "sse~sst" sse[d2]~sst;
output off; /* file end */
screen on;
end;
/*********/
/* Procs */
proc ssec(bet);
local e;
e=(x[.,1]-0.5.*bet[1]/(bet[2]+bet[3]*x[.,4]+bet[4]*x[.,5]-x[.,3]));
//e=x[.,1]-bet[1].*abs(bet[2]+bet[3]*x[.,4]+bet[4]*x[.,5]-x[.,3])^(bet[1]-1)./(abs(bet[2]+bet[3]*x[.,4]+bet[4]*x[.,5]-x[.,3])^bet[1]);
retp(sumc(e^2));
endp;
/*mbet*/
proc mbet(bet,w);
local r, zx,bmloc;
zx=w[.,4];
bmloc=w[.,3];
r=0.5.*bet[1]/(bet[2]+bet[3].*zx+bet[4].*(zx.*zx)-bmloc);
//r=bet[1].*abs(bet[2]+bet[3]*x[.,4]+bet[4]*x[.,5]-bmloc)^(bet[1]-1)./(abs(bet[2]+bet[3]*x[.,4]+bet[4]*x[.,5]-bmloc)^bet[1]);
retp(r);
endp;
proc er(bet);
local e;
e=x[.,1]-0.5.*bet[1]./(bet[2]+bet[3]*x[.,4]+bet[4].*x[.,5]-x[.,3]);
//e=x[.,1]-bet[1].*abs(bet[2]+bet[3]*x[.,4]+bet[4]*x[.,5]-x[.,3])^(bet[1]-1)./(abs(bet[2]+bet[3]*x[.,4]+bet[4]*x[.,5]-x[.,3])^bet[1]);
retp(e);
endp;
[此贴子已经被作者于2008-2-15 1:35:39编辑过]
你的东西又臭又长，谁有耐心看完呢？！ 哈哈！
不过建议您可以试试 Gauss-Netwon Regression，
因为它具有 One-step efficient estimator 的特性；
Gauss code 您就斟酌自己修改一下，只要 GNR 前一步是 any consistent estimator 即可。
扫码加好友，拉您进群
全部版块
我的主页
收藏
