library(ggplot2)
library(parallel)
source("common/functions.r", chdir = TRUE)
data = ProcessRawData()
sp500 = data$sp50026 Running backtests
We are no longer estimating risk one day at a time. We are building a forecast engine that runs across thousands of trading days and produces the dataset the backtests need.
That means fixing one set of parameters for all four models, running them in parallel and saving the results in a format both Unix and Windows can read. The next chapter, Chapter 27, takes that dataset and tests it.
26.1 GARCH estimation speed
Backtesting requires fitting GARCH models thousands of times across rolling windows. The choice of estimation package determines total runtime. We benchmark GARCH(1,1) estimation across R, Python and Julia using a 2000-observation window from S&P 500 returns, the same window size used for the backtest estimation below.
| Platform | Package | Normal (ms) | Student-t (ms) |
|---|---|---|---|
| R | tseries | 0.7 | — |
| Julia | ARCHModels | 0.8 | 1.6 |
| Julia | Manual | 1.5 | — |
| Python | arch | 4.7 | 6.9 |
| R | fGarch | 13 | 41 |
| R | rugarch | 42 | 66 |
| R | tsgarch | 180 | 296 |
These timings are illustrative and relative rather than a reproducible benchmark, since hardware and package versions were not recorded. R tseries and Julia ARCHModels are the fastest implementations for normal errors. The manual Julia implementation using Nelder-Mead is roughly twice as slow despite being compiled code.
For Student-t distributions, Julia ARCHModels was the quickest in this run at 1.6 ms. Python arch at 7 ms is next, followed by R fGarch at 41 ms. The R frameworks rugarch and tsgarch, which support a wider range of GARCH variants and distributions, are slower, reflecting their more extensive modelling infrastructure.
For rolling window backtests with thousands of re-estimations the differences above compound, and on this hardware ARCHModels was the quickest for both distributions. In R, tseries is fastest for normal errors but only supports that distribution. For Student-t, fGarch is the fastest R option.
26.2 Data and libraries
Backtesting adds parallel processing and serialisation to the libraries used so far.
tsgarchprovides GARCH estimation (loaded viafunctions.r);parallelenables multicore processing viamclapply().
import numpy as np
import pandas as pd
import sys
sys.path.insert(0, 'common')
from functions import (ProcessRawData,
Risk_HS, Risk_EWMA, Risk_nGARCH, Risk_tGARCH, RunOneDay)
import multiprocessing
data = ProcessRawData()
sp500 = data['sp500']
y = sp500['y'].values- Risk functions (
Risk_HS,Risk_EWMA,Risk_nGARCH,Risk_tGARCH) andRunOneDayare imported fromcommon/functions.py; multiprocessingis loaded here for core detection, andjoblibis introduced later for the parallel run.
using Statistics, Distributions, DataFrames, Dates
using ARCHModels
using Serialization
include("common/functions.jl");
data = ProcessRawData();
sp500 = data["sp500"];
y = sp500[!, :y];- Risk functions (
Risk_HS,Risk_EWMA,Risk_nGARCH,Risk_tGARCH) andRunOneDayare loaded fromcommon/functions.jl; ARCHModelsprovides fast GARCH model estimation (see Table 26.1);- Julia’s built-in threading via
Threads.@threadsprovides parallel processing.
26.3 Backtesting framework setup
All four models use the same estimation window, probability level and portfolio value, so their forecasts can be compared directly.
par = list()
par$asset = "S&P 500" # asset being tested
par$probability = 0.01 # probability level (1% VaR)
par$value = 1000 # portfolio value in currency units
par$WE = 2000 # estimation window size
par$T = length(sp500$y) # total sample size
par$WT = par$T-par$WE # testing window size
par$lambda = 0.94 # decay parameter for EWMA
par$Methods = c("HS","EWMA","nGARCH","tGARCH") # models we runpar = {
'asset': 'S&P 500',
'probability': 0.01,
'value': 1000,
'WE': 2000,
'T': len(y),
'lambda': 0.94,
'Methods': ['HS', 'EWMA', 'nGARCH', 'tGARCH']
}
par['WT'] = par['T'] - par['WE']par = Dict(
"asset" => "S&P 500",
"probability" => 0.01,
"value" => 1000,
"WE" => 2000,
"T" => length(y),
"lambda" => 0.94,
"Methods" => ["HS", "EWMA", "nGARCH", "tGARCH"]
);
par["WT"] = par["T"] - par["WE"];26.4 Backtesting functions
Systematic backtesting requires consistent methodology — the same estimation windows and the same set of models, applied to every trading day in the sample. The risk functions and RunOneDay() are defined in the common function files:
- R:
common/functions.r—Risk_HS,Risk_EWMA,Risk_nGARCH,Risk_tGARCH,RunOneDay - Python:
common/functions.py—Risk_HS,Risk_EWMA,Risk_nGARCH,Risk_tGARCH,RunOneDay - Julia:
common/functions.jl—Risk_HS,Risk_EWMA,Risk_nGARCH,Risk_tGARCH,RunOneDay
These functions implement Historical Simulation, EWMA, GARCH with normal errors and GARCH with Student-t errors. The RunOneDay() function calculates all risk measures for a single trading day, returning VaR and ES for each method.
26.5 Saving results
We argued in Section 16.3.1 that intermediate files should normally be avoided. Backtesting provides an exception. Re-estimating GARCH models across thousands of days takes time, and the subsequent statistical analysis in the next chapter is separate from the computation here. We save the backtest results to files so the validation chapter can load them directly without re-running the estimation.
We use language-specific binary formats (.RData, .pkl, .jls) for speed and precision, as discussed in Section 10.5. For workflows involving multiple languages, use Parquet. All three languages can read it.
26.5.1 Metadata for reproducibility
Save metadata with the results. It complements reproducible environments by recording the versions that generated the output. Knowing which package versions were used helps diagnose why results might differ when re-running code months or years later, or when results differ across machines. For regulatory purposes, being able to demonstrate exactly how and when results were generated matters, and sharing results with colleagues is easier when they can verify they have compatible software versions.
We save the backtest results along with the parameters used and metadata including the run date, machine name and package versions. This creates a self-documenting output that can be understood and reproduced later.
26.6 Parallel processing
Each day’s forecast is independent of others within the rolling window framework, making backtesting well-suited for multi-core implementation. Most modern computers have at least 4 cores, and some have 16 or more. Most languages use only one core by default. Parallel backtesting can therefore deliver large speedups.
26.6.1 Available cores
cores = detectCores()
cat("Available cores:", cores, "\n")Available cores: 10cores = multiprocessing.cpu_count()
print(f"Available cores: {cores}")Available cores: 10cores = Sys.CPU_THREADS;
println("Available cores: $cores");Available cores: 426.6.2 Parallel implementation
Each language provides different mechanisms for parallel processing. The examples below show complete workflows that run the backtest in parallel and save the results with metadata.
26.6.2.1 R: mclapply (Mac/Linux)
The mclapply() function provides efficient parallel processing on Unix systems through forking. It distributes work across cores with minimal overhead.
# Run parallel backtest
backtest = mclapply(
(par$WE+1):par$T,
RunOneDay,
y = sp500$y,
par = par,
Methods = par$Methods,
mc.cores = cores,
mc.preschedule = TRUE
)
backtest = data.frame(do.call(rbind, backtest))
names(backtest) = RunOneDay(par = par, Methods = par$Methods, HeaderOnly = TRUE)
backtest$date = sp500$date[backtest$index]
backtest$date.t = ymd(backtest$date)
# Create metadata
metadata = list(
run_date = Sys.time(),
machine = Sys.info()["nodename"],
start_date = as.character(min(backtest$date)),
end_date = as.character(max(backtest$date)),
R_version = R.version.string,
tsgarch_version = packageVersion("tsgarch")
)
# Save with metadata
results = list(backtest = backtest, par = par, metadata = metadata)
save(results, file = "Risk/backtest.RData")26.6.2.2 R: foreach with doParallel (all platforms)
For Windows compatibility, or when forking is unavailable, the foreach package with doParallel provides a cluster-based alternative that works on all operating systems.
library(foreach)
library(doParallel)
cl = makeCluster(cores)
registerDoParallel(cl)
clusterEvalQ(cl, {
library(reshape2)
source("common/functions.r", chdir = TRUE)
})
clusterExport(cl, c("par", "sp500"))
backtest = foreach(i = (par$WE+1):par$T) %dopar% {
RunOneDay(i, y = sp500$y, par = par, Methods = par$Methods)
}
stopCluster(cl)
backtest = data.frame(do.call(rbind, backtest))
names(backtest) = RunOneDay(par = par, Methods = par$Methods, HeaderOnly = TRUE)
backtest$date = sp500$date[backtest$index]
backtest$date.t = ymd(backtest$date)
# Create metadata (created in the parent process after the loop, once
# results have been collected back from the workers)
metadata = list(
run_date = Sys.time(),
machine = Sys.info()["nodename"],
start_date = as.character(min(backtest$date)),
end_date = as.character(max(backtest$date)),
R_version = R.version.string,
tsgarch_version = packageVersion("tsgarch")
)
results = list(backtest = backtest, par = par, metadata = metadata)
save(results, file = "Risk/backtest.RData")The cluster approach requires explicitly loading libraries and exporting variables to worker processes, unlike mclapply() which inherits the parent environment through forking.
26.6.2.3 Python: joblib
Python’s joblib library provides straightforward parallel processing with Parallel and delayed. The n_jobs=-1 parameter uses all available cores.
from joblib import Parallel, delayed
import pickle
import platform
from datetime import datetime
import arch
# Run parallel backtest
results = Parallel(n_jobs=-1)(
delayed(RunOneDay)(t, y, par, par['Methods'])
for t in range(par['WE'], par['T'])
)
backtest = pd.DataFrame(results)
backtest['date'] = sp500['date'].values[par['WE']:par['T']]
backtest['date.t'] = pd.to_datetime(backtest['date'].astype(str), format='%Y%m%d')
# Create metadata
metadata = {
'run_date': datetime.now().isoformat(),
'machine': platform.node(),
'start_date': str(backtest['date'].min()),
'end_date': str(backtest['date'].max()),
'python_version': platform.python_version(),
'arch_version': arch.__version__,
'numpy_version': np.__version__,
'pandas_version': pd.__version__
}
# Save with metadata
with open("Risk/backtest.pkl", "wb") as f:
pickle.dump({'backtest': backtest, 'par': par, 'metadata': metadata}, f)26.6.2.4 Julia: Threads.@threads
Julia provides native multi-threading. The @threads macro distributes loop iterations across available threads. Julia threads share memory, avoiding the serialisation overhead of process-based parallelism.
using Pkg, Serialization
# Run parallel backtest
results = Vector{Dict{String, Any}}(undef, par["WT"])
Threads.@threads for i in 1:par["WT"]
t = par["WE"] + i
results[i] = RunOneDay(t, y, par; Methods=par["Methods"])
end
backtest = DataFrame(results)
backtest[!, :date] = [sp500[Int(idx), :date] for idx in backtest[!, :index]]
# Create metadata
metadata = Dict(
"run_date" => string(Dates.now()),
"machine" => gethostname(),
"start_date" => string(minimum(backtest[!, :date])),
"end_date" => string(maximum(backtest[!, :date])),
"julia_version" => string(VERSION),
"ARCHModels_version" => let
deps = Pkg.dependencies()
uuid = first(k for (k, v) in deps if v.name == "ARCHModels")
string(deps[uuid].version)
end
)
# Save with metadata
serialize("Risk/backtest.jls", Dict("backtest" => backtest, "par" => par, "metadata" => metadata))26.7 Sanity check
A quick plot of VaR against actual returns serves as a sanity check. The VaR line should track volatility clusters, and returns should occasionally breach the threshold but not too often.
26.8 Next steps
We now have a complete backtest dataset with VaR and ES forecasts from four models (HS, EWMA, GARCH, GARCH-t) across thousands of trading days. The data is saved in native formats for each language:
- R:
Risk/backtest.RData - Python:
Risk/backtest.pkl - Julia:
Risk/backtest.jls
The next chapter, Chapter 27, loads these files and implements statistical tests for backtesting validation, including:
- Violation frequency analysis to test if violations occur at the expected rate;
- Independence testing to check if violations are clustered;
- Comparative analysis across HS, EWMA, GARCH and GARCH-t models;
- ES evaluation using the Acerbi-Székely approach.