# Contradictory results when estimating GJR-GARCH(1,1) with rugarch package

I am using financial stock data (1588 observations, daily returns) to estimate
a GARCH(1,1) and a GJR-GARCH(1,1) model.

A GARCH(1,1) takes the form: $$\sigma_t^2 = \omega+\alpha_1 a_{t-1}^2+\beta_1\sigma_{t-1}^2$$,

while a GJR-GARCH(1,1) takes the form:

$$\sigma_t^2 = \omega+(\alpha_1+\gamma_1N_{t-1})a_{t-1}^2+\beta_1\sigma_{t-1}^2$$, where $$N_{t-1}$$ is an indicator function with

$$N_{t-1} = \begin{cases} 1 & \text{if} \; a_{t-i} < 0, \\ 0 & \text{if} \; a_{t-i} \geq 0. \\ \end{cases}$$, and $$\alpha_i$$,$$\gamma_i$$ and $$\beta_j$$ are non-negative parameters.

I am using the rugarch package in R to estimate the models (I assume constant mean):

spec1 = ugarchspec(variance.model = list(model = "sGARCH", garchOrder = c(1,1)),
mean.model = list(armaOrder = c(0,0)))
garch = ugarchfit(spec = spec1, data = data)
spec2 = ugarchspec(variance.model = list(model = "fGARCH", garchOrder = c(1,1),
submodel = "GJRGARCH"),
mean.model = list(armaOrder = c(0,0)))
gjr = ugarchfit(spec = spec2, data = data)


For the GARCH(1,1) model, I get following output:

Optimal Parameters
------------------------------------
Estimate  Std. Error  t value Pr(>|t|)
mu      0.000309    0.000305   1.0135  0.31083
omega   0.000010    0.000001   7.7836  0.00000
alpha1  0.142502    0.013095  10.8818  0.00000
beta1   0.820859    0.016215  50.6237  0.00000


Everything seems to work fine. The covariance-stationarity condition is fulfilled as $$\alpha_1+\beta_1=0.9633609<1$$
Naturally, this value coincides with the line of code: garch@fit$persistence. For the GJR-GARCH(1,1), I get following output: Optimal Parameters ------------------------------------ Estimate Std. Error t value Pr(>|t|) mu 0.000065 0.000301 0.21742 0.827883 omega 0.000011 0.000001 14.45405 0.000000 alpha1 0.140940 0.013092 10.76573 0.000000 beta1 0.810593 0.015614 51.91367 0.000000 eta11 0.206881 0.052294 3.95610 0.000076  The covariance-stationarity condition is not fulfilled for this model as $$\alpha_1+\beta_1+0.5\gamma_1=1.0549735>1$$ (Note that this condition follows from the fact that I am assuming a symmetric distribution, the standard normal distribution). My first "conflicting" result and question is: How is it possible that the fitted GARCH(1,1) fulfills the covariance-stationarity condition, but the GJR-GARCH(1,1) does not using the exact same data? My second conflicting result is that if a run the line of code gjr@fit$persistence I get the value 0.9575659, which does not coincide with the "manually" calculated value of 1.0549735. The formula for the persistence of a GJR in the vignette of the rugarch package is exactly as I stated it above: https://cran.r-project.org/web/packages/rugarch/vignettes/Introduction_to_the_rugarch_package.pdf

Does someone have a clue, where the inconsistency could come from?

There are two different parametrizations of the GJR-GARCH model in rugarch, and you're applying the formula for the persistence from one parametrization to the other.

For GJR-GARCH(1,1), the first one is the one you've shown, which in the documentation is written like this:

$$\sigma_t^2 = w + (\alpha + \gamma I_{t-1}) \varepsilon_{t-1}^2 + \beta \sigma_{t-1}^2$$

where $$I_t = 1$$ when $$\varepsilon_t \leq 0$$, and 0 otherwise.

As you say, when the standardized residuals are drawn from a symmetric distribution, the persistence is given by:

$$\alpha + \beta + \frac{1}{2}\gamma$$

There is an alternate, equivalent parametrization, which arises as a special case of the family-GARCH model, which goes like this:

$$\sigma_t^2 = w + \alpha' \sigma^2_{t-1} (|z_{t-1}|-\eta z_{t-1})^2 + \beta \sigma_{t-1}^2$$

Since the sign of $$z_t$$ and $$\varepsilon_t$$ is the same and since $$\varepsilon_t = \sigma_t z_t$$, you can show that this is the same as:

$$\sigma_t^2 = w + \left(\alpha'(1-\eta)^2 + \alpha'\left[(-1-\eta)^2-(1-\eta)^2\right]I_{t-1} \right)\varepsilon^2_{t-1} + \beta \sigma_{t-1}^2$$

Reading off the coefficients, we have the relationship between the two parametrizations:

$$\alpha = \alpha'(1-\eta)^2$$ $$\gamma = \alpha'\left[(-1-\eta)^2-(1-\eta)^2\right]$$

Therefore, in the second parametrization the persistence is given by:

$$\alpha'(1-\eta)^2 + \beta + \frac{1}{2}\alpha'\left[(-1-\eta)^2-(1-\eta)^2\right]$$

The more general version of this appears in the documentation you've linked to as equations 27-28.

In rugarch, the first parametrization is obtained from rugarch::ugarchspec by passing :

variance.model = list(model = "gjrGARCH", garchOrder = c(1,1))


The second, which is the one you have, by passing:

variance.model = list(model = "fGARCH", garchOrder = c(1,1), submodel = "GJRGARCH")


Plugging in your coefficient estimates (to the precision that you've displayed them), I get that the persistence of your GJR-GARCH model is $$0.9575652$$, which is very close to what rugarch::persistence gave you.

• I was also wondering why I got the parameter "eta11" in the output of the GJR-GARCH(1,1), but your answer has clarified this point as well. Using the code that you stated for the first parametrization, I now get consistent results and get the parameter "gamma1" instead of "eta11" in the output of the GJR-GARCH(1,1). Thank you very much! Commented Nov 16, 2021 at 11:16