0) Model assumptions

In order to keep our example simple, we make the following assumptions.

• There is no natural catastrophe events risk NE (this can be easily added by means of an additional distribution to consider in the below lines of codes)
• There are only single model segments (modelling retrocession might require multiple model segments, but this is not the focus of this article)
• We ignore discounting (meaning we assume forward values = present values)
• The FINMA StandRe documentation as at Jan 2021 is still applicable
• The parameters of the distributions (AER, AEP, IE1 and IE2) have been estimated (which is the main purpose of the FINMA StandRe template)

1) Import parameters

We first import the values needed to simulate the non-life risk, define the parameters of the t-copula as prescribed by FINMA and set the number of simulations as well as a seed for reproducibility.

# IF ON A MOBILE, TURN THE SCREEN TO LANDSCAPE FOR READABILITY
# Import and define parameters
df_parameters           <- read.table(
                           file   = "standre.txt", 
                           sep    = "\t", 
                           header = TRUE)
rownames(df_parameters) <- df_parameters[,1]
t_rho                   <- 0.23
t_degreefreedom         <- 4.00
N                       <- as.integer(2.0e6)
set.seed(100)

The following table shows our invented parameters of the distributions which we are going to simulate.

2) Simulate IE1 and IE2

Second, we simulate the compound Poisson-Pareto distributions.

# Simulate IE1 and IE2
draw_ie <- data.frame(matrix(nrow = N, ncol = 2))
for(i in 1:2){
 model       <- paste0("IE", toString(i))
 theta       <- (df_parameters["pareto_tail", model] 
                        * df_parameters["threshold", model] 
                        / df_parameters["pareto_initial", model])
 draw_ie[,i] <- actuar::rcompound(N,
                 model.freq = stats::rpois(df_parameters["freq",model]),
                 model.sev  = actuar::rpareto2(
                      min   = df_parameters["threshold",model], 
                      scale = theta, 
                      shape = df_parameters["pareto_tail",model])
                 )
}

3) Simulate AER and AEP

Third, we draw from a 2-dimensional t-copula on which we apply the inverse lognormal CDFs to obtain draws for AER and AEP.

# Simulate AER and AEP
draw_ae           <- data.frame(matrix(nrow = N, ncol = 2))
colnames(draw_ae) <- c("R", "P")
draw_copula       <- copula::rCopula(
                      N, 
                      copula::tCopula(
                      param   = t_rho,
                      dim     = 2, 
                      df      = t_degreefreedom))
draw_ae[,"R"]     <- stats::qlnorm(
                      draw_copula[,1], 
                      meanlog = df_parameters["mu", "AER"], 
                      sdlog   = df_parameters["sigma", "AER"] 
                      )
draw_ae[,"P"]     <- stats::qlnorm(draw_copula[,2],
                      meanlog = df_parameters["mu", "AEP"],
                      sdlog   = df_parameters["sigma", "AEP"]
                      )

In order to get an idea of the imposed dependency between AER and AEP, let's plot 2000 draws from the t-copula,

4) Aggregate the components and calculate key statistics

The 4th and final step is to aggregate the different distributions and calculate key statistics for the FINMA StandRe template (mean, standard devation, value at risk VaR(99%), expected shortfall ES(99%)). We also show the demeaned expected shortfall (we call it d-ES) as this corresponds to the definition of risk capital according to FINMA.

# Function to compute the results
compute_results <- function(draw_vector) {
    result_EL          <- mean(draw_vector)
    result_SD          <- sd(draw_vector)
    draw_sort          <- sort(draw_vector, decreasing = TRUE)[1:
                            as.integer(0.01 * length(draw_vector))] 
    result_VaR         <- draw_sort[length(draw_sort)]
    result_ES          <- mean(draw_sort)
    result_d_ES        <- result_ES - result_EL
    return(c(result_EL, result_SD, result_VaR, result_ES, result_d_ES))
}
# Set up output data.frame
df_results <- data.frame(matrix(nrow = 5, ncol = 7))
colnames(df_results)   <- c("Total","AEP_IE1","AER_IE2","AEP","AER","IE1","IE2")
rownames(df_results)   <- c("MEAN","SD","VaR","ES","d-ES")
# Aggregate
aep_ie1                <- draw_ae[,"P"] + draw_ie[,1]
aer_ie2                <- draw_ae[,"R"] + draw_ie[,2]
total                  <- aep_ie1 + aer_ie2
# Compute results
df_results[,"IE1"]     <- compute_results(draw_ie[,1])
df_results[,"IE2"]     <- compute_results(draw_ie[,2])
df_results[,"AEP"]     <- compute_results(draw_ae[,"P"])
df_results[,"AER"]     <- compute_results(draw_ae[,"R"])
df_results[,"AER_IE2"] <- compute_results(aer_ie2)
df_results[,"AEP_IE1"] <- compute_results(aep_ie1)
df_results[,"Total"]   <- compute_results(total)

5) Dispay the results and conclude

We conclude by presenting a table with the aggregations 'Total' (= AEP + AER + IE1 + IE2), AEP + IE1 and AER + IE2. If we were to compute the entire SST, the next step would be to transfer the aggregated distribution 'Total' into the FINMA SST template, but this is out of scope of this article. Finally, we recommend to test the sufficiency of the number of simulations; for our example we used 2 million simulations and running it for the seeds 1, 100 and 200 resulted in variations in d-ES of the 'Total' distribution of less than 0.5%.