## ----setup, include=FALSE-----------------------------------------------------
knitr::opts_chunk$set(collapse = TRUE, comment = "#>")
library(actuarialfitdist)
set.seed(20261001)

## ----auto-limit-data----------------------------------------------------------
n <- 1200
true_loss <- rlnorm(n, meanlog = log(9000), sdlog = 0.95)
policy_limit <- sample(c(10000, 25000, 50000, 100000), n, replace = TRUE,
                       prob = c(0.15, 0.35, 0.35, 0.15))

pd <- data.frame(
  loss = pmin(true_loss, policy_limit),
  deductible = 0,
  limit = policy_limit
)

mean(pd$loss >= pd$limit)

## ----auto-limit-fit-----------------------------------------------------------
pd_candidates <- fit_severity_candidates(
  pd,
  loss = "loss",
  limit = "limit",
  distributions = c("lognormal", "weibull", "gamma", "burr", "loglogistic")
)

pd_candidates$comparison

## ----auto-limit-factors-------------------------------------------------------
pd_factors <- severity_factors(
  pd_candidates,
  deductible = 0,
  limit = c(10000, 25000, 50000, 100000, 250000),
  base_deductible = 0,
  base_limit = 25000
)

head(pd_factors)

## ----auto-limit-side-by-side--------------------------------------------------
compare_factors(
  pd_candidates,
  metric = "rating_factor",
  deductible = 0,
  limit = c(10000, 25000, 50000, 100000),
  base_deductible = 0,
  base_limit = 25000
)

## ----collision-data-----------------------------------------------------------
n <- 1500
collision_ground_up <- rgamma(n, shape = 2.2, scale = 3200)
collision_ded <- sample(c(250, 500, 1000, 2000), n, replace = TRUE,
                        prob = c(0.10, 0.35, 0.40, 0.15))
keep <- collision_ground_up > collision_ded

collision <- data.frame(
  loss = collision_ground_up[keep],
  deductible = collision_ded[keep],
  limit = Inf
)

## ----collision-fit------------------------------------------------------------
collision_fit <- fit_severity(
  collision,
  loss = "loss",
  deductible = "deductible",
  limit = "limit",
  distribution = "gamma"
)

summary(collision_fit)

## ----collision-factors--------------------------------------------------------
severity_factors(
  collision_fit,
  deductible = c(250, 500, 1000, 2000, 5000),
  limit = Inf,
  base_deductible = 500,
  base_limit = Inf
)

## ----homeowners-data----------------------------------------------------------
n <- 1600
coverage_a <- sample(seq(150000, 900000, by = 25000), n, replace = TRUE)
home_ded <- sample(c(1000, 2500, 5000, 10000), n, replace = TRUE)
shape <- 1.35
scale_coef <- 0.075
latent_scale <- scale_coef * coverage_a
home_ground_up <- rweibull(n, shape = shape, scale = latent_scale)
home_limit <- coverage_a
keep <- home_ground_up > home_ded

home <- data.frame(
  loss = pmin(home_ground_up[keep], home_limit[keep]),
  deductible = home_ded[keep],
  limit = home_limit[keep],
  coverage_a = coverage_a[keep]
)

## ----homeowners-fit-----------------------------------------------------------
home_fit <- fit_severity(
  home,
  loss = "loss",
  deductible = "deductible",
  limit = "limit",
  distribution = "weibull",
  scale_by = "coverage_a",
  scale_intercept = FALSE
)

summary(home_fit)

## ----homeowners-factors-------------------------------------------------------
severity_factors(
  home_fit,
  scale_value = c(200000, 400000, 600000, 800000),
  deductible = c(1000, 2500, 5000),
  limit = c(100000, 250000, 500000),
  base_deductible = 1000,
  base_limit = 250000
)

## ----homeowners-fit-plot, fig.width=7, fig.height=4---------------------------
plot_fit(
  home_fit,
  type = "density",
  scale_range = c(300000, 450000)
)

## ----homeowners-scale-plot, fig.width=7, fig.height=4-------------------------
plot_scaling_fit(home_fit, bins = 10, nsim = 3, seed = 12)

## ----splice-data--------------------------------------------------------------
n <- 1500
threshold_true <- 30000
is_body <- runif(n) < 0.90

# Draw body values from a lognormal conditional on being below the threshold.
Fb <- plnorm(threshold_true, log(9000), 0.75)
body_draw <- qlnorm(runif(n, 0, Fb), log(9000), 0.75)

# Draw Pareto-II tail values conditional on exceeding the threshold.
shape_tail <- 2.1
scale_tail <- 18000
Ft <- 1 - (1 + threshold_true / scale_tail)^(-shape_tail)
u_tail <- runif(n, Ft, 1)
tail_draw <- scale_tail * ((1 - u_tail)^(-1 / shape_tail) - 1)

splice_ground_up <- ifelse(is_body, body_draw, tail_draw)
splice_ded <- sample(c(0, 1000, 2500), n, replace = TRUE)
splice_lim <- sample(c(50000, 100000, 250000, Inf), n, replace = TRUE)
keep <- splice_ground_up > splice_ded

splice_dat <- data.frame(
  loss = pmin(splice_ground_up[keep], splice_lim[keep]),
  deductible = splice_ded[keep],
  limit = splice_lim[keep]
)

## ----splice-fit---------------------------------------------------------------
splice_fit <- fit_spliced_severity(
  splice_dat,
  loss = "loss",
  deductible = "deductible",
  limit = "limit",
  body = "lognormal",
  tail = "pareto2",
  threshold_method = "search",
  threshold_probs = c(0.60, 0.95),
  threshold_grid = 4
)

splice_fit
splice_fit$threshold_search

## ----splice-factors-----------------------------------------------------------
severity_factors(
  splice_fit,
  deductible = c(0, 1000, 5000),
  limit = c(25000, 50000, 100000, 250000),
  base_deductible = 0,
  base_limit = 50000
)

## ----valuation----------------------------------------------------------------
severity_factors(
  home_fit,
  scale_value = 400000,
  deductible = 2500,
  limit = 25000,
  valuation_factor = c(1.00, 0.80, 0.60, 0.40, 0.20),
  base_deductible = 2500,
  base_limit = 25000,
  base_valuation_factor = 1
)

## ----weights------------------------------------------------------------------
weighted_collision <- collision[1:1000, ]
weighted_collision$w <- sample(1:3, nrow(weighted_collision), replace = TRUE)

weighted_fit <- fit_severity(
  weighted_collision,
  loss = "loss",
  deductible = "deductible",
  limit = "limit",
  distribution = "gamma",
  weights = "w",
  fixed = list(shape = 2.2)
)

coef(weighted_fit)

## ----hessian------------------------------------------------------------------
vcov(collision_fit)
confint(collision_fit)

## ----bootstrap, eval=FALSE----------------------------------------------------
# boot <- bootstrap_fit(collision_fit, R = 5, type = "parametric", seed = 99)
# head(boot$estimates)
# boot$success_rate

## ----custom-data--------------------------------------------------------------
d <- fit_plot_data(home_fit, scale_range = c(300000, 450000), nsim = 3, seed = 5)
head(d)

## ----new-api, fig.width=7, fig.height=4---------------------------------------
set.seed(20261008)
no_terms <- data.frame(loss = rlnorm(250, log(12000), .8))
fit_no_terms <- fit_severity(no_terms, loss = "loss",
  distribution = "lognormal", hessian = FALSE)
print(fit_no_terms)
plot_fit(fit_no_terms, x_scale = "log10",
  xlim = c(1000, 100000))

## ----excluded-rows------------------------------------------------------------
problem <- no_terms
problem$deductible <- 0
problem$deductible[3] <- problem$loss[3] + 10
clean_fit <- suppressWarnings(fit_severity(problem, loss = "loss",
  deductible = "deductible", distribution = "lognormal", hessian = FALSE))
clean_fit$excluded_rows
clean_fit$excluded_details

## ----candidate-comparison-plot, fig.width=8, fig.height=4---------------------
plot_candidate_fits(pd_candidates, type = "density",
  x_scale = "log10", xlim = c(200, 150000))