Insurance embedded options and guarantees are not traded on markets. It’s therefore a bit difficult to define precisely what really is their current price; seen as a best estimate reserve from the insurer’s perspective.

The principle of market consistent actuarial valuation is to price these options and guarantees as we would price traded financial derivatives with similar features. A simple example of an insurance contract with embedded option, is a one offering a minimum guaranteed rate. In such a contract, the insurer will have to pay at a future date \(T\), some interests equal to : \[ max(L(T, S), K) \] where \(K\) is a fixed guaranteed rate defined at the inception of the contract for the maturity \(T\), and \(L(T, S)\) is a stochastic interest rate prevailing at time T, for a future maturity \(S\).

This expression can be rewritten as : \[ K + (L(T, S) - K)^+ \] where \(x \mapsto (x)^+\) is a function giving the positive part of \(x\) (i.e \(x\) if \(x > 0\), and \(0\) otherwise). The contract can thus be seen as a porftolio consisting in a zero-coupon bond and an interest rate caplet (here with nominal value equal to \(1\)).

Where a closed-form expression exists for pricing the option, it could directly be used. But the embedded options values don’t usually have closed-form expressions, and Monte Carlo simulation is required.

Hence, risks factors are projected after being calibrated to derivatives, where a deep and liquid market exists for these products. For fixed income modeling, for example, interest rates models parameters are usually calibrated to caps (sum of caplets) or swaptions (option on swaps) prices. These derivatives offer payoffs that are similar to insurance embedded options.

After the model calibration, one kind of model validation is based on studying the Monte Carlo simulation errors made by the calibrated model on the estimation of zero-coupons prices and caps/swaptions prices. The tests carried out are often referred to as ‘martingale tests’ (but they aren’t) or ‘market consistency tests’, but there’s no new concept here : it’s about studying the Monte Carlo simulation errors.

In this post, we consider the \(G2++\) short rate model (a 2-factor Hull & White model). The simulation of the model is made with R package ESGtoolkit. For more resources on ESGtoolkit, see the package vignette, or these slides. We’ll study the Monte Carlo errors made by the simulation on the estimation of zero-coupons prices.

The G2++ short rate model is defined as : \[ dx_t = -a x_t + \sigma dW^{(x)}_t \] \[ dy_t = -b y_t + \eta dW^{(y)}_t \] \[ dW^{(x)}_tdW^{(y)}_t = \rho dt \] \[ r_t = x_t + y_t + \phi_t \]

where \(x_0 = 0\), \(y_0 = 0\), \((W^{(x)}_t)_{t \geq 0}\) and \((W^{(y)}_t)_{t \geq 0}\) are correlated standard brownian motions with instantaneous correlation equal to \(\rho\). \(\phi_t\) is a function depending on market instantaneous forward rates (for more details on this model and \(\phi_t\), see this book or the example in R below), ensuring theoretically that the zero-coupon prices generated by the model are equal to market prices. That is : \[ P(0, T) = P^M(0, T) \] where \(P(0, T)\) is the price at time \(t=0\) of a zero-coupon with maturity \(T\) derived from the model as : \[ \mathbb{E}\left[ exp \left( -\int_0^T r_s d_s \right) \right] =: P(0, T) \] and \(P^M(0, T)\) is the market (bootstrapped) price at time \(t=0\) of a zero-coupon with maturity \(T\).

Here, the G2++ was previously calibrated to ATM Euro caps prices observed on December 31, 2011 (I used R package mcGlobaloptim), and we’ll investigate if \(P(0, T)\) is significantly close to \(P^M(0, T)\) for each \(T\), at a given \(5\%\) confidence level. That is, constructing confidence intervals for the mean of \[ exp \left( -\int_0^T r_s d_s \right) - P^M(0, T) \]

rm(list=ls())
# Observed maturities
u <- 1:30
# Yield to maturities
txZC <- c(0.01422,0.01309,0.01380,0.01549,0.01747,0.01940,
          0.02104,0.02236,0.02348, 0.02446,0.02535,0.02614,
          0.02679,0.02727,0.02760,0.02779,0.02787,0.02786,
          0.02776,0.02762,0.02745,0.02727,0.02707,0.02686,
          0.02663,0.02640,0.02618,0.02597,0.02578,0.02563)
# Zero-coupon prices
p <- c(0.9859794,0.9744879,0.9602458,0.9416551,0.9196671,
       0.8957363,0.8716268,0.8482628,0.8255457,0.8034710,
       0.7819525,0.7612204,0.7416912,0.7237042,0.7072136
       ,0.6922140,0.6785227,0.6660095,0.6546902,0.6441639,
       0.6343366,0.6250234,0.6162910,0.6080358,0.6003302,
       0.5929791,0.5858711,0.5789852,0.5722068,0.5653231)
# Creating a function for the simulation of G2++
suppressPackageStartupMessages(library(ESGtoolkit))
# Function of the number of scenarios, the method (
# antithetic or not)
simG2plus <- function(n, methodyc = c("NS", "SV", "HCSPL", "SW"))
{
    # Horizon, number of simulations, frequency
    horizon <- 20
    freq <- "semi-annual" 
    delta_t <- 1/2
    # Parameters found for the G2++
    a_opt <- 0.50000000
    b_opt <- 0.35412030
    sigma_opt <- 0.09416266
    rho_opt <- -0.99855687
    eta_opt <- 0.08439934
    # Simulation of gaussian correlated shocks
    eps <- ESGtoolkit::simshocks(n = n, horizon = horizon,
                     frequency = "semi-annual",
                     family = 1, par = rho_opt)
    # Simulation of the factor x
    x <- ESGtoolkit::simdiff(n = n, horizon = horizon, 
                 frequency = freq,  
                 model = "OU", 
                 x0 = 0, theta1 = 0, theta2 = a_opt, theta3 = sigma_opt,
                 eps = eps[[1]])
    # Simulation of the factor y
    y <- ESGtoolkit::simdiff(n = n, horizon = horizon, 
                 frequency = freq,  
                 model = "OU", 
                 x0 = 0, theta1 = 0, theta2 = b_opt, theta3 = eta_opt,
                 eps = eps[[2]])
    # Instantaneous forward rates, with spline interpolation
    methodyc <- match.arg(methodyc)
    fwdrates <- ESGtoolkit::esgfwdrates(n = n, horizon = horizon, 
    out.frequency = freq, in.maturities = u, 
    in.zerorates = txZC, method = methodyc)
    fwdrates <- window(fwdrates, end = horizon)
    # phi
    t.out <- seq(from = 0, to = horizon, 
                 by = delta_t)
    param.phi <- 0.5*(sigma_opt^2)*(1 - exp(-a_opt*t.out))^2/(a_opt^2) + 
    0.5*(eta_opt^2)*(1 - exp(-b_opt*t.out))^2/(b_opt^2) +
      (rho_opt*sigma_opt*eta_opt)*(1 - exp(-a_opt*t.out))*
      (1 - exp(-b_opt*t.out))/(a_opt*b_opt)
    param.phi <- ts(replicate(n, param.phi), 
                    start = start(x), deltat = deltat(x))
    phi <- fwdrates + param.phi
    colnames(phi) <- c(paste0("Series ", 1:n))
    # The short rates
    r <- x + y + phi
    colnames(r) <- c(paste0("Series ", 1:n))
    return(r)
}

Simulation of the short rates :

set.seed(3)
# 4 types of simulations, by changing the number of scenarios
ptm <- proc.time()
r.SW <- simG2plus(n = 10000, methodyc = "SW")
proc.time() - ptm
##    user  system elapsed 
##    1.61    0.08    1.68
par(mfrow=c(1, 2))
ESGtoolkit::esgplotbands(r.SW, xlab = 'time', ylab = 'short rate')
matplot(r.SW, type = 'l', xlab = 'time', ylab = 'short rate')

plot of chunk unnamed-chunk-2

Stochastic discount factors : \[ D(0, T) = exp \left( -\int_0^T r_s d_s \right) \]

# Stochastic deflators :
deltat_r <- deltat(r.SW)
Dt.SW <- ESGtoolkit::esgdiscountfactor(r = r.SW, X = 1)
Dt.SW <- window(Dt.SW, start = deltat_r, deltat = 2*deltat_r)
# Market prices
horizon <- 20
marketprices <- p[1:horizon]
# Monte Carlo prices
montecarloprices.SW <- rowMeans(Dt.SW)

# Confidence intervals
estim.SW <- apply(Dt.SW - marketprices, 1, function(x) t.test(x)$estimate)
conf.int.SW <- apply(Dt.SW - marketprices, 1, function(x) t.test(x)$conf.int)
conf.int.SW
##            [,1]       [,2]       [,3]       [,4]       [,5]       [,6]
## [1,] -1.004e-04 -0.0001715 -3.648e-04 -6.356e-04 -9.319e-04 -0.0011440
## [2,]  1.251e-06  0.0000536  6.396e-06  6.928e-06  6.672e-05  0.0002472
##            [,7]       [,8]       [,9]      [,10]      [,11]      [,12]
## [1,] -0.0012674 -0.0014203 -0.0017221 -0.0020859 -0.0025032 -0.0028276
## [2,]  0.0005139  0.0007302  0.0007682  0.0007136  0.0005705  0.0004902
##           [,13]      [,14]      [,15]      [,16]      [,17]      [,18]
## [1,] -0.0030000 -0.0032192 -0.0034290 -0.0035552 -0.0036973 -0.0037038
## [2,]  0.0005308  0.0005006  0.0004624  0.0004963  0.0005002  0.0006303
##           [,19]      [,20]
## [1,] -0.0036614 -0.0036352
## [2,]  0.0008046  0.0009548
par(mfrow=c(1, 2))
plot(marketprices, col = "blue", type = 'l', 
     xlab = "time", ylab = "prices", main = "Prices")
points(montecarloprices.SW, col = "red")
matplot(t(conf.int.SW), type = "l", xlab = "time", 
        ylab = "", main = "Confidence Interval \n for the price difference")
lines(estim.SW, col = "blue")

plot of chunk unnamed-chunk-4