Introduction to valuer

Payoff of guarantees

To initialize a product class the user needs to pass a guarantee_payoff class object. For example, a payoff_rollup class is a type of guarantee_payoff (a sub class actually) which represents a payoff which is the maximum between the account value at the time of the payment and the rollup of the premium at a given rate from issue time till the time of the payment. Currently the payoff classes implemented are:

• payoff_rollup
  The payoff will be the highest between the account value and the roll-up of the premium at a given rate,
• payoff_ratchet
  The payoff will be the highest account value recorded at some specified times.

So let’s begin by initializing a payoff_rollup object after choosing a 2% yearly roll-up rate and a single premium of $100 to finance the underlying fund of the policy:

library(valuer)

## Loading required package: orthopolynom

## Loading required package: polynom

rate <- constant_parameters$new(0.02)

premium <- 100

Above we have instantiated a constant_parameters (?constant_parameters) class object to hold the rollup rate. At this point we pass these info into the initialize (new) method of the payoff_rollup class.

rollup <- payoff_rollup$new(premium, rate)

We now have an object of class payoff_rollup called rollup.
As an example of using the object rollup, we’re going to calculate the payoff of a contract starting Jan 1st with expiration on Dec 31st 2020 assuming the account value is 108 at expiration. This is done by means of the get_payoff public method of the rollup object. The method takes as argument the account value and a vector with the initial and expiration times as timeDate(?timeDate) objects.

t1 <- timeDate::timeDate("2016-01-01")

t2 <- timeDate::timeDate("2020-12-31")

rollup$get_payoff(108, c(t1,t2))

## [1] 110.5231

Above the initial time t1 and expiration t2 were defined as timeDate objects. The public method get_payoff of the rollup object is passed the account value and a vector with these times to get the payoff. It’s clear in this case the roll up of premium from 2016-01-01 to 2016-12-31 was payed since greater than the account value.

Contract fees

To finance the guarantees a fee is deducted from the policy account on a yearly basis. As said, in valuer an object of a VA product class holds all info about the contract and the policyholder. So in order to initialize an object representing a VA contract we need to pass the fee to the new method of this object. In addition, since valuer implements a simple state-dependent fee structure where the fee is payed only when the account is below a certain barrier, we would need to pass the value of that barrier as well. In case a barrier is not provided, the fee is assumed to be constant (not state-dependent).

# A constant fee of 2% per year
fee <- constant_parameters$new(0.02)

#Barrier for a state-dependent fee. The fee will be applied only if
#the value of the account is below the barrier
barrier <- 200

Surrender guarantee

A VA contract may embed a surrender guarantee in case the policyholder is allowed to surrender the contract. Usually in this case a penalty is applied and the amount withdrawn would be the account value at the time of surrender multiplied by the penalty. Currently valuer supports a surrender penalty which is held constant during the time the insured can surrender or two types of penalties which are decreasing functions of time. The idea behind decreasing penalties is that the penalty should be set very high at contract inception to discourage early surrenders and progressively reduces to zero. A specific class named penalty_class is used to model these penalties.

#Withdrawal penalty applied in case the insured surrenders the contract

penalty <- penalty_class$new(type = 1, 0.01)

This is a constant penalty or type 1.

To specify a penalty of type 2 we must provide the time (in years) by which the insured can surrender. Usually this is the end of the accumulation period.

penalty <- penalty_class$new(type = 2, const = 0.08, T = 5)

penalty$get(0)

## [1] 0.08

penalty$get(5)

## [1] 0

The type 2 penalty is calculated using the formula

\[ const\left (1-\frac{t}{T} \right )^3 \]

so it’s 8% at t = 0 and decreases to 0 at time t = 5.

Finally type 3 is calculated using the formula

\[ 1 - \exp\left ( -(const / T) \left ( T - \min\left (t, T \right ) \right )\right ) \]

penalty <- penalty_class$new(type = 3, const = 0.08, T = 5)

penalty$get(0)

## [1] 0.07688365

penalty$get(5)

## [1] 0

Initializing a VA product object

We now need to set the time-line of the product. In order to do it we first need to specify the issue date, the end of the accumulation period date and the end of the benefit payment date in case it doesn’t coincide with the end of accumulation date, for example with the GMIB or GMWB riders. Then, we need to set the age of the insured when the contract is issued and pass all other contractual features discussed above. For example if we want to set up a VA with a GMAB rider:

#Five years time-line
begin <- timeDate::timeDate("2016-01-01")
end <- timeDate::timeDate("2020-12-31")

#Age of the policyholder.
age <- 60

#Sets up a VA contract with GMAB guarantee.
contract <- GMAB$new(rollup, t0 = begin, t = end, age = age, fee = fee, barrier = barrier, penalty = penalty)

The contract fee, state-dependent fee barrier, initial age of the insured can be changed later on by means of the corresponding set methods. Their values can be checked using the get methods:

contract$set_age(50)
contract$get_age()

## [1] 50

contract$set_barrier(200)
contract$get_barrier()

## [1] 200

contract$set_fee(constant_parameters$new(0.03))

It’s also possible to change the constant of the surrender penalties while it’s not possible to switch between penalty types unless we redefine the penalty_class object.

contract$set_penalty(0.04)

head(contract$get_penalty())

## [1] 0.003992011 0.003989834 0.003987657 0.003985479 0.003983302 0.003981125

penalty2 <- penalty_class$new(type = 2, const = 0.08, T = 5)

contract$set_penalty_object(penalty2)

head(contract$get_penalty())

## [1] 0.08000000 0.07986892 0.07973799 0.07960720 0.07947656 0.07934605

These set methods are useful when we want to get estimates of the contract value corresponding to different features of the contract but using the same simulated paths of the financial and demographic processes. As an example of this usage scenario, the set_fee method would be called repeatedly while estimating the fair fee of a contract.

Pricing with va_bs_engine

In order to set up a va_bs_engine class object we need to pass to its new method the constant interest rate, parameters to specify the geometric Brownian motion and parameters for the Weibull intensity of mortality function, these are:

Financial parameters

the constant volatility, constant dividend rate and initial value of the fund. Volatility and dividend rate must be passed as constant_parameters objects, while the initial value is a scalar and usually set equal to the single premium. It has to be noted that the volatility and dividend parameters should be obtained by calibrating the GBM (risk neutral) model against some option market data relevant to the actual VA we want to price.

Mortality parameters

Two scalar parameters are needed called c1 and c2. These should be calculated by fitting the survival probabilities obtained by a life table of the population our insured belongs to.

Of course, besides the above parameters, we need to pass the VA product object initialized previously.

#Interest rate
r <- constant_parameters$new(0.03)

#Initial value of the underlying fund
spot <- premium

#Volatility
vol <- constant_parameters$new(0.2)

#Dividend rate
div <- constant_parameters$new(0.0)

#Sets up the pricing engine specifying the va_contract, the interest rate
#the parameters of the Weibull intensity of mortality, the initial fund
#value, the volatility and dividends rate
engine <- va_bs_engine$new(contract, r, c1=90.43, c2=10.36, spot,
volatility=vol, dividends=div)

Now that we have an engine we can use its public methods do_static and do_mixed to price the contract. While do_static assumes the policyholder cannot surrender the contract, do_mixed assumes the possibility of surrendering the contract. We need to pass to these methods the number of Monte Carlo simulations we want to run and a mc_gatherer object (see ?mc_gatherer). The mc_gatherer holds the point estimates of each Monte Carlo run and has methods to print the final estimate, Monte Carlo standard error and a convergence table. It is also possible to call its method plot to produce a Monte Carlo convergence graph.

#Number of paths to simulate
no_of_paths <- 1e3

#Gatherer for the MC point estimates
the_gatherer <- mc_gatherer$new()

engine$do_static(the_gatherer, no_of_paths)
the_gatherer$get_results()

##       mean        se
## 1 106.5022 0.8315978

By means of the method do_mixed it is possible to estimate the value of the contract with the assumption the policyholder can surrender. This is similar to estimating the value of an American derivative by means of Least Squares Monte Carlo. In fact, the algorithm of do_mixed is a slightly modified LSMC procedure which takes into account the demographic risk as outlined in BMOP2011. Therefore, besides the number of simulations and the gatherer, we need to pass parameters relevant to LSMC:

• degree - the maximum degree of the Laguerre polynomials used by least squares in LSMC. The default value is three,

• freq - the frequency of the surrender decision. The default is “3m” which corresponds to a surrender decision taken every three months.

In addition, we can specify by means of the simulate boolean parameter if the simulation of the stochastic processes has to start from scratch or not. The default is TRUE but should be changed to FALSE in case we want to run do_mixed on the same paths used by do_static.

engine$do_mixed(the_gatherer, no_of_paths, degree = 3, freq = "3m", simulate = FALSE)

the_gatherer$get_results()

##       mean        se
## 1 108.3584 0.6211143

the_gatherer$plot()

Pricing with va_sde_engine

The class va_sde_engine leverages on yuima to simulate the process of the underlying fund, interest rate and intensity of mortality. The later is assumed to be independent from the previous two processes. By means of yuima it is possible to specify the stochastic differential equations defining these processes by means of an elegant and natural syntax. The generality and flexibility of this approach comes at the expense of a sensible increase of run time of the simulation and requires a fairly good knowledge of yuima. The initial setup of a va_sde_engine object is far more difficult than the setup of a va_bs_engine one. In fact, along with the VA product object we also need to pass two fairly convoluted lists of parameters:

financial_parms

a list with three elements to set up the financial stochastic differential equations in yuima. The first is a list of parameters for the yuima’s simulate function. The second is a list with parameters for the yuima’s setModel function. The third element is a vector with indices indicating the interest rate and log price in the solve.variable argument of setModel.

mortality_parms

a list with three elements to set up the stochastic differential equations of the intensity of mortality process in yuima. The first is a list of parameters for the yuima’s simulate function. The second is a list with the parameters for the yuima’s setModel function. The third is a vector with indices indicating the intensity of mortality in the solve.variable argument of setModel.

Unfortunately it is not straightforward to write these parameter lists and the user will have to carefully review yuima documentation to do so. With regards to that, good references are IAC2011 and BROU2014. To provide some examples, valuer comes with two set of parameter lists. The first set, financials_BMOP2011 and mortality_BMOP2011, initializes yuima with the system of SDEs from BMOP2011 while the second, financials_BBM2010 and mortality_BBM2010, implements the SDEs described in BBM2010. These example parameters were derived by the authors of the cited articles by calibrating the financial and demographic risk models using market data and specific life tables. In general, an user of valuer may want to do the same using market data relevant to his/her problem.

financials_BMOP2011

## [[1]]
## [[1]]$a
## [1] 0.6
## 
## [[1]]$b
## [1] 0.03
## 
## [[1]]$sigma.r
## [1] 0.03
## 
## [[1]]$c
## [1] 1.5
## 
## [[1]]$d
## [1] 0.04
## 
## [[1]]$sigma.K
## [1] 0.4
## 
## [[1]]$rho.SK
## [1] -0.7
## 
## [[1]]$xinit
##      r0      K0      S0 
## 0.03000 0.04000 4.60517 
## 
## 
## [[2]]
## [[2]]$sol
## [1] "r" "K" "S"
## 
## [[2]]$drift
## [1] "a*(b - r)" "c*(d - K)" "r - 0.5*K"
## 
## [[2]]$diffusion
##      [,1]            [,2]            [,3]                        
## [1,] "sigma.r*sq(r)" "0"             "0"                         
## [2,] "0"             "sigma.K*sq(K)" "0"                         
## [3,] "0"             "rho.SK*sq(K)"  "sqrt(1 - (rho.SK^2))*sq(K)"
## 
## 
## $ind
## [1] 1 3

mortality_BMOP2011

## [[1]]
## [[1]]$x
## [1] 60
## 
## [[1]]$e
## [1] 0.5
## 
## [[1]]$sigma.mu
## [1] 0.03
## 
## [[1]]$c1
## [1] 90.43
## 
## [[1]]$c2
## [1] 10.36
## 
## [[1]]$xinit
## [1] 0.002462962
## 
## 
## [[2]]
## [[2]]$drift
## [1] "e*(mu(t, mortality_BMOP2011[[1]]$x, mortality_BMOP2011[[1]]$c1, mortality_BMOP2011[[1]]$c2) - m)"
## 
## [[2]]$diffusion
## [1] "sigma.mu*sq(m)"
## 
## [[2]]$solve.variable
## [1] "m"
## 
## 
## $ind
## [1] 1

Given these parameter lists it is possible to set up the pricer by typing

engine <- va_sde_engine$new(contract, financials_BMOP2011, mortality_BMOP2011)

## Warning in yuima.warn("'delta' (re)defined."): 
## YUIMA: 'delta' (re)defined.

The value of the contract is again estimated by means of the engine do_static and do_mixed methods

engine$do_static(the_gatherer, no_of_paths)
the_gatherer$get_results()

engine$do_mixed(the_gatherer, no_of_paths, degree = 3, freq = "3m", simulate = FALSE)
the_gatherer$get_results()

Parallel simulation

With the va_sde_engine the simulation of the financial and intensity of mortality paths takes quite some time and it’s where most of the time is spent during the execution of do_static and do_mixed methods. Luckily, va_sde_engine supports the parallel execution of the simulations by means of the foreach package. By executing them on multiple processors/cores the overall execution time will decrease.
It is important to highlight that if we want foreach to execute the code in parallel a parallel back-end must be registered first. Otherwise, the simulation will be done sequentially and foreach will issue a warning to inform about that (just the first time it’s used without a back-end). In order to register a parallel back-end you can use the package doParallel. Please check the introductory vignette of doParallel for more details.
The following code will start a cluster with two CPU cores, register it with doParallel so that foreach can use it during the do_mixed method:

library(doParallel)

cl <- makeCluster(2)

registerDoParallel(cl)

engine$do_mixed(the_gatherer, no_of_paths)

In a system with at least two CPU cores, you should see a sensible decrease of the execution time. Of course, foreach should be installed first, otherwise the engine will run the simulations sequentially.