求助：向量BEKK－GARCH（1，1）模型参数极大似然估计BHHH算法的程序

是 否 +2 论坛币

k人 参与回答

经管之家送您一份

应届毕业生专属福利!

求职就业群

赵安豆老师微信：zhaoandou666

经管之家联合CDA

送您一个全额奖学金名额~ !

立即领取

感谢您参与论坛问题回答

经管之家送您两个论坛币！

+2 论坛币

做Bivariate BEKK－GARCH（1，1）模型，做到最后系数估计不出来了，尤其那个搜索步长是怎么编的，急救！谢了！

我不会EViews编程，现学matlab编，希望能得到各位高手指教！

扫码加我 拉你入群

请注明：姓名-公司-职位

以便审核进群资格，未注明则拒绝

分享0 收藏6 回帖

关键词：极大似然估计 GARCH 似然估计 极大似然 BHHH 模型 算法 向量 极大似然估计 BHHH

http://www.anst.uu.se/matieklo/ox_code.htm

/* BHHH PROCEDURES */

/* Procedure for calculating one sided gradients of individual log-likelihood functions. */

fnGradient(const fnF,const vX0, vF0) { decl dEps,vF1,vF2,dN,vE,vEi,mG,dh, bUseTwoSided=FALSE; dEps = 1e-7; dN = rows(vX0); vE = zeros(dN,1); if (vF0==0) fnF(vX0,&vF0,0,0);; /* Compute function only if not supplied in argument */ mG = zeros(rows(vF0),dN); /* Create a matrix of zeros to store */ /* the first derivatives df(x0)/dx_i */

for (decl i=0;i<dN;++i){ /* Loop i from 1 to n=rows(vX) */ vEi = vE; vEi= 1; /* Construct e_i */ dh = max((fabs(vX0)*dEps),dEps); /* Calculate perturbation factor h */ if (!fnF(vX0+vEi*dh,&vF2,0,0)) println("Gradient evaluation failed."),exit(1); if (bUseTwoSided){ /* Calculate two-sided finite difference */ if (!fnF(vX0-vEi*dh,&vF1,0,0)) println("Gradient evaluation failed."),exit(1); mG[] = (vF2-vF1)/(2*dh); } else{ /* Calculate one-sided finite difference */ vF1 = vF0; mG[] = (vF2-vF1)/ dh; }

} /* Next i */ return mG; /* Return matrix of individual gradients */ }

MaxBHHH(const fnF, const avX, const adFunc, const amInvHess, const fNumDer) { decl vX0,vX1, /* Parameter vectors */ vF0,vF1,vF2, /* Function values (n*1)-vectors */ mG0,vG0,mG1,mG2, /* Gradients in matrix and vector form, Hessian, direction */ mH,mHi,vS, dEps1,dEps2, /* Tolerance for convergence tests */ fConv=MAX_NOCONV,fStepOK, /* Flags for indicating convergence etc*/ dLambda; /* Line search parameter, step length */ adFunc[0] = M_NAN; amInvHess[0] = M_NAN*unit(sizerc(avX[0]));

decl iMaxIter,iPrint,bCompact; [iMaxIter,iPrint,bCompact]=GetMaxControl(); // println("/* INITIALIZATION */"); vX0 = avX[0]; [dEps1,dEps2] = GetMaxControlEps();

// println("/* Check initial function values */"); decl dStart = timer(); fnF(vX0,&vF0,0,0); decl sTimeSpan = sprint((timer() - dStart) / 100); // println("\n\t- One function evaluation ",fNumDer ? "":"(excluding analytical gradients) ", // "took ", sTimeSpan," sec.");

dStart = timer(); decl fInitFunc = fnF(vX0,&vF0,fNumDer ? 0 : &mG0,0); // println("\n\t- One function evaluation ",fNumDer ? "":"(including analytical gradients) ", // "took ", sprint((timer() - dStart) / 100)," sec.");

if (!fInitFunc){ format(90); println("\nInitial function evaluation failed!", "\nInitial position:",vX0', "Initial function values:",vF0' ); return MAX_FUNC_FAIL; }

for (decl i=0;i<iMaxIter;++i){ if (isnan(vF0)) print("Function evaluation failed! vX0'=",vX0'); decl dIterStartTime = timer();

// println("\n/* STEP 1: CALCULATE FUNCTION VALUE, GRADIENT AND HESSIAN */"); // println("vF0 = ",double(sumc(vF0))); // mG = fnGradient(fnF,vX0,vF0); if (fNumDer) mG0 = fnGradient(fnF,vX0,vF0); if (isnan(mG0)){ println(vX0',mG0); return MAX_FUNC_FAIL; } vG0 = sumc(mG0); mH = mG0'mG0;

// println("/* STEP 2: CALCULATE DIRECTION */"); // mHi = invertsym(mH); /* This should/could be updated to LU-decomp solution */ // vS = mHi*vG0';

decl vscale = diagonal(mH)', /* determine new scale factors */ ml,vdiag,fsingular; vscale = vscale .> 0 .? 1 ./ sqrt(vscale) .: 1; // Override scale // println(vscale'); /* solve the system q=Qb as Sq=SQSa, Sa=b */ if (decldl(vscale' .* mH .* vscale, &ml, &vdiag) == 1) { /* solve system for scaled score; solution stored in vdelta */ vS = solveldl(ml, vdiag, vG0' .* vscale) .* vscale; fsingular = FALSE; } else { ml = mH; /* remember for return value */ println(mH); vS = -vG0'; /* steepest decent */ // println("vX0' = ",vX0',"vF0' = ",vF0',"mH = ",mH,"vS' = ",vS'); fsingular = TRUE; } if (vG0*vS<0){ println("Changing direction");vS=-vS; }

// println("/* STEP 2B: LINE SEARCH, Contraction or extension of step length */"); // println("vF0 = ","%12.8f",double(sumc(vF0))); dLambda = 1; /* Was 1 initialize step length dLambda */ decl fStepOK=FALSE, bContraction, bF1 = fnF(vX0+dLambda*vS,&vF1,0,0), // We don't need the derivatives in the line search bF2 = fnF(vX0+2*dLambda*vS,&vF2,0,0);

if (!bF1 || !bF2 || sumc(vF2)<sumc(vF1) ) // Full step lengt implies function fail or lower likelihood bContraction = TRUE; // Reduce step lengt else bContraction = FALSE; // otherwise try to increase step lengt

while (dLambda >1e-15 && dLambda < 1e5){ // Min/max step lengths dLambda = bContraction ? 0.5*dLambda : 2*dLambda; // New potential dLambda (step length)

decl bF2 = fnF(vX0+dLambda*vS,&vF2,0,0);

// println( "bF1 = ",bF1," sumc(vF1)=","%12.8f",double(sumc(vF1)), // "\nbF2 = ",bF2," sumc(vF2)=","%12.8f",double(sumc(vF2))," dLambda=",double(dLambda)); // Case 1. dLambda imples missing value => Continue if we are in contraction phase, // other wise stop. if (bF2==FALSE || bF1==FALSE) if (bContraction){ // println("Case 1a. (bF2==FALSE or bF1==FALSE and contraction => Continue to contract) ",dLambda); bF1=bF2; vF1=vF2; } else{ // println("Case 1b. (bF2==FALSE and expansion => Stop expanding) ",dLambda); dLambda = bContraction ? 2*dLambda : 0.5*dLambda; // Go one step back. fStepOK = TRUE; fnF(vX0+dLambda*vS,&vF2,fNumDer ? 0 : &mG1,0); // We need gradient at chosen step length break; } // Case 2. dLambda implies lower likelihood => Stop else if (sumc(vF2)<sumc(vF1)){ // println("Case 2. (sumc(vF2)<sumc(vF1) ",double(sumc(vF2))," ",double(sumc(vF1)), " dLambda=",dLambda); dLambda = bContraction ? 2*dLambda : 0.5*dLambda; // Go one step back. fStepOK = TRUE; fnF(vX0+dLambda*vS,&vF2,fNumDer ? 0 : &mG1,0); // We need gradient at chosen step length break; }

// Case 3. dLambda implies higher likelihood => Continue the search // in the same direction. else if (sumc(vF2)>sumc(vF1)){ // println("Case 3. (sumc(vF2)>sumc(vF1) ",double(sumc(vF2))," ",double(sumc(vF1))," dLambda=",dLambda); bF1=bF2; vF1=vF2; } } // vF1 and mG1 now contains the function values evaluated at vX1=vX0+dLambda*vS

// println("/* STEP 3: UPDATE PARAMETER VECTOR */"); vX1 = vX0 + dLambda*vS;

// println("/* STEP 4: TEST STRONG CONVERGENCE AT vX0, rel grad.<eps1 and rel. diff. in vX<eps2 */"); decl vpp = vX0 .? 1.0 / sqr(sizerc(vX0)) + fabs(vX0) .: 1; fConv = (fabs(vG0' .* vpp) < dEps1 && fabs((vX1-vX0) ./ vpp) < 10 * dEps1) ? MAX_CONV : MAX_NOCONV;

if (iPrint > 0 && imod(i, iPrint) == 0){ MaxMonitor("BHHH", i, vX0, vG0', vS, double(sumc(vF0)), double(sumc(vF0)), double(dLambda), 0, fConv, bCompact); // println("Iteration timer = ",timespan(dIterStartTime)); }

if (fConv==MAX_CONV || !fStepOK){ break; } else if (i==iMaxIter-1){ return MAX_MAXIT; } else{ // Update old parameter and function values vX0 = vX1; vF0 = vF1; if (fNumDer==FALSE) mG0=mG1; } } // println("// Update input arguments"); avX[0] = vX0; adFunc[0] = double(sumc(vF0)); mHi = invertsym(mH); amInvHess[0] = mHi ? mHi : M_NAN*mH;

if (fConv==MAX_CONV) return MAX_CONV; else{ decl vpp = vX0 .? 1.0 / sqr(sizerc(vX0)) + fabs(vX0) .: 1; fConv = fabs(vG0' .* vpp) < dEps2 ? MAX_WEAK_CONV : MAX_NOCONV; return fConv; } }

[此贴子已经被作者于2005-8-14 18:03:54编辑过]

try als

http://econ.ucsd.edu/software/ARCH.html

拜谢两位侠客指点，不胜感激！

可惜我对编程了解的不深，这么长的程序看了半天还是一头雾水，可否请再次赐教？这个程序是用什么语言编的？数据是怎么才能调进去的？可否直接出结果的？

第二位大哥给的网站打不开也不知道是怎么回事，甚为遗憾！

I have emailed the first two downloaded zip files to you.

check ur mail box.

You said that 可惜我对编程了解的不深，这么长的程序看了半天还是一头雾水.

given this, why do you choose this topic? this is not an easy one if you are Not well in programming.

good luck anyway.

用splus的garch module就可以做！！

大哥可真的是好人啦！邮件也收到了，谢谢了，也谢了您的祝愿！做这个题目确实是有点赶鸭子上阵之嫌，但我确实也是有这个兴趣来做这个模型，而且都做有段时间了，不想这么就放弃了。

splus真的可以很容易的实现吗？没玩过哎，不防一试！

谢谢两位！！

hlllkkkkkkkkkk

尴尴尬尬

请问，谁知道在eviews6.0中使用BHHH算法对EC-EGARCH模型进行检验的步骤啊？？？急切的盼望热心人回复啊