A Two Factor Forward Curve Model with Stochastic Volatility for Commodity Prices

Abstract: We describe a model for evolving commodity forward prices that incorporates three important dynamics which appear in many commodity markets: mean reversion in spot prices and the resulting Samuelson effect on volatility term structure, decorrelation of moves in different points on the forward curve, and implied volatility skew and smile. This model is a "forward curve model" - it describes the stochastic evolution of forward prices - rather than a "spot model" that models the evolution of the spot commodity price. Two Brownian motions drive moves across the forward curve, with a third Heston-like stochastic volatility process scaling instantaneous volatilities of all forward prices. In addition to an efficient numerical scheme for calculating European vanilla and early-exercise option prices, we describe an algorithm for Monte Carlo-based pricing of more generic derivative payoffs which involves an efficient approximation for the risk neutral drift that avoids having to simulate drifts for every forward settlement date required for pricing.

• The paper proposes a two-factor forward curve model incorporating Heston-like stochastic volatility to capture the Samuelson effect and implied volatility skew in commodity markets.
• It utilizes numerical integration of characteristic functions and solves ODEs to efficiently price European vanilla and early exercise options.
• The Monte Carlo simulation with drift approximation reduces computational costs while accurately replicating forward prices and volatility under extreme market conditions.

Stochastic Volatility Model for Commodity Forward Curves

This paper introduces a two-factor forward curve model with stochastic volatility for commodity prices, addressing key dynamics observed in commodity markets: mean reversion in spot prices (Samuelson effect), decorrelation of forward price movements, and implied volatility skew and smile. The model extends existing forward curve models by incorporating a Heston-like stochastic volatility process, thereby enabling the modeling of implied volatility skew and smile, which are not captured by models with lognormally distributed forward prices.

Model Definition and Dynamics

The model defines the evolution of a forward price $F(t, T)$ through the following stochastic differential equations:

$$\frac{dF(t,T)}{F(t,T)} = \sqrt{v(t)} \sigma \left( e^{-\beta_1 (T-t)} dz_1(t) + R e^{-\beta_2 (T-t)} dz_2(t) \right)$$

$$dv(t) = \beta (1 - v(t)) dt + \alpha \sqrt{v(t)} dz_3(t)$$

$E[dz_1(t) dz_2(t)] = \rho dt$

$E[dz_1(t) dz_3(t)] = \rho_1 dt$

$E[dz_2(t) dz_3(t)] = \rho_2 dt$

The forward prices are driven by two Brownian motions, $dz_1(t)$ and $dz_2(t)$, with their impact on different points of the forward curve ($T$) modulated by exponential terms incorporating mean reversion parameters $\beta_1$ and $\beta_2$. This captures the Samuelson effect, where the volatility of forward prices decreases with longer settlement times. The correlation between these two Brownian motions is given by $\rho$, allowing for decorrelation of forward price movements. The instantaneous volatility of all forward prices is scaled by a Heston-like stochastic volatility factor $v(t)$, driven by a third Brownian motion $dz_3(t)$. The parameters $\alpha$ and $\beta$ control the volatility of this factor and its mean reversion speed, respectively. The correlations $\rho_1$ and $\rho_2$ between the volatility factor and the forward price factors determine the implied volatility skew.

Vanilla and Early Exercise Option Pricing

The paper outlines a method for pricing European vanilla and early exercise options under the model. The approach involves calculating the characteristic function of the log-forward price, $x(t,T) = \ln \left(F(t,T)/F(0,T) \right)$, at the option's expiration time $t_e$. Unlike the original Heston model, the characteristic function cannot be derived in closed form and requires one numerical integration. The call price $C(K)$ for strike $K$ is given by

$C(K) = D(T) (F(0,T) - \frac{K}{2} - \frac{K}{\pi} \int_{\theta=0}^\infty {\Re[\frac{f(\theta) e^{-i \theta \ln K/F(0,T)}{\theta^2 + i \theta}]})$

where $D(T)$ is the discount factor to settlement time $T$.

Figure 1: At-the-money forward implied volatility for vanilla options as a function of expiration time, demonstrating the decay of implied volatility due to the Samuelson effect.

The method necessitates solving a pair of ordinary differential equations for functions $A(\tau)$ and $B(\tau)$ numerically, where $\tau = t_e - t$:

$\frac{dA}{d\tau} = \beta B$

$$\frac{dB}{d\tau} = -\frac{1}{2} (\theta^2 + i \theta) \sigma^2_F(t,T) - \beta B + \frac{\alpha^2}{2} B^2 + i \theta B \alpha \sigma (e^{-\beta_1 \tau} e^{-\beta_1 (T-t_e)} \rho_1 + R e^{-\beta_2 \tau} e^{-\beta_2 (T-t_e)} \rho_2)$$

These equations are integrated with initial conditions $A(0) = B(0) = 0$. The resulting characteristic function is then used in a second numerical integration to obtain the option price. This approach remains computationally efficient compared to finite difference methods or Monte Carlo simulation. The model can replicate the volatility term structure and skew observed in commodity markets, as demonstrated in the paper.

Figure 2: Vanilla option implied volatility as a function of option strike price, illustrating the model's ability to generate implied volatility skew.

Monte Carlo Simulation with Drift Approximation

For pricing exotic derivatives that depend on multiple forward prices, the paper introduces a Monte Carlo simulation algorithm. The key challenge is the computational cost of simulating the integrated drift term, which depends on the path of the stochastic volatility factor $v(t)$ and the settlement time $T$. To address this, the paper proposes an approximation for the drift term:

$$I(t,T) \approx \int_{s=0}^t {\sigma^2_F(s,T) ds} + k(t,T) \int_{s=0}^t {w(s) ds}$$

where $I(t,T)$ is the integrated drift term, $w(t) = v(t) - 1$, and $k(t,T)$ is a deterministic "drift factor" chosen to match the variance of the approximate drift with the true drift. The drift factor is computed as:

$$k^2(t,T) = \frac{2 \int_{s_2=0}^t \int_{s_1=0}^{s_1} \sigma^2_F(s_1,T) \sigma^2_F(s_2,T) E[w(s_1) w(s_2)] ds_1 ds_2}{2 \int_{s_2=0}^t \int_{s_1=0}^{s_1} E[w(s_1) w(s_2)]}$$

where

$$E[w(s_1) w(s_2)] = \frac{\alpha^2}{2 \beta} (1 - e^{-2 \beta s_1}) e^{-\beta (s_2 - s_1)}$$

This approximation significantly reduces the computational burden, as only the integrated Heston factor $\int_{s=0}^t {w(s) ds}$ needs to be tracked as an additional path variable.

Figure 3: Approximation error in the simulated forward with extreme market and model parameters, showing minimal error in basis points.

Figure 4: Standard error on the expected value of the forward in the Monte Carlo simulation with extreme market and model parameters, indicating higher error than approximation error.

Figure 5: Approximation error in the simulated at-the-money implied volatility with extreme market and model parameters, demonstrating very small error in percentage volatility.

Figure 6: Standard error on the expected value of the at-the-money implied volatility in the Monte Carlo simulation with extreme market and model parameters, showing higher error than approximation error.

Figure 7: Approximation error in the simulated out-of-the-money implied volatility with extreme market and model parameters, for a strike 1.4 times the forward, indicating tiny error in percentage volatility.

Figure 8: Standard error on the expected value of the out-of-the-money implied volatility in the Monte Carlo simulation with extreme market and model parameters, for a strike 1.4 times the forward, showing higher error than approximation error.

Validation and Conclusion

The paper validates the drift approximation by comparing Monte Carlo simulations with and without the approximation under extreme parameter settings. The approximation demonstrates high accuracy, with errors significantly smaller than the simulation standard errors for forward prices, at-the-money implied volatilities, and out-of-the-money implied volatilities. The model provides a framework for pricing a wide range of commodity derivatives, balancing accuracy and computational efficiency.

The paper suggests future research directions, including:

• Adding deterministic time dependency to $\sigma$ for a better fit to at-the-money implied volatility.
• Adding deterministic time dependency to $\alpha$, $\rho_1$, and/or $\rho_2$ to calibrate the model to market levels of implied volatility smile and skew.
• Adding a second volatility factor to control the correlation of implied volatility moves across different expiration tenors.