cmlmt calculation error

Just now, I debugged my program and found that the value of the matrixs become a dot after a number of looping. I dont know why and how to correct it. I post my code here, its about 250th loop after which error occurs.

p.s. Could someone tell me how to attach files in this forum?

Thanks for your patience.


library cmlmt;

#include cmlmt.sdf;

// load data
r1day1 = xlsReadM("SHIBOR.xlsx","b2:b2059",1,0);
r1day2 = xlsReadM("SHIBOR.xlsx","b3:b2060",1,0);

struct DS d0;
// create a default data structure
d0 = reshape(dsCreate,2,1);
d0[1].datamatrix = r1day1;
d0[2].datamatrix = r1day2;

struct PV p0;
// creates a default parameter structure
p0 = pvCreate;
p0 = pvPackm(p0,0|1,"meanR1",1|1);
// mean parameters in regime 1, arranged as a1,b1
p0 = pvPackm(p0,0|1,"meanR2",1|1);
// mean parameters in regime 2, arranged as a2,b2
p0 = pvPackm(p0,0|.3|.4|1|1,"varR1",1|1|1|1|1);
/* variance parameters in regime1, arranged as
omega1,alpha1,beta1,mu1,v1 */
p0 = pvPackm(p0,0|.4|.3|1|1,"varR2",1|1|1|1|1);
/* variance parameters in regime2, arranged as
omega2,alpha2,beta2,mu2,v2 */
p0 = pvPackm(p0,1|1,"TP",1|1);
// parameters for transition matrix, e1,e2
// p11 = 1/(1+exp(e1)), p22 = 1/(1+exp(e2))
// p0 = pvPackm(p0,1|1,"hstart",1|1);
// parameters for starting value of h(1,0)&h(2,0)
p0 = pvPackm(p0,0|1,"hf1",1|1);
// parameters for calculating f(delata) of regime 1
p0 = pvPackm(p0,0|1,"hf2",1|1);
// parameters for calculating f(delata) of regime 2

struct cmlmtControl c0;
// creates a default structure
c0 = cmlmtControlCreate;
c0.Bounds = { -50 50, -50 50,
-50 50, -50 50,
-10 10, .0001 50, .0001 50, .0001 50, .0001 50,
-10 10, .0001 50, .0001 50, .0001 50, .0001 50,
-50 50, -50 50,
// .0001 10, .0001 10,
-50 50, -50 50,
-50 50, -50 50 };
c0.NumObs = 2058;

// prints the results
call cmlmtprt(cmlmt(&rs_garch,p0,d0,c0));

proc rs_garch(struct pv p, struct ds d, ind);
local meanR1,meanR2,varR1,varR2,TP,hstart,hf1,hf2,tm,
p11,p22,p1,rowN,ex1,ex2,r1dayF,r1dayL,
fvector,pvector,pFvector,pIvector,hvector,hIvector,dIvector,tfvector;
struct modelResults mm;

meanR1 = pvUnpack(p,"meanR1");
/* a1 = meanR1[1];
b1 = meanR1[2]; */
meanR2 = pvUnpack(p,"meanR2");
/* a2 = meanR2[1];
b2 = meanR2[2]; */
varR1 = pvUnpack(p,"varR1");
/* omega1 = varR1[1];
alpha1 = varR1[2];
beta1 = varR1[3];
mu1 = varR1[4];
v1 = varR1[5]; */
varR2 = pvUnpack(p,"varR2");
/* omega2 = varR2[1];
alpha2 = varR2[2];
beta2 = varR2[3];
mu2 = varR2[4];
v2 = varR2[5]; */
TP = pvUnpack(p,"TP");
/* e1 = TP[1];
e2 = TP[2]; */
// hstart = pvUnpack(p,"hstart");
/* h1_0 = hstart[1];
h2_0 = hstart[2]; */
hf1 = pvUnpack(p,"hf1");
/* c1 = hf1[1];
d1 = hf1[2]; */
hf2 = pvUnpack(p,"hf2");
/* c2 = hf2[1];
d2 = hf2[2]; */

p11 = 1/(1+exp(TP[1]));
P22 = 1/(1+exp(TP[2]));
p1 = (1-p22)/(2-p11-p22); // invariant probability of being in regime 1
tm = p11 ~ 1-p11 | 1-p22 ~ p22;// transition matrix

r1dayL = d[1].datamatrix;
r1dayF = d[2].datamatrix;

rowN = rows(r1dayL);

if ind[1];

fvector = zeros(rowN,2); // a matrix consists of [f(r(t+1)|r(t~),s(t)=1),f(r(t+1)|r(t~),s(t)=2)]
pvector = zeros(rowN,2); // a vector consists of [p(s(t)=1|r(t~)),p(s(t)=2|r(t~))]
pIvector = zeros(rowN,2); // a matrix consists of [p1(1,t),p2(1,t)]
pFvector = zeros(rowN,2); // a matrix consists of [p(s(t-1)=1|r(t~)),p(s(t-1)=2|r(t~))]
hvector = zeros(rowN,2); // a matrix consists of [h(1,t),h(2,t)]
hIvector = zeros(rowN,2); // a matrix consists of [hi(1,t),hi(2,t)]
dIvector = zeros(rowN,2); // a matrix consists of [delta_i(1,t),delta_i(2,t)]
tfvector = zeros(rowN,1); // a vector consists of f(r(t+1)|r(t~))

ex1 = r1dayF-(meanR1[1]+meanR1[2]*r1dayL);
ex2 = r1dayF-(meanR2[1]+meanR2[2]*r1dayL);

hvector[1,1] = sumc(ex1^2)/rowN;
hvector[1,2] = sumc(ex2^2)/rowN;

fvector[1,1] = exp(-ex1[1]^2/(2*hvector[1,1]))/(2*pi*hvector[1,1])^.5;
fvector[1,2] = exp(-ex2[1]^2/(2*hvector[1,2]))/(2*pi*hvector[1,2])^.5;

pvector[1,.] = p1~(1-p1);
pFvector[1,.] = pvector[1,.].*fvector[1,.]./sumr(pvector[1,.].*fvector[1,.]);
pIvector[1,1] = p11*pvector[1,1]/(p11*pvector[1,1]+(1-p22)*pvector[1,2]);
pIvector[1,2] = (1-p11)*pvector[1,1]/((1-p11)*pvector[1,1]+p22*pvector[1,2]);
// hvector[1,.] = hstart';
hIvector[1,1] = pIvector[1,1]*hvector[1,1]+(1-pIvector[1,1])*hvector[1,2]
+pIvector[1,1]*(1-pIvector[1,1])*(ex1[1]-ex2[1])^2;
hIvector[1,2] = pIvector[1,2]*hvector[1,1]+(1-pIvector[1,2])*hvector[1,2]
+pIvector[1,2]*(1-pIvector[1,2])*(ex1[1]-ex2[1])^2;
dIvector[1,1] = pIvector[1,1]*ex1[1]/hvector[1,1]^.5+(1-pIvector[1,1])*ex2[1]/hvector[1,2]^.5;
dIvector[1,2] = pIvector[1,2]*ex1[1]/hvector[1,1]^.5+(1-pIvector[1,2])*ex2[1]/hvector[1,2]^.5;

tfvector[1] = sumr(fvector[1,.].*pvector[1,.]);

for i (2,rowN,1);
hvector[i,1] = (varR1[1]+varR1[2]*(hIvector[i-1,1])^(varR1[4]/2) * (abs(dIvector[i-1,1]-hf1[2])-hf1[1]*(dIvector[i-1,1]-hf1[2]))^varR1[5]
+varR1[3]*(hIvector[i-1,1])^(varR1[4]/2))^(2/varR1[4]);
hvector[i,2] = (varR2[1]+varR2[2]*(hIvector[i-1,2])^(varR2[4]/2) * (abs(dIvector[i-1,2]-hf2[2])-hf2[1]*(dIvector[i-1,2]-hf2[2]))^varR2[5]
+varR2[3]*(hIvector[i-1,2])^(varR2[4]/2))^(2/varR2[4]);

fvector[i,1] = exp(-ex1[i]^2/(2*hvector[i,1]))/(2*pi*hvector[i,1])^.5;
fvector[i,2] = exp(-ex2[i]^2/(2*hvector[i,2]))/(2*pi*hvector[i,2])^.5;

pvector[i,.] = pFvector[i-1,.]*tm;
// pvector[i,1] = p11*pFvector[i-1,1]+(1-p22)*pFvector[i-1,2];
// pvector[i,2] = p22*pFvector[i-1,2]+(1-p11)*pFvector[i-1,1];

pFvector[i,.] = pvector[i,.].*fvector[i,.]./sumr(pvector[i,.].*fvector[i,.]);
// pFvector[i,1] = pvector[i,1]*fvector[i,1]/(pvector[i,1]*fvector[i,1]+pvector[i,2]*fvector[i,2]);
// pFvector[i,2] = pvector[i,2]*fvector[i,2]/(pvector[i,1]*fvector[i,1]+pvector[i,2]*fvector[i,2]);

pIvector[i,1] = p11*pvector[i,1]/(p11*pvector[i,1]+(1-p22)*pvector[i,2]);
pIvector[i,2] = (1-p11)*pvector[i,1]/((1-p11)*pvector[i,1]+p22*pvector[i,2]);

hIvector[i,1] = pIvector[i,1]*hvector[i,1]+(1-pIvector[i,1])*hvector[i,2]
+pIvector[i,1]*(1-pIvector[i,1])*(ex1[1]-ex2[1])^2;
hIvector[i,2] = pIvector[i,2]*hvector[i,1]+(1-pIvector[i,2])*hvector[i,2]
+pIvector[i,2]*(1-pIvector[i,2])*(ex1[1]-ex2[1])^2;
dIvector[i,1] = pIvector[i,1]*ex1[1]/hvector[i,1]^.5+(1-pIvector[i,1])*ex2[1]/hvector[i,2]^.5;
dIvector[i,2] = pIvector[i,2]*ex1[1]/hvector[i,1]^.5+(1-pIvector[i,2])*ex2[1]/hvector[i,2]^.5;

tfvector[i] = sumr(fvector[i,.].*pvector[i,.]);

endfor;

mm.function = sumc(ln(tfvector));
endif;
retp(mm);
endp;

7 Answers



0



You are saying it is around the 250th iteration of the loop inside of rs_garch that starts with:

for i (2,rowN,1);

Is this correct? Also, do you know which line is the first to return a '.'?

aptech

1,338


0



I debugged it again, the first row after which dot occurs is now 92, the change is possibly due to the change of its starting value. The result for one of the matrix present like follows:


 

88            0.4838164443136511             0.7044154633233479

89                                 0                                                    0

90                                  0                                                    0

91                                  0                                                    0

92                                  .                                                        .


I dont know where the error be.

Thanks!

 



0



When estimating a garch model like:

r(t) = a + b * r(t-1) + u(t), where u(t) = h(t)^0.5*Z(t), Z(t) is N(0,1)

h(t) = alpha + beta * h(t-1) + zeta * u(t-1).

The parameter beta and zeta should be constrained to be positive, but

if I also constrain alpha to be positive, the result would be that the estimator of alpha is at 0.0001 which is the constraint.

If I constrain alpha to be bigger than some negative number (like -0.1), Gauss would possibly return that function evaluation failed.

Could you tell me how to deal with this issue?

Thanks a lot.



0



In a GARCH model you need alpha>0, zeta>0 and beta>=0 for identification. I personally don't use bounds for the parameters, especially when I only need to guarantee positivity. What about estimating x without constraints and set alpha=exp(x)? Then you have more room for alpha being very small but always positive.



0



What does the "INF" in the correlation matrix of the parameters mean?

For example, following is a result given by GAUSS:


1 -0.96203422 -0.1253878 0.18308101 . 0.41191654 +INF . -0.64990466 -0.79824687 -INF +INF -0.55452406 1.3747147
-0.96203422 1 0.093765118 -0.17692098 . -0.44130707 -INF . 0.76530259 0.96352451 +INF -INF 0.6663268 -1.5238216
-0.1253878 0.093765118 1 -0.95337865 . -0.068382173 +INF . -0.32826896 0.10443985 -INF +INF -0.62167148 0.33967216
0.18308101 -0.17692098 -0.95337865 1 . -0.0050880999 -INF . 0.57599187 -0.22418101 +INF -INF 0.89164092 -0.73449623
. . . . . . . . . . . . . .
0.41191654 -0.44130707 -0.068382173 -0.0050880999 . 1 -INF . -0.19454829 -0.45717146 +INF +INF -0.2895939 0.42120494
+INF -INF +INF -INF . -INF -INF . -INF +INF +INF -INF +INF -INF
. . . . . . . . . . . . . .
-0.64990466 0.76530259 -0.32826896 0.57599187 . -0.19454829 -INF . 1 0.18133867 -INF +INF -0.35575276 0.56408893
-0.79824687 0.96352451 0.10443985 -0.22418101 . -0.45717146 +INF . 0.18133867 1 -INF +INF -0.32258855 -0.60829708
-INF +INF -INF +INF . +INF +INF . -INF -INF -INF +INF -INF +INF
+INF -INF +INF -INF . +INF -INF . +INF +INF +INF -INF +INF -INF
-0.55452406 0.6663268 -0.62167148 0.89164092 . -0.2895939 +INF . -0.35575276 -0.32258855 -INF +INF 1 0.38185242
1.3747147 -1.5238216 0.33967216 -0.73449623 . 0.42120494 -INF . 0.56408893 -0.60829708 +INF -INF 0.38185242 1



0



INF means infinity. It is probably caused by a division by 0.

aptech

1,338


0



When using CMLMT library to estimate a complex model which I couldn't provide the analytical gradient and hessian, Gauss would possibly give the result "The covariance of the parameters failed to invert", then, I want to know how to modify so as to get a reasonable result?

Your Answer


You must login to post answers.

Have a Specific Question?

Get a real answer from a real person

Need Support?

Get help from our friendly experts.

Try GAUSS for 30 days for FREE

See what GAUSS can do for your data

© Aptech Systems, Inc. All rights reserved.

Privacy Policy