garch() uncertainty

As part of an on-going paper with Kerrie Mengersen and Pierre Pudlo, we are using a GARCH(1,1) model as a target. Thus, the model is of the form

$y_t=\sigma_t \epsilon_t \qquad \sigma^2_t = \alpha_0 + \alpha_1 y_{t-1}^2 + \beta_1 \sigma_{t-1}^2$

which is a somehow puzzling object: the latent (variance) part is deterministic and can be reconstructed exactly given the series and the parameters. However, estimation is not such an easy task and using the garch() function in the tseries package leads to puzzling results! Indeed, simulating data shows some high variability of the procedure against starting values:

genedata=function(para,nobs){

   pata=epst=sigt=rnorm(nobs)
   sigt[1]=sqrt(para[1])
   pata[1]=epst[1]*sigt[1]
   for (t in 2:nobs){
      sigt[t]=sqrt(para[1]+para[2]*pata[t-1]^2+para[3]*sigt[t-1]^2)
      pata[t]=epst[t]*sigt[t]
      }
   list(pata=pata,sigt=sigt,epst=epst)
}
> x = genedata(c(1, 0.3, 0.2),1000)$pata
> garch(x,trace=FALSE)

Call:
garch(x = x, trace = FALSE)

Coefficient(s):
       a0         a1         b1
4.362e+00  1.976e-01  6.805e-14
> garch(x,trace=FALSE,start=c(1,.3,.2))

Call:
garch(x = x, trace = FALSE, start = c(1, 0.3, 0.2))

Coefficient(s):
    a0      a1      b1
0.8025  0.2592  0.3255
> simgarch=genedata(c(1, 0.2, 0.7),1000)

Call:
garch(x = simgarch$pat, trace = FALSE)

Coefficient(s):
a0         a1         b1
8.044e+00  1.826e-01  4.051e-14
> garch(simgarch$pat,trace=FALSE,star=c(1, 0.2, 0.7))

Call:
garch(x = simgarch$pat, trace = FALSE, star = c(1, 0.2, 0.7))

Coefficient(s):
a0      a1      b1
1.1814  0.2079  0.6590

The above code clearly shows the huge impact of the starting value on the final estimate….

This entry was posted on May 17, 2012 at 12:12 am and is filed under R, Statistics, University life with tags GARCH, R, times series, tseries. You can follow any responses to this entry through the RSS 2.0 feed. You can leave a response, or trackback from your own site.

3 Responses to “garch() uncertainty”

Andres Trujillo-Barrera Says:
May 17, 2012 at 7:43 am

I’ve never used the package tseries, but have you considered using rugarch or fGarch, they seem way more robust for this kind of models

Reply
- xi'an Says:
  May 17, 2012 at 9:31 am
  
  Thanks. I just tried to install both of them but they both failed for failed dependencies…
  
  Reply
William Volterman Says:
May 17, 2012 at 12:54 am

It seems to be a convergence issue, I ran the code and looked at the likelihood. Pity that the function doesn’t return any information on convergence like most optimizers.

> genedata=function(para,nobs){
+
+ pata=epst=sigt=rnorm(nobs)
+ sigt[1]=sqrt(para[1])
+ pata[1]=epst[1]*sigt[1]
+ for (t in 2:nobs){
+ sigt[t]=sqrt(para[1]+para[2]*pata[t-1]^2+para[3]*sigt[t-1]^2)
+ pata[t]=epst[t]*sigt[t]
+ }
+ list(pata=pata,sigt=sigt,epst=epst)
+ }
>
> set.seed(1000)
> x y y$coef;y$n.likeli
a0 a1 b1
1.795077e+00 1.681971e-01 1.422666e-14
[1] 797.3634
>
> y y$coef;y$n.likeli
a0 a1 b1
1.10351778 0.33430225 0.08950214
[1] 780.3891
>
> y y$coef;y$n.likeli
a0 a1 b1
1.10349913 0.33430273 0.08950677
[1] 780.3891
>
> y y$coef;y$n.likeli
a0 a1 b1
2.492743e+00 1.026425e-01 1.242472e-15
[1] 838.3297

Reply

Xi'an's Og