---
title: "Claim Counts Prediction Using Individual Data"
author: "Gabriele Pittarello"
date: "`r Sys.Date()`"
output: rmarkdown::html_vignette
vignette: >
%\VignetteIndexEntry{Claim Counts Prediction Using Individual Data}
%\VignetteEngine{knitr::rmarkdown}
%\VignetteEncoding{UTF-8}
---
```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
```
```{r eval=FALSE, include=TRUE}
library(ReSurv)
library(reticulate)
use_virtualenv("pyresurv")
library(ggplot2)
```
In this vignette we show a pipeline for predicting individual claim counts using the `ReSurv` package.
This vignette was created for the *2024 Reserving Call Paper Program on Technology and the Reserving Actuary* of the Casualty Actuarial Society (CAS) and the CAS Reserves Working Group. We would like to thank the team of reviewers that provided feedback that improved our package.
# Installation
To install the developer version of the `ReSurv` package.
```{r eval=FALSE, include=TRUE}
library(devtools)
devtools::install_github("edhofman/resurv")
library(ReSurv)
packageVersion("ReSurv")
```
To handle the Python dependencies for *the first package installation*, we suggest to create a dedicated virtual environment. The remaining part of this Section can be disregard from users that are not interested in using our models based on Neural Networks.
```{r eval=FALSE, include=TRUE}
ReSurv::install_pyresurv()
```
The default name of the virtual environment is `"pyresurv"`.
We then suggest to refresh the R session and to import the `ReSurv` package in R using
```{r eval=FALSE, include=TRUE}
library(ReSurv)
reticulate::use_virtualenv("pyresurv")
```
## Managing Multiple Package Dependencies
For the case of multiple packages using isolated-package-environments we suggest the following procedure. Below an example for `pysparklyr`.
```{r eval=FALSE, include=TRUE}
envname <- "./venv"
ReSurv::install_pyresurv(envname = envname)
pysparklyr::install_pyspark(envname = envname)
```
# Data Simulation
We simulate a synthetic data set from scenario Alpha.
```{r eval=FALSE, include=TRUE}
input_data <- data_generator(random_seed = 7,
scenario = "alpha",
time_unit = 1 / 360,
years = 4,
period_exposure = 200)
```
We inspect the data set.
```{r eval=FALSE, include=TRUE}
str(input_data)
```
# Data Pre-Processing
The data is pre-processed using the `IndividualDataPP` function.
```{r eval=FALSE, include=TRUE}
individual_data <- IndividualDataPP(data = input_data,
categorical_features = "claim_type",
continuous_features = "AP",
accident_period = "AP",
calendar_period = "RP",
input_time_granularity = "days",
output_time_granularity = "quarters",
years = 4)
```
# Selection of the hyper-parameters
This Sections first shows the usage of the `ReSurvCV` function for cross-validation. Secondly, we show how to combine the `ReSurvCV` function into the `ParBayesianOptimization`package routine for performing the Bayesian Optimization of Hyperparameters.
## XGB: K-Fold cross-validation (CV)
Example of K-Fold cross-validation (CV) for the XGB model.
```{r eval=FALSE, include=TRUE}
resurv_cv_xgboost <- ReSurvCV(IndividualDataPP = individual_data,
model = "XGB",
hparameters_grid = list(booster = "gbtree",
eta = c(.001),
max_depth = c(3),
subsample = c(1),
alpha = c(0),
lambda = c(0),
min_child_weight = c(.5)),
print_every_n = 1L,
nrounds = 1,
verbose = FALSE,
verbose.cv = TRUE,
early_stopping_rounds = 1,
folds = 5,
parallel = TRUE,
ncores = 2,
random_seed = 1)
```
## NN: K-Fold cross-validation
Example of K-Fold cross-validation (CV) for the NN model.
```{r eval=FALSE, include=TRUE}
resurv_cv_nn <- ReSurvCV(IndividualDataPP = individual_data,
model = "NN",
hparameters_grid = list(num_layers = c(1, 2),
num_nodes = c(2, 4),
optim = "Adam",
activation = "ReLU",
lr = .5,
xi = .5,
eps = .5,
tie = "Efron",
batch_size = as.integer(5000),
early_stopping = "TRUE",
patience = 20),
epochs = as.integer(300),
num_workers = 0,
verbose = FALSE,
verbose.cv = TRUE,
folds = 3,
parallel = FALSE,
random_seed = as.integer(Sys.time()))
```
## ReSurv and Bayesian Parameters Optimisation
### XGB: Bayesian Parameters Optimisation
Example of Bayesian Optimization of Hyperparameters for the XGB model.
```{r eval=FALSE, include=TRUE}
bounds <- list(eta = c(0, 1),
max_depth = c(1L, 25L),
min_child_weight = c(0, 50),
subsample = c(0.51, 1),
lambda = c(0, 15),
alpha = c(0, 15))
obj_func <- function(eta,
max_depth,
min_child_weight,
subsample,
lambda,
alpha) {
xgbcv <- ReSurvCV(
IndividualDataPP = individual_data,
model = "XGB",
hparameters_grid = list(
booster = "gbtree",
eta = eta,
max_depth = max_depth,
subsample = subsample,
alpha = lambda,
lambda = alpha,
min_child_weight = min_child_weight
),
print_every_n = 1L,
nrounds = 500,
verbose = FALSE,
verbose.cv = TRUE,
early_stopping_rounds = 30,
folds = 3,
parallel = FALSE,
random_seed = as.integer(Sys.time())
)
lst <- list(Score = -xgbcv$out.cv.best.oos$test.lkh,
train.lkh = xgbcv$out.cv.best.oos$train.lkh)
return(lst)
}
bayes_out <- bayesOpt(FUN = obj_func,
bounds = bounds,
initPoints = 50,
iters.n = 1000,
iters.k = 50,
otherHalting = list(timeLimit = 18000))
```
### NN: Bayesian Parameters Optimisation
Example of Bayesian Optimization of Hyperparameters for the XGB model.
```{r eval=FALSE, include=TRUE}
bounds <- list(num_layers = c(2L, 10L),
num_nodes = c(2L, 10L),
optim = c(1L, 2L),
activation = c(1L, 2L),
lr = c(0.005, 0.5),
xi = c(0, 0.5),
eps = c(0, 0.5))
obj_func <- function(num_layers,
num_nodes,
optim,
activation,
lr,
xi,
eps) {
optim <- switch(optim,
"Adam",
"SGD")
activation <- switch(activation, "LeakyReLU", "SELU")
batch_size <- as.integer(5000L)
number_layers <- as.integer(num_layers)
num_nodes <- as.integer(num_nodes)
deepsurv_cv <- ReSurvCV(IndividualData = individual_data,
model = "NN",
hparameters_grid = list(num_layers = number_layers,
num_nodes = num_nodes,
optim = optim,
activation = activation,
lr = lr,
xi = xi,
eps = eps,
tie = "Efron",
batch_size = batch_size,
early_stopping = "TRUE",
patience = 20),
epochs = as.integer(300),
num_workers = 0,
verbose = FALSE,
verbose.cv = TRUE,
folds = 3,
parallel = FALSE,
random_seed = as.integer(Sys.time()))
lst <- list(Score = -deepsurv_cv$out.cv.best.oos$test.lkh,
train.lkh = deepsurv_cv$out.cv.best.oos$train.lkh)
return(lst)
}
bayes_out <- bayesOpt(FUN = obj_func,
bounds = bounds,
initPoints = 50,
iters.n = 1000,
iters.k = 50,
otherHalting = list(timeLimit = 18000))
```
# Estimation
In this Section we show the usage of the `ReSurv` method for the estimation of our models parameters. The hyperparameters for NN and XGB were derived with the Bayesian Optimisation procedure.
## COX
```{r eval=FALSE, include=TRUE}
resurv_fit_cox <- ReSurv(individual_data,
hazard_model = "COX")
```
## XGB
```{r eval=FALSE, include=TRUE}
hparameters_xgb <- list(
params = list(
booster = "gbtree",
eta = 0.9611239,
subsample = 0.62851,
alpha = 5.836211,
lambda = 15,
min_child_weight = 29.18158,
max_depth = 1
),
print_every_n = 0,
nrounds = 3000,
verbose = FALSE,
early_stopping_rounds = 500
)
resurv_fit_xgb <- ReSurv(individual_data,
hazard_model = "XGB",
hparameters = hparameters_xgb)
```
## NN
```{r eval=FALSE, include=TRUE}
hparameters_nn = list(num_layers= 2,
early_stopping= TRUE,
patience=350,
verbose=FALSE,
network_structure=NULL,
num_nodes= 10,
activation ="LeakyReLU",
optim ="SGD",
lr =0.02226655,
xi=0.4678993,
epsilon= 0,
batch_size= 5000L,
epochs= 5500L,
num_workers= 0,
tie="Efron" )
resurv_fit_nn <- ReSurv(individual_data,
hazard_model = "NN",
hparameters = hparameters_nn)
```
# Prediction
In this Section we show how to use our models for predictions with the `predict` method and for different output time granularities (quarterly, yearly and monthly). The output can be displayed using the methods `summary` and `print`.
```{r eval=FALSE, include=TRUE}
resurv_fit_predict_q <- predict(resurv_fit_cox)
individual_data_y <- IndividualDataPP(input_data,
id = "claim_number",
continuous_features = "AP_i",
categorical_features = "claim_type",
accident_period = "AP",
calendar_period = "RP",
input_time_granularity = "days",
output_time_granularity = "years",
years = 4,
continuous_features_spline = NULL,
calendar_period_extrapolation = FALSE)
resurv_fit_predict_y <- predict(resurv_fit_cox,
newdata = individual_data_y,
grouping_method = "probability")
individual_data_m <- IndividualDataPP(input_data,
id = "claim_number",
continuous_features = "AP_i",
categorical_features = "claim_type",
accident_period = "AP",
calendar_period = "RP",
input_time_granularity = "days",
output_time_granularity = "months",
years = 4,
continuous_features_spline = NULL,
calendar_period_extrapolation = FALSE)
resurv_fit_predict_m <- predict(resurv_fit_cox,
newdata = individual_data_m,
grouping_method = "probability")
model_s <- summary(resurv_fit_predict_y)
print(model_s)
```
Similar functions for the other models
```{r eval=FALSE, include=TRUE}
resurv_predict_xgb <- predict(resurv_fit_xgb)
resurv_predict_nn <- predict(resurv_fit_nn)
```
# Plot the predicted development factors
In this Section we show an example of development factors plot for the COX model.
```{r eval=FALSE, include=TRUE}
resurv_fit_predict_q$long_triangle_format_out$output_granularity %>%
filter(AP_o==15 & claim_type == 1) %>%
filter(row_number()==1) %>%
select(group_o)
plot(resurv_fit_predict_q,
granularity = "output",
title_par = "COX: Accident Quarter 15 Claim Type 1",
group_code = 30)
```
# CRPS
In this Section we show the CRPS calculation using the `survival_crps` method.
```{r eval=FALSE, include=TRUE}
crps <- survival_crps(resurv_fit_cox)
m_crps <- mean(crps$crps)
m_crps
```
# Plot the development factors
We first extract the output group code for an arbitrary group.
```{r eval=FALSE, include=TRUE}
resurv_fit_predict_q$long_triangle_format_out$output_granularity %>%
filter(AP_o==15 & claim_type == 1) %>%
filter(row_number()==1) %>%
select(group_o)
```
We use the group code to plot the feature dependent development factors.
```{r eval=FALSE, include=TRUE}
plot(resurv_fit_predict_q,
granularity = "output",
title_par = "COX: Accident Quarter 15 Claim Type 1",
group_code = 30)
```
# Appendix
## Code to replicate Figure 1
```{r eval=FALSE, include=TRUE}
p1 <- input_data %>%
as.data.frame() %>%
mutate(claim_type = as.factor(claim_type)) %>%
ggplot(aes(x = RT - AT, color = claim_type)) +
stat_ecdf(size = 1) +
labs(title = "", x = "Notification delay (in days)", y = "ECDF") +
xlim(0.01, 1500) +
scale_color_manual(
values = c("royalblue", "#a71429"),
labels = c("Claim type 0", "Claim type 1")
) +
scale_linetype_manual(values = c(1, 3),
labels = c("Claim type 0", "Claim type 1")) +
guides(
color = guide_legend(
title = "Claim type",
override.aes = list(
color = c("royalblue", "#a71429"),
linewidth = 2
)
),
linetype = guide_legend(
title = "Claim type",
override.aes = list(linetype = c(1, 3), linewidth = 0.7)
)
) +
theme_bw() +
theme(
axis.text = element_text(size = 20),
axis.title.y = element_text(size = 20),
axis.title.x = element_text(size = 20),
legend.text = element_text(size = 20)
)
p1
p2 <- input_data %>%
as.data.frame() %>%
mutate(claim_type = as.factor(claim_type)) %>%
ggplot(aes(x = claim_type, fill = claim_type)) +
geom_bar() +
scale_fill_manual(
values = c("royalblue", "#a71429"),
labels = c("Claim type 0", "Claim type 1")
) +
guides(fill = guide_legend(title = "Claim type")) +
theme_bw() +
labs(title = "", x = "Claim Type", y = "") +
theme(
axis.text = element_text(size = 20),
axis.title.y = element_text(size = 20),
axis.title.x = element_text(size = 20),
legend.text = element_text(size = 20)
)
p2
```
## Code for plotting feature-dependent development factors
Chain-Ladder case.
```{r eval=FALSE, include=TRUE}
clmodel <- resurv_fit_cox$IndividualDataPP$training.data %>%
mutate(DP_o = 16 -
DP_rev_o + 1) %>%
group_by(AP_o, DP_o) %>%
summarize(I = sum(I), .groups = "drop") %>%
group_by(AP_o) %>%
arrange(DP_o) %>%
mutate(I_cum = cumsum(I), I_cum_lag = lag(I_cum, default = 0)) %>%
ungroup() %>%
group_by(DP_o) %>%
reframe(df_o = sum(I_cum * (
AP_o <= max(resurv_fit_cox$IndividualDataPP$training.data$AP_o) - DP_o +
1
)) /
sum(I_cum_lag * (
AP_o <= max(resurv_fit_cox$IndividualDataPP$training.data$AP_o) - DP_o +
1
)),
I = sum(I * (
AP_o <= max(resurv_fit_cox$IndividualDataPP$training.data$AP_o) - DP_o
))) %>%
mutate(DP_o_join = DP_o - 1) %>%
as.data.frame()
clmodel %>%
filter(DP_o > 1) %>%
ggplot(aes(x = DP_o, y = df_o)) +
geom_line(linewidth = 2.5,
color = "#454555") +
labs(title = "Chain ladder",
x = "Development quarter",
y = "Development factor") +
ylim(1, 3. + .01) +
theme_bw(base_size = rel(5)) +
theme(plot.title = element_text(size = 20))
##
clmodel_months <- individual_data_m$training.data %>%
mutate(DP_o = 48 -
DP_rev_o + 1) %>%
group_by(AP_o, DP_o) %>%
summarize(I = sum(I), .groups = "drop") %>%
group_by(AP_o) %>%
arrange(DP_o) %>%
mutate(I_cum = cumsum(I), I_cum_lag = lag(I_cum, default = 0)) %>%
ungroup() %>%
group_by(DP_o) %>%
reframe(df_o = sum(I_cum * (
AP_o <= max(individual_data_m$training.data$AP_o) - DP_o + 1
)) /
sum(I_cum_lag * (
AP_o <= max(individual_data_m$training.data$AP_o) - DP_o + 1
)),
I = sum(I * (
AP_o <= max(individual_data_m$training.data$AP_o) - DP_o
))) %>%
mutate(DP_o_join = DP_o - 1) %>%
as.data.frame()
ticks_at <- seq(1, 48, 4)
labels_as <- as.character(ticks_at)
clmodel_months %>%
filter(DP_o > 1) %>%
ggplot(aes(x = DP_o,
y = df_o)) +
geom_line(linewidth = 2.5,
color = "#454555") +
labs(title = "Chain ladder",
x = "Development month",
y = "Development factor") +
ylim(1, 2.5 + .01) +
scale_x_continuous(breaks = ticks_at, labels = labels_as) +
theme_bw(base_size = rel(5)) +
theme(plot.title = element_text(size = 20))
##
clmodel_years <- individual_data_y$training.data %>%
mutate(DP_o = 4 -
DP_rev_o + 1) %>%
group_by(AP_o, DP_o) %>%
summarize(I = sum(I), .groups = "drop") %>%
group_by(AP_o) %>%
arrange(DP_o) %>%
mutate(I_cum = cumsum(I), I_cum_lag = lag(I_cum, default = 0)) %>%
ungroup() %>%
group_by(DP_o) %>%
reframe(df_o = sum(I_cum * (
AP_o <= max(individual_data_m$training.data$AP_o) - DP_o + 1
)) /
sum(I_cum_lag * (
AP_o <= max(individual_data_m$training.data$AP_o) - DP_o + 1
)),
I = sum(I * (
AP_o <= max(individual_data_m$training.data$AP_o) - DP_o
))) %>%
mutate(DP_o_join = DP_o - 1) %>%
as.data.frame()
ticks_at <- seq(1, 4, 1)
labels_as <- as.character(ticks_at)
clmodel_years %>%
filter(DP_o > 1) %>%
ggplot(aes(x = DP_o,
y = df_o)) +
geom_line(linewidth = 2.5,
color = "#454555") +
labs(title = "Chain ladder",
x = "Development year",
y = "Development factor") +
ylim(1, 2.5 + .01) +
scale_x_continuous(breaks = ticks_at, labels = labels_as) +
theme_bw(base_size = rel(5)) +
theme(plot.title = element_text(size = 20))
```
COX, quarterly output.
```{r eval=FALSE, include=TRUE}
ct <- 0
ap <- 12
resurv_fit_predict_q$long_triangle_format_out$output_granularity %>%
filter(AP_o == ap & claim_type == ct) %>%
filter(row_number() == 1) %>%
select(group_o)
plot(
resurv_fit_predict_q,
granularity = "output",
title_par = "COX: Accident Quarter 12 Claim Type 0",
group_code = 23
)
ct <- 0
ap <- 36
resurv_fit_predict_m$long_triangle_format_out$output_granularity %>%
filter(AP_o == ap & claim_type == ct) %>%
filter(row_number() == 1) %>%
select(group_o)
plot(
resurv_fit_predict_m,
granularity = "output",
color_par = "#a71429",
title_par = "COX: Accident Month 36 Claim Type 0",
group_code = 71
)
ct <- 1
ap <- 7
resurv_fit_predict_m$long_triangle_format_out$output_granularity %>%
filter(AP_o == ap & claim_type == ct) %>%
filter(row_number() == 1) %>%
select(group_o)
plot(
resurv_fit_predict_m,
granularity = "output",
color_par = "#a71429",
title_par = "COX: Accident Month 7 Claim Type 1",
group_code = 14
)
ct <- 0
ap <- 2
resurv_fit_predict_y$long_triangle_format_out$output_granularity %>%
filter(AP_o == ap & claim_type == ct) %>%
filter(row_number() == 1) %>%
select(group_o)
plot(
resurv_fit_predict_y,
granularity = "output",
color_par = "#FF6A7A",
title_par = "COX: Accident Year 2 Claim Type 0",
group_code = 3
)
ct <- 1
ap <- 3
resurv_fit_predict_y$long_triangle_format_out$output_granularity %>%
filter(AP_o == ap & claim_type == ct) %>%
filter(row_number() == 1) %>%
select(group_o)
plot(
resurv_fit_predict_y,
granularity = "output",
color_par = "#FF6A7A",
title_par = "COX: Accident Year 3 Claim Type 1",
group_code = 6
)
```
## Code for computing ARE TOT and ARE CAL
The code below can be used to compute the performance metrics in Table 6 of the manuscript.
```{r eval=FALSE, include=TRUE}
conversion_factor <- individual_data$conversion_factor
max_dp_i <- 1440
true_output <- resurv_fit_cox$IndividualDataPP$full.data %>%
mutate(
DP_rev_o = floor(max_dp_i * conversion_factor) -
ceiling(
DP_i * conversion_factor +
((AP_i - 1) %% (
1 / conversion_factor
)) * conversion_factor) + 1,
AP_o = ceiling(AP_i * conversion_factor),
TR_o = AP_o - 1
) %>%
filter(DP_rev_o <= TR_o) %>%
group_by(claim_type, AP_o, DP_rev_o) %>%
mutate(claim_type = as.character(claim_type)) %>%
summarize(I = sum(I), .groups = "drop") %>%
filter(DP_rev_o > 0)
out_list <- resurv_fit_predict_q$long_triangle_format_out
out <- out_list$output_granularity
#Total output
score_total <- out[, c("claim_type", "AP_o", "DP_o", "expected_counts")] %>%
mutate(DP_rev_o = 16 - DP_o + 1) %>%
inner_join(true_output, by = c("claim_type", "AP_o", "DP_rev_o")) %>%
mutate(ave = I - expected_counts, abs_ave = abs(ave)) %>%
# from here it is reformulated for the are tot
ungroup() %>%
group_by(AP_o, DP_rev_o) %>%
reframe(abs_ave = abs(sum(ave)), I = sum(I))
are_tot <- sum(score_total$abs_ave) / sum(score_total$I)
dfs_output <- out %>%
mutate(DP_rev_o = 16 - DP_o + 1) %>%
select(AP_o, claim_type, DP_rev_o, f_o) %>%
mutate(DP_rev_o = DP_rev_o) %>%
distinct()
#Cashflow on output scale.Etc quarterly cashflow development
score_diagonal <- resurv_fit_cox$IndividualDataPP$full.data %>%
mutate(
DP_rev_o = floor(max_dp_i * conversion_factor) -
ceiling(
DP_i * conversion_factor +
((AP_i - 1) %% (
1 / conversion_factor
)) * conversion_factor) + 1,
AP_o = ceiling(AP_i * conversion_factor)
) %>%
group_by(claim_type, AP_o, DP_rev_o) %>%
mutate(claim_type = as.character(claim_type)) %>%
summarize(I = sum(I), .groups = "drop") %>%
group_by(claim_type, AP_o) %>%
arrange(desc(DP_rev_o)) %>%
mutate(I_cum = cumsum(I)) %>%
mutate(I_cum_lag = lag(I_cum, default = 0)) %>%
left_join(dfs_output, by = c("AP_o", "claim_type", "DP_rev_o")) %>%
mutate(I_cum_hat = I_cum_lag * f_o,
RP_o = max(DP_rev_o) - DP_rev_o + AP_o) %>%
inner_join(true_output[, c("AP_o", "DP_rev_o")] %>% distinct()
, by = c("AP_o", "DP_rev_o")) %>%
group_by(AP_o, DP_rev_o) %>%
reframe(abs_ave2_diag = abs(sum(I_cum_hat) - sum(I_cum)), I = sum(I))
are_cal_q <- sum(score_diagonal$abs_ave2_diag) / sum(score_diagonal$I)
# scoring XGB ----
out_xgb <- resurv_predict_xgb$long_triangle_format_out$output_granularity
score_total_xgb <- out_xgb[, c("claim_type",
"AP_o",
"DP_o",
"expected_counts")] %>%
mutate(DP_rev_o = 16 - DP_o + 1) %>%
inner_join(true_output, by = c("claim_type", "AP_o", "DP_rev_o")) %>%
mutate(ave = I - expected_counts, abs_ave = abs(ave)) %>%
# from here it is reformulated for the are tot
ungroup() %>%
group_by(AP_o, DP_rev_o) %>%
reframe(abs_ave = abs(sum(ave)), I = sum(I))
are_tot_xgb <- sum(score_total_xgb$abs_ave) / sum(score_total$I)
dfs_output_xgb <- out_xgb %>%
mutate(DP_rev_o = 16 - DP_o + 1) %>%
select(AP_o, claim_type, DP_rev_o, f_o) %>%
mutate(DP_rev_o = DP_rev_o) %>%
distinct()
#Cashflow on output scale.Etc quarterly cashflow development
score_diagonal_xgb <- resurv_fit_cox$IndividualDataPP$full.data %>%
mutate(
DP_rev_o = floor(max_dp_i * conversion_factor) -
ceiling(
DP_i * conversion_factor +
((AP_i - 1) %% (
1 / conversion_factor
)) * conversion_factor) + 1,
AP_o = ceiling(AP_i * conversion_factor)
) %>%
group_by(claim_type, AP_o, DP_rev_o) %>%
mutate(claim_type = as.character(claim_type)) %>%
summarize(I = sum(I), .groups = "drop") %>%
group_by(claim_type, AP_o) %>%
arrange(desc(DP_rev_o)) %>%
mutate(I_cum = cumsum(I)) %>%
mutate(I_cum_lag = lag(I_cum, default = 0)) %>%
left_join(dfs_output_xgb, by = c("AP_o", "claim_type", "DP_rev_o")) %>%
mutate(I_cum_hat = I_cum_lag * f_o,
RP_o = max(DP_rev_o) - DP_rev_o + AP_o) %>%
inner_join(true_output[, c("AP_o", "DP_rev_o")] %>% distinct()
, by = c("AP_o", "DP_rev_o")) %>%
group_by(AP_o, DP_rev_o) %>%
reframe(abs_ave2_diag = abs(sum(I_cum_hat) - sum(I_cum)), I = sum(I))
are_cal_q_xgb <- sum(score_diagonal_xgb$abs_ave2_diag) / sum(score_diagonal$I)
# scoring NN ----
out_nn <- resurv_predict_nn$long_triangle_format_out$output_granularity
score_total_nn <- out_nn[, c("claim_type",
"AP_o",
"DP_o",
"expected_counts")] %>%
mutate(DP_rev_o = 16 - DP_o + 1) %>%
inner_join(true_output, by = c("claim_type", "AP_o", "DP_rev_o")) %>%
mutate(ave = I - expected_counts, abs_ave = abs(ave)) %>%
# from here it is reformulated for the are tot
ungroup() %>%
group_by(AP_o, DP_rev_o) %>%
reframe(abs_ave = abs(sum(ave)), I = sum(I))
are_tot_nn <- sum(score_total_nn$abs_ave) / sum(score_total$I)
dfs_output_nn <- out_nn %>%
mutate(DP_rev_o = 16 - DP_o + 1) %>%
select(AP_o, claim_type, DP_rev_o, f_o) %>%
mutate(DP_rev_o = DP_rev_o) %>%
distinct()
#Cashflow on output scale.Etc quarterly cashflow development
score_diagonal_nn <- resurv_fit_cox$IndividualDataPP$full.data %>%
mutate(
DP_rev_o = floor(max_dp_i * conversion_factor) -
ceiling(
DP_i * conversion_factor +
((AP_i - 1) %% (
1 / conversion_factor
)) * conversion_factor) + 1,
AP_o = ceiling(AP_i * conversion_factor)
) %>%
group_by(claim_type, AP_o, DP_rev_o) %>%
mutate(claim_type = as.character(claim_type)) %>%
summarize(I = sum(I), .groups = "drop") %>%
group_by(claim_type, AP_o) %>%
arrange(desc(DP_rev_o)) %>%
mutate(I_cum = cumsum(I)) %>%
mutate(I_cum_lag = lag(I_cum, default = 0)) %>%
left_join(dfs_output_nn, by = c("AP_o", "claim_type", "DP_rev_o")) %>%
mutate(I_cum_hat = I_cum_lag * f_o,
RP_o = max(DP_rev_o) - DP_rev_o + AP_o) %>%
inner_join(true_output[, c("AP_o", "DP_rev_o")] %>% distinct()
, by = c("AP_o", "DP_rev_o")) %>%
group_by(AP_o, DP_rev_o) %>%
reframe(abs_ave2_diag = abs(sum(I_cum_hat) - sum(I_cum)), I = sum(I))
are_cal_q_nn <- sum(score_diagonal_nn$abs_ave2_diag) / sum(score_diagonal$I)
# Scoring Chain-Ladder ----
true_output_cl <- individual_data$full.data %>%
mutate(
DP_rev_o = floor(max_dp_i * conversion_factor) -
ceiling(
DP_i * conversion_factor +
((AP_i - 1) %% (
1 / conversion_factor
)) * conversion_factor) + 1,
AP_o = ceiling(AP_i * conversion_factor)
) %>%
filter(DP_rev_o <= TR_o) %>%
mutate(DP_o = max(individual_data$training.data$DP_rev_o) - DP_rev_o + 1) %>%
group_by(AP_o, DP_o, DP_rev_o) %>%
summarize(I = sum(I), .groups = "drop") %>%
filter(DP_rev_o > 0)
latest_observed <- individual_data$training.data %>%
filter(DP_rev_o >= TR_o) %>%
mutate(DP_o = max(individual_data$training.data$DP_rev_o) - DP_rev_o + 1) %>%
group_by(AP_o) %>%
mutate(DP_max = max(DP_o)) %>%
group_by(AP_o, DP_max) %>%
summarize(I = sum(I), .groups = "drop")
clmodel <- individual_data$training.data %>%
mutate(DP_o = max(individual_data$training.data$DP_rev_o) - DP_rev_o + 1) %>%
group_by(AP_o, DP_o) %>%
summarize(I = sum(I), .groups = "drop") %>%
group_by(AP_o) %>%
arrange(DP_o) %>%
mutate(I_cum = cumsum(I), I_cum_lag = lag(I_cum, default = 0)) %>%
ungroup() %>%
group_by(DP_o) %>%
reframe(df = sum(I_cum * (
AP_o <= max(individual_data$training.data$AP_o) - DP_o + 1
)) /
sum(I_cum_lag * (
AP_o <= max(individual_data$training.data$AP_o) - DP_o + 1
)), I = sum(I * (
AP_o <= max(individual_data$training.data$AP_o) - DP_o
))) %>%
mutate(DP_o_join = DP_o) %>%
mutate(DP_rev_o = max(DP_o) - DP_o + 1)
predictions <- expand.grid(AP_o = latest_observed$AP_o,
DP_o = clmodel$DP_o_join) %>%
left_join(clmodel[, c("DP_o_join", "df")], by = c("DP_o" = "DP_o_join")) %>%
left_join(latest_observed, by = "AP_o") %>%
rowwise() %>%
filter(DP_o > DP_max) %>%
ungroup() %>%
group_by(AP_o) %>%
arrange(DP_o) %>%
mutate(df_cum = cumprod(df)) %>%
mutate(I_expected = I * df_cum - I * lag(df_cum, default = 1)) %>%
select(DP_o, AP_o, I_expected)
conversion_factor <- individual_data$conversion_factor
max_dp_i <- 1440
score_total <- predictions %>%
inner_join(true_output_cl, by = c("AP_o", "DP_o")) %>%
mutate(ave = I - I_expected, abs_ave = abs(ave)) %>%
# from here it is reformulated for the are tot
ungroup() %>%
group_by(AP_o, DP_rev_o) %>%
reframe(abs_ave = abs(sum(ave)), I = sum(I))
are_tot <- sum(score_total$abs_ave) / sum(score_total$I)
score_diagonal <- individual_data$full.data %>%
mutate(
DP_rev_o = floor(max_dp_i * conversion_factor) -
ceiling(
DP_i * conversion_factor +
((AP_i - 1) %% (
1 / conversion_factor
)) * conversion_factor) + 1,
AP_o = ceiling(AP_i * conversion_factor)
) %>%
group_by(claim_type, AP_o, DP_rev_o) %>%
mutate(claim_type = as.character(claim_type)) %>%
summarize(I = sum(I), .groups = "drop") %>%
group_by(claim_type, AP_o) %>%
arrange(desc(DP_rev_o)) %>%
mutate(I_cum = cumsum(I)) %>%
mutate(I_cum_lag = lag(I_cum, default = 0)) %>%
mutate(DP_o = max(individual_data$training.data$DP_rev_o) - DP_rev_o + 1) %>%
left_join(CL[, c("DP_o", "df")], by = c("DP_o")) %>%
mutate(I_cum_hat = I_cum_lag * df,
RP_o = max(DP_rev_o) - DP_rev_o + AP_o) %>%
inner_join(true_output_cl[, c("AP_o", "DP_rev_o")] %>% distinct()
, by = c("AP_o", "DP_rev_o")) %>%
group_by(AP_o, DP_rev_o) %>%
reframe(abs_ave2_diag = abs(sum(I_cum_hat) - sum(I_cum)), I = sum(I))
are_cal_q <- sum(score_diagonal$abs_ave2_diag) / sum(score_diagonal$I)
```