SECTION 1: LEARNING OBJECTIVES

By the end of this lesson, you will be able to:

  • Understand the concept of volatility clustering – the tendency of large price changes to be followed by large changes (of either sign) – and why it matters for risk.

  • Define the ARCH and GARCH models and their parameters.

  • Derive the GARCH(1,1) model and interpret its coefficients (persistence, mean reversion, and the unconditional variance).

  • Estimate GARCH models using Maximum Likelihood Estimation (MLE) in Python.

  • Use GARCH to forecast volatility over multiple horizons.

  • Combine GARCH volatility forecasts with VaR/ES to produce dynamic risk measures that adapt to changing market conditions.

  • Evaluate GARCH forecasts using loss functions (e.g., MSE, QLIKE) and compare to simpler models (e.g., Historical Volatility, EWMA).

  • Apply GARCH models to financial data for risk management, option pricing, and portfolio optimisation.

  • Understand extensions such as GARCH-M, EGARCH (for asymmetric volatility), and GJR-GARCH.


SECTION 2: VOLATILITY CLUSTERING – THE STYLIZED FACT

Financial returns exhibit volatility clustering: periods of high volatility are followed by high volatility, and periods of low volatility by low volatility. This violates the i.i.d. assumption of constant volatility.

Examples:

  • Stock market crashes (e.g., 1987, 2008) show extreme volatility followed by elevated volatility for months.

  • Calm periods (e.g., 2017) show persistently low volatility.

This persistence implies that volatility is time-varying and can be forecasted using past information.

The challenge: Volatility is latent (unobservable). We need models to capture its dynamics.


SECTION 3: ARCH AND GARCH MODELS – MATHEMATICAL FOUNDATION

3.1 The ARCH(q) Model (Engle, 1982)

The ARCH (AutoRegressive Conditional Heteroskedasticity) model assumes that the conditional variance depends on past squared errors.

Let rt be the return, with zero mean (or demeaned), and assume:

rt=σtεt,εt∼N(0,1)

where the conditional variance is:

σt2=α0+α1rt−12+α2rt−22+⋯+αqrt−q2

with α0>0αi≥0, and ∑i=1qαi<1 for stationarity.

Limitation: ARCH(q) often requires a large q to capture persistence.

3.2 The GARCH(p,q) Model (Bollerslev, 1986)

GARCH (Generalised ARCH) extends ARCH by allowing past conditional variances to enter the equation.

The GARCH(p,q) model:

σt2=ω+∑i=1qαirt−i2+∑j=1pβjσt−j2

  • ω>0

  • αi≥0 (ARCH coefficients)

  • βj≥0 (GARCH coefficients)

  • ∑i=1qαi+∑j=1pβj<1 (stationarity condition)

Interpretation:

  • The unconditional (long-run) variance is: σ2=ω1−∑αi−∑βj

  • The persistence is measured by α+β. If it is close to 1, volatility shocks decay slowly.

3.3 The GARCH(1,1) Model – Most Common

The workhorse model:

σt2=ω+αrt−12+βσt−12

  • ω : constant term

  • α : reaction to past squared returns (news impact)

  • β : persistence of volatility

  • Condition: α+β<1, with α,β,ω>0.

Forecasting:

  • 1-step ahead: σt+1∣t2=ω+αrt2+βσt2

  • h-step ahead (for h>1): the forecast converges to the unconditional variance:

σt+h∣t2=σ2+(α+β)h−1(σt+1∣t2−σ2)


SECTION 4: ESTIMATION – MAXIMUM LIKELIHOOD (MLE)

Given a sample of returns r1,r2,…,rT, the likelihood for GARCH(1,1) with normal errors is:

L=∏t=1T12πσt2exp⁡(−rt22σt2)

The log-likelihood (ignoring constants) is:

ℓ=−12∑t=1T(log⁡(σt2)+rt2σt2)

We maximise this with respect to ω,α,β. Because the model is recursive (volatility depends on past volatilities), we initialise σ12 with the unconditional variance or the sample variance.

Implementation: The arch package in Python provides a robust implementation.


SECTION 5: EXTENSIONS – EGARCH AND GJR-GARCH

EGARCH (Exponential GARCH) – Nelson (1991):
Models the log of variance, allowing for asymmetric effects (leverage effect): negative shocks increase volatility more than positive shocks of the same magnitude.

ln⁡(σt2)=ω+∑i=1qαi(∣rt−i∣σt−i−2π)+∑k=1pβkln⁡(σt−k2)+∑i=1qγirt−iσt−i

  • The term with γ captures the asymmetric effect.

GJR-GARCH (Glosten-Jagannathan-Runkle):
Adds a leverage term:

σt2=ω+αrt−12+βσt−12+γIt−1rt−12

where It−1=1 if rt−1<0, else 0.

Both models are widely used in finance to capture the leverage effect (negative returns increase future volatility).


SECTION 6: IMPLEMENTATION IN PYTHON – GARCH(1,1) WITH arch

We will fit a GARCH(1,1) model to the portfolio returns from Lesson 1 and use it for volatility forecasting and dynamic VaR.

python
# ===================================================================
# MODULE 5, LESSON 2: GARCH MODELS FOR VOLATILITY FORECASTING
# ===================================================================

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
from arch import arch_model
from scipy.stats import norm, t
import warnings
warnings.filterwarnings('ignore')

# Set style
sns.set_style("whitegrid")
np.random.seed(42)

print("="*70)
print("GARCH(1,1) VOLATILITY FORECASTING")
print("="*70)

# ----------------------------------------------------------------
# PART A: GENERATE DATA WITH VOLATILITY CLUSTERING (using GARCH process)
# ----------------------------------------------------------------

# Simulate a GARCH(1,1) process for demonstration
n = 1000
omega = 0.00001
alpha = 0.1
beta = 0.85

# Initialize
sigma2 = np.zeros(n)
r = np.zeros(n)
sigma2[0] = omega / (1 - alpha - beta)  # unconditional variance
r[0] = np.sqrt(sigma2[0]) * np.random.randn()

# Simulate
for t in range(1, n):
    sigma2[t] = omega + alpha * r[t-1]**2 + beta * sigma2[t-1]
    r[t] = np.sqrt(sigma2[t]) * np.random.randn()

# Convert to DataFrame
df_garch = pd.DataFrame({'return': r})
print("Simulated GARCH(1,1) data generated.")
print(f"True parameters: omega={omega}, alpha={alpha}, beta={beta}")

# Plot returns and volatility
fig, axes = plt.subplots(2, 1, figsize=(14, 8))

ax = axes[0]
ax.plot(df_garch['return'], linewidth=0.8)
ax.set_title('Simulated Returns with Volatility Clustering', fontsize=12)
ax.set_ylabel('Return')
ax.grid(True, alpha=0.3)

ax = axes[1]
ax.plot(np.sqrt(sigma2), color='red', linewidth=1.5)
ax.set_title('True Conditional Volatility', fontsize=12)
ax.set_ylabel('Volatility')
ax.grid(True, alpha=0.3)

plt.tight_layout()
plt.savefig('garch_simulation.png', dpi=300)
plt.show()

# ----------------------------------------------------------------
# PART B: ESTIMATE GARCH(1,1) MODEL
# ----------------------------------------------------------------

print("\n" + "-"*60)
print("PART B: GARCH(1,1) Estimation")
print("-"*60)

# Fit GARCH(1,1) with normal errors
model = arch_model(df_garch['return'], mean='zero', vol='GARCH', p=1, q=1)
res = model.fit(update_freq=5, disp='off')

print(res.summary())

# Extract parameters
omega_est = res.params['omega']
alpha_est = res.params['alpha[1]']
beta_est = res.params['beta[1]']

print(f"\nEstimated parameters:")
print(f"  omega = {omega_est:.6f}")
print(f"  alpha = {alpha_est:.6f}")
print(f"  beta  = {beta_est:.6f}")
print(f"  persistence (alpha+beta) = {alpha_est + beta_est:.6f}")

# Conditional volatility (fitted)
cond_vol = res.conditional_volatility
df_garch['fitted_vol'] = cond_vol

# Unconditional variance
uncond_var = omega_est / (1 - alpha_est - beta_est)
uncond_vol = np.sqrt(uncond_var)
print(f"Unconditional volatility (long-run): {uncond_vol:.6f}")

# ----------------------------------------------------------------
# PART C: FORECASTING VOLATILITY
# ----------------------------------------------------------------

print("\n" + "-"*60)
print("PART C: Volatility Forecasting")
print("-"*60)

# Forecast next 10 days
forecast = res.forecast(horizon=10)
forecast_variance = forecast.variance.iloc[-1].values  # last row, horizon

print("Volatility forecasts (next 10 days):")
for h in range(1, 11):
    vol_forecast = np.sqrt(forecast_variance[h-1])
    print(f"  Day {h}: {vol_forecast:.6f}")

# Plot forecast
fig, ax = plt.subplots(figsize=(12, 5))
ax.plot(np.sqrt(sigma2[-100:]), label='True Volatility', linewidth=2, color='blue')
ax.plot(cond_vol[-100:], 'r--', label='Estimated Volatility', linewidth=2, alpha=0.7)
# Add forecast point
ax.axvline(x=len(sigma2)-1, color='black', linestyle=':', alpha=0.7)
ax.plot(range(len(sigma2)-1, len(sigma2)+10-1), np.concatenate([[cond_vol[-1]], np.sqrt(forecast_variance)]), 'go-', label='Forecast', markersize=8)
ax.set_title('Volatility Forecast – GARCH(1,1)', fontsize=12)
ax.set_xlabel('Time')
ax.set_ylabel('Volatility')
ax.legend()
ax.grid(True, alpha=0.3)
plt.tight_layout()
plt.savefig('garch_forecast.png', dpi=300)
plt.show()

# ----------------------------------------------------------------
# PART D: DYNAMIC VaR USING GARCH FORECASTS
# ----------------------------------------------------------------

print("\n" + "-"*60)
print("PART D: Dynamic VaR with GARCH")
print("-"*60)

# For each day, compute 1-day VaR using estimated conditional volatility
confidence = 0.95
alpha_level = 1 - confidence
z_alpha = norm.ppf(alpha_level)

# We assume mean return is zero (or use estimated mean if any)
# VaR = - (mu + z_alpha * sigma_t)  (loss positive)
df_garch['VaR_95'] = -(0 + z_alpha * cond_vol)  # positive loss

# Compare actual returns to VaR
exceptions = (df_garch['return'] < -df_garch['VaR_95']).sum()
print(f"Number of exceptions (returns < -VaR): {exceptions} out of {len(df_garch)}")
print(f"Expected exceptions: {alpha_level * len(df_garch):.0f}")

# Plot returns and VaR
fig, ax = plt.subplots(figsize=(14, 6))
ax.plot(df_garch['return'], label='Returns', linewidth=0.8, alpha=0.7)
ax.plot(df_garch['VaR_95'], color='red', linewidth=1.5, label='95% VaR (GARCH)')
ax.fill_between(df_garch.index, -df_garch['VaR_95'], df_garch['VaR_95'], alpha=0.1, color='red')
ax.set_title('Dynamic VaR from GARCH(1,1)', fontsize=12)
ax.set_ylabel('Return / VaR')
ax.legend()
ax.grid(True, alpha=0.3)
plt.tight_layout()
plt.savefig('dynamic_var_garch.png', dpi=300)
plt.show()

# ----------------------------------------------------------------
# PART E: COMPARISON WITH HISTORICAL VOLATILITY (EWMA)
# ----------------------------------------------------------------

print("\n" + "-"*60)
print("PART E: Comparison with EWMA Volatility")
print("-"*60)

# EWMA (Exponentially Weighted Moving Average) – lambda = 0.94 (RiskMetrics)
lambda_ewma = 0.94
sigma2_ewma = np.zeros(n)
sigma2_ewma[0] = df_garch['return'].var()
for t in range(1, n):
    sigma2_ewma[t] = lambda_ewma * sigma2_ewma[t-1] + (1 - lambda_ewma) * df_garch['return'].iloc[t-1]**2
vol_ewma = np.sqrt(sigma2_ewma)

# Plot comparison
fig, ax = plt.subplots(figsize=(12, 5))
ax.plot(np.sqrt(sigma2), label='True Volatility', linewidth=1.5, color='black')
ax.plot(cond_vol, label='GARCH(1,1)', linewidth=1.5, color='red')
ax.plot(vol_ewma, label='EWMA (λ=0.94)', linewidth=1.5, color='blue', alpha=0.7)
ax.set_title('Volatility Comparison: GARCH vs EWMA', fontsize=12)
ax.set_xlabel('Time')
ax.set_ylabel('Volatility')
ax.legend()
ax.grid(True, alpha=0.3)
plt.tight_layout()
plt.savefig('volatility_comparison.png', dpi=300)
plt.show()

# ----------------------------------------------------------------
# PART F: APPLYING GARCH TO REAL DATA (SIMULATED PORTFOLIO)
# ----------------------------------------------------------------

print("\n" + "-"*60)
print("PART F: GARCH on Realistic Portfolio Returns (from Lesson 1)")
print("-"*60)

# Use the portfolio returns from Lesson 1
# We already have df_returns['Portfolio'] as port_returns

port_series = pd.Series(port_returns, name='return')
model_real = arch_model(port_series, mean='zero', vol='GARCH', p=1, q=1)
res_real = model_real.fit(update_freq=5, disp='off')

print("GARCH(1,1) on Real Portfolio Returns:")
print(res_real.summary())

# Extract and interpret
omega_r = res_real.params['omega']
alpha_r = res_real.params['alpha[1]']
beta_r = res_real.params['beta[1]']
persistence = alpha_r + beta_r
uncond_var_r = omega_r / (1 - persistence)

print(f"\nPersistence (alpha+beta) = {persistence:.4f}")
print(f"Unconditional volatility = {np.sqrt(uncond_var_r):.4f}")

# Forecast 1-day VaR using GARCH for the next day
last_vol = res_real.conditional_volatility[-1]
forecast_var = res_real.forecast(horizon=1).variance.iloc[-1, 0]
next_vol = np.sqrt(forecast_var)
var_95_garch = - (0 + z_alpha * next_vol)
print(f"Next-day 95% VaR (GARCH) = {var_95_garch*100:.4f}%")

# Compare with historical VaR from Lesson 1
hist_var_95 = hist_results[0.95]['VaR']
print(f"Historical 95% VaR (from previous lesson) = {hist_var_95*100:.4f}%")
print(f"Difference: GARCH VaR is {'higher' if var_95_garch > hist_var_95 else 'lower'} by {abs(var_95_garch - hist_var_95)*100:.4f}%")

# ----------------------------------------------------------------
# PART G: MODEL DIAGNOSTICS – STANDARDISED RESIDUALS
# ----------------------------------------------------------------

print("\n" + "-"*60)
print("PART G: Diagnostics – Standardised Residuals")
print("-"*60)

# Standardised residuals = return / conditional volatility
std_resid = res_real.resid / res_real.conditional_volatility

# Ljung-Box test for autocorrelation of residuals and squared residuals
from statsmodels.stats.diagnostic import acorr_ljungbox

lb_resid = acorr_ljungbox(std_resid, lags=10, return_df=True)
lb_sq_resid = acorr_ljungbox(std_resid**2, lags=10, return_df=True)

print("Ljung-Box Test on Standardised Residuals (p-values):")
print(lb_resid['lb_pvalue'].tolist()[:5])
print("Ljung-Box Test on Squared Standardised Residuals (p-values):")
print(lb_sq_resid['lb_pvalue'].tolist()[:5])

print("\nInterpretation:")
print("  • If p-values are > 0.05, we fail to reject the null of no autocorrelation.")
print("  • This indicates the GARCH model has adequately captured volatility dynamics.")

# QQ plot of standardised residuals
from scipy import stats
fig, ax = plt.subplots(figsize=(8, 6))
stats.probplot(std_resid, dist="norm", plot=ax)
ax.set_title('Q-Q Plot of Standardised Residuals', fontsize=12)
plt.tight_layout()
plt.savefig('garch_qqplot.png', dpi=300)
plt.show()

SECTION 7: BUSINESS APPLICATIONS IN BANKING

 
 
Application Role of GARCH Benefit
Market Risk (VaR/ES) Dynamic volatility forecasts for VaR/ES. More accurate risk capital than constant volatility.
Option Pricing Volatility input for Black-Scholes (or local volatility models). Better pricing and hedging.
Risk Management Volatility forecasts for stress testing and scenario analysis. Improved risk assessment.
Portfolio Optimisation Dynamic covariances for mean-variance optimisation. Enhanced asset allocation.
Hedge Ratio Calculation Time-varying volatilities for optimal hedging. More effective hedging strategies.

SECTION 8: REGULATORY CONSIDERATIONS

  • Basel III allows banks to use internal models for market risk, which require accurate volatility forecasts. GARCH is a common choice.

  • Model validation (SR 11-7) requires that GARCH models are backtested, and their forecasts are compared with realised volatility.

  • IFRS 9 / CECL may use GARCH to forecast volatility for ECL calculations (though typically simpler models are used).

  • FRTB requires that internal models (including GARCH) be subject to rigorous backtesting at the trading desk level.


SECTION 9: SUMMARY FOR THE DATA PRACTITIONER

  • Volatility is time-varying and predictable – GARCH models capture this persistence.

  • GARCH(1,1) is the industry standard, with parameters ω,α,β.

  • Persistence ( α+β ) near 1 indicates slow mean reversion.

  • Volatility forecasts from GARCH can be used to compute dynamic VaR/ES.

  • Extensions (EGARCH, GJR-GARCH) capture asymmetry (leverage effect).

  • Model validation involves checking standardised residuals for autocorrelation and normality.

  • In practice, GARCH is widely used for risk management, option pricing, and portfolio construction.


SECTION 10: RECOMMENDED NEXT STEPS

  1. Fit GARCH models to actual asset returns (e.g., S&P 500, FX rates).

  2. Compare GARCH(1,1) with other volatility models (EWMA, HAR-RV).

  3. Implement rolling-window forecasts and backtest VaR.

  4. Explore multivariate GARCH (DCC-GARCH) for portfolio volatility.

  5. Study Extreme Value Theory (EVT) for tail risk modelling.

  6. Prepare for the next lesson on Model Validation and Stress Testing.


[END OF LESSON 2 – MODULE 5]