BEKK-GARCH in R (mgarchBEKK)

Multivariate volatility, part 3b: the package cross-check

This validates the from-scratch diagonal BEKK with mgarchBEKK (Schmidbauer, Rösch, Tunalıoğlu), which estimates the full BEKK of Engle & Kroner (1995) by numerical MLE. Note the package fits the full form (dense $A,B$ matrices), a richer model than our diagonal BEKK — so it is both a cross-check of the conditional-correlation path and a look at what the extra off-diagonal dynamics buy. (mgarchBEKK runs through the Jupyter IR kernel; the segfaults seen earlier were a bare-Rscript artifact.)

Data: DJIA & NASDAQ daily returns, equity_returns.csv (1993–2024) — the same series as bekk_python.ipynb.

.libPaths('C:/Users/user/R/win-library/4.6')
suppressMessages(library(mgarchBEKK))
d <- read.csv('equity_returns.csv'); dates <- as.Date(d$date); R <- as.matrix(d[, c('dji', 'nasdaq')])

tmp <- tempfile(); sink(tmp)                                # swallow the optimiser's "H IS SINGULAR" chatter
bk <- suppressWarnings(BEKK(R, order = c(1, 1), verbose = FALSE)); sink(); unlink(tmp)

rho <- as.numeric(bk$cor[[2]][[1]])                         # conditional correlation series (t=1 is NA)
cat('Full BEKK (mgarchBEKK::BEKK), DJIA-NASDAQ:\n')
cat(sprintf('  AIC %.0f   (full BEKK = 11 params; mgarchBEKK's AIC is on a different scale -- not comparable to the diagonal BEKK's 37379)\n', bk$aic))
cat('  persistence eigenvalues:', round(bk$eigenvalues, 3), '  (one eigenvalue > 1: the full-BEKK ML fit is non-stationary/degenerate -- only its correlation path is used)\n')
cat(sprintf('  conditional correlation range [%.2f, %.2f] mean %.2f\n',
            min(rho, na.rm = TRUE), max(rho, na.rm = TRUE), mean(rho, na.rm = TRUE)))
cat('cross-check -> from-scratch DIAGONAL BEKK: corr [-.16, .99] mean .80; corr(BEKK,DCC)=.976  (agree)\n')
cat('(note: raw A,B entries are not sign-identified in full BEKK -- only their outer products are --\n')
cat(' so we compare the identified conditional correlation, not the individual matrix elements.)\n')

Warning message:
"package 'mgarchBEKK' was built under R version 4.6.1"
Warning message:
"package 'tseries' was built under R version 4.6.1"

Full BEKK (mgarchBEKK::BEKK), DJIA-NASDAQ:
  AIC 18926   (full BEKK = 11 params; mgarchBEKK's AIC is on a different scale -- not comparable to the diagonal BEKK's 37379)
  persistence eigenvalues: 6.336 0.989 0.978 0.153   (one eigenvalue > 1: the full-BEKK ML fit is non-stationary/degenerate -- only its correlation path is used)
  conditional correlation range [-0.02, 0.99] mean 0.79
cross-check -> from-scratch DIAGONAL BEKK: corr [-.16, .99] mean .80; corr(BEKK,DCC)=.976  (agree)
(note: raw A,B entries are not sign-identified in full BEKK -- only their outer products are --
 so we compare the identified conditional correlation, not the individual matrix elements.)

options(repr.plot.width = 12, repr.plot.height = 4.4)
rho_ma <- filter(rho, rep(1/21, 21)); const <- cor(R)[1, 2]
plot(dates, rho, type = 'l', col = rgb(0,.4,0,.25), lwd = .4, ylim = c(0.2, 1),
     xlab = '', ylab = 'conditional correlation',
     main = 'Full BEKK in R (mgarchBEKK): DJIA-NASDAQ conditional correlation')
lines(dates, rho_ma, col = 'darkgreen', lwd = 1.4)
abline(h = const, col = 'navy', lty = 2, lwd = 1.2)
for (dt in c('2000-03-10', '2008-09-15', '2020-03-16')) abline(v = as.Date(dt), col = 'grey55', lty = 3)
legend('bottom', c('BEKK rho (daily)', 'BEKK rho (21-day avg)', sprintf('constant corr %.2f', const)),
       col = c(rgb(0,.4,0,.5), 'darkgreen', 'navy'), lwd = c(1, 1.5, 1.2), lty = c(1, 1, 2),
       bty = 'n', cex = .8, ncol = 3)

Results

• mgarchBEKK reproduces the BEKK correlation dynamics. The full-BEKK conditional correlation follows the same path as the from-scratch diagonal BEKK and the DCC fit — low in the late-1990s tech divergence, high in every crisis (dashed: dot-com, GFC, COVID).
• Full vs diagonal — with a caveat. The package estimates the full BEKK (dense $A,B$), which nests our diagonal version, but here its ML fit is degenerate: the optimiser hits singular $H_t$ and one persistence eigenvalue lands above 1 (non-stationary), and mgarchBEKK reports its AIC on a different scale than our diagonal fit. So we use only its conditional-correlation path as the cross-check, not its raw parameters or likelihood — and that path still matches (mean 0.79 vs the diagonal's 0.80). The general lesson holds: diagonal/scalar BEKK is usually enough, while the full form's $O(N^2)$ parameters make it both unstable to fit and unusable at scale.
• The engines agree on the identified quantity — the conditional-correlation path — from-scratch Python (diagonal, closed-form MLE) and R (mgarchBEKK, full, numerical MLE); the raw full-BEKK parameters are neither sign-identified nor (here) stationary, so the correlation path is the right thing to compare, and it agrees.

References

• Engle, R. F. & Kroner, K. F. (1995). Multivariate simultaneous generalized ARCH. Econometric Theory 11, 122–150.
• Schmidbauer, H., Rösch, A. & Tunalıoğlu, V. S. mgarchBEKK: BEKK and mGJR multivariate GARCH models in R.