I am supposed to solve the below question using Thieles differential equation, and I am having a hard time formulating Thieles differential equation for this specific insurance case and would appreciate some help in doing this. This is a question in a course called Life Insurance Mathematics.
Determine the single premium to be paid from the time point when a healthy individual turns 50, in order for this to be a fair premium, if the last time point when benefits can be received corresponds to the individuals 65 birthday. Do these calculations using a discretization corresponding to the time units 1 day, 1 month, and 1 year. Report how the results are affected by the discretization, and illustrate how the reserves evolve as a function of time.
The information I have available in order to solve this is the following:
We consider a health insurance with three states, healthy ($2$), sick ($1$) and dead ($0$) with the following transitions being allowed: $2 \leftrightarrow 1$, $1 \to 0$, $2 \to 0$. At time $0$, an annual sickness benefit of $15\%$ of a monthly salary of $30 000$ is paid for a person being sick for one year. These payments are assumed to occur continuously in time, and there is an indexation corresponding to an increase of $1.5\%$ per annum (deterministic). That is, $$b_1(t) = 12 · 30 000 · 0.15 · 1.015^t.$$
We also assume that the mortality is described by the cohort mortality given by the Makeham-function $$\mu(t)=a+bexp\{ct\}$$
with parameters $a=3.5 \cdot 10^{-4}$, $b=7 \cdot 10^{-8}$, and $c=0.157$. Lastly, I also have the parameters $\mu_{12}=50$, $\mu_{21}=2$ and the constant interest rate $r(t)=0.04$. Now, let's continue with what I have done so far.
I know that a fair premium (Actuarially fair insurance) means that the price of this specific health insurance should equal the expected payouts, and my guess is that you find this fair premium by adjusting an initial guess until the Thiele differential equation is satisfied using an R code. So, now it remains to determine the Thiele differential equation so that I can write this code, and this is where I get stuck.
So, the Thiele differential equation for reserves is given by $$\frac{d}{dt}V^J(t)=r(t)V^J(t)-b^J(t)-\sum_{k \neq j} \mu_{jk}^{t}(b^{J^k}(t)+V^k(t)-V^J(t)),$$ so I wrote the following equation for this in R:
# Calculate fair premium using Thiele differential equation for reserves while (t <= t_max) {
mu <- a + b * exp(c * (50 + t))
b_sick <- 12 * 30000 * 0.15 * 1.015^t
# Calculate reserves using Thiele differential equation
dV_total <- r * V_total - b_1 - mu_12 * (V_total - fair_premium_guess) + mu_21 * V_total
V_total <- V_total + dV_total * dt
reserves <- c(reserves, V_total)
time_points <- c(time_points, t)
t <- t + dt }
where V_total is serving as a cumulative measure of reserves, including both benefits and interest, and it is being updated over time, r is the constant interest rate $r(t)$, and I use fair_premium_guess since the fair premium is determined by adjusting an initial guess until the Thiele differential equation is satisfied. Now I wonder, is dV_total <- r * V_total - b_1 - mu_12 * (V_total - fair_premium_guess) + mu_21 * V_total the correct way to express the Thiele differential function for reserves in this specific case of health insurance?
I know that we also can find the mortality via the Makeham-function, and maybe this is the one I should use? So, in summary, I am having a hard time expressing the Thiele differential equation (and not the R code) and would appreciate help with this!
Thank you.