求助garch模型参数估计的代码程序

求助garch模型参数估计的代码程序

See also
========================
library cml;
#include cml.ext;
cmlset;
dsn = "garch";
depvar = { INDEX };
indvar = { };

p = 2;
q = 2;
output file = garch.out reset;

z = loadd(dsn);
if not scalmiss(indvar);
{ depvar,dv,indvar,iv } = indices2(dsn,depvar,indvar);
else;
{ depvar,dv } = indices(dsn,depvar);
endif;
if scalmiss(dv);
errorlog "error: variable not found in dataset "$+dsn;
end;
endif;
if not scalmiss(indvar);
z = z[.,dv]~ones(rows(z),1)~z[.,iv];
numx = rows(indvar) + 1;
_cml_ParNames = "CONST" | indvar | "kappa" |
(0$+"delta"$+ftocv(seqa(1,1,p),1,0)) |
(0$+"alpha"$+ftocv(seqa(1,1,q),1,0));
else;
z = z[.,dv]~ones(rows(z),1);
numx = 1;
_cml_ParNames = "CONST" | "kappa" |
(0$+"delta"$+ftocv(seqa(1,1,p),1,0)) |
(0$+"alpha"$+ftocv(seqa(1,1,q),1,0));
endif;

_ww_ = { -1e250 1e250 };
_cml_Bounds = ones(numx+1+p+q,2).*_ww_;
_cml_Bounds[numx+1,1] = .01;
_cml_Bounds[numx+2:numx+p+q+1,1] = zeros(p+q,1);
_cml_c = zeros(1,rows(_cml_Parnames));
_cml_c[1,numx+2:numx+p+q+1] = -ones(1,p+q);
_cml_d = -.99;


r = maxc(p|q);
__weight = (rows(z)/(rows(z)-r))*ones(rows(z),1);
__weight[1:r] = zeros(r,1);

if not scalmiss(indvar);
start = z[.,1] / z[.,2:cols(z)];
else;
start = meanc(z[.,1]);
endif;
start = start | 1 | .01 * ones(p+q,1);
_cml_GradProc = &grd;
{ b,f,g,h,ret } = cml( z, 0, &garch, start );
__title = "Wald confidence limits";
cl1 = cmlclimits(b,h);
call cmlclprt(b,f,g,cl1,ret);
output off;

proc garch( b, x );
local u2,kappa,delta,alpha,h;
u2 = (x[.,1] - x[.,2:cols(x)] * b[1:numx])^2;
kappa = b[numx+1];
delta = b[numx+2:numx+p+1];
alpha = b[numx+p+2:numx+p+q+1];
h = recserar(kappa+_shape(u2,q,1)*alpha,meanc(u2)*ones(p,1),delta);
retp( -0.5*( (u2 ./ h) + ln(2 * pi) + ln(h) ) );
endp;
proc grd(b,x);
local j,u,u2,kappa,delta,alpha,r,h,s,g;
u = x[.,1] - x[.,2:cols(x)] * b[1:numx];
u2 = u .* u;
kappa = b[numx+1];
delta = b[numx+2:numx+p+1];
alpha = b[numx+p+2:numx+p+q+1];
g = ones(rows(x),rows(b));
s = _shape(u2,q,1);
h = recserar(kappa+s*alpha,meanc(u2)*ones(p,1),delta);
g[.,1:numx] = -2*_dshape(u.*x[.,2:cols(x)],q,alpha);
g[.,numx+2:numx+p+1] = _shape(h,p,-1);
g[.,numx+p+2:numx+p+q+1] = s;
g = (0.5*(u2./h - 1)./h) .*
recserar(g,zeros(p,cols(g)),delta.*ones(1,cols(g)));
g[.,1:numx] = g[.,1:numx] + (u./h).*x[.,2:cols(x)];
retp(g);
endp;
proc _shape(z,m,r);
local y,n,v;
n = rows(z) - m - 1;
y = zeros(rows(z),m);
if r == 1;
v = seqa(1,1,m)' + seqa(0,1,n+1);
else;
v = seqa(m,-1,m)' + seqa(0,1,n+1);
endif;
y[m+1:rows(z),.] = reshape(z[v],n+1,m);
retp(y);
endp;
proc _dshape(z,m,a);
local y,j,n,v;
n = rows(z) - m - 1;
y = zeros(rows(z),cols(z));
v = seqa(1,1,m)' + seqa(0,1,n+1);
j = 1;
do until j > cols(z);
y[m+1:rows(z),j] = reshape(z[v,j],n+1,m) * a;
j = j + 1;
endo;
retp(y);
endp;

如以上回答内容为低俗、色情、不良、暴力、侵权、涉及违法等信息，可以点下面链接进行举报！

点此我要举报以上问答信息