Projet A.V.

Part 2 : Dynamic mortality modeling

1. Fit the model (range of ages: 20 to 100 years old) with the log and logit links and compare results. Analyze the goodness of t and based on this comparison, determine which link log function is the most appropriate.

After fitting the model for the two links, we obtain the Log-Likelihood, the Akaike information criteria and the Bayesian information criteria, which are summarized in the following table :

Model	Loglikelihood	AIC	BIC
Log	-10027.30	20644.40	22344.32
Logit	-10075.3	20740.04	22439.76

The fitting of the LC model has been realized using two different links : log and logit. After evaluation of the quality of the fitting through the aforementioned criteria’s, we conclude that the Log-model is the most appropriate.

Figure 1 : Parameters for the Lee-Carter model fitted with log link to the Belgium Female population for ages 20-100 from 1990 to 2018.

Figure 2 : Parameters for the Lee-Carter model fitted with logit link to the Belgium Female population for ages 20-100 from 1990 to 2018.

2. For the log model, analyze the features of ( \(Kt^1\)), \((Kt^2)\) and select an appropriate model to forecast the time series.

In the LC model extended to two factors :

• \(Kt^1\) represents the overall mortality level in time. Therefore, we expect to observe in almost every country a quasi linear decrease. As for the model selection, a random walk with drift has been shown to provide a reasonable fit, that is,

• \(Kt^2\) , captures an effect of shape by age, modifying the slope or convexity of mortality.

With the simple code provided in the paper “StMoMo: An R Package for Stochastic Mortality Modeling”, we obtain the appropriates ARIMA processes used to model the factors:

Series: k1  ARIMA(2,1,0) with drift   Coefficients:           ar1      ar2    drift       -0.9309  -0.4369  -1.2883 s.e.   0.1717   0.1704   0.0993  sigma^2 = 1.634:  log likelihood = -45.5 AIC=99.01   AICc=100.75   BIC=104.34 Series: k2  ARIMA(0,1,0)   sigma^2 = 0.04366:  log likelihood = 4.11 AIC=-6.22   AICc=-6.06   BIC=-4.89

The factor \(Kt^1\) is modeled as an ARIMA (2,1,0) process with drift reflecting a non-stationary dynamic and a steady decline in the general level of mortality.

The factor \(Kt^2\) is modeled as an ARIMA(0,1,0) process, indicating erratic dynamics without drift or significant AR structure.

We forecast the period indexes \(Kt^1\) and \(Kt^2\) using the ARIMA models selected. The fan charts below show the expected future trajectory and the associated uncertainty.

\(Kt^1\) continues its long-term decreasing trend, indicating ongoing mortality improvements.

\(Kt^2\) remains relatively stable, with increasing uncertainty over time.

This confirms that the selected ARIMA specifications provide reasonable long-term forecasts.

Figure 5 : Forecast of the two period indexes over 50 years

Figure 6 : Temporal evolution of the period factors

Figure 7 : Fanchart for mortality rates qxt at ages x=20,x=40, x=60 andx=80 from the LC model fitted to the Belgium female population 1990-2018

3. For the log model, simulate 1000 mortality curves (or more) for the year 2023 and analyze them (average, standard deviation, confidence interval by ages). Using the average forecast death probabilities in 2023, calculate the life expectancy at ages; 20, 30, 40, 50, 60, and 70 years old. Compare your results with the life expectancy (same ages) computed with the real raw 2023 mortality rates and legal mortality tables.

The package StMoMo provides the function simulate for simulating trajectories from GAPC stochastic mortality models. We use this function to simulate mortality curves for the next five years starting in 2018.

From the simulated 2023 mortality rates \(\{m_x\}\), we construct the corresponding life table using standard actuarial formulas.

Conversion from mortality rates to one-year death probabilities

For each age \(x\), the mortality rate from the model is denoted \(m_x\).
We convert it into a one-year death probability using:

\[ q_x = 1 - e^{-m_x}. \]

This transformation is consistent with the assumption of a constant hazard over each age interval.

Survival function

We fix the radix of the life table at:

\[ l_{x_0} = 1, \]

and compute the number of survivors at age \(x+1\) recursively:

\[ l_{x+1} = l_x(1 - q_x). \]

This yields the survival curve \(\{l_x\}\) for the simulated year 2023.

Person-years lived

For all ages except the last open interval:

\[ L_x = \frac{l_x + l_{x+1}}{2}. \]

For the last age \(\omega\), we use:

\[ L_\omega = \frac{l_\omega}{q_\omega}. \]

This gives the full sequence \(\{L_x\}\) of person-years lived.

Total future lifetime

The cumulated future lifetime is:

\[ T_x = \sum_{k=x}^{\omega} L_k. \]

We compute it using a backward cumulative sum.

Life expectancy

The remaining life expectancy at age \(x\) is finally obtained as:

\[ e_x = \frac{T_x}{l_x}. \]

After collecting the necessary data for the year 2023, we are now able to compare it with the other data in our disposition, summarized in the following table :

Age	model	raw	legal
20	64.69263	64.48	61.30540
30	54.60	54.60	51.66329
40	45.034	44.82	42.11134
50	35.48638	35.29	32.79595
60	26.39998	26.15	23.98544
70	18.02005	17.71	16.10289

Now we can compare our results using the following graphs.

Figure 8 : Graphs of the remaining expectancy according to our simulations.

First, we can see that life expectancy decreases in a quasi-linear trend as age increase. Furthermore, we can also see that our model prediction is consistent with actual data and, unsurprisingly, we can visualize the prudence of legal mortality tables.

Figure 9 : Mortality rates mx(x, 2023) across ages predicted by the model, with 95% PI. Comparison with the legal and raw mortality tables

Here, we can see that :

• the mortality rate increases exponentially with age

• the forecast is more precise at younger ages, more uncertainty at older ages

• the legal table mortality rates constantly lie above the raw and simulated data except from age 90+

  

Figure 10 : histograms of mortality rates simulated by the model for ages 20, 40 and 60

Figure 11 : Boxplots of the simulated death probabilities at given ages

Figure 12: Graph of the ratio mx(model/mx(raw)

Figure 13 : Graph of the standard deviation of the probabilities of death by ages.

In conclusion, the stochastic simulations produced by the LC model for 2023 generate realistic mortality patterns, demonstrating increasing uncertainty with advancing age. The distributions of the simulated mortality rates (mx) and probabilities of death (qx) are narrower at younger ages and widen progressively to reflect the natural increase in variability with age.

Life expectancies computed from the forecasted 2023 mortality rates consistently fall between the raw HMD values and the legal mortality table, which is intentionally conservative, confirming that the model provides reliable projections. The ratio of mx(model)/mx(raw) varies around 1, with greater deviation at younger ages.

Overall, the LC model provides coherent and realistic mortality forecasts for 2023, with uncertainty behaving as expected across the age range.

4. Under the assumption that i = 1.25%, perform at least 1000 simulations in order to price a temporary death insurance: age x = 60 years old in 2018, duration n = 20 years and capital C = 100000. Calculate the standard deviation and the 1%, 5%, 95% and 99% percentile of this price. Explain clearly the method! Compare its value with the one computed using a legal mortality table.

To answer this question, we must our log-model to simulate a 1000 mortality curves and, for each curve, calculate the price of a temporary death insurance.

Price based on the statutory mortality table

We use the previously calibrated log-link model. (fit_log) to simulate future mortality trajectories:

• horizon : \(n = 20\) years (2019–2038)
  
• numbre of simulations : \(N = 1000\)

For each simulation \(k\), we obtain a mortality rate trajectory \(\mu_{x,t}^{(k)}\) for ages 60 to 79.
We convert to annual death probabilities
\[ {}_1q_{60+t-1}^{(k)} = 1 - \exp\big(-\mu_{60+t-1}^{(k)}\big), \quad t = 1,\dots,n. \]

We then create the survival probaibilities
\[ {}_0p_{60}^{(k)} = 1,\qquad {}_{t}p_{60}^{(k)} = \prod_{j=0}^{t-1} \big(1 - {}_1q_{60+j}^{(k)}\big). \]

The pure premium for simulation \(k\) is then: \[ PV^{(k)} = C \sum_{t=1}^{n} v^t\, {}_{t-1}p_{60}^{(k)}\, {}_1q_{60+t-1}^{(k)}. \]

By repeating this for \(k = 1,\dots,N\), we obtain a price distribution from which we calculate:

• the mean,

• the standard deviation,

• the 1%, 5%, 95% and 99% quantiles.

Price based on legal life tables

The statutory mortality table is now used to calculate the pure premium for term life insurance.

Let’s note \(q_x^{\text{legal}}\) the probability of death at age \(x\) according to the statutory table.
For an inssured aged \(60\) years old in 2018 on a span of \(n = 20\) years, we need the probabilities of death at ages \(60,61,\dots,79\).

The annual probabilities of death used in the contract are therefore:\[ {}_1q_{60+t-1}^{\text{legal}} = q_{60+t-1}^{\text{legal}}, \qquad t = 1,\dots,20. \]

Based on these probabilities, survival probabilities are constructed: \[ {}_0p_{60}^{\text{legal}} = 1, \] and, for \(t \ge 1\), \[ {}_{t}p_{60}^{\text{legal}} = \prod_{j=0}^{t-1} \left(1 - {}_1q_{60+j}^{\text{legal}}\right). \]

The actuarial present value of term life insurance based on the statutory table is therefore: \[ PV_{\text{legal}} = C \sum_{t=1}^{20} v^t \, {}_{t-1}p_{60}^{\text{legal}} \, {}_1q_{60+t-1}^{\text{legal}}, \]

Table :

Mean	Sd	Q1%	Q5%	Q95%	Q99%	Legal Price
16535.672	1105.608	14032.404	14698.502	18388.197	19308.307	28185.248