Time Series - HW 6

Lecture Note of Bus 41202, Spring 2012:

More Volatility Models. Mr. Ruey Tsay

The GARCH-M model

rt = μ + cσ 2 t + at, at = σt�t, σ

2 t = α0 + α1a

2 t−1 + β1σ

2 t−1

where c is referred to as risk premium, which is expected to be posi-

tive.

Example: A GARCH(1,1)-M model for the monthly excess returns

of S&P 500 index from January 1926 to December 1991. For numer-

ical stability, I use percentage returns. The fitted model is

rt = 0.743 + 0.048σ 2 t + at, σ

2 t = 0.812 + 0.123a

2 t−1 + .854σ

2 t−1.

Std err of risk premium is 0.141 so that the estimate is not statistically

significant at the usual 5% level. R demonstration > source("garchM.R")

> mn=garchM(sp5*100)

[1] 2399.822

0: 2380.0229: 0.422452 0.00561297 0.806149 0.121976 0.854361

3: 2379.9525: 0.427596 0.00652832 0.806107 0.121466 0.855756

6: 2378.4838: 0.606807 0.0375055 0.801889 0.126193 0.846833

9: 2377.9587: 0.673997 0.00166937 0.798191 0.125517 0.851915

12: 2377.8474: 0.692209 0.0444464 0.802574 0.122278 0.854141

15: 2377.7922: 0.742796 0.0480959 0.812187 0.122530 0.853507

Maximized log-likehood: 2377.792

Coefficient(s):

Estimate Std. Error t value Pr(>|t|)

mu 0.7427956 0.1540336 4.82229 1.4192e-06 ***

gamma 0.0480959 0.1408765 0.34140 0.7327991

omega 0.8121873 0.2858128 2.84168 0.0044877 **

alpha 0.1225297 0.0220596 5.55449 2.7843e-08 ***

beta 0.8535072 0.0219079 38.95885 < 2.22e-16 ***

---

1

Remarks: This R script is relatively slow. It takes longer time due

to its use of recursive loop in evaluating likelihood function.

The EGARCH model

Asymmetry in responses to past positive and negative returns:

g(�t) = θ�t + γ[|�t|− E(|�t|)],

with E[g(�t)] = 0.

To see asymmetry of g(�t), rewrite it as

g(�t) =

⎧⎪⎪⎨ ⎪⎪⎩ (θ + γ)�t − γE(|�t|) if �t ≥ 0, (θ − γ)�t − γE(|�t|) if �t < 0.

An EGARCH(m, s) model:

at = σt�t, ln(σ 2 t ) = α0 +

1 + β1B + · · · + βs−1Bs−1 1 − α1B −···− αmBm

g(�t−1).

Some features of EGARCH models:

• uses log trans. to relax the positiveness constraint • asymmetric responses

Consider an EGARCH(1,1) model

at = σt�t, (1 − αB) ln(σ2t ) = (1 − α)α0 + g(�t−1),

Under normality, E(|�t|) = √ 2/π and the model becomes

(1 − αB) ln(σ2t ) = ⎧⎪⎪⎨ ⎪⎪⎩ α∗ + (θ + γ)�t−1 if �t−1 ≥ 0, α∗ + (θ − γ)�t−1 if �t−1 < 0

2

where α∗ = (1 − α)α0 − √ 2 π γ.

This is a nonlinear fun. similar to that of the threshold AR model of

Tong (1978, 1990).

Specifically, we have

σ2t = σ 2α t−1 exp(α∗)

⎧⎪⎪⎪⎪⎪⎨ ⎪⎪⎪⎪⎪⎩

exp[(θ + γ) at−1√ σ2t−1

] if at−1 ≥ 0, exp[(θ − γ) at−1√

σ2t−1 ] if at−1 < 0.

The coefs (θ + γ) & (θ − γ) show the asymmetry in response to positive and negative at−1. The model is, therefore, nonlinear if

θ �= 0. Thus, θ is referred to as the leverage parameter. Focus on the function g(�t−1). The leverage parameter θ shows the

effect of the sign of at−1 whereas γ denotes the magnitude effect.

See Nelson (1991) for an exmaple of EGARCH model.

Another example: Monthly log returns of IBM stock from Jan-

uary 1926 to December 1997 for 864 observations.

For textbook, an AR(1)-EGARCH(1,1) is obtained (RATS program):

rt = 0.0105 + 0.092rt−1 + at, at = σt�t

ln(σ2t ) = −5.496 + g(�t−1)

1 − .856B, g(�t−1) = −.0795�t−1 + .2647[|�t−1|−

√ 2/π],

Model checking:

For ãt: Q(10) = 6.31(0.71) and Q(20) = 21.4(0.32)

For ã2t: Q(10) = 4.13(0.90) and Q(20) = 15.93(0.66)

3

Discussion:

Using √ 2/π ≈ 0.7979 ≈ 0.8, we obtain

ln(σ2t ) = −1.0 + 0.856 ln(σ2t−1) + ⎧⎪⎪⎨ ⎪⎪⎩ 0.1852�t−1 if �t−1 ≥ 0 −0.3442�t−1 if �t−1 < 0.

Taking anti-log transformation, we have

σ2t = σ 2×0.856 t−1 e

−1.001 × ⎧⎪⎪⎨ ⎪⎪⎩ e0.1852�t−1 if �t−1 ≥ 0 e−0.3442�t−1 if �t−1 < 0.

For a standardized shock with magnitude 2, (i.e. two standard devi-

ations), we have

σ2t (�t−1 = −2) σ2t (�t−1 = 2)

= exp[−0.3442 × (−2)]

exp(0.1852 × 2) = e 0.318 = 1.374.

Therefore, the impact of a negative shock of size two-standard de-

viations is about 37.4% higher than that of a positive shock of the

same size.

Forecasting: some recursive formula available

Another parameterization of EGARCH models

ln(σ2t ) = α0 + α1 |at−1| + γ1at−1

σt−1 + β1 ln(σ

2 t−1),

where γ1 denotes the leverage effect.

Below, I re-analyze the IBM log returns by extending the data to

December 2009. The sample size is 1008.

The fitted model is

rt = 0.012 + at, at = σt�t

4

ln(σ2t ) = −0.611 + 0.231|at−1|− 0.250at−1

σt−1 + 0.92 ln(σ2t−1).

Since EGARCH and TGARCH (below) share similar objective and

the latter is easier to estimate. We shall use TGARCH model.

The Threshold GARCH (TGARCH) or GJR Model A

TGARCH(s, m) or GJR(s, m) model is defined as

rt = μt + at, at = σt�tσ 2 t = α0 +

s∑ i=1

(αi + γiNt−i)a 2 t−i +

m∑ j=1

βjσ 2 t−j,

where Nt−i is an indicator variable such that

Nt−i =

⎧⎪⎪⎨ ⎪⎪⎩ 1 if at−i < 0,

0 otherwise.

One expects γi to be positive so that prior negative returns have

higher impact on the volatility.

The Asymmetric Power ARCH (APARCH) Model.

This model was introduced by Ding, Engle and Granger (1993) as a

general class of volatility models. The basic form is

rt = μt + at, at = σt�t, �t ∼ D(0, 1) σδt = ω +

s∑ i=1

αi(|at−i|− γiat−i)δ + m∑ j=1

βjσ δ t−j

where δ is a non-negative real number. In particular, δ = 2 gives rise

to the TGARCH model and δ = 0 corresponds to using log(σt).

Theoretically, one can use any power δ to obtain a model. In practice,

two things deserve further consideration. First, δ will also affect the

5

specification of the mean equation, i.e., model for μt. Second, it is

hard to interpret δ, except for some special values such as 0, 1, 2.

In R, one can fix the value of δ a priori using the subcommand

include.delta=F, delta = 2.

Here I pre-fix δ = 2. Thus, we can use APARCH model to estimate

TGARCH model. Consider the percentage log returns of monthly

IBM stock from 1926 to 2009.

R demonstration

> da=read.table("m-ibm2609.txt",header=T)

> head(da)

date ibm

1 19260130 -0.010381

.....

6 19260630 0.068493

> ibm=log(da$ibm+1)*100

> m1=garchFit(~aparch(1,1),data=ibm,trace=F,delta=2,include.delta=F)

> summary(m1)

Title:

GARCH Modelling

Call:

garchFit(formula = ~aparch(1, 1), data = ibm, delta = 2, include.delta = F,

trace = F)

Mean and Variance Equation:

data ~ aparch(1, 1)

[data = ibm]

Conditional Distribution: norm

Coefficient(s):

mu omega alpha1 gamma1 beta1

1.18659 4.33663 0.10767 0.22732 0.79468

Std. Errors: based on Hessian

Error Analysis:

Estimate Std. Error t value Pr(>|t|)

mu 1.18659 0.20019 5.927 3.08e-09 ***

omega 4.33663 1.34161 3.232 0.00123 **

6

alpha1 0.10767 0.02548 4.225 2.39e-05 ***

gamma1 0.22732 0.10018 2.269 0.02326 *

beta1 0.79468 0.04554 17.449 < 2e-16 ***

---

Log Likelihood:

-3329.177 normalized: -3.302755

Standardised Residuals Tests:

Statistic p-Value

Jarque-Bera Test R Chi^2 67.07416 2.775558e-15

Shapiro-Wilk Test R W 0.9870142 8.597234e-08

Ljung-Box Test R Q(10) 16.90603 0.07646942

Ljung-Box Test R Q(15) 24.19033 0.06193099

Ljung-Box Test R Q(20) 31.89097 0.04447407

Ljung-Box Test R^2 Q(10) 4.591691 0.9167342

Ljung-Box Test R^2 Q(15) 11.98464 0.6801912

Ljung-Box Test R^2 Q(20) 14.79531 0.7879979

LM Arch Test R TR^2 7.162971 0.8466584

Information Criterion Statistics:

AIC BIC SIC HQIC

6.615430 6.639814 6.615381 6.624694

> plot(m1) <= shows normal distribution is not a good fit.

>

> m1=garchFit(~aparch(1,1),data=ibm,trace=F,delta=2,include.delta=F,cond.dist="std")

> summary(m1)

Title:

GARCH Modelling

Call:

garchFit(formula = ~aparch(1, 1), data = ibm, delta = 2, cond.dist = "std",

include.delta = F, trace = F)

Mean and Variance Equation:

data ~ aparch(1, 1)

[data = ibm]

Conditional Distribution: std

Coefficient(s):

mu omega alpha1 gamma1 beta1 shape

1.20476 3.98975 0.10468 0.22366 0.80711 6.67329

Std. Errors: based on Hessian

Error Analysis:

Estimate Std. Error t value Pr(>|t|)

7

mu 1.20476 0.18715 6.437 1.22e-10 ***

omega 3.98975 1.45331 2.745 0.006046 **

alpha1 0.10468 0.02793 3.747 0.000179 ***

gamma1 0.22366 0.11595 1.929 0.053738 .

beta1 0.80711 0.04825 16.727 < 2e-16 ***

shape 6.67329 1.32779 5.026 5.01e-07 ***

---

Log Likelihood:

-3310.21 normalized: -3.283938

Standardised Residuals Tests:

Statistic p-Value

Jarque-Bera Test R Chi^2 67.82336 1.887379e-15

Shapiro-Wilk Test R W 0.9869698 8.212564e-08

Ljung-Box Test R Q(10) 16.91352 0.07629962

Ljung-Box Test R Q(15) 24.08691 0.06363224

Ljung-Box Test R Q(20) 31.75305 0.04600187

Ljung-Box Test R^2 Q(10) 4.553248 0.9189583

Ljung-Box Test R^2 Q(15) 11.66891 0.7038973

Ljung-Box Test R^2 Q(20) 14.18533 0.8209764

LM Arch Test R TR^2 6.771675 0.872326

Information Criterion Statistics:

AIC BIC SIC HQIC

6.579782 6.609042 6.579711 6.590898

> plot(m1)

Make a plot selection (or 0 to exit):

1: Time Series

2: Conditional SD

3: Series with 2 Conditional SD Superimposed

4: ACF of Observations

5: ACF of Squared Observations

6: Cross Correlation

7: Residuals

8: Conditional SDs

9: Standardized Residuals

10: ACF of Standardized Residuals

11: ACF of Squared Standardized Residuals

12: Cross Correlation between r^2 and r

13: QQ-Plot of Standardized Residuals

Selection: 13

For the percentage log returns of IBM stock from 1926 to 2009, the

8

fitted GJR model is

rt = 1.20 + at, at = σt�t, �t ∼ t∗6.67 σ2t = 3.99 + 0.105(|at−1|− 0.224at−1)2 + .807σ2t−1,

where all estimates are significant, and model checking indicates that

the fitted model is adequate.

Note that, we can obtain the model for the log returns as

rt = 0.012 + at, at = σt�t

σ2t = 3.99 × 10−4 + 0.105(|at−1|− 0.224at−1)2 + .807σ2t−1.

The sample variance of the IBM log returns is about 0.005 and the

empirical 2.5% percentile of the data is about −0.130. If we use these two quantities for σ2t−1 and at−1, respectively, then we have

σ2t (−) σ2t (+)

= 0.0004 + 0.105(0.130 + 0.224 × 0.130)2 + 0.807 × 0.005 0.0004 + 0.105(0.130 − 0.224 × 0.130)2 + 0.807 × 0.005

= 1.849.

In this particular case, the negative prior return has about 85% higher

impact on the conditional variance.

Stochastic volatility model

A (simple) SV model is

at = σt�t, (1 − α1B −···− αmBm) ln(σ2t ) = α0 + vt where �t’s are iid N(0, 1), vt’s are iid N(0, σ

2 v), {�t} and {vt} are

independent.

9

�3 �2 �1 0 1 2 3

� 2

0 2

4

qnorm � QQ Plot

Theoretical Quantiles

S a m

p le

Q u a n til

e s

Figure 1: Normal probability plot for TGARCH(1,1) model fitted to monthly percentage log returns of IBM stock from 1926 to 2009

10

�4 �2 0 2 4

� 2

0 2

4

qstd � QQ Plot

Theoretical Quantiles

S a m

p le

Q u a n til

e s

Figure 2: QQ plot for TGARCH(1,1) model fitted to monthly percentage log returns of IBM stock from 1926 to 2009.

11

Long-memory SV model

A simple LMSV is

at = σt�t, σt = σ exp(ut/2), (1 − B)dut = ηt where σ > 0, �t’s are iid N(0, 1), ηt’s are iid N(0, σ

2 η) and indepen-

dent of �t, and 0 < d < 0.5.

The model says

ln(a2t) = ln(σ 2) + ut + ln(�

2 t)

= [ln(σ2) + E(ln �2t)] + ut + [ln(� 2 t) − E(ln �2t)]

≡ μ + ut + et.

Thus, the ln(a2t) series is a Gaussian long-memory signal plus a non-

Gaussian white noise; see Breidt, Crato and de Lima (1998).

Application

see Examples 3.4 & 3.5

12