Dynamics of a fractal-fractional mathematical model for the management of waste plastic in the ocean with four different numerical approaches

1.   Introduction

In the preface of the book of ^[1], it is stated that reinsurance is a fascinating field, where several of the challenges of classical insurance are amplified, particularly when it comes to dealing with extreme situations like large claims and rare events. However, leaving this enthusiastic phrase aside, the reinsurance field is also exciting because of its practicality and ability to solve many financial problems for insurance companies. As is well known by the actuarial community, a reinsurance contract establishes an agreement or contract between the insurer and the reinsurer, whereby the former (also called the ceding company) indemnifies the latter a part of the risk subscribed by the insurer, with a third party (usually the customer) paying, in turn, a reinsurance premium. As such, it is a contract similar to the one that an insurance company signs with a client to cover damages of any kind. The difference, however, is that the reinsurer covers only a part of the compensation, called priority. This type of reinsurance is called excess-of-loss reinsurance. For a detailed study of reinsurance, in all its forms, readers are referred to the works of ^[1,2,3,4,5] and the references therein. Recent works on reinsurance include contributions such as ^[6,7,8,9]. In these recent works, the reinsurance technique is also used in non-life insurance, but its scope is broadened, including, for example, non-life multiline insurers under the Solvency Ⅱ umbrella, risk transfers in the context of the 2007-08 global crisis, reinsurance using two layers and multiple reinsurers, and reinsurance in a game-theoretic setting, considering one insurer and two reinsurers.

As part of the collective risk model, in excess-of-loss reinsurance, it is usual to assume that the risk parameters associated with the random variables, the number of claims, and the claim size are independent. This is the main idea in ^[10], inspired by the ideas developed in ^[11]. For example, we can describe this situation by taking the Poisson frequency and the Pareto severity as appropriate distributional assumptions. The use of the Bayesian methodology facilitates the computation of predictive and posterior moments of the accumulated claims and the interpretation of the parameters involved. ^[12] also studied the Bayesian estimators in Paretian excess-of-loss reinsurance. Although these derivations have traditionally been built on the assumption of independence (see ^[13]), there has been great interest in modeling dependent risks due to the increasing complexity of the current insurance and reinsurance products.

In this paper, we will pay attention to the calculation of the prior and posterior moments of the size of the total claims, assuming some degree of dependence derived from the prior distributions of the risk parameters (profiles) involved in the Bayesian framework. This degree of dependence will be modeled using the Sarmanov-Lee family of distributions; specific details of this type of copula can be seen in ^[14,15]. This model allows us to study the impact of dependence on the next year's total claims costs, showing that a low degree of correlation produces highly sensitive results. Furthermore, the case where the moments are calculated separately without assuming dependence is also considered. In recent years, some actuarial works have considered the Sarmanov family of bivariate distributions: see, for example, ^[16,17,18]. Also, ^[19] considered, from a Bayesian point of view, the calculation of credibility premiums for the compound loss under a bivariate prior for the risk parameters of the frequency and severity distributions.

The impact of assuming the dependence feature is a matter of concern. At first sight, it would seem that the effect of including dependence has a low impact in the first and second moments; in fact, the reinsurer should be prudent, since ignoring the hypothesis of dependence might lead to budget imbalance, especially when there is a change of scale in the monetary units used. This work assumes the dependence hypothesis between the risk profiles instead of the claims number and the size of the corresponding claims. The former approach not only retains the main features of the latter case but also reduces the mathematical burden of the independence assumption between the number and the severity.

The rest of the paper has the following structure: The basic collective risk model considered in this paper is introduced in Section 2. Here, the bivariate prior distribution assumed for the risk profile (dependent on two parameters) is also presented and studied. In the original paper of ^[10], a classical gamma distribution was assumed, but this seems inadequate considering this risk parameter's support. To make the paper self-contained, we recall the main results of ^[10] in Section 3. Here, we also give closed-form expressions for the risk, collective, and Bayes net premiums in the excess-of-loss reinsurance model. An extension of the basic collective risk model and the reinsurance model, which includes some degree of dependence, is introduced in Section 4. Here, we also provide properties of this extended model together with the expressions of the collective and Bayesian net premiums, which are written in closed-form expressions. The cases with and without a layer are both studied. Numerical applications are provided in Section 5, and final comments are given in the last section.

2.   The specific basic collective risk model

Let $N(t)$ be the random variable denoting the number of claims reported to the insurer during $[0, t]$ and let $Y_1, Y_2, \dots, Y_{N(t)}$ be the corresponding individual claim amounts up to time $t$. Assume the following: $(i)$ $\left\{N(t), \, t\geq 0\right\}$ is a time homogeneous Poisson process with intensity $\lambda > 0$, i.e., $p_n(t) = \exp(-\lambda t)(\lambda t)^n/n!$ and $(ii)$ $\left\{Y_i\right\}_{i = 1, 2, \dots}$ are stochastically independent and independent of $N(t)$ with common survival function $\bar F(y|c, \psi) = (c/y)^{\psi}$, $y\geq c$, i.e., the classical Pareto distribution with scale parameter $c > 0$ and shape parameter $\psi > 0$. Its mean is given by $\mathbb{E}(Y) = c\psi/(\psi-1)$, $\psi > 1$. The compound process described above has a cumulative distribution function (see, for example, ^[20]) given by

where $F^{n\ast}(\cdot)$ represents the $n$th-fold convolution of $F(\cdot|c, \psi)$. In practice, it is difficult to derive a closed-form expression for the distribution function given in (2.1). Therefore, numerical approximations have been proposed based on the normal and gamma distribution and fast Fourier transform. Other approximations that have been considered are based on Edgeworth expansions and the Gram-Charlier series. Furthermore, an exact recursive expression for (2.1) can be found by using Panjer's recursion formula (see ^[21] and [22, Chapter 4], among others). A closed-form expression for the Pareto type Ⅲ (Lomax) distribution has been derived for the convolution above in ^[23], and the asymptotic distribution of an infinite sum of i.i.d. classical Pareto type Ⅰ random variables is provided in ^[24].

Nevertheless, in the expression (2.1), it is not necessary to calculate premiums based on the moments of the distribution. We will consider that the parameter $c$ is known, that the risk profiles $\lambda$ and $\psi$ are unknown, and that they take values from the random variables $\Lambda$ and $\Psi$, respectively. As it is well-known, a premium calculation principle (^[25,26]; among others) assigns to each vector of risk profiles $\Theta = (\lambda, \psi)$ a premium within the set ${P}\in\mathbb{R^+}$, the action space. Using the net premium principle, it is easy to see that, in this case, it gives the unknown risk premium

We have followed the idea of ^[10] in order to define the net premium. That is, the net premium is defined as the sum of losses instead of the sum of discounted losses, which leads to the exclusion of the time value of the monetary losses.

In the absence of experience, the actuary computes the collective (prior) premium, $\mathbb{E}_{\pi}[P(\Theta)]$, where $\pi$ is the prior distribution assigned to the risk profile $\Theta$. On the other hand, if experience is available, the actuary takes a sample

and uses this information to estimate the unknown risk premium $P(\Theta)$ through the Bayes premium $P^{\ast}$, obtained by computing the posterior expectation $\mathbb{E}_{\pi^{\ast}}[P(\Theta)]$. Here, $\pi^{\ast}$ is the posterior distribution of the risk parameter, $\Theta$, given the sample information provided in (2.2). For that purpose, we need also the likelihood function, which is given by

where $z = \sum_{i = 1}^{n}\log(y_t/c)$.

We will assume that the distribution of $\Lambda$ is gamma with shape parameter $\nu > 0$ and rate parameter $\tau > 0$ and that the distribution of $\Psi$ is Erlang with shape parameter $\gamma\in\mathbb{N}_+ = \{1, 2, \dots\}$ and rate parameter $\xi > 0$. That is, it is assumed that $\Lambda\sim{\cal G}(\nu, \tau)$ and $\Psi\sim{\cal E}(\gamma, \xi, k)$, a $k > 0$ translated Erlang variate, with probability density functions (pdfs) given by, ^*

^*[10] assumes a gamma prior for the parameter $\psi$ that seems inappropriate since from Proposition 3, we have that $\psi > k$. A translated gamma or Erlang distribution should be more accurate in the case considered here.

where $\Gamma(\cdot)$ is the Euler gamma function and $\Gamma(\gamma) = (\gamma-1)!$. The choice of the prior distributions, in addition to mathematical convenience (calculations are provided), is done in this way following Hesselager's work ^[10]. The fundamental objective of this work is to study the variation of the premiums assuming a certain degree of dependence between the risk profiles. Mean and variance of these two distributions are $\nu/\tau$, $\nu/\tau^2$, $k+\gamma/\xi$, and $\gamma/\xi^2$, respectively.

Finally, the joint pdf of $\Theta = (\Lambda, \Psi)$ is assumed to be the product of these densities. Thus,

i.e., both random variables are assumed to be independent. The next result provides the posterior distribution of (2.6) given the sample information provided in (2.3). As we will see, the posterior distribution is conjugate for $\lambda$ and almost conjugate (quasi-conjugate) for the $\psi$ risk profile.

Proposition 1. Let $\left\{X(t)\right\}_{t\geq 0}$ be the compound Poisson process defined above. Suppose that $(\Lambda, \Psi)$ follows the prior distribution with pdf given in (2.6). Then, the posterior distribution of $(\Lambda, \Psi)$ given the sample information $\aleph = \left\{ N_c(T) = n, \; \; Y_{c, t} = y_{ti}, \; \; i = 1, 2, \dots, n\right\}$, assumed to be independent and identically distributed, is given by,

where

and the updated parameters are given by

being $\kappa\equiv\kappa(\gamma, \xi^{\ast}, k)$ the normalization constant given by

Proof. The posterior pdf $\pi^{\ast}_1(\lambda|\nu^{\ast}, \tau^{\ast})$ is straightforwardly obtained from (2.3), (2.4), and Bayes' theorem. Now, $\pi_2^{\ast}(\psi|\gamma^{\ast}, \xi^{\ast}, k)$ is proportional to

where the constant of proportionality is given by

Hence the result. □

As we can see in the next result, bivariate, prior, and posterior distributions can now be used to get the collective (prior) and posterior (Bayes) premiums of the aggregated model described above based on the first moment (the net premium).

Proposition 2. The collective (prior) and Bayes (posterior) net premiums, i.e., $k = 1$, for the compound aggregate model described above are given by,

respectively.

Proof. Expression (2.12) is obtained as follows,

from which we get the result. Observe that we have used the fact that the mean of the random variable $\Psi$ following the pdf $\pi_2(\psi|\gamma, \xi)$ is $1+(\gamma-1)/\xi$.

To get (2.13), we need to compute

Now, taking into account (2.7) and (2.8), the result is obtained directly. □

3.   The reinsurance model under independence

^[10] proposed the aggregate model described in the previous section to calculate reinsurance premiums in the collective risk model. For that, let $N_u(t) = \sum_{i = 1}^{N(t)}I(Y_i > u)$ be the number of claims exceeding the level $u$ during $[0, t]$. The latter expression is again a Poisson process with intensity $\lambda_u = \lambda \Pr(Y_i > u) = \lambda (\sigma/u)^{\psi}$, for $u\geq \sigma$.

Consider an excess-of-loss cover for the layer $b-a$ in excess of $a$. For a claim $Y_a$, exceeding the priority or deductible $a$, the reinsurer will cover the amount $Z_i = \min\left\{Y, b\right\}-a$. Our interest now lies on the random variable, the reinsurer's total claim amount during period $[0, t]$, $X(t) = \sum_{i = 1}^{N_{{a}}(t)}Z_i$. It is clear that $\left\{X(t)\right\}_{t\geq 0}$ is a compound Poisson process with intensity $\lambda_a$ and distribution function given by

Now, at time $T$, the practitioner wishes to predict the reinsurer's total claims cost for the next year, i.e. $X = X(T+1)-X(T)$. Therefore, we must first compute the so-called risk premium of order $k$, $k = 1, 2, \dots$, i.e., premiums based on non-central moments. Observe that the first moment is the so-called net risk premium (see, for instance, ^[25]). Then, we need

where $\Theta = (\Lambda, \Psi)$ is the vector of risk parameters and $\mu_k(\psi) = \mathbb{E}(Z_i^k|\psi)$.

Then, it is necessary to compute $\mu_k(\psi) = \mathbb{E}(Z_i^k|\psi)$. A closed-form expression for this $k$th moment is provided in ^[10] giving

The latter result can be written as a closed-form expression regarding the incomplete Beta function, as shown in the following proposition.

Proposition 3. It is verified that the non-central moments $\mu_k(\psi) = \mathbb{E}(Z_i^k|\psi)$ can be written as

with $\psi > k$, where $I_z(r, s)$ is the incomplete Beta function, given by

Proof. We have that

Let now the change of variable be $z = 1-a/y$, and we get that the last integral can be rewritten as

from which we have the result. □

The main idea of Hesselager's model (see ^[10]) is to consider that the parameters $\lambda$ and $\psi$ (often termed as risk profiles) are unknown and will be regarded as outcomes of random parameters $\Lambda$ and $\Psi$, respectively, and to treat the problem in a Bayesian framework. For that reason, the prior distributions of $\Lambda$ and $\Psi$ are required, and the premium provided in (3.1) is called risk premium, net risk premium when $k = 1$.

Proposition 4. The risk, collective (prior), and Bayes (posterior) net premiums, i.e., $k = 1$, are given by,

respectively and where $a\left\{\begin{array}{c} < \\ = \\ > \end{array}\right\}c$.

Proof. By taking $t = 1$, we have that $P(\lambda, \psi) = \mathbb{E}(N(t))\mathbb{E}(Y) = \lambda \mu_1(\psi)$. Now, $\mu_1(\psi)$ is obtained directly from (3.2) or (3.3) by taking $k = 1$, from which we easily get (3.1). To get (3.5), we have to compute ${P}(\nu, \tau, \gamma, \xi) = \mathbb{E}_p[P(\Lambda, \Psi)]$, where the expectation is taken with respect to $p(\lambda, \psi|\nu, \tau, \gamma, \xi, 1)$. Now, due to the hypothesis of independence between the risk profiles, we have that

where $\mathbb{E}_{\pi_1}(\Lambda) = \nu/\tau$. On the other hand, we have that

from which we get (3.5) after some simple algebra. Expression (3.6), which is the posterior expectation $P^{\ast}(\nu^{\ast}, \tau^{\ast}, \gamma^{\ast}, \xi^{\ast}, k) = \mathbb{E}_{\pi^{\ast}}[P(\Lambda, \Psi)]$, is obtained from (3.1) using the updated parameters provided in (2.7)–(2.10) and arranging the constant of normalization. □

As a special case of the model studied before, we now consider that $b\to \infty$. This is the basic excess-of-loss reinsurance model considered in [27, p. 11] and [1, p. 49]. In this case, an unlimited excess-of-loss treaty exists with retention $a > 0$, and the mean excess amount gives the expected amount to be paid by the reinsurer. Thus, $z > 0$ and the risk premium can be obtained from (3.1) as

Furthermore, from (3.5) and (3.6) and assuming that $b\to\infty$, the collective and Bayes premiums are given by,

An anonymous reviewer has pointed out this important detail: suppose the deductible $a$ is sufficiently large. The excess loss should be roughly distributed as the generalized Pareto distribution (Pickands-Balkema-de Haan theorem), which might be of more interest because the underlying distribution of individual severity can be assumed arbitrarily. See, for instance ^[28].

4.   Breaking the independence hypotheses

We focus here on the collective and Bayes net premiums for the aggregate amount of claims under a compound model, assuming some degree of dependence between the random variables $\Lambda$ and $\Psi$. In our case, it is challenging to implement this assumption in the first stage of the model, involving Poisson and Pareto random variables. As pointed out by ^[29], an attractive way to model the dependence between risks is via the dependence between the risk profiles. Thus, we translate the non-dependence to the second stage (the prior distributions of the risk profiles).

In the statistical literature, there exists a lot of research papers dealing with the construction of bivariate distributions (see, for example, ^[30,31,32]). Most of them are based on the idea of building a bivariate distribution with given marginal distributions, and the Sarmanov-Lee family of distributions ^[14,15] follows this idea. Specifically, assume $f_{1}(x_{1})$ and $f_{2}(x_{2})$ are univariate probability density functions for the continuous case or probability mass functions for the discrete one, with supports defined on subsets of $\mathbb{R}$ (which can be the total real line). Let $\phi_{i}(t)$, $i = 1, 2$, be bounded nonconstant functions such that ${\int}_{-\infty}^{\infty}\phi_{i}(t)f_{i}(t)\ dt = 0$ (they are usually named mixing functions). Then, provided that $\omega$ is a real number that satisfies the condition $1+\omega\ \phi_{1}(x_{1}) \phi_{2}(x_{2})\geq 0$ for all $x_{1}$ and $x_{2}$, the function given by

which is a genuine bivariate joint density (or a probability mass function in the discrete case) with specified marginals $f_{1}(x_{1})$ and $f_{2}(x_{2})$. It is well-known that the family given in (4.1) is a special case of the construction provided in ^[33] and includes some of the well-known Farlie-Gumbel-Morgenstern distributions as special cases. For details, see ^[34,35,36]. Observe that the dependence is captured through the mixing functions and the $\omega$ parameter.

Let us assume that a certain degree of dependence between the random variables $\Lambda$ and $\Psi$ exist, and specify a marginal prior gamma distribution for the first random variable and a translated Erlang distribution for the latter one, which, in both cases, are the conjugate a priori distribution of the Poisson and classical Pareto distributions. See the work of ^[37], where conjugate analysis for the Pareto distribution is studied.

We first need the following Lemma to obtain the prior distribution after assuming some degree of dependence.

Lemma 1. Let us consider the following functions

Then, it is verified that

if $\omega_1\leq \omega\leq \omega_2$, where

being

Proof. To see this, observe that because $d(\varphi_1(\lambda|\nu, \tau))/d\lambda < 0$, we have that $\varphi_1(\lambda|\nu, \tau)$ is a decreasing function on $\lambda$; also, since the support of $\lambda$ is in the interval $(0, \infty)$, we have that the range of variation of $\varphi_1(\lambda|\nu, \tau)$ is $(-\delta_1(\nu, \tau), 1-\delta_1(\nu, \tau))$. Using the same argument, we get that the range of variation of $\varphi_2(\psi|\gamma, \xi, k)$ is given by $(-\delta_2(\gamma, \xi, k), \exp(-k)-\delta_2(\gamma, \xi, k))$. Now we have that the expression $H(\nu, \tau, \gamma, \xi, \omega)\geq 0$ if

from which, taking into account that $0 < \delta_1(\nu, \tau) < 1$ and $0 < \delta_2(\gamma, \xi, k) < 1$, it is a simple exercise to see that the range of $\omega$ is given by $(\omega_1, \omega_2)$ provided above. □

The following result provides the prior distribution under the Sarmanov-Lee family of distributions.

Proposition 5. The expression given by

where $H(\nu, \tau, \gamma, \xi, k, \omega)$ is given in (4.2), defines a genuine prior bivariate distribution of $(\Lambda, \Psi)$ with marginal distributions given by $\pi_1(\lambda|\nu, \tau)$ and $\pi_2(\psi|\gamma, \xi, k)$, which always satisfies that $\omega_1\leq\omega\leq \omega_2$ and $\omega_1$ and $\omega_2$ are given by (4.3) and (4.4), respectively.

Proof. When $\omega_1 < \omega < \omega_2$, we have that $\pi(\lambda, \psi|\nu, \tau, \gamma, \xi, \omega) > 0$. Now, taking into account that

we get in a straightforward way that

Hence the result. □

As can be seen, the prior bivariate distribution (4.5) can be written as a linear combination of the product of gamma and Erlang distributions. Additionally, the special case $\omega = 0$ provides the prior bivariate distribution given by $\pi(\lambda, \psi|\nu, \tau, \gamma, \xi, k) = \pi_1(\lambda|\nu, \tau)\, \pi_2(\psi|\gamma, \xi, k)$, i.e. it represents the case in which we assume independence between $\Lambda_{c}$ and $\Psi$. Figure 1 shows the pdf given in (4.5) for parameters values given by $(\tau, \nu, \gamma, \xi) = (5, 10, 10, 4)$ and different values of $\omega$ moving in the range $[-1.24, 30.19]$ obtained from (4.3)-(4.4). The discrepancies between the different graphs are best appreciated in the contour plots given in Figure 1 below.

4.1.   Properties of the bivariate distribution

Some properties of the bivariate distribution proposed in this section are now studied.

Lemma 2. The following relations are verified,

Proof. See the Appendix. □

The next result shows the cross moment and correlation coefficient between the random variables $\Lambda_{c}$ and $\Psi$ when the prior distribution (4.5) is used.

Proposition 6. For the prior bivariate distribution (4.5), the cross moment and correlation coefficient between $\Lambda_{c}$ and $\Psi$ are given by

Proof. By using Theorems 1 and 2 in ^[15], we have that

where

$\mu_1$, $\mu_2$, $\sigma_1$ and $\sigma_2$, the means and standard deviations of the gamma pdf (2.4) and Erlang pdf (2.5), respectively. These are given by $\mu_1 = \nu/\tau$, $\sigma_1 = \sqrt{\nu}/\tau$, $\mu_2 = k+\gamma/\xi$, and $\sigma_2 = \gamma/\xi^2$.

Now, by using the expressions (4.6) and (4.7) provided in Lemma 2 and some straightforward algebra, we get the result. □

Observe that the correlation coefficient in (4.8) is directly proportional to $\omega$. The proportionality factor, given by $\sqrt{\nu\gamma}\delta_1(\nu+1, \tau)\delta_2(\gamma, \xi, k)/(\tau(1+\xi))$, is always positive, and the sign of the correlations will depend on the sign of $\omega$, which allows values in the range provided by (4.3)-(4.4). A study of the linear correlation coefficient as a function of the values of the hyperparameters is carried out in the next proposition.

Proposition 7. Assume that $\omega > 0$. Then, the following is verified:

(1) $\rho(\Lambda, \Psi)$ is an increasing (decreasing) function on $\nu$ if $\tau > (<)\left[\exp(1/(2\nu)-1\right]^{-1}$, for all $\zeta > 0$ and $\gamma\in\mathbb{N}$.

(2) $\rho(\Lambda, \Psi)$ is an increasing (decreasing) function on $\tau$ if $\nu > (<)\tau$, for all $\zeta > 0$ and $\gamma\in\mathbb{N}$.

(3) $\rho(\Lambda, \Psi)$ is an increasing (decreasing) function on $\gamma$ if $\zeta > (<)\left[\exp(1/(2\gamma)-1\right]^{-1}$, for all $\tau > 0$ and $\nu > 0$.

(4) $\rho(\Lambda, \Psi)$ is an increasing (decreasing) function on $\zeta$ if $\gamma > (<)\zeta$, for all $\tau > 0$ and $\nu > 0$.

Proof. See the Appendix. □

Figure 2 plots the correlation given as a function of $\omega$. The correlation coefficients seem limited to values between -0.15 and 0.15, coherent with the correlation in most empirical cases found in actuarial statistics.

The posterior distribution can be computed easily, and the next result shows that this posterior distribution can be written as a linear combination of the posterior distribution studied in Section 2.

Proposition 8. The posterior distribution of $(\Lambda_{c}, \Psi)$ given the sample information $\aleph$ is a linear combination of products of posterior univariate gamma and Erlang densities, given by

where $\kappa^{\ast}\equiv\kappa(\gamma^{\ast}, \xi^{\ast}+1, k)$ is given in (2.11) and the normalization constant, $\kappa_p\equiv \kappa_p(\gamma^{\ast}, \xi^{\ast}, k)$, satisfies

Proof. The result follows after some intensive algebra by using (2.3), (4.5), Bayes' Theorem and arranging parameters. □

4.2.   Premiums under the dependence model

Expressions given in (4.5) and (4.9) can be used now to get the collective and Bayes net premiums under the non-independent model. These are shown in the following result.

Proposition 9. Under the bivariate model with dependence, the prior (collective) and Bayesian (posterior) premiums are given by,

where $P(\nu, \tau, \gamma, \xi)$ is the collective premium under the different models studied previously, the basic specific model (expression given in (2.12)), the reinsurance model with a layer (expression given in (3.5)), and the reinsurance model without a layer (expression given in (3.7)). Alternatively, $P^{\ast}(\nu^{\ast}, \tau^{\ast}, \gamma^{\ast}, \xi^{\ast}, 1)$ is the net Bayes premium under the models studied, expressions (2.13), (3.6), and (3.8).

Proof. See the Appendix. □

In the following corollary, we give an alternative expression for (4.10).

Corollary 1. The collective premium given in (4.10) can be rewritten as

Proof. See the Appendix. □

We get the corresponding collective and Bayes premium under the independent model by taking $\omega = 0$ in (4.10) and (4.11).

Note that the Bayesian premium has been calculated sequentially rather than in a batch, which means that the count of the number of losses may not be a homogeneous Poisson over time.

5.   Numerical applications

In this section, two numerical examples will be considered. We first used the same dataset for comparison in ^[10,11]. The second dataset corresponds to catastrophe data regarding the total damage caused by a storm during 1956–1977, which is displayed in Table 1. The amounts and the year of occurrence are exhibited in this table. We have only considered losses below 6 billion dollars. A comprehensive version of this dataset is shown in [38, p. 128]. For this reason, in our analysis, we only consider 10 entries.

The elicited hyperparameter values of the prior distributions are shown in Table 2, together with the bounds of the $\omega$ parameter and other statistics needed to obtain the collective and Bayesian net premiums. The maximum likelihood estimation method requires knowledge of the distribution obtained from (2.1), which is not simple; hence, the hyperparameters were estimated by the method of moments, where for the parameters $(\nu, \tau)$, the sample mean and index of dispersion (ratio of variance to the mean assumed to be approximately 1.25) were used to calculate the estimates. Furthermore, for the parameters $(\xi, \gamma)$, we used the sample mean with the coefficient of variation and the approximation of the value of $\gamma$ to the next integer. The sample sizes, the values of the layers, and the corresponding time periods for each data set are also given.

The collective and Bayesian net premiums computed using (4.10) and (4.11) are shown in Table 3 for the two examples considered and the basic collective model. Observe that the case $\omega = 0$ reduces to the model where independence between the risk profiles is assumed. The correlation coefficient obtained for each value of $\omega$ obtained via (4.8) is also shown in this table. It can be observed that incorporating dependence has a low impact on both premiums, as judged by the changes in $\omega$ when this dependence parameter increases or decreases (i.e., more separated from the hypotheses of independence).

The collective and Bayesian net premiums, including a layer of 1.5 million in excess of the priority $a$ for Example 1 and a layer of 2 million in excess of the priority $a$ for Example 2, are shown in Tables 4 and 5, respectively. Again, the case $\omega = 0$ reduces to the model where independence between the risk profiles is assumed, and the correlation coefficient is also included. It can be noticed that increasing the value of the deductible $a$ decreases the value of both premiums. In cases where the deductible is 1.5 million, the Bayes premium is higher than the collective premiums in both examples.

In addition, the collective and Bayesian net premiums without including a layer in excess of the priority $a$ for both examples are shown in Tables 6 and 7, respectively. Again, the case $\omega = 0$ reduces to the model where independence between the risk profiles is assumed, and the correlation coefficient is also shown in these tables. It can be noticed that increasing the value of the deductible $a$ decreases the value of both premiums. In cases where the deductible is 1.5 million, the Bayes premium is higher than the collective premiums in both examples. In Example 1, it can be observed that increasing the value of the deductible $a$ decreases the value of both premiums. However, this is not verified in the second example, where the Bayes premiums increase for the deductible value $a = 1.5$ million and decrease again for $a = 2.2$. Also, the collective premiums increase for the deductible of $a = 2.2$ million.

Even though, for all the tables above, the predicted correlation coefficients seem to be minor, this is not a trivial issue since it could have a significant impact when changes of scales in the monetary units are considered; therefore, the analysis we have carried out in this work shows that although the correlation is usually low, the impact of assuming independence might have potential financial consequences. It is important to note that the dependence hypothesis is established between the risk profiles and not directly between the number of claims and the corresponding severity. However, given the mathematical difficulty that the latter scenario would entail, this can be replaced by transferring this hypothesis to the risk profiles, which retain characteristics of the former case. Even though when calculating bonus-malus premiums in the automobile insurance portfolio, independence between the number of claims and the amount of the claim is usually assumed, some recent works such as ^[39] have discussed the calculation of bonus-malus premiums for the compound loss by incorporating a certain degree of dependence between the frequency and severity. However, in practice, the correlation between both variables is usually low. For example, by considering the automobile insurance portfolio used in ^[40] and taken from ^[41] that corresponds to the years 2004-05 and contains data of 67,856 insureds, of whom 4,624 filed a claim, i.e., a positive claims amount, the correlation between the number of claims and the amount for the whole portfolio is 0.4818. Nevertheless, when it is restricted to the 4,624 policyholders that have declared at least one claim, the correlation drops to 0.0762, consistent with the values shown in Table 3.

6.   Final comments and future works

We have presented an excess-of-loss reinsurance model derived from the compound Poisson distribution. The classical Pareto distribution was assumed for the severity of claims, and dependence between the number of claims and claims size was placed on the risk profiles. In this work, initially, we considered independence between the risk profiles whose parameters follow gamma and shifted Erlang distributions. Next, the hypothesis of independence was broken by assuming some degree of dependence between the random variables associated with these risk profiles via a prior bivariate distribution based on the Sarmanov-Lee probabilistic family. This model was used to study the effect of incorporating some degree of dependence of these variables on the collective and Bayesian net premiums. The results obtained lead us to conclude that even at moderate levels of correlation between the risk profiles, these aforementioned reinsurance premiums are more sensitive than those computed, assuming the independence of the risk profiles.

The methodology proposed in this paper is suitable for application in the collective model, in which the two relevant variables used in calculating premiums are used simultaneously: the number of claims and the amount associated with them. The difficulty of working simultaneously with these two variables is well-known, especially when obtaining the total claim amount distribution. Two scenarios that should be considered are, on the one hand, the study of the sensitivity of the obtained premiums to slight modifications in the choice of the prior distribution. To that end, we are discussing carrying out a Bayesian robustness analysis (local or global). On the other hand, choosing the Sarmanov-Lee family to construct the prior distribution that allows for a certain degree of dependence between the two parameters can be very rigid (it has the advantage of being almost conjugate) and could be replaced by a bivariate distribution based on another type of copula.

Use of Generative-AI tools declaration

The author declares he has not used Artificial Intelligence (AI) tools in the creation of this article.

Acknowledgments

EGD was partially funded by grant PID2021-127989OB-I00 (Ministerio de Economía y Competitividad, Spain).

Conflict of interest

Prof. Emilio Gómez-Déniz is the Guest Editor of special issue "Advances in Probability Distributions and Social Science Statistics" for AIMS Mathematics. Prof. Emilio Gómez-Déniz was not involved in the editorial review and the decision to publish this article.

The author declares no conflict of interest.