using Dates
using MortalityTables
using DifferentialEquations
using Plots
using ColorSchemes # for Turbo colors, which emphasize readability
using ActuaryUtilities
using LaTeXStrings
Introduction
This demonstrates an example of defining an universal life policy value roll like a discrete differential equation.
It then uses the SciML DifferentialEquations package to “solve” the policy projection given a single point, but also to see how the policy projection behaves under different premium and interest rate conditions.
The first time this is run, it will download and precompile some large packages.
Let’s use the 2001 CSO table as the basis for cost of insurance charges:
cso = MortalityTables.table("2001 CSO Super Preferred Select and Ultimate - Male Nonsmoker, ANB")
MortalityTable (CSO/CET):
Name:
2001 CSO Super Preferred Select and Ultimate - Male Nonsmoker, ANB
Fields:
(:select, :ultimate, :metadata)
Provider:
American Academy of Actuaries
mort.SOA.org ID:
1076
mort.SOA.org link:
https://mort.soa.org/ViewTable.aspx?&TableIdentity=1076
Description:
2001 Commissioners Standard Ordinary (CSO) Super Preferred Select and Ultimate Table – Male Nonsmoker. Basis: Age Nearest Birthday. Minimum Select Age: 0. Maximum Select Age: 99. Minimum Ultimate Age: 16. Maximum Ultimate Age: 120
Next, policy mechanics are coded. It’s essentially a discrete differential equation, so it leverages DifferentialEquations.jl
The projection is coded in the Discrete DifferentialEquation format:
\[u_{n+1} = f(u,p,t_{n+1})\]
In the code below, this translates to:
u
is the state of the system. We will track three variables to represent thestate
:state[1]
is the account valuestate[2]
is the premium paidstate[3]
is the policy duration
p
are the parameters of the system.t
is the time, which will represent days since policy issuance
function policy_projection(state, p, t)
# grab the state from the inputs
av = state[1]
# calculated variables
cur_date = p.issue_date + Day(t)
dur = duration(p.issue_date, cur_date)
att_age = p.issue_age + dur - 1
# lapse if AV <= 0
lapsed = (av <= 0.0) & (t > 1)
if !lapsed
monthly_coi_rate = (1 - (1 - p.mort_assump[att_age])^(1 / 12))
## Periodic Policy elements
# annual events
if Dates.monthday(cur_date) == Dates.monthday(p.issue_date) ||
cur_date == p.issue_date + Day(1) # OR first issue date
premium = p.annual_prem
else
premium = 0.0
end
# monthly_events
if Dates.day(cur_date) == Dates.day(p.issue_date)
coi = max((p.face - av) * monthly_coi_rate, 0.0)
else
coi = 0.0
end
# daily events
int(av) = av * ((1 + p.int_rate)^(1 / 360) - 1.0)
# av
new_av = max(0.0, av - coi + premium + int(av - coi))
# new state
return [new_av, premium, dur] # AV, Prem, Dur
else
# new state
return [0.0, 0.0, dur] # AV, Prem, Dur
end
end
policy_projection (generic function with 1 method)
The following function will create a named tuple of parameters given a varying prem
(premium) and int
(credit rate).
params(prem, int) = (
int_rate=int,
issue_date=Date(2010, 1, 1),
face=1e6,
issue_age=25,
mort_assump=MortalityTables.table("2001 CSO Super Preferred Select and Ultimate - Male Nonsmoker, ANB").ultimate,
projection_years=75,
annual_prem=prem,
)
params (generic function with 1 method)
Runing the system
This results in the following plot. The tracked output variables u1 and u2 represent the two vars that we tracked above: account value and cumulative premium.
begin
p = params(
8000.0, # 8,000 annual premium
0.08 # 8% interest
)
# calculate the number of days to project
projection_end_date = p.issue_date + Year(p.projection_years)
days_to_project = Dates.value(projection_end_date - p.issue_date)
# the [0.0,..] are the initial conditions for the tracked variables
prob = DiscreteProblem(policy_projection, [0.0, 0.0, 0], (0, days_to_project), p)
proj = solve(prob, FunctionMap())
plot(proj)
end