SAS中文论坛

标题: 向高手请教极大似然估计程序问题 [打印本页]

作者: shiyiming 时间: 2008-3-4 09:28
标题: 向高手请教极大似然估计程序问题
各位高手，这是朱世武老师《金融计算与建模》中第23章里的程序，要解决的是CIR模型利率期限结构拟合的问题。其中B2W是7日回购利率两周加权指数平均数，是时间序列数据。我用的数据和朱老师用的时期不一样，但是C2、C3、C4、C5、C6和书上算出来的数差距不大，Ｃ1差距非常大，书上是零点几，我算出来是九点几。另一个问题是下面的“ｆ＝”这段程序提交后总是说语句错误或顺序不对，可是也没有发现什么错误啊。还看不懂con={0 0 0 0 -5,1 50 1 1 5};　x={0.2 5 0.02 0.2 0.8};这两个语句。请各位高手帮帮忙啊，急切想弄清这个问题！
data Resdat.B2W;
set Resdat.Bchmkir;
where code='B2W';
keep date ir;
run;
data B2W;
set Resdat.B2W;
where '05jan2001'd<=date<='13aug2003'd;
run;
data B2W;
set B2W;
dir=dif(ir);
lagir=lag(ir);
ivir=1/ lagir;
if dir=. then delete;
run;
proc iml;
reset deflib=work;
use B2W;
list;
list all;
read all var {lagir} into ir;
read all var {dir} into dir;
read all var {ivir} into ivir;
print ir dir ivir;
store ir dir ivir;
run;
quit;
/* 23.1.5 CIR模型利率期限结构拟合 */
proc iml;
reset deflib=work;
load ir dir ivir;
[color=#FF0040:1vyn82iv]c1[/color:1vyn82iv]=sum(dir##2#365#ivir);
c2=-2#sum(dir#ivir);
c3=2*sum(dir);
c4=sum((1/365)#ivir);
c5=-2*656/365;
c6=sum(ir#(1/365));
print c1 c2 c3 c4 c5 c6;
quit;
proc iml;
reset deflib=work;
start F_BETTS(x);
[color=#FF0040:1vyn82iv]f[/color:1vyn82iv]=-(656/2)*log(x[1]**2)-1/2*(x[1]**(-2))*(0.0352866416+(0.3776176352)*x[2]*x[3]+(-0.008852)*x[2]+(82.768352682)*x[2]*x[2]*x[3]*x[3]+(-3.594520548)*x[2]*x[2]*x[3]+(0.0392739233)*x[2]*x[2])-(4/2)*log(x[4]**2) /* = */
-(1/2)*(x[4]**(-2))*((log(0.9465)-log(2*((x[2]+x[5])**2+2*x[1]**2)**(1/2)*exp((x[2]+x[5]+((x[2]+x[5])**2+2*x[1]**2)**(1/2))*2.202740/2)/((x[2]+x[5]+((x[2]+x[5])**2+2*x[1]**2)**(1/2))*(exp(((x[2]+x[5])**2+2*x[1]**2)**(1/2)*2.202740)-1)+2*((x[2]+x[5])**2+2*x[1]**2)**(1/2)))*2*x[2]*x[3]/(x[1]**2)-2*(exp(((x[2]+x[5])**2+2*x[1]**2)**(1/2)*2.202740)-1)/((x[2]+x[5]+((x[2]+x[5])**2+2*x[1]**2)**(1/2))*(exp(((x[2]+x[5])**2+2*x[1]**2)**(1/2)*2.202740)-1)+2*((x[2]+x[5])**2+2*x[1]**2)**(1/2))*0.021027)**2)
-(1/2)*(x[4]**(-2))*((log(0.9802)-log(2*((x[2]+x[5])**2+2*x[1]**2)**(1/2)*exp((x[2]+x[5]+((x[2]+x[5])**2+2*x[1]**2)**(1/2))*0.880837/2)/((x[2]+x[5]+((x[2]+x[5])**2+2*x[1]**2)**(1/2))*(exp(((x[2]+x[5])**2+2*x[1]**2)**(1/2)*0.880837)-1)+2*((x[2]+x[5])**2+2*x[1]**2)**(1/2)))*2*x[2]*x[3]/(x[1]**2)-2*(exp(((x[2]+x[5])**2+2*x[1]**2)**(1/2)*0.880837)-1)/((x[2]+x[5]+((x[2]+x[5])**2+2*x[1]**2)**(1/2))*(exp(((x[2]+x[5])**2+2*x[1]**2)**(1/2)*0.880837)-1)+2*((x[2]+x[5])**2+2*x[1]**2)**(1/2))*0.021027)**2)
-(1/2)*(x[4]**(-2))*((log(0.9800)-log(2*((x[2]+x[5])**2+2*x[1]**2)**(1/2)*exp((x[2]+x[5]+((x[2]+x[5])**2+2*x[1]**2)**(1/2))*0.88/2)/((x[2]+x[5]+((x[2]+x[5])**2+2*x[1]**2)**(1/2))*(exp(((x[2]+x[5])**2+2*x[1]**2)**(1/2)*0.88)-1)+2*((x[2]+x[5])**2+2*x[1]**2)**(1/2)))*2*x[2]*x[3]/(x[1]**2)-2*(exp(((x[2]+x[5])**2+2*x[1]**2)**(1/2)*0.88)-1)/((x[2]+x[5]+((x[2]+x[5])**2+2*x[1]**2)**(1/2))*(exp(((x[2]+x[5])**2+2*x[1]**2)**(1/2)*0.88)-1)+2*((x[2]+x[5])**2+2*x[1]**2)**(1/2))*0.021027)**2)
-(1/2)*(x[4]**(-2))*((log(0.9905)-log(2*((x[2]+x[5])**2+2*x[1]**2)**(1/2)*exp((x[2]+x[5]+((x[2]+x[5])**2+2*x[1]**2)**(1/2))*0.405427/2)/((x[2]+x[5]+((x[2]+x[5])**2+2*x[1]**2)**(1/2))*(exp(((x[2]+x[5])**2+2*x[1]**2)**(1/2)*0.405427)-1)+2*((x[2]+x[5])**2+2*x[1]**2)**(1/2)))*2*x[2]*x[3]/(x[1]**2)-2*(exp(((x[2]+x[5])**2+2*x[1]**2)**(1/2)*0.405427)-1)/((x[2]+x[5]+((x[2]+x[5])**2+2*x[1]**2)**(1/2))*(exp(((x[2]+x[5])**2+2*x[1]**2)**(1/2)*0.405427)-1)+2*((x[2]+x[5])**2+2*x[1]**2)**(1/2))*0.021027)**2)
-(1/2)*(x[4]**(-2))*((log(0.984)-log(2*((x[2]+x[5])**2+2*x[1]**2)**(1/2)*exp((x[2]+x[5]+((x[2]+x[5])**2+2*x[1]**2)**(1/2))*0.681383/2)/((x[2]+x[5]+((x[2]+x[5])**2+2*x[1]**2)**(1/2))*(exp(((x[2]+x[5])**2+2*x[1]**2)**(1/2)*0.681383)-1)+2*((x[2]+x[5])**2+2*x[1]**2)**(1/2)))*2*x[2]*x[3]/(x[1]**2)-2*(exp(((x[2]+x[5])**2+2*x[1]**2)**(1/2)*0.681383)-1)/((x[2]+x[5]+((x[2]+x[5])**2+2*x[1]**2)**(1/2))*(exp(((x[2]+x[5])**2+2*x[1]**2)**(1/2)*0.681383)-1)+2*((x[2]+x[5])**2+2*x[1]**2)**(1/2))*0.021027)**2)
;
return(f);
finish F_BETTS;
con={0 0 0 0 -5,1 50 1 1 5};
x={0.2 5 0.02 0.2 0.8};
optn={1 3};
call nlpnra(rc,xres,"F_BETTS",x,optn,con);
store xres;
quit;

欢迎光临 SAS中文论坛 (https://mysas.net/forum/)