Although ARIMA, ADL, and ARCH/GARCH models appear different,
they all belong to the same family of linear dynamic models —
each describes how a current value depends on past information.
The ARIMA(p,d,q) model captures persistence in the mean of a time series.
\[ y_t = c + \phi_1 y_{t-1} + \dots + \phi_p y_{t-p} + \theta_1 \varepsilon_{t-1} + \dots + \theta_q \varepsilon_{t-q} + \varepsilon_t \]
The Autoregressive Distributed Lag (ADL) model extends ARIMA to include other predictors.
\[ y_t = \beta_0 + \sum_{i=1}^{p}\beta_i y_{t-i} + \sum_{j=0}^{q}\gamma_j x_{t-j} + \varepsilon_t \]
ADL is basically an ARIMA on \(y_t\), but with exogenous regressors \(x_t\).
The ARCH(p) and GARCH(p,q) models apply the same dynamic logic to the variance of errors, not the mean.
\[ \sigma_t^2 = \alpha_0 + \sum_{i=1}^{p}\alpha_i u_{t-i}^2 + \sum_{j=1}^{q}\phi_j \sigma_{t-j}^2 \]
In spirit, a GARCH model is an ADL applied to the error variance.
| Model | Dependent variable | What it models | Dynamic structure |
|---|---|---|---|
| ARIMA | \(y_t\) | Mean of a single series | Past \(y_t\)’s and shocks |
| ADL | \(y_t\) | Mean with external drivers | Past \(y_t\)’s and past \(x_t\)’s |
| ARCH/GARCH | \(\sigma_t^2\) | Conditional variance | Past \(u_t^2\)’s and \(\sigma_t^2\)’s |
All three share the same recursive idea:
> “Today depends on yesterday.”
Only the object of interest (mean vs variance) changes.
| Level | Model | Focus |
|---|---|---|
| 1 | ADL / ARIMA | Predicting the expected level of a series |
| 2 | ARCH / GARCH | Predicting the expected variability of that series |
So, GARCH ≈ ARIMA logic applied to volatility,
and ADL ≈ ARIMA logic extended to multiple variables.
ARIMA → dynamics of the mean
ADL → dynamics between variables
GARCH → dynamics of the variance
All are linear-in-parameters models that explain persistence —
in levels (ARIMA), across variables (ADL), or in volatility (GARCH).
In ordinary regression or ARIMA models, we assume that the variance of the errors is constant:
\[u_t \sim N (0, \sigma^2)\]
That’s called homoskedasticity. But in real-world financial or economic data (returns, inflation, exchange rates, etc.), you often see volatility clustering — periods of calm and periods of turbulence.
Example:
So the variance (volatility) is changing over time. That’s called heteroskedasticity.
Engle (1982) introduced the ARCH model (Autoregressive Conditional Heteroskedasticity).
Idea: Let’s model the variance itself as something that depends on past squared shocks.
\[\sigma_t^2 = \alpha_0 + \alpha_1 u_{t-1}^2 + \alpha_2 u_{t-2}^2 + \dots + \alpha_p u_{t-p}^2\]
Here:
So volatility today depends on how big shocks were yesterday, the day before, etc.
In the ARCH model, we use squared shocks (\(u_{t-i}^2\)) instead of the raw shocks (\(u_{t-i}\)).
There are two main reasons:
Variance must be positive
Variance (\(\sigma_t^2\)) represents volatility — it cannot be negative.
Squaring ensures that each term contributes positively: \[
\sigma_t^2 = \alpha_0 + \alpha_1 u_{t-1}^2 + \alpha_2 u_{t-2}^2 + \dots + \alpha_p u_{t-p}^2,
\quad \text{where } \alpha_0 > 0, \; \alpha_i \ge 0.
\]
We care about the magnitude, not the direction, of shocks
Both large positive and large negative shocks increase volatility.
Squaring \(u_t\) removes the sign and keeps only the size of the movement.
# Simulate 100 random shocks
set.seed(123)
u_t <- rnorm(100, mean = 0, sd = 1)
#Compute squared shocks
u_t2 <- u_t^2
#Plot
par(mfrow = c(2,1), mar = c(3,4,2,1))
plot(u_t,
type = "h",
col = "steelblue",
main = "Raw shocks (u_t)",
ylab = "u_t")
plot(u_t2,
type = "h",
col = "firebrick",
main = "Squared shocks (u_t^2)",
ylab = "u_t^2")
If volatility persists for a long time (as in financial data), ARCH needs many lags (\(p\)) to capture that persistence.
That means too many parameters → not efficient.
GARCH solves that by adding lagged variances themselves:
\[\sigma_t^2 = \alpha_0 + \sum_{i=1}^{p}\alpha_i u_{t-i}^2 + \sum_{j=1}^{q}\phi_j \sigma_{t-j}^2\]
This means volatility depends both on:
The GARCH(1,1) model is written as:
\[ R_t = \beta_0 + u_t, \quad u_t \sim N(0, \sigma_t^2) \]
\[ \sigma_t^2 = \alpha_0 + \alpha_1 u_{t-1}^2 + \phi_1 \sigma_{t-1}^2 \]
where:
| Symbol | Meaning | Expected Sign |
|---|---|---|
| \(\beta_0\) | Constant (mean of returns) | Any |
| \(\alpha_0\) | Long-run average variance | \(> 0\) |
| \(\alpha_1\) | Reaction to last period’s shock (ARCH effect) | \(\ge 0\) |
| \(\phi_1\) | Volatility persistence (GARCH effect) | \(\ge 0\) |
Describes the average return level — how returns behave on average over time.
Captures the time-varying volatility — how uncertainty evolves over time.
This discussion focuses on modeling volatility clustering with ARCH and GARCH. These topics equip you with advanced tools for working with structured data and financial or economic series that exhibit time-varying variance.
Required Readings:
Exploring Volatility Clustering:
Read Chapter 16.4 of Econometrics with R for an introduction to volatility clustering and the motivation behind ARCH and GARCH models.
Select or simulate a time series where you suspect time-varying volatility, such as financial returns, commodity prices, or any series with fluctuating variance.
I select MSFT monthly returns and fit a GARCH model from fGarch library.
Registered S3 method overwritten by 'quantmod':
method from
as.zoo.data.frame zoo ── Attaching core tidyquant packages ─────────────────────── tidyquant 1.0.11 ──
✔ PerformanceAnalytics 2.0.8 ✔ TTR 0.24.4
✔ quantmod 0.4.28 ✔ xts 0.14.1
── Conflicts ────────────────────────────────────────── tidyquant_conflicts() ──
✖ zoo::as.Date() masks base::as.Date()
✖ zoo::as.Date.numeric() masks base::as.Date.numeric()
✖ PerformanceAnalytics::legend() masks graphics::legend()
✖ quantmod::summary() masks base::summary()
ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
✔ dplyr 1.1.4 ✔ readr 2.1.5
✔ forcats 1.0.0 ✔ stringr 1.5.1
✔ ggplot2 4.0.0 ✔ tibble 3.3.0
✔ lubridate 1.9.4 ✔ tidyr 1.3.1
✔ purrr 1.1.0
── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
✖ dplyr::filter() masks stats::filter()
✖ dplyr::first() masks xts::first()
✖ dplyr::lag() masks stats::lag()
✖ dplyr::last() masks xts::last()
ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errorsconflicted::conflicts_prefer(dplyr::lag)[conflicted] Will prefer dplyr::lag over any other package.conflicted::conflicts_prefer(dplyr::last)[conflicted] Will prefer dplyr::last over any other package.conflicted::conflicts_prefer(dplyr::filter)[conflicted] Will prefer dplyr::filter over any other package.?tq_get
# Get Netflix stock prices
df_ts <- tq_get("MSFT",
get = "stock.prices",
from = "1990-01-01"
)
str(df_ts)tibble [9,031 × 8] (S3: tbl_df/tbl/data.frame)
$ symbol : chr [1:9031] "MSFT" "MSFT" "MSFT" "MSFT" ...
$ date : Date[1:9031], format: "1990-01-02" "1990-01-03" ...
$ open : num [1:9031] 0.606 0.622 0.62 0.635 0.622 ...
$ high : num [1:9031] 0.616 0.627 0.639 0.639 0.632 ...
$ low : num [1:9031] 0.598 0.615 0.616 0.622 0.615 ...
$ close : num [1:9031] 0.616 0.62 0.638 0.622 0.632 ...
$ volume : num [1:9031] 5.30e+07 1.14e+08 1.26e+08 6.96e+07 5.90e+07 ...
$ adjusted: num [1:9031] 0.377 0.379 0.39 0.381 0.387 ...# Calculate daily percentage change
df_ts <- df_ts %>%
mutate(
pct_change = (close - lag(close)) / lag(close) * 100
)
# Preview
head(df_ts, 10)# A tibble: 10 × 9
symbol date open high low close volume adjusted pct_change
<chr> <date> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
1 MSFT 1990-01-02 0.606 0.616 0.598 0.616 53035200 0.377 NA
2 MSFT 1990-01-03 0.622 0.627 0.615 0.620 113774400 0.379 0.564
3 MSFT 1990-01-04 0.620 0.639 0.616 0.638 125740800 0.390 2.94
4 MSFT 1990-01-05 0.635 0.639 0.622 0.622 69566400 0.381 -2.45
5 MSFT 1990-01-08 0.622 0.632 0.615 0.632 58982400 0.387 1.53
6 MSFT 1990-01-09 0.632 0.639 0.627 0.630 70300800 0.386 -0.275
7 MSFT 1990-01-10 0.625 0.634 0.612 0.613 103766400 0.375 -2.75
8 MSFT 1990-01-11 0.616 0.622 0.590 0.601 95774400 0.368 -1.98
9 MSFT 1990-01-12 0.591 0.606 0.583 0.598 148910400 0.366 -0.433
10 MSFT 1990-01-15 0.595 0.604 0.592 0.598 62467200 0.366 0 # Plot closing prices over time
ggplot(data = df_ts,
mapping = aes(x = date,
y = pct_change)
) +
geom_line(color = "blue") +
labs(
title = "MSFT Daily Price Pct Change Over Time",
x = "Date",
y = "Closing Price (USD)"
) +
theme_minimal()Warning: Removed 1 row containing missing values or values outside the scale range
(`geom_line()`).
# Convert to monthly and calculate % change
df_monthly <- df_ts %>%
group_by(year_month = floor_date(x = date, unit = "month")) %>%
summarise(close = last(close), .groups = "drop") %>%
mutate(
pct_change = (close - lag(close)) / lag(close) * 100
) %>% filter(!is.na(pct_change))
# Plot closing prices over time
ggplot(df_monthly, aes(x = year_month, y = pct_change)) +
geom_line(color = "blue") +
labs(
title = "MSFT Monthly Price Pct Change Over Time",
x = "Date",
y = "Closing Price (USD)"
)
# plot sample autocorrelation of daily percentage changes
acf(df_monthly$pct_change,
main = "ACF"
)
We see that autocorrelations are weak - so that it is difficult to predict future outcomes using, e.g., an AR model.
The tall bar at lag 0 = 1 is expected — every ACF plot has it.
The other bars (lags 1–25) hover near zero, meaning autocorrelation is weak. → The series behaves like white noise — difficult to predict using an AR model.
Modeling Conditional Variance:
Fit an ARCH or GARCH model to your series using R (for example, the fGarch or rugarch package).
Plot the estimated conditional variance and compare it to the unconditional variance from simpler models.
Interpret the significance of the parameters and how volatility evolves over time.
fGarch
Univariate or multivariate GARCH time series fitting. Estimates the parameters of a univariate ARMA-GARCH/APARCH process.
NOTE: Packages 'fBasics', 'timeDate', and 'timeSeries' are no longer
attached to the search() path when 'fGarch' is attached.
If needed attach them yourself in your R script by e.g.,
require("timeSeries")??garchFit
# estimate GARCH(1,1) model of daily percentage changes
GARCH_MSFT <- garchFit(data = df_monthly$pct_change,
trace = F # optimization process of fitting the model parameters should be printed ?
)
summary(GARCH_MSFT)
Title:
GARCH Modelling
Call:
garchFit(data = df_monthly$pct_change, trace = F)
Mean and Variance Equation:
data ~ garch(1, 1)
<environment: 0x1682c5008>
[data = df_monthly$pct_change]
Conditional Distribution:
norm
Coefficient(s):
mu omega alpha1 beta1
1.944520 2.867759 0.091883 0.864059
Std. Errors:
based on Hessian
Error Analysis:
Estimate Std. Error t value Pr(>|t|)
mu 1.94452 0.35302 5.508 3.62e-08 ***
omega 2.86776 1.24339 2.306 0.021088 *
alpha1 0.09188 0.02693 3.412 0.000644 ***
beta1 0.86406 0.03596 24.028 < 2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Log Likelihood:
-1497.807 normalized: -3.483272
Description:
Sat Nov 8 22:40:35 2025 by user:
Standardised Residuals Tests:
Statistic p-Value
Jarque-Bera Test R Chi^2 33.4237353 5.522403e-08
Shapiro-Wilk Test R W 0.9833088 7.367492e-05
Ljung-Box Test R Q(10) 11.3320110 3.322430e-01
Ljung-Box Test R Q(15) 14.0513585 5.216376e-01
Ljung-Box Test R Q(20) 16.7053539 6.720150e-01
Ljung-Box Test R^2 Q(10) 8.2238795 6.069789e-01
Ljung-Box Test R^2 Q(15) 10.5940497 7.808050e-01
Ljung-Box Test R^2 Q(20) 16.5345516 6.829572e-01
LM Arch Test R TR^2 11.5729227 4.805587e-01
Information Criterion Statistics:
AIC BIC SIC HQIC
6.985149 7.022952 6.984978 7.000076 After running the GARCH model, you need to interpret two linked equations — the mean returns equation and the conditional heteroskedasticity equation (the two-equation system). Part 2, the heteroskedasticity component, can become more complex depending on the chosen GARCH(p, q) specification.
First, I write both equations using Greek symbols, and then replace the Greeks with the estimated coefficients from the model output. Be aware that different software packages use different notation, so it is important to interpret the parameters carefully.
I have adjusted some coefficient symbols from the textbook to match the notation used by the garchFit() command in R.
\[ R_t = 1.97218 + \mu_t, \quad \mu_t \sim N(0, \sigma_t^2) \]
Where:
garcgFit, but called \(\beta_0\) in textbook system of equations 16.3 (see screenshot above).\[ \sigma_t^2 = 2.91377 + 0.09266 \epsilon_{t-1}^2 + 0.86307 \sigma_{t-1}^2 \]
Where:
\(\omega\) = 2.91377: Long-run average variance (baseline risk level).
\(\alpha_1\) = 0.09266: Impact of the most recent shock. A large \(\mu_{t-1}^2\) will increase \(\sigma_t^2\) in the next period.
\(\beta_1\) = 0.86307: Persistence of volatility from the previous period. A high value indicates slow decay of volatility after shocks.
This formulation allows the volatility to change over time, capturing periods of high and low market turbulence.
The persistence parameter:
\[ \alpha_1 + \beta_1 = 0.95573 \]
is close to 1, indicating high volatility persistence — shocks to volatility tend to last for a long time.
# compute deviations of the percentage changes from their mean
dev_mean_MSFT <- df_monthly$pct_change - GARCH_MSFT@fit$coef[1]
# plot deviation of percentage changes from mean
plot(x = df_monthly$year_month, y = dev_mean_MSFT,
type = "l",
col = "steelblue",
ylab = "Percent",
xlab = "Date",
main = "Estimated Bands of +- One Conditional Standard Deviation",
cex.main=0.8,
lwd = 0.2)
# add horizontal line at y = 0
abline(0, 0)
# add GARCH(1,1) confidence bands (one standard deviation) to the plot
lines(x = df_monthly$year_month,
GARCH_MSFT@fit$coef[1] + GARCH_MSFT@sigma.t,
col = "darkred",
lwd = 0.5)
lines(x = df_monthly$year_month,
GARCH_MSFT@fit$coef[1] - GARCH_MSFT@sigma.t, ,
col = "darkred",
lwd = 0.5)
The light blue line appears to be the actual returns or residuals from a financial time series (likely % deviations from a mean).
The two red lines represent the estimated bands of ± 1 conditional standard deviation from a volatility model (likely a GARCH-type model).
Conditional variance changes over time based on recent information. In GARCH models, this means periods of high volatility widen the bands, while quiet periods narrow them.
Unconditional variance assumes volatility is constant over time — from a simple model like historical variance. This would be a flat horizontal band (same width everywhere).
Comparison:
If you plotted the unconditional variance, it would likely underestimate volatility during high-stress periods (e.g., late 1990s, 2008) and overestimate during calm periods.
The GARCH conditional variance adapts to spikes in volatility and periods of calm, making it better for risk estimation.
Discussion Questions:
In your dataset, how did the ARCH/GARCH model improve the understanding of volatility compared to traditional constant-variance models?
How would you incorporate volatility modeling into forecasting workflows alongside ETS or ARIMA models?
How ARCH/GARCH Improved Understanding of Volatility
The GARCH(1,1) model:
Captured volatility clustering
The plot shows that conditional standard deviation bands widen during turbulent market periods (e.g., dot-com bubble, 2008 financial crisis, COVID-19 shock) and narrow during calm periods.
A constant-variance model would have used one average variance level for all periods, missing this dynamic.
Quantified persistence of volatility
α₁ (0.09266) → About 9% of last period’s shock size carries into today’s volatility.
β₁ (0.86307) → About 86% of yesterday’s volatility persists today.
Persistence α₁ + β₁ ≈ 0.956 → Shocks have long-lasting effects.
Constant-variance models cannot measure or use this persistence in risk estimates.
Improved risk measurement
Prediction intervals adapt to market conditions — wide when volatility is high, narrow when low — unlike ARIMA/ETS with fixed-width intervals.
This improves Value-at-Risk (VaR) and stress-testing accuracy.
Better model diagnostics
Ljung–Box and ARCH–LM tests show no remaining autocorrelation in squared residuals → volatility structure is well captured.
Constant-variance models leave significant autocorrelation in squared residuals, meaning they systematically miss volatility patterns.
Incorporating Volatility Modeling with ETS/ARIMA
ARIMA–GARCH Hybrid (Finance Standard)
Step 1: Fit ARIMA to model the conditional mean of returns.
Step 2: Fit GARCH to ARIMA residuals to model conditional variance.
Step 3: Forecast mean & volatility jointly.
ARIMA forecast → central path.
GARCH forecast → time-varying prediction intervals or VAR.
rugarch
Another popular way to implement the GARCH.
The rugarch package aims to provide a flexible and rich univariate GARCH modelling and testing environment. Modelling is a simple process of defining a specification and fitting the data.
Loading required package: parallel?rugarch
??ugarchspec
?fGarch::garchSpec
# Specify and fit a GARCH(1,1) model
spec <- ugarchspec(
variance.model = list(model = "sGARCH", garchOrder = c(1, 1)),
mean.model = list(armaOrder = c(0, 0), include.mean = TRUE),
distribution.model = "norm" # standard normal distribution
)
fit <- ugarchfit(spec = spec,
data = df_monthly$pct_change
)
show(fit)
*---------------------------------*
* GARCH Model Fit *
*---------------------------------*
Conditional Variance Dynamics
-----------------------------------
GARCH Model : sGARCH(1,1)
Mean Model : ARFIMA(0,0,0)
Distribution : norm
Optimal Parameters
------------------------------------
Estimate Std. Error t value Pr(>|t|)
mu 1.944303 0.353009 5.5078 0.000000
omega 2.870657 1.246834 2.3024 0.021315
alpha1 0.091512 0.026836 3.4101 0.000649
beta1 0.864090 0.036068 23.9574 0.000000
Robust Standard Errors:
Estimate Std. Error t value Pr(>|t|)
mu 1.944303 0.333575 5.8287 0.000000
omega 2.870657 1.297646 2.2122 0.026953
alpha1 0.091512 0.024405 3.7498 0.000177
beta1 0.864090 0.031507 27.4255 0.000000
LogLikelihood : -1497.793
Information Criteria
------------------------------------
Akaike 6.9851
Bayes 7.0229
Shibata 6.9849
Hannan-Quinn 7.0000
Weighted Ljung-Box Test on Standardized Residuals
------------------------------------
statistic p-value
Lag[1] 2.460 0.1168
Lag[2*(p+q)+(p+q)-1][2] 2.475 0.1949
Lag[4*(p+q)+(p+q)-1][5] 3.923 0.2638
d.o.f=0
H0 : No serial correlation
Weighted Ljung-Box Test on Standardized Squared Residuals
------------------------------------
statistic p-value
Lag[1] 0.8411 0.3591
Lag[2*(p+q)+(p+q)-1][5] 3.6034 0.3078
Lag[4*(p+q)+(p+q)-1][9] 5.0811 0.4176
d.o.f=2
Weighted ARCH LM Tests
------------------------------------
Statistic Shape Scale P-Value
ARCH Lag[3] 4.031 0.500 2.000 0.04466
ARCH Lag[5] 4.834 1.440 1.667 0.11226
ARCH Lag[7] 5.239 2.315 1.543 0.20144
Nyblom stability test
------------------------------------
Joint Statistic: 0.821
Individual Statistics:
mu 0.2670
omega 0.2536
alpha1 0.4205
beta1 0.3795
Asymptotic Critical Values (10% 5% 1%)
Joint Statistic: 1.07 1.24 1.6
Individual Statistic: 0.35 0.47 0.75
Sign Bias Test
------------------------------------
t-value prob sig
Sign Bias 0.9214 0.3574
Negative Sign Bias 0.9660 0.3346
Positive Sign Bias 0.9571 0.3391
Joint Effect 1.8637 0.6012
Adjusted Pearson Goodness-of-Fit Test:
------------------------------------
group statistic p-value(g-1)
1 20 19.95 0.3974
2 30 25.72 0.6404
3 40 40.51 0.4035
4 50 57.67 0.1852
Elapsed time : 0.1495972 Note sGARCH is standard GARCH, not seasonal GARCH.
There are lots of GARCH models - EGARCH, APARCH, TGARCH, QGARCH, GJRGARCH, IGARCH, AVGARCH, PGARCH, GARCH-M and the standard GARCH.
You can add refinement/structure to the mean model and variance model by changing the parameters.
Try including an autoregressive (AR) term in the mean equation and observe how the estimated mean equation changes.
Because the residuals (\(u_t\)) are derived from the mean model, this change will also affect the variance equation indirectly.
# Specify and fit a GARCH(1,1) model with an AR(1) mean equation
spec2 <- ugarchspec(
variance.model = list(model = "sGARCH", garchOrder = c(1, 1)),
mean.model = list(armaOrder = c(1, 0), include.mean = TRUE),
distribution.model = "norm" # standard normal innovations
)
fit2 <- ugarchfit(
spec = spec2,
data = df_monthly$pct_change
)
show(fit2)
*---------------------------------*
* GARCH Model Fit *
*---------------------------------*
Conditional Variance Dynamics
-----------------------------------
GARCH Model : sGARCH(1,1)
Mean Model : ARFIMA(1,0,0)
Distribution : norm
Optimal Parameters
------------------------------------
Estimate Std. Error t value Pr(>|t|)
mu 1.934127 0.328464 5.8884 0.000000
ar1 -0.081618 0.051847 -1.5742 0.115443
omega 2.762622 1.223137 2.2586 0.023906
alpha1 0.089354 0.026514 3.3701 0.000751
beta1 0.867601 0.035604 24.3680 0.000000
Robust Standard Errors:
Estimate Std. Error t value Pr(>|t|)
mu 1.934127 0.329930 5.8622 0.000000
ar1 -0.081618 0.045735 -1.7846 0.074333
omega 2.762622 1.227563 2.2505 0.024418
alpha1 0.089354 0.023757 3.7612 0.000169
beta1 0.867601 0.029859 29.0570 0.000000
LogLikelihood : -1496.564
Information Criteria
------------------------------------
Akaike 6.9840
Bayes 7.0313
Shibata 6.9838
Hannan-Quinn 7.0027
Weighted Ljung-Box Test on Standardized Residuals
------------------------------------
statistic p-value
Lag[1] 0.001980 0.9645
Lag[2*(p+q)+(p+q)-1][2] 0.007669 1.0000
Lag[4*(p+q)+(p+q)-1][5] 1.446899 0.8594
d.o.f=1
H0 : No serial correlation
Weighted Ljung-Box Test on Standardized Squared Residuals
------------------------------------
statistic p-value
Lag[1] 0.7965 0.3721
Lag[2*(p+q)+(p+q)-1][5] 3.0707 0.3943
Lag[4*(p+q)+(p+q)-1][9] 4.4756 0.5102
d.o.f=2
Weighted ARCH LM Tests
------------------------------------
Statistic Shape Scale P-Value
ARCH Lag[3] 3.036 0.500 2.000 0.08145
ARCH Lag[5] 4.025 1.440 1.667 0.17168
ARCH Lag[7] 4.443 2.315 1.543 0.28715
Nyblom stability test
------------------------------------
Joint Statistic: 0.9213
Individual Statistics:
mu 0.30060
ar1 0.09275
omega 0.25608
alpha1 0.40934
beta1 0.36894
Asymptotic Critical Values (10% 5% 1%)
Joint Statistic: 1.28 1.47 1.88
Individual Statistic: 0.35 0.47 0.75
Sign Bias Test
------------------------------------
t-value prob sig
Sign Bias 1.059 0.2900
Negative Sign Bias 1.234 0.2177
Positive Sign Bias 1.037 0.3003
Joint Effect 2.643 0.4500
Adjusted Pearson Goodness-of-Fit Test:
------------------------------------
group statistic p-value(g-1)
1 20 19.86 0.4030
2 30 30.88 0.3709
3 40 40.33 0.4115
4 50 43.72 0.6864
Elapsed time : 0.0209651 The estimated equations based on the optimal parameters are as follows:
Mean equation:
\[ r_t = 0.3006 + 0.0928\,r_{t-1} + u_t \]
Conditional variance equation:
\[ \sigma_t^2 = 0.2561 + 0.4093\,u_{t-1}^2 + 0.3689\,\sigma_{t-1}^2 \]
The mean equation includes both a constant term (\(\mu = 0.3006\)) and an autoregressive term (\(\text{AR}(1) = 0.0928\)).
This implies a small positive autocorrelation in returns — today’s return depends slightly on yesterday’s return.
The ARCH term (\(\alpha_1 = 0.4093\)) shows how strongly volatility reacts to recent shocks (\(u_{t-1}^2\)).
Larger shocks yesterday lead to higher variance today.
The GARCH term (\(\beta_1 = 0.3689\)) captures volatility persistence — how much of yesterday’s volatility carries forward to today.
The total persistence is: \[ \alpha_1 + \beta_1 = 0.7782 \] which is less than 1, meaning volatility is mean-reverting (it gradually returns to normal after shocks).
The long-run (unconditional) variance is: \[ \frac{\omega}{1 - \alpha_1 - \beta_1} = \frac{0.2561}{1 - 0.7782} \approx 1.15 \] giving a long-run volatility of \(\sqrt{1.15} \approx 1.07\).
| Term | Symbol | Value | Interpretation |
|---|---|---|---|
| Mean (intercept) | \(\mu\) | 0.3006 | Average return |
| AR(1) term | \(\phi_1\) | 0.0928 | Dependence on \(r_{t-1}\) |
| Constant (variance intercept) | \(\omega\) | 0.2561 | Baseline variance |
| ARCH(1) | \(\alpha_1\) | 0.4093 | Reaction to new shocks (\(u_{t-1}^2\)) |
| GARCH(1) | \(\beta_1\) | 0.3689 | Persistence of volatility (\(\sigma_{t-1}^2\)) |
In short:
The model suggests modest return autocorrelation and strong but decaying volatility persistence —
periods of high volatility subside gradually toward a long-run variance of about 1.15.
Increase the lag orders in the variance equation to estimate a more complex GARCH(\(p\), \(q\)) model. This allows volatility to depend on multiple past shocks and lagged variances.
# Specify and fit a GARCH(3,2) model
spec2 <- ugarchspec(
variance.model = list(model = "sGARCH", garchOrder = c(3, 2)),
mean.model = list(armaOrder = c(0, 0), include.mean = TRUE),
distribution.model = "norm" # standard normal distribution
)
fit2 <- ugarchfit(spec = spec2,
data = df_monthly$pct_change
)
show(fit2)
*---------------------------------*
* GARCH Model Fit *
*---------------------------------*
Conditional Variance Dynamics
-----------------------------------
GARCH Model : sGARCH(3,2)
Mean Model : ARFIMA(0,0,0)
Distribution : norm
Optimal Parameters
------------------------------------
Estimate Std. Error t value Pr(>|t|)
mu 1.894470 0.354161 5.34918 0.000000
omega 4.514070 2.448657 1.84349 0.065258
alpha1 0.022001 0.042639 0.51599 0.605861
alpha2 0.059165 0.056502 1.04715 0.295032
alpha3 0.066420 0.066628 0.99688 0.318823
beta1 0.485763 0.483566 1.00454 0.315116
beta2 0.296979 0.415483 0.71478 0.474745
Robust Standard Errors:
Estimate Std. Error t value Pr(>|t|)
mu 1.894470 0.323101 5.86340 0.000000
omega 4.514070 2.593365 1.74062 0.081750
alpha1 0.022001 0.035733 0.61572 0.538080
alpha2 0.059165 0.056121 1.05425 0.291767
alpha3 0.066420 0.055307 1.20092 0.229782
beta1 0.485763 0.240991 2.01569 0.043832
beta2 0.296979 0.191163 1.55354 0.120294
LogLikelihood : -1496.411
Information Criteria
------------------------------------
Akaike 6.9926
Bayes 7.0588
Shibata 6.9921
Hannan-Quinn 7.0187
Weighted Ljung-Box Test on Standardized Residuals
------------------------------------
statistic p-value
Lag[1] 2.956 0.08557
Lag[2*(p+q)+(p+q)-1][2] 2.964 0.14353
Lag[4*(p+q)+(p+q)-1][5] 4.436 0.20446
d.o.f=0
H0 : No serial correlation
Weighted Ljung-Box Test on Standardized Squared Residuals
------------------------------------
statistic p-value
Lag[1] 0.0009432 0.9755
Lag[2*(p+q)+(p+q)-1][14] 5.0244426 0.7665
Lag[4*(p+q)+(p+q)-1][24] 8.9404252 0.8020
d.o.f=5
Weighted ARCH LM Tests
------------------------------------
Statistic Shape Scale P-Value
ARCH Lag[6] 0.3535 0.500 2.000 0.5521
ARCH Lag[8] 1.0593 1.480 1.774 0.7478
ARCH Lag[10] 1.6861 2.424 1.650 0.8286
Nyblom stability test
------------------------------------
Joint Statistic: 1.1627
Individual Statistics:
mu 0.2614
omega 0.2218
alpha1 0.3566
alpha2 0.4509
alpha3 0.3843
beta1 0.3686
beta2 0.3692
Asymptotic Critical Values (10% 5% 1%)
Joint Statistic: 1.69 1.9 2.35
Individual Statistic: 0.35 0.47 0.75
Sign Bias Test
------------------------------------
t-value prob sig
Sign Bias 1.0209 0.3079
Negative Sign Bias 1.4782 0.1401
Positive Sign Bias 0.4694 0.6390
Joint Effect 2.4096 0.4919
Adjusted Pearson Goodness-of-Fit Test:
------------------------------------
group statistic p-value(g-1)
1 20 18.09 0.5162
2 30 32.00 0.3199
3 40 44.79 0.2418
4 50 45.81 0.6031
Elapsed time : 0.03688502 The estimated equations based on the optimal parameters are as follows:
Mean equation:
\[ r_t = 1.894470 + u_t \]
Conditional variance equation:
\[ \begin{aligned} \sigma_t^2 &= 4.514070 + 0.022001\,u_{t-1}^2 + 0.059165\,u_{t-2}^2 + 0.066420\,u_{t-3}^2 \\ &\quad + 0.485763\,\sigma_{t-1}^2 + 0.296979\,\sigma_{t-2}^2 \end{aligned} \]
The mean equation indicates an average return level of 1.894, with shocks captured by \(u_t\).
The ARCH terms (\(\alpha_1 = 0.022\), \(\alpha_2 = 0.059\), \(\alpha_3 = 0.066\)) measure how recent shocks affect current volatility.
Each squared residual \(u_{t-i}^2\) represents the impact of past return shocks on today’s variance.
The GARCH terms (\(\beta_1 = 0.486\), \(\beta_2 = 0.297\)) measure the persistence of volatility — how much past volatility carries over into the current period.
The total persistence is: \[ \alpha_1 + \alpha_2 + \alpha_3 + \beta_1 + \beta_2 = 0.9303 \] which is less than 1, implying that volatility is mean-reverting but decays slowly.
The long-run (unconditional) variance is: \[ \frac{\omega}{1 - \sum \alpha_i - \sum \beta_j} = \frac{4.514070}{1 - 0.930328} \approx 64.79 \] giving a long-run volatility of \(\sqrt{64.79} \approx 8.05\).
| Term | Symbol | Value | Interpretation |
|---|---|---|---|
| Mean (intercept) | \(\mu\) | 1.894 | Average return |
| Constant (variance intercept) | \(\omega\) | 4.514 | Baseline variance |
| ARCH(1) | \(\alpha_1\) | 0.022 | Response to \(u_{t-1}^2\) |
| ARCH(2) | \(\alpha_2\) | 0.059 | Response to \(u_{t-2}^2\) |
| ARCH(3) | \(\alpha_3\) | 0.066 | Response to \(u_{t-3}^2\) |
| GARCH(1) | \(\beta_1\) | 0.486 | Persistence of \(\sigma_{t-1}^2\) |
| GARCH(2) | \(\beta_2\) | 0.297 | Persistence of \(\sigma_{t-2}^2\) |
In short:
The model captures moderate short-term shock effects and strong volatility persistence,
suggesting that periods of high volatility decay gradually toward a long-run average variance of about 64.8.
In GARCH models, we don’t assume that shocks are always normally distributed — financial data often show fat tails (extreme events more likely than under normality) and sometimes asymmetry (negative shocks having different effects from positive ones).
So we can specify how the conditional errors are distributed —
“norm”: simplest; if residuals’ histogram looks bell-shaped.
“std”: use if you see fat tails in residuals (high kurtosis, large outliers).
“sstd”: use if both fat tails and asymmetry appear (negative shocks produce larger volatility).
“ged”/“sged”: slightly less common but can fine-tune tail behavior.
Variance depends on past squared shocks:
\[ \sigma_t^2 = \alpha_0 + \alpha_1 u_{t-1}^2 + \alpha_2 u_{t-2}^2 + \dots \]
Limitation: Needs many lags to capture persistent volatility
Adds lagged variance terms to make it efficient:
\[ \sigma_t^2 = \alpha_0 + \alpha_1 u_{t-1}^2 + \phi_1 \sigma_{t-1}^2 \]
Mean equation:
\[
r_t = \mu + u_t
\]
Shock structure:
\[
u_t = \sigma_t \cdot \varepsilon_t, \quad \varepsilon_t \sim N(0,1)
\]
Variance equation (GARCH):
\[
\sigma_t^2 = \alpha_0 + \alpha_1 u_{t-1}^2 + \phi_1 \sigma_{t-1}^2
\]
Therefore:
\[
u_t \sim N(0, \sigma_t)
\]
where \(\sigma_t\) is time-varying and predicted by GARCH.
Two components:
\[ \sigma_t^2 = \alpha_0 + \alpha_1 u_{t-1}^2 + \phi_1 \sigma_{t-1}^2 \]
| Term | What it does |
|---|---|
| \(\alpha_0\) | Baseline variance / minimum level of volatility that exists even when there are no shocks |
| \(\alpha_1\) | Reaction to new shocks Big shock yesterday → higher variance today |
| \(\phi_1\) | Persistence of volatility High variance yesterday → stays high today |
Condition for stability: \(\alpha_1 + \phi_1 < 1\).
Volatility Clustering: Large changes follow large changes; small follow small
Persistence: \(\alpha_1 + \phi_1\) measures how long volatility stays elevated
Long-run variance:
\[
\frac{\alpha_0}{1 - \alpha_1 - \phi_1}
\]
gives the average variance over time. It accounts for both the baseline level and the persistence of volatility.
Mean reversion: After a shock, variance gradually returns to the long-run average
Mean equation (ARIMA, AR, etc.): Models expected returns — gives residuals \(u_t\)
Variance equation (GARCH): Models expected volatility — uses those \(u_t\) to predict \(\sigma_t^2\)
Historical note:
The parameter nameomegaoriginates from Bollerslev’s (1986) original formulation of the GARCH model.
In his notation, the conditional variance equation was written as
\[ \sigma_t^2 = \omega + \alpha_1 u_{t-1}^2 + \beta_1 \sigma_{t-1}^2, \]
where \(\omega\) represents the constant term — the baseline or “intercept” in the variance equation.Later textbooks often replaced \(\omega\) with \(\alpha_0\) for consistency with the ARCH model notation,
but most software packages (includinggarchFit()andrugarch) retainedomegato stay true to Bollerslev’s original paper.