Options on bonds, swaptions and caps

rl.ZeroCouponBondOption, rl.CouponBondOption, rl.Swaption and rl.CapFloor under SwitchingVasicek, SwitchingHullWhite or SwitchingG2.

The problem

Under the Vasicek model, a call expiring at \(T\) with strike \(K\) on the zero-coupon bond maturing at \(S\) is

\[C = \mathbb{E}\Big[e^{-\int_0^Tr_s\,ds}\,\big(P(T, S) - K\big)^+\Big].\]

The bond at expiry depends on the regime then: on \(\{y_T = j\}\), \(P(T, S) = A_j\,e^{-b\,r_T}\) with \(b = B(S - T)\) and \(A_j = a_j(S - T)\) from the bond’s reduced system. So the option pays when \(r_T < r_j^* = \log(A_j/K)/b\): there is one exercise boundary per regime.

Reduction

The price is a sum over the regime at expiry of terms of the form

\[\mathbb{E}\big[e^{-\int_0^Tr - c\,r_T}\,\mathbf 1\{r_T < r_j^*\}\,\mathbf 1\{y_T = j\}\big], \qquad c = b \text{ and } c = 0 .\]

Without the indicator on \(r_T\) this is a bond-like quantity with two changes. The payoff contains \(e^{-c\,r}\), so the coefficient of the rate starts at \(c\) rather than zero,

\[B_c(t) = c\,e^{-at} + \frac{1 - e^{-at}}{a},\]

and it is paid only in regime \(j\), so the terminal vector is the unit vector \(e_j\) rather than \(\mathbf 1\):

\[\Psi_j(c) = \mathbb{E}\big[e^{-\int_0^Tr - c\,r_T}\,\mathbf 1\{y_T = j\} \mid y_0 = i\big] = e^{-B_c(T)\,r_0}\,a_i(T), \qquad a' = (Q + \operatorname{diag} g)\,a, \quad a(0) = e_j,\]

\[g_i(t) = -a\,b_i\,B_c + \tfrac12\sigma_i^2\,B_c^2 .\]

The indicator on \(r_T\) is restored by Gil–Pelaez inversion in the complex exponent \(c - iu\):

\[\mathbb{E}\big[e^{-\int_0^Tr - c\,r_T}\,\mathbf 1\{r_T < r_j^*\}\,\mathbf 1\{y_T = j\}\big] = \tfrac12\Psi_j(c) - \frac1\pi\int_0^\infty\operatorname{Im}\big(e^{-iu\,r_j^*}\,\Psi_j(c - iu)\big)\frac{du}{u}.\]

The regime at expiry costs one solve per regime and nothing else.

Closed form, first order

With the terminal vector \(e_j = \tfrac12\mathbf 1 + \tfrac12s_j(1, -1)^\top\), \(s_1 = 1\), \(s_2 = -1\), the stationary half is carried as for a bond and the other half relaxes, leaving a first-order imprint through \(\tilde g(0)\), which is not zero here because \(B_c(0) = c\):

\[\Psi_j(c) = \tfrac12\,e^{-B_c(T)\,r_0}\exp\Big(\int_0^T\bar g + \frac\varepsilon2\int_0^T\tilde g^{\,2}\Big) \Big(1 \pm \frac\varepsilon2\,\tilde g(T) + s_j\,\frac\varepsilon2\,\tilde g(0)\Big) + O(\varepsilon^2).\]

The term in \(\tilde g(T)\) is the memory of the regime at the start; the term in \(\tilde g(0)\) is its mirror image, the memory of the regime at expiry. Because \(B_c = c\,e^{-at} + B\), the first-order factor is a polynomial in \(c\), \(\Pi_j(c) = \sum_{k=0}^4\pi_{jk}c^k\), and a factor \(c^k\) on the transform is the \(k\)-th derivative of the Gaussian density of \(r_T\) under the averaged model’s \(T\)-forward measure, \(r_T \sim N(\mu, v)\). The price is then

\[C = \bar P(0, T)\sum_{j=1}^2\tfrac12\Big(R_0(A_j) + \varepsilon\sum_{k=0}^4(-1)^k\pi_{jk}\,R_k(A_j)\Big) + O(\varepsilon^2),\]

\[R_0(A) = A\,e^{-b\mu + \frac12b^2v}\,\Phi(z^* + b\sqrt v) - K\,\Phi(z^*), \qquad z^* = \frac{\log(A/K)/b - \mu}{\sqrt v},\]

\[R_k(A) = v^{-k/2}\Big(A\,e^{-b\mu + \frac12b^2v}\sum_{i=0}^k\binom ki(-b\sqrt v)^{k-i}H_i(z^* + b\sqrt v) - K\,H_k(z^*)\Big), \qquad H_0 = \Phi, \quad H_i(w) = -\mathrm{He}_{i-1}(w)\,\varphi(w).\]

\(R_0\) is Jamshidian’s formula, so at \(\varepsilon = 0\) this is the averaged model’s option price. The ingredients are

\[\mu = e^{-aT}r_0 + a\bar b\,B(T) - \bar s\,M_{1,1}, \qquad v = \bar s\,M_{2,0}, \qquad M_{k,m} = \int_0^Te^{-kat}B(t)^m\,dt,\]

\[\Pi_j(c) = \tfrac12\int_0^T\tilde g_c^{\,2} \pm \tfrac12\,\tilde g_c(T) + \tfrac12\,s_j\,\tilde g_c(0), \qquad \tilde g_c = -a\tilde b\,B_c + \tfrac12\tilde s\,B_c^2 .\]

Python

The closed form is regimelib._engine.bond_option_explicit.call; the library’s engines solve the Gil–Pelaez integrals without expanding.

from regimelib._engine.bond_option_explicit import call as closed_form_call
r0, a, b, sigma, T, S = 0.04, 0.5, [0.05, 0.03], [0.015, 0.010], 1.0, 4.0

print(" lam     K   Jamshidian, averaged   closed form, first order   library, numerical")
for lam in (10.0, 20.0, 40.0):
    model = rl.SwitchingVasicek(rl.RegimeChain.twoState(lam, lam), r0, a, b, sigma)
    for K in (0.86, 0.90):
        option = rl.ZeroCouponBondOption("call", K, T, S)
        option.setPricingEngine(rl.NumericalSwitchingEngine(model, regime=0))
        averaged = closed_form_call(a, b, sigma, r0, T, S, K, lam, start=0, order=0)
        first = closed_form_call(a, b, sigma, r0, T, S, K, lam, start=0, order=1)
        print(f"{lam:4.0f}  {K:.2f}   {averaged:.8f}             {first:.8f}                 {option.NPV():.8f}")

lam     K   Jamshidian, averaged   closed form, first order   library, numerical
0.86   0.02638714             0.02622112                 0.02621545
0.90   0.00135934             0.00134597                 0.00134605
0.86   0.02638714             0.02630380                 0.02630239
0.90   0.00135934             0.00135215                 0.00135217
0.86   0.02638714             0.02634539                 0.02634504
0.90   0.00135934             0.00135562                 0.00135562

Each doubling of the switching rate divides the first-order error by about four.

Coupon bonds, swaptions and caps

A coupon bond \(\sum_kc_kA_{kj}e^{-b_kr}\) is decreasing in \(r\) in every regime, so Jamshidian’s decomposition holds regime by regime. In regime \(j\) there is one \(r_j^*\) at which the coupon bond equals the strike, and

\[\Big(\sum_kc_kA_{kj}e^{-b_kr_T} - K\Big)^+ = \sum_kc_k\big(A_{kj}e^{-b_kr_T} - K_{kj}\big)^+ \quad\text{on } \{y_T = j\}, \qquad K_{kj} = A_{kj}e^{-b_kr_j^*}.\]

The strikes differ by regime, which is the only change from the model without regimes. A receiver swaption is the call on the coupon bond at par and a payer swaption the put. A caplet on \([t_{k-1}, t_k]\) with strike \(K\) is \(1 + \tau_kK\) puts on the bond maturing at \(t_k\), struck at \(1/(1 + \tau_kK)\) and expiring at \(t_{k-1}\). Under Hull–White the deterministic shift scales each cash flow; under G2++ the single Gaussian variable \(B_ax_T + B_bz_T\) takes the place of \(r_T\).

lam = 10.0
model = rl.SwitchingVasicek(rl.RegimeChain.twoState(lam, lam), 0.03, 0.5, [0.06, 0.02], [0.015, 0.008])
swaption = rl.Swaption("payer", 2.0, [3.0, 4.0, 5.0, 6.0, 7.0], 0.035, notional=100.0)
swaption.setPricingEngine(rl.NumericalSwitchingEngine(model, regime=0))
reference = swaption.NPV()
print(f"2x5 payer swaption, numerical   {reference:.6f}")
for order in (0, 1, 2, 4):
    swaption.setPricingEngine(rl.FastSwitchingEngine(model, order=order, regime=0))
    print(f"  expansion, order {order}            {swaption.NPV():.6f}   error {abs(swaption.NPV() - reference):.1e}")

cap = rl.CapFloor("cap", [0.5 * k for k in range(1, 11)], 0.04, notional=100.0)
cap.setPricingEngine(rl.NumericalSwitchingEngine(model, regime=0))
print(f"semiannual cap at 4%, 0.5y to 5y  {cap.NPV():.6f}")

2x5 payer swaption, numerical   1.924124
  expansion, order 0            1.895914   error 2.8e-02
  expansion, order 1            1.923805   error 3.2e-04
  expansion, order 2            1.924123   error 1.9e-06
  expansion, order 4            1.924124   error 3.7e-08
semiannual cap at 4%, 0.5y to 5y  1.171073

Bermudan swaptions depend on the path of the rate and are priced on the coupled equations on a grid (Other instruments).