Incomplete market evaluation
of pension liabilities
University of Amsterdam
Julius Linssen
December 13, 2018
Verklaring eigen werk
Hierbij verklaar ik, Julius Linssen, dat ik deze scriptie zelf geschreven heb en dat ik de volledige verantwoordelijkheid op me neem voor de inhoud ervan.
Ik bevestig dat de tekst en het werk dat in deze scriptie gepresenteerd wordt origineel is en dat ik geen gebruik heb gemaakt van andere bronnen dan die welke in de tekst en in de referenties worden genoemd.
De Faculteit Economie en Bedrijfskunde is alleen verantwoordelijk voor de begelei-ding tot het inleveren van de scriptie, niet voor de inhoud.
Contents
1 Introduction 4
2 Forecasting Mortality 6
2.1 Lee-Carter Model . . . 6
2.2 Estimating the Lee-Carter Model . . . 7
2.3 Li-Lee Model . . . 9
2.4 Estimating the Li-Lee Model . . . 10
3 Applying Model Ambiguity 12 3.1 Model Dynamics And Solution . . . 12
3.2 Bounding Model Ambiguity . . . 15
4 The life insurance contract 18 4.1 Specifying the contract . . . 18
4.2 Applying the Li-Lee and Lee-Carter model . . . 18
4.3 Applying model ambiguity . . . 19
4.4 Results . . . 21
5 Conclusion 27 6 References 29 References 29 A Lee-Carter and Li-Lee results 30 A.1 The Lee-Carter Model . . . 30
B Log-normal deaths assumption 36
C Solution of the robust objective 41
D Proof of solution 46
1
Introduction
Since the 1970’s life expectancies have steadily been rising, which has caused pen-sion funds to encounter increasing costs in their payments to the penpen-sioners. The uncertainty around these increases are of a large concern to the pension industry. Several mortality models are used in practice to capture this longevity risk. Among these models, the Lee-Carter model is regarded as the benchmark method since it appeared in 1992 (R. D. Lee and L. Carter [7]). Since then many expansions of the model have been made, one of which intends to avoid diverging projected mortality rates in similar populations. This extension to the Lee-Carter model is called the Li-Lee model (N. Li and R. Lee [8]) and this is used by the Actuarial Society in the Netherlands to calculates dutch mortality figures while taking the mortality figures of other western European countries into account.
Although these models provide us with great methods of calculating future mor-tality rates, they become less useful at older ages and less certain when faced with very long-term liabilities. Besides the uncertainty in mortality calculations pension funds face another problem, pension liabilities are not products that one can easily trade in making the market for pension liabilities incomplete. When markets are complete there exists a unique way of pricing such that all contracts have a unique price. This is known as the First Fundamental Theorem of Asset Pricing. The Second Fundamental Theorem of Asset Pricing [11] states that the lack of this completeness implies no unique pricing kernels exist, but that there are multiple correct ways of pricing. This means that, even when all mortality figures are known perfectly, there still exists a chance of mismatch between optimal prices based on a complete market assumption and the consequence of the accompanying asset management strategy in practice. Balter and Pelsser [1] have studied a method of pricing liabilities in such a way that this mismatch is minimized. They show that under a few assumptions a unique analytical solution can be found, even when an insurer or pension funds is unsure about the model of mortality.
Optimizing trading strategies for profit maximization in an incomplete market gen-erally leads to an unbounded problem. A concave utility function can solve the unboundedness which has been considered in different settings. Bordigoni et al. [4] introduce uncertainty as an additional penalty term measured by its entropy. Ma-toussi et al. [10] focus on the uncertainty of the volatility. Hern´andez-Hern´andez and Schied [5] consider logarithmic and power utility and Mania and Tevzadze [9] solve the utility maximisation via duality. Balter and Pelsser [1] show that a utility
function is not needed in a setting where the goal is to maximize profit under model uncertainty, if the uncertainty is bounded.
In this thesis we explain the Lee-Carter and Li-Lee model and show how they can be used to calculate future mortality models from data in the Human Mortality Database. Anyone who subscribes to the organisation that manages this database can have access to its data and it is used by many instances, including the Actu-arial Society of the Netherlands, to calculate mortality figures. We also attempt to explain the theory of Balter and Pelsser and how to use their apply their model to Lee-Carter and Li-Lee parameters. Unfortunately it is not able to find analytical results with these parameters, so we also show a way to use the model of Balter and Pelsser to get incomplete market prices that are consistent with the Lee-Carter and Li-Lee results, but use a less desirable model to define the incomplete markets. This thesis aims to show the effect of an incomplete market to value the price of a life insurance contract. Furthermore, we want to compare the use of the Lee-Carter mdoel and the Li-Lee model and show how they impact the expected price of future liabilities given a complete or incomplete market. In section 2 we discuss the Lee-Carter and Li-Lee model and how to use deaths and exposures from the Human Mortality Database to calculate their parameters. Section 3 discusses the theory of the model that Balter and Pelsser introduced to calculate liabilities robustly using model uncertainty and defines the alternative models that capture our model uncer-tainty. We derive all the results we obtain using the previously discussed theory and finally, section 5 concludes.
2
Forecasting Mortality
In order to be able to calculate the liabilities that depend on mortality, we must model mortality. Modeling mortality has been an active area of research in disci-plines such as actuarial science and demology for many years. Due to the increasing proportion of aged people in Western countries, this forecasting of mortality became even more important, since aged people are more expensive with respect to insurance or pension funds. Where first models would treat the mortality distribution as static and provide only point estimates of future mortality rates, nowadays mortality is expressed as a function of time while taking into account uncertainty in the mortal-ity rates. The focus hereby lies on modeling the yearly central death rate µ(x(t), t) which is allowed to change over time t and depends on the age x of the individual or multiple individuals at hand. Consider T the time of death for an individual. For this individual we define their chance of dying before a certain time t as
P (T ≤ t) = 1 − e−R0tµ(x(s),s)ds.
Note that the age x(t) of individuals also change over time so x depends on t. This rate parameter µ(x, t) is referred to as the central death rate and has been the focus of mortality models ever since Lee and Carter based their model around it in 1992 [7]. In the next subsections we explain how this model is defined and how it can be applied to forecast mortality figures for future times, based on publically available data from the Netherlands.
In the Netherlands the Dutch Actuarial Society uses another model called the Li-Lee model to forecast mortality. This model allows other datasets to be included in the forecasts, which is used to include figures from other countries in the EU with similar living conditions as in the Netherlands. The last subsections will show how this model is defined and also how this model can be used to forecast mortality based on publically available data from the Netherlands and similar countries.
2.1
Lee-Carter Model
First we will use the Lee-Carter model to find estimates for µ(x, t). The Lee-Carter model assumes that the central rate is constant in an age-period, which describes all the individuals of a certain age in a certain calendar year. For example, all the people of age 30 in 1970 form an age-period, and are thus assumed to have the same central death rate over that period. It then estimates the logarithm of µx,t of a person with
age x in year t by:
where αx and βx are time-invariant parameters and κt is often called the mortality
index and is independent of age. The parameter αx describes the age-specific
mortal-ity on average. The parameter κt is a stochastic process and determines the change
in mortality rates over time, and this change differs when βx differs, meaning that
an average improvement in the mortality index can have a different effect on people aged 30 than children of 1 year old.
Note that if we model log(µx,t) = ˆαx+ ˆβxκˆtwith estimated values for αx, βx and
κt, that we find the same model if we multiply ˆβx with a constant c and divide ˆκt
by the same constant c. The same is true if we add c to all values ˆκt and subtract
ˆ
βx∗ c from all the values ˆαx. To specify a unique model we need to impose additional
restrictions and we choose the restrictionsP
tˆκt= 0 and
P
xβˆx = 1.
In order to calculate any future mortality rates the Lee-Carter model imposes a stochastic structure on κt, usually a random walk with drift
κt+1= κt+ δ + t (2)
where t are independently normally distributed variables. Now we explain how the
values for αx, βx, κt, δ and the variance of t can be estimated.
2.2
Estimating the Lee-Carter Model
As input for the modelling process of the Dutch mortailty figures of the Lee-Carter model we take discrete observations in the following forms:
• Exposures Ex,t, the amount of time that people of ages in [x, x + 1[ have lived
in the period [t, t + 1[
• Deaths Dx,t, the number of deaths among people of ages in [x, x + 1[ in the
period [t, t + 1[
We do this for the years 1970 until 2014 for all ages up to 90.This information is available for the Netherlands in the Human Mortality Database (HMDB) which anyone can access for free after registering. There is more data available in the database for the Netherlands, it has data up to 2016 at time of writing. However, for some countries that will later be used to estimate the Li-Lee mortality figures there is only data up to 2014 so for all the countries we will use the data up to 2014. How this data can be used to find values for the Lee-Carter model is explained next.
As mentioned, we assume the force of mortality µx,t in an age-period is
e−R0tµ(x(s),s)ds we then find that, for an individual, the likelihood of surviving for
ti < 1 years in an age-period and then dying is dtdF (t)|t=ti = e
µtiµ. Also we can find
that the likelihood of surviving for ti = 1 years in this age-period (so not dying in
this ageperiod) is 1 − F (ti) = e−µx,tti = e−µx,t. If we now define Ti to be the random
variable indicating years spent alive in the cell by person i and Di equal to 1 if the
individual dies, 0 if not we can write both probabilities as (µ)Die−µTi.
Now consider all the individuals in a certain age period (x, t). If we sum all the variables Ti of these individuals we retrieve our observation Ex,t =
P
iTi and if we
sum up the variables Di of these individuals we retrieve the observed Dx,t =PiDi.
Assuming that the deaths of individuals are all independently distributed we can now consider multiple individuals and use these identities to find the likelihood that Dx,t
deaths occured while a group of people spend Ex,t years alive during a certain
age-period. This likelihood can be expressed as a function of Dx,t, Ex,tand µx,t = eαx+βxκt
by L(Dx,t, Ex,t, µx,t(αx, βx, κt)) = Y i (µ)Die−µx,tTi = (µx,t) P iDieµx,tPiTi = (µ x,t)Dx,te−µx,tEx,t.
From now on, we will just refer to this likelihood as L. The log likelihood for all x and t is now given by
log(L) =X
x
X
t
(Dx,tlog(µx,t) − µx,tEx,t)
which, if we replace µx,t by writing it in terms of αx, βx and κt becomes
X
x
X
t
(Dx,t(αx+ βxκt) − Ex,teαx+βxκt).
We use this expression to find the estimators for αx, βx, and κt by maximizing the
likelihood. This is achieved by taking the partial derivatives to αx, βx, and κt and
setting them to 0: 0 = ∂ ∂αx logL( ˆαx, ˆβx, ˆκt) → 0 = X t (Dx,t− Ex,teαˆx+ ˆβxˆκt), 0 = ∂ ∂βx logL( ˆαx, ˆβx, ˆκt) → 0 = X t ˆ κt(Dx,t− Ex,teαˆx+ ˆBxκˆt), 0 = ∂ ∂κt logL( ˆαx, ˆβx, ˆκt) → 0 = X x ˆ βx(Dx,t− Ex,teαˆx+ ˆβxκˆt).
This system cannot be solved direcly, but the Newton/Raphson process [6] can be used to solve it iteratively. After this process gives the estimates ( ˆαx, ˆβx,ˆκx) we
would like to fit the stochastic structure we imposed on ˆκt with (2) where t are
independently normally distributed variables, the so called random walk with drift. Values for δ and the variance of tcan quickly be found by the rwf () function in R,
which concludes this section.
2.3
Li-Lee Model
In the Netherlands, whenever longevity calculations need to be made, the results of the Actuarial Society, the AG in Dutch for short, are often used as a reference point. As this thesis tries to do relevant calculations, we will also look into the same model they use for our longevity calculations: the Li-Lee model [8]. The Li-Lee model uses a similar approach as the Lee-Carter model , but considers a larger dataset and computes Lee-Carter values for this data. In our case, we will be looking at all the European countries with a GDP larger than the European average. For this data we estimate values Ax, Bx and Kt in the same way as is done for the Lee-Carter model:
log(µEurx,t ) = Ax+ BxKt+ εx,t
where we apply the restrictions P
tKt = 0 and
P
xBx = 1. To obtain the
mortality figures for Dutch individuals we must add the deviation of the Dutch data to the European data to the model as well:
log(µx,t) = Ax+ BxKt+ αx+ βxκt+ ε0x,t
where P
tκt = 0 and
P
xβx = 1. At this point it is important to note that this
is not exactly how Li and Lee introduced their model in 2005. They assume that the variables Ax found when modeling mortality the larger dataset are sufficient for
estimating mortality for the smaller dataset as well. Therefore they do not include αxin their model. As Li and Lee assume that in the long run the mortality figures for
the smaller and larger dataset converge, that is a very valid assumption. However, the AG chooses to include the αx as it does prove to be a significant difference, thus
giving a better within-sample fit, and we choose to resemble the AG-model as much as possible. Adding the αx term also allows us to fit the mortality model to the
larger and smaller dataset in very similar ways, as we will see later.
Here the future values for Kt are found by imposing a random walk with drift
structure on it. For the values κt we impose a First-order autoregressive structure.
European and Dutch datasets indeed converge. They are defined by the regressive equations
Kt+1 = Kt+ θ + t,
κt+1 = aκt+ δt,
where both t and δt are normally distributed with mean 0 and variances σ2 and
σ2
δ respectively. This once again gets you a model that fits the Dutch mortality data,
but this time we have taken two stochastic processes into account: The development of the European mortality and the development of the deviation of the Netherlands. This approach is taken because from 1970 onwards the mortality probablities have all similarly decreased. For this reason the AG has chosen to consider all the Eu-ropean countries with a GDP above the EuEu-ropean average. This ensures that the projection will not depend exclusively on data relating to the Netherlands, in which specific fluctuations may have occurred in the past, which do not necessarily say anything about future developments. The estimate is that the long-term increase in life expectancy in the Netherlands can be predicted more precisely by including a broader European population. The successive projections are also expected to be more stable than they would be if only data pertaining to the Netherlands were used. Besides the statistical benefits, it makes sense to model mortality this way, as one can assume that the decline in mortality rates comes partly from improvements in technology and medicine. Improvements in these fields and other relevant factors aren’t restricted to single countries. A higher longevity in one country therefore probably guarantees a higher longevity in other similar countries as well. We also know that welfare is positively correlated with longevity, which is why the choice for countries with similar welfare makes sense.
2.4
Estimating the Li-Lee Model
As input for our modeling process of the European figures we once again take discrete observations in the following forms:
• Exposures Ex,t, the amount of time that people of ages in [x, x + 1[ have lived
in the period [t, t + 1[
• Deaths Dx,t, the number of deaths among people of ages in [x, x + 1[ in the
We do this for the years 1970 until 2014 for all ages up to 90 and for all of the following countries: Belgium, Denmark, Germany, Finland, France, Ireland, Ice-land, Luxembourg, The Netherlands, Norway, Austria, United Kingdom, Sweden and Switzerland. It is worth mentioning that the figures that were found by the AG for estimating the values of the Li-Lee model were only based partially [12] on the data in the mortality database. At the time they made this choice because not all data up to 2014 was in the HMDB yet. At this time that data is in the database, so we choose to base our figures only on data available in the HMDB.
Once again we try to find estimates for the central death rate µx,t which for ease
of reference we will this time call µLL
x,t. From this point onward, when we refer to
the estimate of µx,t based on the Lee-Carter model we will use µLCx,t instead. As
mentioned in the previous subsection, the values for Ax, Bx and Kt are found in the
same way as the values αx, βx and κt were found in the Lee-Carter model, except
this time they are based of European data, not just Dutch data. The same goes for the value of θ and the variance of the error term t, that are part of the stochastic
structure of Kt. After these values are found we can estimate the parameters that
define the deviation of the Dutch model: αx, βx, κt, a and the variance of δt. To find
these values recall from 2 subsections before that we find estimates for the Lee-Carter values by maximizing the likelihood
logLEU =X
x
X
t
(Dx,tEUlog(µLLx,t) − µLLx,tEx,tEU).
This likelihood is for general µx,t and thus is also found for the Li-Lee model. By
using European values for Dx,tand Ex,tthe previously discussed estimates are found.
If we use these estimates and use Dutch values for Dx,t and Ex,twe can find estimates
for the parameters αx, βx, κt, a and δt: logLN L =
P x P t(D N L x,t log(µLLx,t) − µLLx,tEx,tN L).
If we now fill in log(µLLx,t) = Ax+ BxKt+ αx+ βxκt and use the estimated values we
already have we get the following likelihood: logLN L =X x X t (DN Lx,t ( ˆAx+ ˆBxKˆt+ αx+ βxκt) − e( ˆAx+ ˆBx ˆ Kt+αx+βxκt)EN L x,t ).
When we take the partial derivatives with respect to the values αx, βx and κx we see
that the term DN L
x,t ( ˆAx + ˆBxKˆt) falls away. Therefore we see that when we use the
modified Dutch exposure ˜EN L
x,t = Ex,t ∗ e( ˆAx+ ˆBx ˆ
Kt we once again retrieve the same
partial derivatives for the estimations of αx, βx and κt. For the values of a and the
3
Applying Model Ambiguity
The Lee-Carter and the Li-Lee model apply powerful ways to model and forecast mortality. Unfortunately, these models can still be wrong and random deviations from forecasts can have a huge impact on liabilities. We therefore want to be prudent and price liabilities robustly, such that even in a worst-case scenario we are able to pay off our liabilities. To be able to do this we must pick a sensible set of allowable deviations from the standard model. Balter and Pelsser [1] have studied a pricing method for liabilities in incomplete market settings that take the uncertainty of the base model into account. Since liabilities depending on future mortality are not often traded in, the market for them can certainly be considered incomplete, making the model of Balter and Pelsser very useful for our problem. It is still standard practice for pension funds and insurance companies to model liabilities by computing the expected value of their contracts using mortality figures found with a Lee-Carter or Li-Lee model. We would like to compare the differences between these expected values and the robust value found by the model Balter and Pelsser defined. In this section we discuss the theory of Balter and Pelsser. Their model allows for multiple traded assets and multiple non-traded assets to be part of the evaluation of liabilities. We try to explain their theory as general as possible so that any reader who would like to make interest rates, stock or multiple non-traded assets part of their model will still find this section very useful. However, in our results we focus only on the effect model ambiguity has on mortality estimates, taking only 1 non-traded asset and no stocks or interest rate into account.
3.1
Model Dynamics And Solution
Consider an economic agent who has a liability L(·) he has to pay at time T . The agent seeks to hedge his liability that depends on both hedgeable and unhedgeable risk factors as best as possible. The hedging portfolio A(·) is represented by the positions in the hedgeable risk factors, the traded assets. We assume there is at least 1 non-traded asset, otherwise the market wouldn’t be incomplete. Let there be n tradeable assets X = (X1, . . . , Xn) and l untradeable assets called risk factors
Y = (Y1, . . . , Yl) and let there be a bank account X0 on which the agent can go short
or long for the interest rate r(t, X0(t), X(t)). These three processes are modeled by
the following stochastic differential equation:
d X0(t) X(t) Y (t) = r(t, X0(t), X(t)) µx(t, X(t), Y (t)) µy(t, X(t), Y (t)) dt + 0n 0l ΣXX(t, X(t), Y (t)) ΣXY(t, X(t), Y (t)) ΣY X(t, X(t), Y (t)) ΣY Y(t, X(t), Y (t)) 1/2 dW x(t) Wy(t) . (3)
Here µX(t, X(t), Y (t)) and µY(t, X(t), Y (t)) are respectively Rn and Rl-valued. The
dimension of covariance matrices are Rn×n, Rn×l, Rl×nand Rl×lfor ΣXX(t, X(t), Y (t)), ΣXY(t, X(t), Y (t)), ΣY X(t, X(t), Y (t)) and ΣY Y(t, X(t), Y (t)) respectively. Both
drift vectors are allowed to depend on the tradeable and untradeable assets. WX
and WY are the n-dimensional and l-dimensional Brownian motions under the base model P, respectively. It is assumed that the covariance matrix Σ(t, X(t), Y (t)) is positive definite and thus invertible. The subscript attached to 0n and 0l indicate
the dimension of these zero vectors. The interest rate that drives the bank account is the integral over the short rate (r(t, X0(t), X(t)) that is assumed to be progressively
measurable and which maps to R, so negative interest rates are not excluded. We will the refer to the combined tradeable and non-tradeable assets as Z(t) = (X(t), Y (t))0. Without the bank account, the above SDE reduces to the matrix-vector form
dZ(t) = µ(t, Z(t))dt + Σ(t, Z(t))1/2dWZ(t).
The liabilities L are assumed to be a function of time T and Z(T ) at time T . Now the dynamics of the theory have been introduced we can introduce the objective. Often we find in literature that the agent finds an optimal investment strategy, leading to the (nonrobust) objective
max
θX(t)∈ΘE
P[A(T, X(T )) − L(T, X(T ), Y (T ))]. (4)
where Θ are all admissible trading strategies and P denotes the baseline model. The surplus at terminal time T is a known function. Of course, risk measure P is derived from the baseline model and we allow for other models by changing this risk measure to the risk measure L. We identify an alternative risk measure L by considering the Radon-Nikodym derivative R(t) = dLdP that follows the PDE dR(t) = λ(t, X, Y )R(t)dWZ(t). The change of measure is driven by the parameter λ(t, X, Y ), which is allowed to be stochastic. We use parameter λ to define a limit of the amount of ambiguity, by only allowing values for λ s.t. |λ(t, X, Y )| ≤ k. We call the set off all alternative risk measuresL and by Girsanov’s theorem the set L can be identified by
L ={L ∼ P : R(t) = dL
dP|t, dR(t) = λ(t, X, Y )R(t)dW
Z(t)
with |λ(t, X, Y )| ≤ k}
for some constant k. The width of the set L of alternatives indicates the amount of ambiguity and is represented by this scalar. The larger k, the larger the set of
alternative models and the more uncertain the agent is.
Note that when the worst-case drift of the risky assets is larger than the interest rate, then it is optimal to always invest as much as possible in the risky asset. In this theoretical setting this would lead to an infinite amount being invested in the risky asset so the worst-case drift must from the baseline model must be smaller than the interest rate. To assure this we put a restriction on k with assumption 1, which incorporates the interest rate in the uncertainty in such a way that the problem always stays bounded.
Assumption 1 We assume that the radius on the uncertainty setL is large enough such that it includes interest rate
k2 − (µX(t, Z) − r(, X0, X))0(ΣXX(t, Z))−1(µX(t, Z)) − r(t, X0, X)) ≥ 0
This assumption ensures that k2 is larger than the market price of risk, which means
that the worst-case scenario always gives a negative market price of risk. In other words, the investor would have been better of by not investing in the risky assets. This way factor k still is allowed to vary, which means that different agents can have different degrees of uncertainty, i.e. a different value of k, as long as assumption 1 is not violated. Now we have a set of measures L that lead to different prices of the liability. We now want to price our liability in such a way that we are always certain we end up with a profit, no matter what measure turns out to be the correct measure. To ensure this we assume the worst measure turns out to be the correct one and optimize for that scenario. This gives us the robust objective:
max
θX(t)∈ΘminL∈LE
L[A(T, X(T )) − L(T, X(T ), Y (T ))]. (5)
Here θ(t) is the investment strategy at time t. This strategy equals θ(t) := (θ1(t), . . . , θn(t),
01, 02, . . . , 0l) which defines the amount invested in the traded assets X1 up to Xn.
The restriction on the unhedgeable part is often represented as Ξθ(t) = 0l where
Ξ = [0l×n|IPl × l] is an [l × (n + l)] matrix with to the left the [l × n] dimensional
zero-matrix and on the right the l-dimensional identity matrix.
The robust objective can be interpreted as the objective that makes the agent look for an investment strategy that is least effected under plausible scenarios that emerge from the baseline model. To robustify the investment strategy, the inner part of the optimisation is play by what Balter and Pelsser [1] refer to as a ”malevolent mother nature” who minimises the surplus by choosing the worst-case measure L∗. The agent then searches for the strategy that maximises his surpluss at time T given the choice of mother nature.
Balter and Pelsser prove in their paper [1] that, given the model introduced above, the indifference price π(t, X, Y ) that solves the objective (5) is given by the partial differential equation ∂tπ + ∂X0 π · r · X + ∂ 0 Yπ · (µ Y − ΣY X(ΣXX)−1 qX)+ c q ∂0 Yπ · (ΣY Y − ΣY X(ΣXX)−1)ΣXY) · ∂Yπ+ (6) 1 2tr(∂XXπΣ XX + 2∂ XYπΣXY + ∂Y YπΣY Y) − r · π = 0)
with terminal value π(T, X, Y ) = L(T, X, Y ), where qX = µ(t, Z(t))−r(t, X0(t), X(t))Z(t)
and c =pk2− (qX)0(ΣXX)−1qX. Proving that this is the solution to the robust
ob-jective is done by first proving that the solution to the robust obob-jective exists and is unique. Once this has been established a single solution must be found, which then must be the only solution. Balter and Pelsser [1] provide both proofs in their paper but on a level that can still be too high for an inexperienced reader. Proving that the solution exists and is unique requires the use of mathematics that are beyond the scope of this thesis. However, we do try to provide a more insightful proof that, if there exists a unique solution to the robust objective, it is found by solving partial differential equation (6). This proof is provided in Appendix C.
3.2
Bounding Model Ambiguity
What is still left for us to do is setting a bound on the model ambiguity by setting a value for k. Balter and Pelsser discuss this problem as well in their paper on am-biguity bounds [2]. One of their approaches considers a value of k that allows for all alternative models that, when tested against the base model, would have a type I error of 5% or less and a type error of 20% or less. We will use this approach as well and show how to retrieve it in this section.
Assume once again that we use a model that can be described by diffusion pro-cesses. This means that we are considering stochastic processes X that are described by stochastic differential equations of the form
dX(t, ω) = µ(t, ω)dt + σ(t, ω)dW (t, ω). (7) Let (Ω, F , P) be a probability space where P defines the baseline model. For the specification of possible alternative models, we consider the Brownian motion with a stochastic drift dW (t, ω) + λ(t, ω)dt. Denote both the alternative model and the alternative probability measure L. The likelihood ratio dLdP|T, which is based on the
information over the interval [0, T ] is given by the random value LR(T , ω) at time T and can be characterised by the stochastic differential equation
dL(t, ω) = λ(t, ω)L(t, ω)dWP(t, ω) over [0, T ] (8)
with initial condition L(0) = 1. This random variable LR(T , ω) can be represented by LR(T , ω) = exp{−1 2 Z T 0 λ(t, ω)2dt + Z T 0 λ(t, ω)dWP(t, ω)}. (9)
Hence, the value at time T of the likelihood ratio LR(T , ω) is completely determined by the realistion ω of a path of the Brownian motion {WP(t, ω)}
o≤t≤T and the process
{λ(t, ω)}0≤t≤T that characterises the alternative model L along this path.
Based on the realised path of the Brownian motion we could test at time T if model P should be rejected in favour of model L for the specific process λt, ω. We are testing two simple hypotheses, H0 : P versus HA : L and the Neyman-Pearson
Lemma tells us that the most powerful test between two simple hypotheses is a likelihood ratio test. This means we reject the model P if LR(T , ω) is larger than the critical value ξ. It is crucial to realise that before time T the test statistic is a random variable and we can impose probabilities on the erros.
To implement the test procedure, we derive the critical value ξ from the equation P[LR(T , ω) ≥ ξ] = α (10) We set the critical value such that the probability of incorrectly rejecting model P is equal to α. This is known as the Type I error. Usually, α is set at 5%.
Another error that can occur is the Type II error: the error of incorrectly rejecting model L. This probability is denoted by β and can be computed as
L[LR(T , ω) < ξ] = β (11) A typical value for β is 20%. The complement of the Type II error is the probability of accepting model L when L is the true model and is equal to 1 − β.
Now consider once again (9). Here we can see that log(LR(T , ω)) is normally dis-tributed with a variance ofR0T λ(t)2dt. The likelihood ratio test procedure LR(T , ω) > ξ is therefor equivalent to performing a test on the statistic
lr(T , ω) = Z T
0
λ(t)dWP(t, ω) (12)
which is more appealing to work with. The test statistic lr(T , ω) has a normal distribution with mean 0 and variance R0T λ(t)2dt under the baseline model P. The
hypothesis P is rejected if lr(T , ω) ≥ ζ. Under the alternative model L the test statistic has a normal distribution with mean R0T λ(t)2dt and variance RT
0 λ(t) 2dt.
The power of the test statistic can be computed explicitly as
L[lr(T , ω) > ζ] = Φ(Φ−1(α) + ( Z T
0
λ(t)2dt)12) (13)
leading to the following restriction for λ:
( Z T
0
λ(t)2dt)12 ≤ Φ−1(1 − β) − Φ−1(α) (14)
If we now take α = 5% and β = 20% we find (R0T λ(t)2dt ≤ 2.48). Now we can look back at (1) and see that we required |λ(t, X, Y )| to be smaller that k at all times. We can see that choosing the value 2.48√
T for k will give us that (
RT
0 λ(t)
2dt)12 ≤
(R0T k2dt)12 = (T · (2.48√
T )
2)12 = 2.48 which is in line with our boundaries for the
drift term λ given the desired Type I and II errors for hypothetical testing of the alternative models. We therefor choose this value of k = 2.48√
4
The life insurance contract
4.1
Specifying the contract
In this section we will try to show the effect of using an incomplete market model instead of standard mortality calculations to calculate liabilities. We will first show the estimated liabilities based on the Lee-Carter and Li-Lee model. Then we show how to base the incomplete market calculations on the Li-Lee model and see how the liabilities increase. Lastly, we will investigate what difference using either the Li-Lee or Lee-Carter model makes in this regard.
To show this we will analyze the liabilities of a life insurer who insures 100 male participants and an insurer who insures 100 female participants, all of the same age Astart. We assume all contracts start on January 1st of 2015 so that all estimates
depend solely on estimated mortality figures, not realized mortality figures. The contract pays out a value of 1 to all the surviving participants at a later age Aend.
We denote the time that has passed since 2015 in years as t and we denote the time all participants reach the age Aend as time t = T .
We keep track of the number of survivors over time, denoted for respectively female and male participants as Nf emale(t) and Nmale(t). Using this notation we
have Nf emale(0) = Nmale(0) = 100. We will try and calculate the expected value of
the liabilities using the mortality figures given by the Lee-Carter and Li-Lee model and the valuation based on model ambiguity.
4.2
Applying the Li-Lee and Lee-Carter model
Recall that the Li-Lee and Lee-Carter model estimate the central death rates µx,t
for female and male individuals. This death rate relates to the survival chance S(t) up to time t by the relation S(t) = eR0tµ(t)dt. We calculate this chance of survival
for every year, which gives us the expected number of survivors when the contract expires. The algorithm we use for this is given on the following page. The expected remaining number of survivors give us our liability.
Algorithm 1 Li-Lee expectations startage ← Astart
pensionage ← Aend
Nf emale ← 100
Nmale ← 100
for i = 0 to i = pensionage - startage - 1 do x = startage + i
t = 2015 + i
Sf emale = exp(µx,t,f emale)
Smale = exp(µx,t,male)
Nf emale = Nf emale∗ Sf emale
Nmale= Nmale∗ Smale
Liabilityf emale ← Nf emale
Liabilitymale ← Nmale
4.3
Applying model ambiguity
Before we can apply the theory of Balter and Pelsser we must first define our market dynamics. Since we are solely interested in the mortality dynamics and not in market dynamics we will ignore any tradeable assets and set the interest rate at 0. This leaves us with a decision to make on the dynamics of the untradeable assets. We could choose the Li-Lee model for this purpose, which has static variables Ax, Bx, αx
and βx and stochastic variables Kt and κt. κt models the gradual decrease of the
Dutch deviation from the European trend and goes to zero over time. The more interesting stochastic variable is Kt which measures the decrease in mortality over
time. This K(t) variable is the most interesting part in the Li-Lee model we would preferably find a solution to the incomplete market defined by
dK(t) = νdt + βdW (t) π(T ) =
T −1
Y
i=0
exp(−µxstart+i,tstart+i)
where µx,t = exp(A(x) + B(x) ∗ K(t) + α(x) + β(x) ∗ κ(t))
Where all the remaining values of A(x), B(x), α(x), β(x) and κ(t) are based on the estimated values of the Li-Lee model. A similar model can be made for the Lee-Carter model using only the values α(x), β(x) and κ(t) and κ(t) would be modeled as a
random walk with drift . Unfortunately these systems are not solvable directly, which makes it a poor choice for a model in our case. Therefor we will try to find a simpler model to define the incomplete market with. To do this, note that the error term in K(t) is normally distributed and see that log(µx,t) depends linearly on K(t). This
means that the central death rates are approximately log-normally distributed. Since the chance of dying relates to the central death rate as Pdeath = 1 − exp(−µx,t) ≈ µ
for small values of µ we see the the chances of dying in a single year, and therefor the deaths in a single year, are approximately log-normally distributed. Since the product of log-normal variable is also a log-normal variable we can say that approximately all deaths over time are log-normally distributed. This gives us an easily solvable model which we will use, however we do have to check how good our approximation is. This and more is done in Appendix B. Denote D(t) as the amount of deaths that have occurred up to time t. We then define our model as follows:
dD(t) = νD(t)dt + βD(t)dW (t) (15) π(T ) = (N (0) − D(T ))
This system is solvable, but has some drawbacks. One large drawback is that now theoretically the amount of deaths can decrease over time due to the Brownian motion. To assure this does not happen β needs to be a small enough number. We will come back on this drawback when discussing the combined results, for now we focus on the solution of this model.
We now apply the main result of section 3, we can solve the sytem by solving the PDE (6). This PDE is now given by
∂tπ + ∂D0 π ∗ (νD(t)) + c q ∂0 Dπ ∗ (β2D(t)2)∂Dπ +1 2(∂DDπβ 2D(t)2) = 0 where c = k.
This PDE is solved by π(t) = N (0) − D(t) ∗ e((ν−β∗k)∗(T −t)), which is shown in
Appendix D. We could make the constant N (0) the number of individuals at the start N (0) = 100 and have no initial deaths D(0) = 0 for both female and male participants. This leads to the the problem that at time 0 now no deaths have occured, which gives us π(t) = π = 100 for both the male and female results. Of course this is not a very useful approach, some deaths need to have occured at time t = 0. To assure this we use the following logic. Using our mortality figures we can calculate the expected number of survivors at a later time. Now we flip this
idea around. Knowing the number of survivors at time t = 0, we can calculate the expected population at an earlier time! We now must decide on what earlier time we want to pick. Note that our model depends on the assumption that (µx,t)
increases log-linearly. This assumption is true from about the age of 20 as is shown in Appendix A. We therefor choose for this earlier time the time that the participants were aged 20. Using this theoretical population as our constant C and setting D(0) as the expected amount of deaths between now and the time when our participants were aged 20, we have defined our model. The effect of this model choice on our log-linearity assumption is also shown in Appendix B.
To find estimates for our values of ν and β we once again model the log(Deaths) as a random walk with drift. Notice that we defined the amount of deaths that occured over time as a log-normal process dD(t) = νD(t)dt + βD(t)dW (t), using Ito’s lemma we can then find that the log(Deaths) follow the process dlog(D(t)) = (ν −12β2)dt + βdW (t). We recognize this a random walk with drift and use the rwf () function in R to find estimators for these values νlog = (ν − 12β2) and βlog = β and
find the values used for the model ambiguity calcultions by defining ˆν = νlog+12βlog2
and ˆβ = βlog. This procedure is summarized by Algorithm 2.
4.4
Results
Now we have discussed how the results are obtained we can analyze more thoroughly the effect of incomplete market dynamics on our liabilities. To do this we focus on a four situations: (Astart, Aend) ∈ {30, 50} × {65, 71}. This shows us the effects of
different starting and pension ages, but also the effect of the uncertainty factor k =
2.48√
T , as this depended on the time horizon. We show the density distribution of the
amount of deaths at time t = T and, based on the Li-Lee model and the Lee-Carter model and show where the incomplete market liabilities lie on this distribution.
Algorithm 2 Incomplete market expectation startage ← Astart pensionage ← Aend Nf emale ← 100 Nmale ← 100 Nf emalestart ← 100 Nstart male ← 100
survivorsf emale ← rep(100, pensionage − startage + 2)
survivorsmale ← rep(100, pensionage − startage + 2)
k ← 2.48/(Aend− Astart)
for i = 20 to i = startage - 1 do x = i
t = 2015 − startage + i Sf emale = exp(µx,t,f emale)
Smale = exp(µx,t,male)
Nstart
f emale = Nf emalestart /Sf emale
Nstart
male = Nmalestart/Smale
for i = 0 to i = pensionage - startage - 1 do x = startage + i
t = 2015 + i
Sf emale = exp(µx,t,f emale)
Smale = exp(µx,t,male)
Nf emale = Nf emale∗ Sf emale
Nmale= Nmale∗ Smale
survivorsf emale[i + 1] ← Nf emalestart − Nf emale
survivorsmale[i + 1] ← Nmalestart− Nmale
rwf.deaths.f emale ← rwf ((log(survivorsf emale)), drif t = T RU E)
rwf.deaths.male ← rwf ((log(survivorsmale)), drif t = T RU E)
νf emale,log ← as.numeric(rwf.deaths.f emale$model$drif t)
βf emale ← rwf.deaths.f emale$model$sd
νf emale ← νf emale,log+ β
2
2
liabilityf emale ← Nf emalestart −(Nf emalestart −100)∗exp(νf emale+βf emale∗k)∗(Astart−Aend)
νmale,log ← as.numeric(rwf.deaths.male$model$drif t)
βmale ← rwf.deaths.male$model$sd
νmale← νmale,log +β
2
2
The liabilities derived from model ambiguity, based on the figures from the Li-Lee model, are as follows:
Astart Aend νf emale βf emale Lf emale νmale βmale Lmale
30 65 0.0867 0.0118 95.6 0.0763 0.00614 93.6 30 71 0.0833 0.0138 93.4 0.0753 0.00618 90.1 50 65 0.0749 0.00690 95.2 0.0747 0.00188 93.0 50 71 0.0718 0.00777 91.9 0.0735 0.00270 87.6
The results based on the Lee-Carter model are as follows: Astart Aend νf emale βf emale Lf emale νmale βmale Lmale
30 65 0.0913 0.00898 94.6 0.0737 0.00658 94.5 30 71 0.0869 0.0138 92.2 0.0742 0.00627 91.0 50 65 0.0794 0.00591 94.7 0.0707 0.00114 93.7 50 71 0.0751 0.00889 91.2 0.0715 0.00192 88.5
Note that females tend to live longer than males according to our estimates. This is a result that is often found in mortality anlysis, so we are happy that we reproduce it here. Furthermore we see that the results where Aend = 71 are substantially lower
that the results of Aend = 65, given that Astart is the same. An obvious result, but it
shows that the later ages are very important for the calculations of liabilities. What is less obvious is that the liabilities where the Astart = 30 are about equal or even
higher than the liabilities where Astart = 50, where Aend is the same. This means
that the people that reach older ages 20 years later will die substantially less. Even so much less that the all deaths from age 30 to 71 combined are still less than the deaths that occur for people who are now of age 50 until they reach the same age of 71. This once again shows the large impact decreasing mortality rates have on the liabilities of pension funds.
Where the models are concerned we see that the male estimates of the Lee-Carter model show a lower mortality and the female estimates a higher mortality than the estimates of the Li-Lee model, but these differences are not very large. The observations made in the previous paragraph are true for both the Lee-Carter and the Li-Lee model. As far as model ambiguity is concerned the Lee-Carter model and the Li-Lee model do not seem to differ too much.
Next we plot the density distributions for the number of survivors at time t = T based on the Li-Lee (black) and the Lee-Carter (brown). This nicely shows the dif-ference between the two models for our case. We also plot both the values based
on model ambiguity from the incomplete market model from Balter and Pelsser. To show the effect of ambiguity we take two different values for k, klow = 0 and
khigh = 2 ∗ k. We see clearly that using a factor k of zero retrieves the mean value
which is what we expect if we have no uncertainty. A factor k finds an even more robust value that is further away from the mean. The results are given on the next page, female results left and male results right. We once again notice that the Lee-Carter estimates show a lower mortality for the males (and thus a higher liability), but a higher mortality for the female participants as became clear from the tables presented before.
We mention as well that these values of ν and β are based on the expected figures for future mortality, which do not take the uncertainty of K(t) or κt into account.
However, this uncertainty does not have a large effect on the estimates of β and ν, the standard deviation β is more influenced by how well the expected deaths fit the geometic brownian motion (15). To show this we have included appendix E that show the results if simulations were used to find the values of ν and β.
Figure 1: Upper: liabilities of female and male participants together with the density distribution of the amount of deaths at time t = T with a startage of 30 and pensionage of 65. Lower: startage of 30 and pensionage of 71
IME lines give the value of the liabilities at time t = 0 based on the incom-plete market model using model ambiguity.
Figure 2: Upper: liabilities of female and male participants together with the density distribution of the amount of deaths at time t = T with a startage of 50 and pensionage of 65. Lower: startage of and pensionage of 71
IME lines give the value of the liabilities at time t = 0 based on the incom-plete market model using model ambiguity.
5
Conclusion
The Lee-Carter and Li-Lee model both provide powerful ways to forecast mortality. Using publicly available data we have been able to fit these models and find results in very similar ways. This method retrieved the well known results that males have a slightly higher mortality rate than females and that there is a slight increase at young adolescent ages in the central death rate. The assumptions that mortality tends to decline log-linearly and that the deviation of Dutch data compared to European data tends to converge over time have also been checked and seem to be a good fit.
To apply model ambiguity in calculations where these mortality figures are con-cerned, we have looked at the model Balter and Pelsser [1] introduce for calculations in an imcomplete market based on model ambiguity. This theory allows multiple traded assets and an interest rate to be part of the calculations, but we have only concerned ourselves with the effect model ambiguity has on mortality figures. This can be done by modeling the time dependent factor κtand Kt of the Lee-Carter and
Li-Lee model and having your liabilities depend on the central death rate µx,t that
these models estimate. Multiple non-traded assets can be included in the model of Balter and Pelsser, so in the future the factor κt of the Li-Lee model, that estimates
the Dutch deviation of the decline in mortality compared to other European coun-tries, could be included as well. However, these settings do not provide analytical solutions. To be able to find analytical solutions we have tried to model the total amount of deaths that occur over time as a geometric Brownian Motion, which gives us a solution to the for the liabilities at time t based on the theory of Balter and Pelsser [1]: π(t) = N (0) − D(t) ∗ e(ν−β∗k)∗(T −t). This method has the downside that
theoretically the number of deaths over time can decrease, but the standard devia-tion is small enough to ensure that the chance of that happening is negligible. We have used the expected number of deaths based on the Lee-Carter and the Li-Lee model to find values for these parameters as well as simulations, and have concluded that the difference is not too big. The uncertainty coming from the variance β2 depends mostly on how well the log-normality assumption fits the expected number of deaths than on the uncertainty factor of K(t) or κt. This shows that modeling the
deaths that occur over time as a geometric brownian motion does come with some downsides, but it does allow for analytical results, showing the power of the theory that Balter and Pelsser [1] derive.
The define bounds for the uncertainty Balter and Pelsser [2] find a set of plausible alternatives to the base model, which is taken into account by the factor k. We have
chosen the value of k that allows for all alternatives models that, if tested against the base model, have a significance level α of 5% and a probablity of type II error of β = 20%. This way k gets the value of k = 2.48T where T is the time horizon. This is just one possible choice for k and we try to show the significance of this choice by showing values based on different values of k as well. Setting the value of k to 0 gets you the expected liability based on the Lee-Carter or Li-Lee model and having a k twice as big finds a more robust value, about twice as far away from the mean value based on the mortality figures as the standard value of k. We hope this thesis has make the effect of model ambiguity and its theory more clear. Perhaps in the future the effect of model ambiguity can be studied in a more realistic study with stock and interest rates. Perhaps an approximate solution to the PDE (6) can be found with the incomplete market model that truly resembles the Li-Lee or Lee-Carter model, showing the effect model ambiguity has on trading strategies, which were left out in this thesis but can definitely be included in the theory Balter and Pelsser provide.
6
References
References
[1] Anne G. Balter and Antoon Pelsser [2015] ’Pricing and hedging in incomplete markets with model uncertainty’ Netspar
[2] Anne G. Balter and Antoon Pelsser [2015] ’Quantifying ambiguity bounds through hypothetical statistical testing’ Netspar
[3] Biagini S, Pinar MC¸ [2017] ’The robust Merton problem of an ambiguity averse investor’ Mathematics and Financial Economics 11(1):1-24
[4] Bordigone G, Maoussi A, Schweizer M [2007] ’A stochastic control approach to a robust utility maximization problem’ Stochastic Analysis and Applications, 125-151 (Springer)
[5] Hern´andez-Hern´andez D, Schied A [2007] ’Robust utility maximization in a stochastic factor model’ Statistics & Decisions 24(1):109-125
[6] Kendall E. Atkinson [1989] ’An Introduction to Numerical Analysis’ John Wiley & Sons, Inc, ISBN 0-471-62489-6
[7] Lee, R.D., and L. Carter [1992], ’Modeling and Forecasting the Time Series of U.S. Mortality,’ Journal of the American Statistical Association, 87, 659-671 [8] Nan Li and Ronald Lee [2005] ’Coherent mortality forecasts for a group of
pop-ulations: an extension of the lee-carter method’ Demography, 42(3): 575-594 [9] Mania M, Tevzadze R [2008] ’Backward stochastic differential equations and their
applications’ Journal of Mathematical Sciences 153(3):291-380
[10] Matoussi A, Possama¨ı D, Zhou C [2015] ’Robust utility maximization in non-dominated models with 2BSDE: The uncertain volatility model’ Mathematical Finance 25(2):258-287
[11] Y. Miyahara [2011] ’Option pricing in incomplete markets, modeling based on geometric L´evy processes and minimal entropy martingale measures’ World Sci-entific
[12] Royal Dutch Actuarial Association [2016] ’Projection Table AG 2016’ Royal Dutch Actuarial Association
A
Lee-Carter and Li-Lee results
In Section 2.2 and 2.4 we discuss how mortality figures can be found using data from the Human Mortality Database and fitting the Lee-Carter and Li-Lee model to this data. The results of this process are given in this appendix.
A.1
The Lee-Carter Model
We find the values of αx, βx and κt using the following algorithm that summarizes
the theory that was discussed in Section 2.2. Algorithm 3 Lee-Carter fit
αx ← rep(0, 91) βx ← rep(1/91, 91) κt ← (22 : −22) while true do µx,t ← exp( ˆαx+ ˆβx∗ ˆκt) LL ←P x,t(Dx,t∗ log(µx,t) − Ex,t∗ µx,t)
fit new values: ˆ αx ← ˆαx− ( P t(Dx,t− Ex,t∗ µx,t))/ P t(−Ex,t∗ µx,t) µx,t ← exp( ˆαx+ ˆβx∗ ˆκt) applyconstraints() ˆ βx← ˆβx− (Pt(ˆκt∗ (Dx,t− Ex,t∗ µx,t))/Pt(ˆκ2t ∗ (−Ex,t∗ µx,t)) µx,t ← exp( ˆαx+ ˆβx∗ ˆκt) applyconstraints() ˆ κt← ˆκt− (Px( ˆβx∗ (Dx,t− Ex,t∗ µx,t))/Px( ˆβx2∗ (−Ex,t∗ µx,t)) µx,t ← exp( ˆαx+ ˆβx∗ ˆκt) applyconstraints() check improvements: LL2 ←P x,t(Dx,t∗ log(µx,t) − Ex,t∗ µx,t) DLL ← LL − LL2 if DLL < 10−12 then Break else loop
This algorithm gives the following estimates for αx, βx and κt:
Figure 3: Upper: The female estimates for αx, βx and κt in the Lee-Carter model.
Lower: The male estimates for αx, βx and κt in the Lee-Carter model.
The Lee-Carter model expects a linear downward trend in the estimates for κt
which we retrieve in our results. Combining these estimates finds the central death rates µx,t which we show on the next page. We note that the bump in the central
death rate around adolescent ages is higher for males and that the male central death rates reach higher values at the older ages.
Figure 4: Female central death rate progression over time for the years 1970, 1980, 1990, 2000 and 2010
Figure 5: Male central death rate progression over time for the years 1970, 1980, 1990, 2000 and 2010
The Li-Lee model is fitted in two parts. The first part of the fitting process is fitting the Ax, Bxand Ktvalues. This is done by taking all the european data and iteratively
improving the fitted values ˆAx, ˆBxand ˆKtin the same way as the estimates for αx, βx
and κtwere fitted in the Lee-Carter algorithm (3). This procedure gives the following
estimates for Ax, Bx and Kt:
Figure 6: Upper: The female estimates for Ax, Bx and Kt in the Li-Lee model.
Lower: The male estimates for Ax, Bx and Kt in the Li-Lee model.
Again we find a linear downward trend in the estimates for Ktwhich is in accordance
Next we multiply the Dutch exposures with µx,t and use these new exposures to fit
ˆ
αx, ˆβx and ˆκt in the same way. Now the results are as follows:
Figure 7: Upper: The female estimates for αx, βx and κtin the Li-Lee model. Lower:
The male estimates for αx, βx and κt in the Li-Lee model.
For these estimates of κt we wanted to fit a autoregressive model that converges to
0 over time. We indeed confirm in these results that the values of κt tend towards 0
as time moves forward. Combining these results we find the following central death rates for the Li-Lee model that are displayed on the following page. We again note that the bump in the central death rates at adolescent ages are higher for male estimates than female estimates and that the male estimates reach a higher central death rate than the female estimates.
Figure 8: Female central death rate progression over time for the years 1970, 1980, 1990, 2000 and 2010
Figure 9: Male central death rate progression over time for the years 1970, 1980, 1990, 2000 and 2010
B
Log-normal deaths assumption
As discussed in section 4, we assume the deaths occured up to time t are log-normally distributed. We can assume the amount of deaths are 0 at time t = 0. We could also assume the amount of deaths are equal to a number that is consistent with a group of people whose deaths were tracked from when they were aged 20. These assumption lead to different estimates for ν and β. Here we show the progress of the log(deaths) over time and see that the second assumption leads to the best fit. The log(deaths) under the both assumptions are given on the next pages, first for the Lee-Carter model and then for the Li-Lee model. We see in the plots where 0 deaths are assumed to have occured at time t = 0 that the amount of deaths seem to grow relatively quicker than the model would explain. This is caused by the initial amount of deaths being equal to zero, in which case no growth should be expected at all and thus the growth that we do find seems to be relatively faster than a log-normal model explains and then start increasing linearly. We therefore imagine a different amount of clients having entered the contract some time earlier which would allow for a small number of deaths already being present. For this number we choose the expected amount of deaths, had the clients entered at the age of twenty, that is in line with an exact amount of 100 clients left at the age the contract starts for female and male participants. The procedure that gets this number of clients is part of Algorithm 2. We choose the age of 20 because from this age onwards we see in appendix A that the central death rates increase log-linearly and this is in line with the model choice for the incomplete market dynamics. Under this assumption we would expect the log(deaths) to increase linearly from the start, which is indeed the case. The female results do seem to be slighly underestimated by the model, but definitely show an increased fit as well.
Using the theoretical deaths at time t = 0 we use log of the expected deaths to fit the random walk with drift. We plot the fit of the random walk with drift we obtain in these figures as well.
Figure 10: Left: female log(deaths) over time for periods 30-65, 30-71, 50-65 and 50-71, starting at 0 deaths based on the Lee-Carter model. Right: The male results under the same assumptions.
Figure 11: Left: female log(deaths) over time for periods 30-65, 30-71, 50-65 and 50-71, starting with a positive amount of deaths based on the Lee-Carter model. Right: The male results under the same assumptions.
Figure 12: Left: female log(deaths) over time for periods 30-65, 30-71, 50-65 and 50-71, starting at 0 deaths based on the Li-Lee model. Right: The male results under the same assumptions.
Figure 13: Left: female log(deaths) over time for periods 30-65, 30-71, 50-65 and 50-71, starting with a positive amount of deaths based on the Li-Lee model. Right: The male results under the same assumptions.
C
Solution of the robust objective
In Section 3.2 we stated that the robust objective max
θX(t)∈ΘminL∈LE
L[A(T, X(T )) − L(T, X(T ), Y (T ))]
is solved by the solution of the partial differential equation ∂tπ + ∂X0 π · r · X + ∂ 0 Yπ · (µ Y − ΣY X (ΣXX)−1qX)+ c q ∂Y0 π · (ΣY Y − ΣY X(ΣXX)−1)ΣXY) · ∂ Yπ+ 1 2tr(∂XXπΣ XX+ 2∂ XYπΣXY + ∂Y YπΣY Y) − r · π = 0).
In this appendix we try to provide some insight into how this result can be achieved. To be clear, we do not provide a full proof of the statement, but show only that if the robust objective has a unique solution, it is the solution to the PDE.
As a first step we link an expectation to a solution to a partial differential equation via the Feynman-Ka¸c equation (FK). This equation states that a partial differential equation of N dimensions of the form
∂u ∂t(x, t)+ N X i=1 µi(x, t) ∂u ∂x(x, t)+ 1 2 N X i=1 N X j=1 γij(x, t) ∂2u ∂xi∂xj (x, t)−r(x, t)u(x, t)+f (x, t) = 0
with terminal condition u(x, T ) = ψ(x) and γij(x, t) = PNk=1σik(x, t)σjk(x, t) is
solved by the conditional expectation
u(x, t) = EQ[ Z T t e−RtTr(Xτ,τ )dτf (X r, r)dr + e− RT t r(Xτ,τ )dτψ(X T)|Xt = x]
To apply this to our situation we define ˜A(t, X(t)) = A(t,X(t))X
0(t) and ˜L(t, X(t)) =
L(t,X(t))
X0(t) to ignore any terms depending on r and find that the conditional expectation
EL[ ˜A(T, X(T )) − ˜L(T, Z(T ))| ˜A(t) = ˜A, X(t) = X, Y (t) = Y ] = ¯v(t, ˜A, X, Y ), where the measure L represents the measure (dWX(t) + X(t, Z)dt, dWY(t, Z) +
Y(t, Z)dt) satisfies the PDE
∂t¯v + ˜θX(t)0(qX(t, Z) + X(t, Z))∂A˜v + (µ¯ X(T, Z) + X(t, Z))0∂X¯v+ (µY(T, Z) + Y(t, Z))0∂Yv +¯ 1 2tr(∆( ˜A, Z) · ¯v(t, ˜A, Z) · ΦΣ(t, Z)Φ 0 ) = 0
where the operator ∆ givess the partial derivatives and the operator Φ extends the covariance matrix Σ(t, Z) with the covariance matrices with respect to ˜A(t, X)
∆( ˜A, Z) = ∂A ˜˜A ∂AX˜ ∂AY˜ ∂X ˜A ∂XX ∂XY ∂Y ˜A ∂Y X ∂Y Y , Φ = [˜θX 0l]0 I[n+l] (16)
and ∂ is used for the partial derivative operator. Now see that we have defined ¯
v(t, ˜A, X, Y ) as our expected surplus which gives us the new robust objective max
˜ θ(t)∈Θ
min
(t,Z)∈Lv(t, ˜¯ A, X, Y ) = v(t, ˜A, X, Y ). (17)
This robust optimised value is given by the new PDE equation, called the Hamilton-Jacobi-Bellman (HJB) equation, ∂tv + max ˜ θ(t) min (t,Z)[˜θ X(t)0 (qX(t, Z) + X(t, Z))∂A˜v + ∂Z0 v(µ(t, Z) + (t, Z))+ 1 2∂A ˜˜Av ˜θ X (t)0ΣXX(t, Z)˜θX(t) + ∂AX0˜ vΣXX(t, Z)˜θX(t) + ∂AY0˜ vΣY X(t, Z)˜θX(t)]+ 1 2tr(∆(Z) · ν(t, ˜A, Z)Σ(t, Z)) = 0) (18) s.t. (x, Z)0Σ(t, Z)−1(t, Z) ≤ k2 Ξ˜θ = 0l with v(T, ˜A, X, Y ) = ˜A(T, X) − ˜L(T, Y )
where we have written out all elements of the trace that depend on ˜θ in the max-min part of the equation. Now imagine that at some time t we have two situations, with one situation having an asset value of 100 and another one having 101. In the first situation we would of course value our position 1 lower than in the second situation as we simply have 1 value less of assets. This works for any deviation of the assets, meaning that the value of our position is always linear in ˜A. Therefore we can propose that
v(t, ˜A, X, Y ) = ˜A(t, X) − ˜w(t, X, Y )
Now we can fill this into the HJB equation (18) and find that the function ˜w(·) satisfies
− ∂tw + max˜ ˜ θ(t) min (t,Z)[˜θ X(t)0 (qX(t, Z) + X(t, Z)) ∂Z0 w(µ(t, Z) + (t, Z))] −˜ 1 2tr(∆(Z) · ˜w(t, Z)Σ(t, Z)) = 0 s.t. (t, Z)0Σ−1(t, Z) ≤ k2 (19) Ξ˜θ = 0l with ˜w(T, X, Y ) = ˜L(T, Z).
So instead of the initial problem we can solve the problem for an agent who maximizes − ˜w(t, X, Y ) which once again leads to a new robust optimization problem
m(t, X, Y ) = maxθ˜X(t)min(t,Z)[˜θX(t)0(qX(t, Z) + X(t, Z)) − ∂Z0w(µ(t, Z) + (t, Z))].˜ (20)
We will now find the optimal values of ˜θ and and plug these into the HJB equation (19) to get the final result. To do this we solve two the minimisation
min
(t,Z)
˜
θ0 · q + 0(˜θ − ∂Zw)˜ (21)
s.t. 0Σ−1 ≤ k2
this is a linear objective in and can immediately be solved by using a Lagrangian. This results in ? = −k√ Σ(˜θ−∂Zw)˜
(˜θ−∂Zw)˜ 0Σ(˜θ−∂Zw)˜
. However, the full Lagrangian of (21) we find L(˜θ, , λ0, λ) = ˜θ0q + 0(˜θ − ∂Zw) − λ˜ 0 1 2( 0 Σ−1 − k2) − λ0(Ξ˜θ − 0l), (22)
where λ has dimension [l × 1] and λ0 is a scalar. Now the first order conditions are
∂L ∂ ˜θ = q + − Ξ 0 λ = 0n+l ∂L ∂ = −∂Zw + ˜˜ θ − λ0Σ −1 = 0n+l (23) ∂L ∂λ = −Ξ˜θ = 0l ∂L ∂λ0 = −1 2( 0 Σ−1 − k2) = 0. The first of the equations gives us
Plugging this into the second equation gives us the optimal value of ˜θ ˜
θ? = ∂Zw + λ˜ 0Σ−1(Ξ0λ?− q)
Now note that our definition of Ξ gives us that 0l = Ξ(∂Zw + λ˜ 0Σ−1)
and we can use this in the first equation to retrieve ΞΣ−1q + ΞΣ−1 − ΞΣ−1Ξ0λ = 0l
and thus
λ? = (ΞΣ−1Ξ0)−1Ξ(Σ−1q − λ−10 ∂Zw)˜
Now λ? is an expression of λ
0 and ? is and expression of λ? and thus we can solve
the constraint from the fourth equation 0Σ−1 − k2 = 0 which gives us the result λ−10 = ±
s
q0Σ−1Ξ0(ΞΣ−1Ξ0)−1ΞΣ−1q − q0Σ−1q + k2
∂ZwΞ˜ 0(ΞΣ−1Ξ0)−1Ξ∂Zw˜
Since is chosen s.t. the objective is minimised, the second order derivative ∂∂2L2 is
positive so we find that λ0 is the positive square root in the above equation. In the
root we can see the term Ψ = Ξ0(ΞΣ−1Ξ0)−1Ξ twice which we name Ψ. Working this term out first gives the result Ψ = 0[n×n] 0[n×l]
0[l×n] ΣY Y − ΣY X(ΣXX)−1ΣXY
. Call the right-lower matrix S and we get the result for ? in terms of q, Σ and ˜w:
? = Ψ(− s q0Σ−1XΣ−1q − q0Σ−1q + k2 ∂Z0 wX∂˜ Zw˜ ∂Zw + Σ˜ −1q) − q = " −qX S q k2−(qX)0(ΣXX)−1qX ∂Zw(Σ˜ Y Y−ΣY X(ΣXX)−1)ΣXY∂Zw˜∂Zw − Σ˜ Y X(ΣXX)−1qX # (24)
which is found by using the identities ΨΣ−1 = 0[n×n] 0[n×l] −ΣY X(ΣXX)−1 I [l×l] and q0Σ−1ΨΣ−1q − q0Σ−1q = (−(qX)0(ΣXX)−1qX, 0)0. Now we have the θ and for the
max-min problem (19) and if we will them in we find that the solution to (20) satisfies the PDE ∂tw − ∂˜ X0 w · r · X − ∂˜ 0 Yw · (µ˜ Y − ΣY X(ΣXX)−1qX) − 1 2tr(∆(Z) · ˜w · Σ) − c q ∂Y0 w · (Σ˜ Y Y − ΣY X(ΣXX)−1ΣXY) · ∂ Yw = 0.˜ (25)
Recall that ˜w(t, X, Y ) = w(t, X, Y )/X0(t) = e− Rt
0r(s,X0(s),X(s))dsw(t, X, Y ) and since
this is our final result we name this π(t, X, Y ) = w(t, X, Y ). Rephrasing the PDE (25) in terms of π gives ∂tπ − ∂X0 π · r · X − ∂ 0 Yπ · (µ Y − ΣY X(ΣXX)−1 qX) −1 2tr(∂XXπΣ XX+ 2∂ XYπΣXY + ∂Y YπΣY Y) − c q ∂Y0 π · (ΣY Y − ΣY X(ΣXX)−1ΣXY) · ∂ Yπ = 0. (26)
D
Proof of solution
The main result of the paper of Balter and Pelsser [1] is that at time t the value of liablities π(t) can be found by solving the PDE
∂tπ + ∂X0 π · r · X + ∂ 0 Yπ · (µ Y − ΣY X(ΣXX)−1 qX)+ c q ∂Y0 π · (ΣY Y − ΣY X(ΣXX)−1)ΣXY) · ∂ Yπ+ (27) 1 2tr(∂XXπΣ XX+ 2∂ XYπΣXY + ∂Y YπΣY Y) − r · π = 0).
In section 4.3 we wanted to use this results to find a solution to our model
dD(t) = νD(t)dt + βD(t)dW (t) (28) π(T ) = (N (0) − D(T )).
This model does not include any interest rate or tradeable assets. We can identify the drift and variance of the untradeable asset as µY = νD(t) and ΣY Y = β2D(t)2.
This reduces (27) to ∂tπ + ∂D0 π · (νD(t)) + k q ∂D0 π · (β2D(t)2)∂ Dπ + 1 2(∂DDπβ 2D(t)2) = 0
and that we proposed that the solution to this PDE was given by
π(t) = N (0) − D(t) · e((ν − β · k)(T − t)). (29) We first take a look at all the partial derivatives separately:
∂tπ(t) = (ν − β · k)D(t)e((ν−β·k)(T −t)) = (ν − β · k)(N (0) − π(t))
∂Dπ(t) = −e((ν−β·k)(T −t)) = (π(t) − N (0))/D(t)
∂DDπ(t) = 0
Putting this all together in the PDE gives us
(ν−β·k)(N (0)−π(t))+((π(t)−N (0))/D(t))·ν·D(t)+kp((π(t) − N (0))/D(t))2· β2D(t)2
= (ν−β·k)(N (0)−π(t))+((π(t)−N (0))/D(t))·ν·D(t)+k((π(t)−N (0))/D(t))·β·D(t) = (ν − β · k)(N (0) − π(t)) + ((π(t) − N (0))) · ν + k((π(t) − N (0))) · β
= (ν − β · k − ν + β · k)(N (0) − π(t)) = 0
which proves that (29) is indeed a solution to (27), which gives the unique solution of the liabilty at time t.
E
estimating with simulations
In section 4, which discussed the application of model ambiguity on life insurance contracts, we have mentioned that the uncertainty of K(t) and κt does not play a
huge role in estimating the values β and ν that define the model ambiguity setting. To show this we have run a 100 simulations for all the discussed start- and pensionages. These results differ only slightly from the results mentioned in section 4, but they do but do show shortcomings of modeling the deaths over time as a geometric brownian motion.
We start the simulations by simulating values for K(t) and κ for respectively the Li-Lee and Lee-Carter model. The results thereof are given below.
Figure 14: The simulations of the future estimates Kt of the Li-Lee model and the
values κt (ranged yellow to red) and the 95% confidance interval based on the model
