Introduction
ETF options on SPY and QQQ are among the most liquid derivatives globally, making accurate pricing critical. The Black–Scholes model, the long-standing benchmark, assumes constant volatility and continuous price paths, failing to capture key market features: volatility skew, term-structure effects, and sudden jumps during stress events. Merton (1976) added Poisson jumps; Heston (1993) introduced stochastic volatility. Yet a core limitation persists: Poisson jumps are memoryless, unable to replicate the clustering of large moves observed in crises (2008, COVID-19), where one shock triggers others. This misprices tail-sensitive OTM puts. The Hawkes process resolves this via self-exciting dynamics: each jump temporarily raises the probability of further jumps before excitation decays, capturing feedback loops like stop-losses and margin calls. This paper introduces a Hawkes Jump-Diffusion (Hawkes-JD) model for ETF options, calibrates it to SPY data, and benchmarks it against Black–Scholes, Merton-JD, and Heston-SV.
ETF options on SPY and QQQ are among the most liquid derivatives globally, making accurate pricing critical. The Black–Scholes model, the long-standing benchmark, assumes constant volatility and continuous price paths, failing to capture key market features: volatility skew, term-structure effects, and sudden jumps during stress events. Merton (1976) added Poisson jumps; Heston (1993) introduced stochastic volatility. Yet a core limitation persists: Poisson jumps are memoryless, unable to replicate the clustering of large moves observed in crises (2008, COVID-19), where one shock triggers others. This misprices tail-sensitive OTM puts. The Hawkes process resolves this via self-exciting dynamics: each jump temporarily raises the probability of further jumps before excitation decays, capturing feedback loops like stop-losses and margin calls. This paper introduces a Hawkes Jump-Diffusion (Hawkes-JD) model for ETF options, calibrates it to SPY data, and benchmarks it against Black–Scholes, Merton-JD, and Heston-SV.
Model
Hawkes Conditional Intensity
The jump arrival rate λ(t) depends on the full history of past jumps:
λ(t) = μ +
Σ
ti < t
ϕ(t − ti)
where μ > 0 is the baseline (exogenous) jump rate, and ϕ(t − tᵢ) is an excitation kernel measuring the residual influence of past jump tᵢ. Using an exponential kernel
ϕ(t) = αe−βt
gives the closed-form intensity:
λ(t) = μ +
Σ
ti < t
α · e−β(t − ti)
α > 0 is the immediate intensity spike after each jump; β > 0 is the decay speed (high β = short memory, low β = long memory).
Asset Price Dynamics
The ETF price Sₜ follows a jump-diffusion SDE under the physical measure:
The ETF price Sₜ follows a jump-diffusion SDE under the physical measure:
dSt / St− = μS dt + σ dWt + (Yt − 1) dNt
where σ is constant volatility, Wₜ is a Brownian motion, Nₜ is the Hawkes counting process, and Yₜ > 0 is the lognormal jump multiplier (e.g., Yₜ = 0.90 implies a 10% instantaneous price drop).
Risk-Neutral Pricing
To price derivatives, we shift to the risk-neutral measure Q. Because jumps carry a directional return, a compensator
To price derivatives, we shift to the risk-neutral measure Q. Because jumps carry a directional return, a compensator
κ = EQ[Yt − 1]
must be subtracted from the drift. The risk-neutral SDE is:
dSt / St− = (r − λ(t) · κ) dt + σ dW̃t + (Yt − 1) dNt
The Hawkes parameters also shift under Q, with risk-neutral excitation α̃ typically larger than physical α, reflecting the fear premium embedded in option prices. Under Q, the intensity becomes:
λ̃(t) = μ̃ +
Σ
ti < t
α̃ · e−β̃(t − ti)
Monte Carlo Pricing (Ogata Thinning)
Because the Hawkes process is path-dependent, options are priced via Monte Carlo. A European call price is:
Because the Hawkes process is path-dependent, options are priced via Monte Carlo. A European call price is:
C0 = e−rT · EQ[max(ST − K, 0)] ≈ e−rT · (1/M)
M
Σ
j=1
max(ST,j − K, 0)
Jump arrival times are simulated using Ogata's thinning algorithm: propose arrivals from a bounding Poisson process, then accept with probability
λ(tcandidate) / λ̄
Empirical Analysis
Statistical Motivation
Daily SPY log-returns (2000–2024, n = 6,288) strongly reject both Gaussian and Poisson benchmarks. Jarque–Bera rejects normality (stat = 32,765, p ≈ 0; excess kurtosis = 11.17, skewness = −0.26). Ljung–Box on squared returns gives Q(20) = 7,891, confirming volatility clustering. Kolmogorov–Smirnov rejects Poisson inter-arrival times (p < 0.001). A rolling 63-day volatility filter identifies 209 jumps (~8.4/year). MLE gives:
μ̂ = 0.00985,
α̂ = 0.06881,
β̂ = 0.09754
Branching ratio η = α/β = 0.705, implying average cluster size of 3.4 jumps and a half-life of ~7 days. ΔAIC = 342 favors Hawkes over Poisson.
Figure 1: Hawkes intensity λ(t), SPY 2000–2024 - spikes at 2008 GFC, 2011 EU debt crisis, 2020 COVID-19, 2022 rate shock
Calibration Results
Four models are calibrated to WRDS OptionMetrics SPY surfaces (82 monthly dates, Jan 2018–Dec 2024). On a single calm-regime date (VIX ≈ 17–18, Aug 29 2025), in-sample IV-RMSE results are:
Four models are calibrated to WRDS OptionMetrics SPY surfaces (82 monthly dates, Jan 2018–Dec 2024). On a single calm-regime date (VIX ≈ 17–18, Aug 29 2025), in-sample IV-RMSE results are:
| Model | IV-RMSE |
|---|---|
| Hawkes–JD | 0.0183 |
| Merton–JD | 0.0245 |
| Black–Scholes | 0.0426 |
| Heston–SV | 0.1084 |
Table 1: In-sample IV-RMSE, Aug 29 2025 (calm regime, VIX ≈ 17–18). Lower is better
Figure 2: IV-RMSE bar chart - Hawkes-JD lowest, followed by Merton-JD, Black-Scholes, Heston-SV
| Model | Mean IV-RMSE | Median IV-RMSE | % Dates Best |
|---|---|---|---|
| Merton–JD | 0.0242 | 0.0192 | 89.0% |
| Hawkes–JD | 0.0313 | 0.0294 | 7.3% |
| Black–Scholes | 0.0454 | -- | -- |
| Heston–SV | 0.1321 | 0.1232 | 0.0% |
Table 2: Panel IV-RMSE summary across 82 monthly dates (Jan 2018–Dec 2024).
| Regime | Black–Scholes | Hawkes–JD | Merton–JD | Heston–SV |
|---|---|---|---|---|
| Calm | 0.040 | 0.017 | 0.018 | 0.120 |
| Normal | 0.043 | 0.030 | 0.023 | 0.131 |
| Stress | 0.059 | 0.052 | 0.033 | 0.146 |
Table 3: Mean IV-RMSE by market regime. Hawkes-JD narrows the gap vs Merton-JD under stress.
Figure 4: IV smile - post-crash λ₀ = 7λ̄Q produces materially steeper OTM put skew without recalibration
Conclusions
Merton-JD achieves the best overall IV fit (mean RMSE 0.0245, best on 89% of dates), but relies on economically implausible parameters disconnected from historical data. Hawkes-JD accepts a modest fitting penalty in exchange for structural discipline: physical-measure parameters are anchored to historical MLE estimates, and the risk-neutral baseline is derived, not freely optimized. This yields a quantifiable jump-clustering premium coherently tracks regime shifts, a property Merton-JD cannot provide. Under stress, Hawkes-JD matches Merton-JD most closely, validating the self-exciting mechanism exactly when tail risk is greatest. No model prices inside the bid-ask spread consistently, so all results represent relative rather than absolute accuracy.
Merton-JD achieves the best overall IV fit (mean RMSE 0.0245, best on 89% of dates), but relies on economically implausible parameters disconnected from historical data. Hawkes-JD accepts a modest fitting penalty in exchange for structural discipline: physical-measure parameters are anchored to historical MLE estimates, and the risk-neutral baseline is derived, not freely optimized. This yields a quantifiable jump-clustering premium coherently tracks regime shifts, a property Merton-JD cannot provide. Under stress, Hawkes-JD matches Merton-JD most closely, validating the self-exciting mechanism exactly when tail risk is greatest. No model prices inside the bid-ask spread consistently, so all results represent relative rather than absolute accuracy.
Written by Tommaso Domenis, Ales Langr, Vincenzo Volpe, Zikuan Wang