G2++

rl.SwitchingG2(chain, termStructure, a, sigma, b, eta, rho) is QuantLib’s G2(termStructure, a, sigma, b, eta, rho) with both volatilities and the correlation allowed to switch.

The model

\[r_t = x_t + z_t + \varphi(t), \qquad dx_t = -a\,x_t\,dt + \sigma_{y_t}\,dW^1_t, \qquad dz_t = -b\,z_t\,dt + \eta_{y_t}\,dW^2_t, \qquad d\langle W^1, W^2\rangle = \rho_{y_t}\,dt,\]

with \(x_0 = z_0 = 0\) and \(\varphi\) fitted to the market curve.

Reduction

Try \(u_i = e^{-B_a(t)x - B_b(t)z}a_i(t)\) in the pricing equation for \(\mathbb{E}[e^{-\int(x + z)}]\). The terms in \(x\) and in \(z\) cancel separately when \(B_a' = 1 - aB_a\) and \(B_b' = 1 - bB_b\), and

\[g_i(t) = \tfrac12\sigma_i^2\,B_a^2 + \tfrac12\eta_i^2\,B_b^2 + \rho_i\sigma_i\eta_i\,B_aB_b, \qquad B_a = \frac{1 - e^{-at}}{a}, \quad B_b = \frac{1 - e^{-bt}}{b},\]

half the variance rate of \(\int r\) in regime \(i\). As for Hull–White, the shift is fitted with the averaged covariance and cancels the averaged part: \(P_i(0, T) = P^M(0, T)\,e^{-\int\bar g}\,a_i(T)\).

Closed form

The switched quantities are the two variances and the covariance. With half-differences

\[\tilde s = \tfrac12(\sigma_1^2 - \sigma_2^2), \qquad \tilde h = \tfrac12(\eta_1^2 - \eta_2^2), \qquad \tilde c = \tfrac12(\rho_1\sigma_1\eta_1 - \rho_2\sigma_2\eta_2),\]

the square of \(\tilde g\) is a quartic in \(B_a\) and \(B_b\), which needs the integrals

\[M_{pq} = \int_0^TB_a^p\,B_b^q\,dt = \frac{1}{a^p\,b^q}\sum_{j=0}^{p}\sum_{k=0}^{q}(-1)^{j+k}\binom pj\binom qk\,\Phi(ja + kb), \qquad \Phi(0) = T, \quad \Phi(\kappa) = \frac{1 - e^{-\kappa T}}{\kappa}.\]

Then

\[\int_0^T\tilde g^{\,2} = \tfrac14\tilde s^2M_{40} + \tfrac14\tilde h^2M_{04} + \big(\tilde c^2 + \tfrac12\tilde s\tilde h\big)M_{22} + \tilde s\,\tilde c\,M_{31} + \tilde h\,\tilde c\,M_{13},\]

\[\tilde g(T) = \tfrac12\tilde s\,B_a^2 + \tfrac12\tilde h\,B_b^2 + \tilde c\,B_aB_b, \qquad \tilde g(0) = 0, \qquad \tilde g'(T) = \tilde s\,B_aE_a + \tilde h\,B_bE_b + \tilde c\,(E_aB_b + B_aE_b),\]

with \(E_a = e^{-aT}\), \(E_b = e^{-bT}\), and the bond relative to the market is the common formula with \(\int\bar g\) removed.

Python

from math import comb
a, b, sigma, eta, rho, flat, lam, T = 0.5, 0.1, [0.015, 0.006], [0.012, 0.008], [-0.4, -0.7], 0.03, 5.0, 5.0

def M(p, q):
    Phi = lambda kappa: T if kappa == 0 else (1 - np.exp(-kappa * T)) / kappa
    return sum((-1) ** (j + k) * comb(p, j) * comb(q, k) * Phi(j * a + k * b)
               for j in range(p + 1) for k in range(q + 1)) / (a ** p * b ** q)

Ea, Eb = np.exp(-a * T), np.exp(-b * T)
Ba, Bb = (1 - Ea) / a, (1 - Eb) / b
_, st = half([s * s for s in sigma])
_, ht = half([e * e for e in eta])
_, ct = half([rho[i] * sigma[i] * eta[i] for i in range(2)])
ratio = second_order(
    int_gbar=0.0,
    int_gt2=(st * st * M(4, 0) / 4 + ht * ht * M(0, 4) / 4 + (ct * ct + st * ht / 2) * M(2, 2)
             + st * ct * M(3, 1) + ht * ct * M(1, 3)),
    gt_T=st * Ba * Ba / 2 + ht * Bb * Bb / 2 + ct * Ba * Bb,
    gt_0=0.0,
    dgt_T=st * Ba * Ea + ht * Bb * Eb + ct * (Ea * Bb + Ba * Eb),
    lam=lam)

model = rl.SwitchingG2(rl.RegimeChain.twoState(lam, lam), flat, a, sigma, b, eta, rho)
print(f"closed form, second order  {ratio:.10f}")
print(f"library, numerical         {library_bond(model, T) / np.exp(-flat * T):.10f}")

closed form, second order  1.0000322126
library, numerical         1.0000322147

Instruments

Bonds, options on bonds and caps by the characteristic-function engines, where a zero-coupon bond depends on the factors through the single Gaussian variable \(B_ax_T + B_bz_T\); swaptions and Bermudans on the two-factor grid.