Diffusion Approximations of the Ruin Probability for the Insurer–Reinsurer Model Driven by a Renewal Process

Abstract

We introduce here a diffusion-type approximation of the ruin probability both in finite and infinite time for a two-dimensional risk process, where claims and premiums are shared with a predetermined proportion. This type of process is often called the insurer–reinsurer model. We assume that the flow of claims is governed by a general renewal process. A simple ruin probability formula for the model is known only in infinite time for the special case of the Poisson process and exponentially distributed claims. Therefore, there is a need for simple analytical approximations. In the literature, in the infinite-time case, for the Poisson process, a De Vylder-type approximation has already been introduced. The idea of the diffusion approximation presented here is based on the weak convergence of stochastic processes, which enables one to replace the original risk process with a Brownian motion with drift. By applying this idea to the insurer–reinsurer model, we obtain simple ruin probability approximations for both finite and infinite time. We check the usefulness of the approximations by studying several claim amount distributions and comparing the results with the De Vylder-type approximation and Monte Carlo simulations. All the results show that the proposed approximations are promising and often yield small relative errors.

1. Introduction

Risk theory in general and ruin probabilities in particular have been an active area of research since the classical Cramér–Lundberg model introduced in 1903 by the Swedish actuary Filip Lundberg (Asmussen and Albrecher 2010). The Cramér–Lundberg model describes the situation of an insurance company that experiences two opposing cash flows: incoming premiums and outgoing claims. The traditional approach in risk theory is to study the probability of ruin, i.e., the probability that the risk process will ever go below zero. The model is described in terms of the classical risk process U, which is defined by

$$U\left(t\right)=u+ct-\sum _{k=1}^{N\left(t\right)}{X}_k,$$

(1)

where $u\ge 0$ denotes the initial capital, c is a positive premium rate, $N={\left(N\left(t\right)\right)}_{t\ge 0}$ describing a claim flow is a Poisson process with intensity $\lambda$ independent of $\left(X_k\right)$ and claim amounts ${\left(X_k\right)}_{k=1}^{\infty }$ form a sequence of positive independent identically distributed (i.i.d.) random variables, with mean value $\mu$ and variance ${\sigma }^2$. The constant c can be written as $c=(1+\theta )\mu \lambda$, where $\theta >0$ is called the relative safety loading (Grandell 1991; Rolski et al. 1999).

One of the key issues of the collective risk theory concerns calculating the ruin probability, i.e., the probability that the risk process becomes negative. We distinguish between two types of ruin probability. The ruin probability $\psi (u,t)$ in finite time t (or within a finite horizon) of a company described by the risk process (1) is given by

$$\psi (u,t)=P\left(U\right(s)<0\mathrm{for}\mathrm{some}s\le t),0<t<\infty ,u\ge 0.$$

The ruin probability $\psi \left(u\right)$ in infinite time can be defined as $\psi \left(u\right)=\psi (u,\infty )$.

Ruin theory is believed to be important for modern risk management. For solvency purposes, the ruin probability can be used as a rough approximation of the insolvency. Moreover, fixing it at an acceptably low level, the needed capital and the rate of the flow of premiums can be estimated. It can also serve as a useful tool in long-range planning for the use of insurer’s funds. In addition, ruin theory has fruitful methodological links and applications to other fields of applied probability, such as queueing theory and mathematical finance (Asmussen and Albrecher 2010).

The ruin probabilities in infinite and finite time, even for the classical risk process, can only be calculated for a few special cases of the claim amount distribution. For the infinite horizon case, there are well-known elementary results for zero initial capital, and the exponential and mixture of two exponential claim amount distributions, see Grandell (1991); Panjer and Willmot (1992). For the results for general phase-type distributions, in particular for mixture of n exponential distributions, see Asmussen and Albrecher (2010). For the finite horizon case, the only convenient “semi-elementary” formula (involving only a simple integral) exists for the exponential distribution (Grandell 1991; Rolski et al. 1999). However, this case can always be approximated by the Monte Carlo method. Therefore, finding a reliable approximation, especially in the ultimate case, is really important from a practical point of view. Various approximations both in finite and infinite time have already been discussed in the literature, let us just mention a few: the classical Cramér–Lundberg approximation (with time-dependent version being Segerdahl’s), saddlepoint approximations (for example Arfwedson’s), De Vylder (often regarded as the best among “simple” approximations), diffusion approximations, and large claim approximations, see Arfwedson (1955); Asmussen and Albrecher (2010); Barndorff-Nielsen and Schmidli (1995); De Vylder (1978); Grandell (2000). For a review and numerical comparison of approximations, see, e.g., Burnecki et al. (2005); Grandell (1991).

For the finite horizon ruin probabilities, crude Monte Carlo method can be always applied. However, if the considered ruin probabilities values are very small, such as less than $1/1000$ (a rare type event), it is expensive in terms of computer time to obtain reasonably precise estimates of the ruin probability. Variance reduction techniques can improve this situation. In this context, we mention importance sampling techniques, where the idea is to sample after an exponential tilt of the probability measure Blanchet and Lam (2011); Collamore (2002).

Andersen (1957) in a paper to the 1957 International Congress of actuaries in New York proposed a generalization of the classical risk process. Instead of assuming just exponentially distributed independent interoccurrence (waiting) times, he introduced a more general distribution function but retained the assumption of independence. The Sparre Andersen model assumes that the interoccurrence times are i.i.d. random variables, hence the counting process becomes a renewal process. The ruin problem for an infinite time was already considered by Sparre Andersen himself. For the classical finite time ruin probability treatment within the model, see, e.g., Grandell (1991); Thorin (1971). The renewal process model brings much more flexibility in modelling the flow of claims than the Poisson process, see, e.g., Burnecki and Weron (2005), where the renewal process with log-normal waiting times was found superior to the Poisson process in modelling losses from the Property Claim Services dataset.

The first modern treatment of diffusion approximation in risk theory, based on weak convergence, was presented by Iglehart (1969), who already discussed it in the renewal claim counting process framework; see also Bohman (1972); Gluckman (1970); Grandell (1972). The idea is to let the number of claims grow in a unit time interval and at the same time make the claim sizes smaller in such a way that the risk process converges weakly to a Brownian diffusion. It was later extended to the weak convergence of compound renewal processes subject to discounting, with deterministic or stochastic interest rate, see Braun (1986); Garrido (1988); Harrison (1977). Apart from the weak convergence, another common alternate motivation for diffusion claim models was to see them as solutions of stochastic differential equations; see the results of Abikhalil (1986); Ruohonen (1980); Garrido (1988, 1989).

This idea of the Brownian diffusion approximation was later extended to the stable diffusion approximation which accounts for heavy-tailed claims belonging to the domain of attraction of the $\alpha$-stable law ($1<\alpha <2$), see Furrer et al. (1997). In Burnecki (2000) it was proved that only self-similar processes with stationary increments appear as weak limits of the risk processes and, conversely, every finite mean H-self-similar process with stationary increments, can result as a weak approximation. We also note that in queuing theory the diffusion approximation is known as the “heavy traffic approximation”, see Asmussen and Albrecher (2010).

Recently, multidimensional risk processes have been introduced in the literature to account for multiple lines of business of an insurance company and collaborating insurance companies. The ruin probability can now be defined in several ways, e.g., when all lines or all companies are ruined or at least one. Collamore (1996) was among the first to introduce the multidimensional ruin problem for light-tailed claims and general ruin sets. Multidimensional heavy-tailed processes were first studied by Hult et al. (2005). They mainly concentrated on multivariate regularly varying random walks and calculated sharp boundaries for the asymptotic ruin probability.

Since different risks usually have an effect on a few lines of business at the same time, the statistical dependence among claims in these lines should be taken into account. The multidimensional risk process was specialized to the two-dimensional case with claims shared with a predetermined proportion in (Avram et al. 2008a, 2008b). This case is usually referred to as the insurer–reinsurer model, as it well describes the quota share proportional treaty. It can also be used to model two branches of the same insurance company. The ruin occurs here if one or both companies go bankrupt. The former case, which is more interesting from a practical point of view, is usually analyzed, and the latter can be obtained from the former in a straightforward way. To the best of the authors’ knowledge, the only simple ruin probability formulas for the insurer–reinsurer model were provided for exponentially distributed claims in Avram et al. (2008a) (by explicitly inverting the Laplace transform) and later in Burnecki et al. (2021) (by means of a change-of-measure technique). Based on the latter result, a De Vylder-type approximation of the ruin probability for the insurer–reinsurer model was introduced in Burnecki et al. (2019) for general claim amount distributions with a finite third moment. Another type of dependence was studied in Behme et al. (2020), where the link was established by a random bipartite network. Badescu et al. (2011) introduced an extension to a system of two insurers, where the first insurer is experiences claims arising from two independent compound Poisson processes and the second insurer covers a proportion of the claims. In Michna (2020), a model driven by a general spectrally positive or negative Lévy process was investigated, see also Avram et al. (2008b).

In this paper, we introduce a diffusion-type of approximation of the ruin probability for the insurer–reinsurer model, where the flow of claims is controlled by natural generalization of the Poisson process, namely for a renewal process. We calculate the ruin probability formulas for both finite and infinite time and study the performance of the approximations. The article is organized as follows. In Section 2, the model is presented and ruin probabilities are defined. In Section 3 we recall the idea of the diffusion approximation, which is based on the weak convergence of the risk process to a Brownian diffusion. The main results are presented in Section 4. For the insurer–reinsurer model driven by the renewal process, we derive ruin probability formulas for both finite- and infinite-time horizons. The special case of the Poisson process is also discussed. In order to check the quality of the approximations, in Section 5 we perform a Monte Carlo simulation study by considering mixture of two exponential, gamma, Weibull (light-tailed cases), lognormal, Pareto and Burr (heavy-tailed cases) claim amount distributions. We compute relative errors of the approximations with respect to the ruin probability values calculated by means of Monte Carlo simulations for the infinite-time case taking sufficiently long time horizons. For the infinite-time case, we compare also the efficiency of the approximation with the De Vylder’s-type approximation recently introduced in the literature. Section 6 summarizes our results.

2. Insurer–Reinsurer Model

In this section, we analyse a two-dimensional model that describes capitals of insurance and reinsurance companies or two business lines of the same insurance company with proportional claim sizes in the renewal process framework. The model can be considered as a system of two processes ($U_1\left(t\right)$, $U_2\left(t\right)$) defined as follows:

$$\begin{array}{c}\hfill \left(\begin{array}{c}{U}_1\left(t\right)\\ U_2\left(t\right)\end{array}\right)=\left(\begin{array}{c}{u}_1\\ u_2\end{array}\right)+\left(\begin{array}{c}\delta (1+{\theta }_1)\\ (1-\delta )(1+{\theta }_2)\end{array}\right)\alpha \mu t-\left(\begin{array}{c}\delta \\ 1-\delta \end{array}\right)\sum _{k=1}^{N\left(t\right)}{X}_k,\end{array}$$

(2)

where $N\left(t\right)$ is a renewal process with inter-arrival times with mean $1/\alpha$ and claim amount sequence ${\left(X_k\right)}_{k\ge 1}$ consists of i.i.d. random variables with mean $\mu$ and variance ${\sigma }^2$.

In that system, the claim arrival process is common for both lines. We call it the insurer–reinsurer model since its primary application is to describe the evolution of capital for an insurer and a reinsurer under the quota share reinsurance, where the premiums and claims are divided between the insurer and the reinsurer with constant proportions $\delta$ and $1-\delta$ with $\delta \in (0,1)$, respectively; see, e.g., Avram et al. (2008b); Foss et al. (2017); Michna (2020).

For the insurer–reinsurer model, we define the probability of ruin if at least one of the two companies is ruined, which is also called ’or’ probability in the literature Avram et al. (2008a). The time of ruin is as follows:

$$\begin{array}{c}\hfill \tau (u_1,u_2)=inf\{t\ge 0:U_1\left(t\right)<0\mathrm{or}{U}_2\left(t\right)<0\}.\end{array}$$

(3)

The ruin probability in finite time t can be expressed as:

$$\begin{array}{c}\hfill \psi (u_1,u_2,t)=\mathbf{P}(\tau (u_1,u_2)<t)\end{array}$$

(4)

and the ruin probability in infinite time is given by:

$$\begin{array}{c}\hfill \psi (u_1,u_2)=\mathbf{P}(\tau (u_1,u_2)<\infty ).\end{array}$$

(5)

In this work, we are interested in analytical approximations of the ruin probabilities given by Equations (4) and (5) for the insurer–reinsurer model. Let us observe now that by rescaling the two stochastic processes $U_1\left(t\right)$ and $U_2\left(t\right)$ of the insurer–reinsurer model (2) by factors ${\delta }^{-1}$ and ${(1-\delta )}^{-1}$, respectively, we arrive at a new system of risk processes denoted by $(X_1\left(t\right),X_2\left(t\right))$:

$$\begin{array}{c}\hfill \left(\begin{array}{c}{X}_1\left(t\right)\\ X_2\left(t\right)\end{array}\right)=\left(\begin{array}{c}{U}_1\left(t\right)/\delta \\ U_2\left(t\right)/(1-\delta )\end{array}\right)=\left(\begin{array}{c}{x}_1\\ x_2\end{array}\right)+\left(\begin{array}{c}1+{\theta }_1\\ 1+{\theta }_2\end{array}\right)\alpha \mu t-\sum _{k=1}^{N\left(t\right)}{X}_k.\end{array}$$

(6)

Please note that, due to scaling invariance, the ruin probabilities in finite and infinite time for the normalized insurer–reinsurer risk process $(X_1\left(t\right),X_2\left(t\right))$ are exactly the same as for the original process $(U_1\left(t\right),U_2\left(t\right))$. Moreover, due to the higher acquisition and administration costs of the insurer, it is natural to expect that the premium rate for the insurer is higher than for the reinsurer and therefore from now on we assume that ${\theta }_1>{\theta }_2$, see Burnecki et al. (2021).

In this paper we will also rely on the fact that due to the specific construction of the considered two-dimensional risk process, in which we assume that premiums and claims are split between the insurer and the reinsurer with a fixed proportion, the two-dimensional problem can be reduced to the one-dimensional one, see Avram et al. (2008b).

3. Classical Diffusion Approximation

In this section we present the idea of replacing the risk process with a Brownian motion with drift. First, let us recall the definition of weak convergence of stochastic processes, see, e.g., Grandell (1991).

Definition 1.

Let $D=D[0,\infty )$ be the space of càdlàg functions, i.e., everywhere right-continuous and with left limits everywhere on $[0,\infty )$. A stochastic process $X={\left(X\left(t\right)\right)}_{t\ge 0}$ is said to be in D if all its realizations are in D. A sequence ${\left(X^{\left(n\right)}\right)}_{n\in \mathbf{N}}$ of stochastic processes in D is said to converge weakly in the Skorokhod $J_1$ topology to a stochastic process X if for every bounded continuous functional f on D it follows that

$$\underset{n\to \infty }{lim}Ef\left(X^{\left(n\right)}\right)=Ef\left(X\right).$$

In this case, we write $X^{\left(n\right)}\Rightarrow X$. The weak convergence of $X^{\left(n\right)}$ to X implies, for example, convergence of the finite-dimensional distributions provided that the limit process X is continuous in probability, and that ${inf}_{0\le t\le t_0}{X}^{\left(n\right)}\left(t\right)\stackrel{d}{\to }{inf}_{0\le t\le t_0}X\left(t\right)$ for any $t_0<\infty$.

Let us now formally define a renewal process.

Definition 2.

Let us assume that the inter-arrival (waiting) times ${\left(M_k\right)}_{k\in \mathbf{N}}$ are assumed to be independent positive random variables with mean $1/\alpha$ and variance ${\sigma }_M^2$. The renewal process N is defined as:

$$N\left(t\right)=max\left\{n:\sum _{k=1}^n{M}_k\le t\right\}.$$

It is well-known that for such a process,

$$\underset{t\to \infty }{lim}\frac{Var\left(N\right(t\left)\right)}{t}={\sigma }_M^2{\alpha }^3$$

and that

$$\frac{N\left(nt\right)-\alpha nt}{\sqrt{n}}\Rightarrow {\sigma }_M{\alpha }^{3/2}W\left(t\right),\mathrm{as}n\to \infty ,$$

where W is a standard Brownian motion, see Furrer et al. (1997); Grandell (1991).

If

$$\underset{n\to \infty }{lim}\left(c_n-\alpha n\frac{\mu }{\sqrt{n}}\right)=c-\alpha \mu ,$$

then

$$u+c_{n}t-\frac{1}{\sqrt{n}}\sum _{k=1}^{N\left(nt\right)}{X}_k\Rightarrow u+(c-\alpha \mu )t-\sqrt{\alpha {\sigma }^2+{\mu }^2{\sigma }_M^2{\alpha }^3}W\left(t\right).$$

For the special case of the Poisson process the diffusion approximation is equivalent to replacing the risk process with the Brownian diffusion with drift by matching its first two moments.

As a result, we obtain that:

$$\begin{array}{c}\psi (u,t)\approx \hfill \\ \mathbf{P}\left\{\underset{0\le s\le t}{inf}\left(u+\theta \alpha \mu s-\sqrt{\alpha {\sigma }^2+{\mu }^2{\sigma }_M^2{\alpha }^3}W\left(s\right)\right)<0\right\}=\overline{\mathsf{\Phi }}\left(\frac{u+\theta \alpha \mu t}{\sqrt{(\alpha {\sigma }^2+{\mu }^2{\sigma }_M^2{\alpha }^3)t}}\right)\hfill \\ +exp\left\{-\frac{2\theta \mu u}{{\sigma }^2+{\mu }^2{\sigma }_M^2{\alpha }^2}\right\}\mathsf{\Phi }\left(\frac{-u+\theta \alpha \mu t}{\sqrt{(\alpha {\sigma }^2+{\mu }^2{\sigma }_M^2{\alpha }^3)t}}\right)\hfill \end{array}$$

(7)

and

$$\psi \left(u\right)\approx \mathbf{P}\left\{\underset{s\ge 0}{inf}\left(u+\theta \alpha \mu s-\sqrt{\alpha {\sigma }^2+{\mu }^2{\sigma }_M^2{\alpha }^3}W\left(s\right)\right)<0\right\}=exp\left\{-\frac{2\theta \mu u}{{\sigma }^2+{\mu }^2{\sigma }_M^2{\alpha }^2}\right\},$$

(8)

where $\mathsf{\Phi }$ is the standard normal distribution and $\overline{\mathsf{\Phi }}=1-\mathsf{\Phi }$, e.g., Asmussen and Albrecher (2010); Grandell (1991). To conclude, the approximation is obtained as a limit of a specifically constructed sequence of risk processes out of the original one. The limit is in the weak convergence of the processes sense.

4. Diffusion Approximations for the Insurer–Reinsurer Model

We derive here diffusion approximations for the two-dimensional insurer–reinsurer model, with claims belonging to the domain of attraction of Gaussian law. Let us recall that since both the insurer and the reinsurer process are driven by the same counting process and claims are split by a fixed proportion between the two parties, the time horizon can be divided into intervals where one process is above the other. This allows one to reduce the problem of two dimensions to a single dimension. In particular, if $x_1>x_2$, the ruin is fully controlled by the reinsurer process (it is always below the insurer since we assumed that ${\theta }_1>{\theta }_2$), hence the ruin probability for the insurer–reinsurer model is simply the ruin probability of the one-dimensional process $X_2\left(t\right)$ describing the capital of the reinsurer. Therefore, we concentrate here on the (non-degenerate) case $x_1<x_2$ and by T we denote the transition time point before which the ruin is controlled by the insurer and after which the ruin is triggered by the reinsurer. The transition time $T=T(x_1,x_2)=(x_2-x_1)/\left(\mu \alpha ({\theta }_1-{\theta }_2)\right)$, for details, see Avram et al. (2008b); Burnecki et al. (2021). When considering the ruin probability in finite time t we will also assume that $t>T$. Otherwise, the two-dimensional ruin probability problem is simply reduced to the finite-time ruin probability of the insurer for t being the time horizon.

Theorem 1.

Let $(X_1\left(t\right),X_2\left(t\right))$ be the insurer–reinsurer model given by Equation (6) and $T=T(x_1,x_2)=(x_2-x_1)/\left(\mu \alpha ({\theta }_1-{\theta }_2)\right)$, where $x_1<x_2$. For the probability of ruin in finite time t, we also assume $t>T$. Then,

$$\begin{array}{ccc}\hfill & \psi (x_1,x_2,t)\approx \hfill & \\ \hfill & =1-{\int }_0^{\infty }\left\{\mathsf{\Phi }\left(\frac{z+{\gamma }_2(t-T)}{\sqrt{\zeta (t-T)}}\right)-exp\left\{-\frac{2{\gamma }_{2}z}{\zeta }\right\}\mathsf{\Phi }\left(\frac{-z+{\gamma }_2(t-T)}{\sqrt{{\gamma }_2(t-T)}}\right)\right\}·\hfill & \hfill \\ \hfill & ·\left\{\frac{1}{\sqrt{\zeta T}}{\mathsf{\Phi }}^{\prime }\left(\frac{z-(x_1+{\gamma }_{1}T)}{\sqrt{\zeta T}}\right)-g\left(x_1,0,z,{\gamma }_1/\zeta ,\zeta T\right)\right\}dz\hfill \end{array}$$

(9)

and

$$\begin{array}{ccc}\hfill & \psi (x_1,x_2)\approx \hfill & \\ \hfill & =1-{\int }_0^{\infty }\left(1-exp\left\{-\frac{2{\gamma }_{2}z}{\zeta }\right\}\right)\left\{\frac{1}{\sqrt{\zeta T}}{\mathsf{\Phi }}^{\prime }\left(\frac{z-(x_1+{\gamma }_{1}T)}{\sqrt{\zeta T}}\right)-g\left(x_1,0,z,{\gamma }_1/\zeta ,\zeta T\right)\right\}dz,\hfill \end{array}$$

(10)

where ${\gamma }_i={\theta }_i\alpha \mu$ ($i=1,2$) and $\zeta =\alpha {\sigma }^2+{\mu }^2{\sigma }_M^2{\alpha }^3$.

Proof. 

First, let us now focus on the ruin probability in finite time $t>T$. By rewriting the general formula for the infinite-time ruin probability $\psi (x_1,x_2)$ for the insurer–reinsurer model given in Proposition 1 in Avram et al. (2008b), we obtain for the finite horizon case that:

$$\psi (x_1,x_2,t)=1-{\int }_0^{\infty }\mathbf{P}\left(\underset{0\le s\le t-T}{inf}{X}_2\left(s\right)-x_2+z>0\right)\mathbf{P}\left(\underset{0\le s\le T}{inf}{X}_1\left(s\right)>0,X_1\left(T\right)\in dz\right).$$

(11)

The first probability under the integral can be approximated by using the diffusion approximation presented in Section 3:

$$\begin{array}{cc}\hfill & \mathbf{P}\left(\underset{0\le s\le t-T}{inf}{X}_2\left(s\right)-x_2+z>0\right)\approx \mathsf{\Phi }\left(\frac{z+{\theta }_2\alpha \mu (t-T)}{\sqrt{(\alpha {\sigma }^2+{\mu }^2{\sigma }_M^2{\alpha }^3)(t-T)}}\right)\hfill \\ \hfill & -exp\left\{-\frac{2{\theta }_2\mu z}{{\sigma }^2+{\mu }^2{\sigma }_M^2{\alpha }^2}\right\}\mathsf{\Phi }\left(\frac{-z+{\theta }_2\alpha \mu (t-T)}{\sqrt{(\alpha {\sigma }^2+{\mu }^2{\sigma }_M^2{\alpha }^3)(t-T)}}\right).\hfill \end{array}$$

(12)

To approximate the second probability under the integral, we will apply the following result for the Brownian motion with drift presented in Borodin and Salminen (2002) (Equation (1.1.8), p. 257)):

$$\begin{array}{cc}\hfill & \mathbf{P}\left(\underset{0\le s\le T}{inf}x+\gamma s+W\left(s\right)\le y,x+\gamma T+W\left(T\right)\in dz\right)\hfill \\ \hfill & =\frac{1}{\sqrt{2\pi T}}exp\left\{\gamma (z-x)-\frac{{\gamma }^{2}T}{2}-\frac{{\left(\right|z-y|+x-y)}^2}{2T}\right\}dz=:g(x,y,z,\gamma ,T)dz.\hfill \end{array}$$

(13)

To simplify the presentation of subsequent calculations, we introduce the following notation: ${\gamma }_i={\theta }_i\alpha \mu$ ($i=1,2$) and $\zeta =\alpha {\sigma }^2+{\mu }^2{\sigma }_M^2{\alpha }^3$. It is easy to check that:

$$\begin{array}{cc}\hfill & \mathbf{P}\left(\underset{0\le s\le T}{inf}{X}_1\left(s\right)\le 0,X_1\left(T\right)\in dz\right)\approx \hfill \\ \hfill & \approx \mathbf{P}\left(\underset{0\le s\le T}{inf}x+{\gamma }_{1}s-\sqrt{\zeta }W\left(s\right)\le y,x+{\gamma }_{1}T-\sqrt{\zeta }W\left(T\right)\in dz\right)=\hfill \\ \hfill & =\mathbf{P}\left(\underset{0\le s\le T}{inf}x+{\gamma }_{1}s+\sqrt{\zeta }W\left(s\right)\le y,x+{\gamma }_{1}T+\sqrt{\zeta }W\left(T\right)\in dz\right)\stackrel{|\zeta s=s^{\prime }|}{=}\hfill \\ \hfill & =\mathbf{P}\left(\underset{0\le s^{\prime }/\zeta \le T}{inf}x+\frac{{\gamma }_1{s}^{\prime }}{\zeta }+W\left(s^{\prime }\right)\le y,x+\frac{{\gamma }_1\zeta T}{\zeta }+W\left(\zeta T\right)\in dz\right)=\hfill \\ \hfill & =\mathbf{P}\left(\underset{0\le s^{\prime }\le \zeta T}{inf}x+\frac{{\gamma }_1{s}^{\prime }}{\zeta }+W\left(s^{\prime }\right)\le y,x+\frac{{\gamma }_1\zeta T}{\zeta }+W\left(\zeta T\right)\in dz\right)\stackrel{|\zeta T=T^{\prime }|}{=}\hfill \\ \hfill & =\mathbf{P}\left(\underset{0\le s^{\prime }\le T^{\prime }}{inf}x+\frac{{\gamma }_1{s}^{\prime }}{\zeta }+W\left(s^{\prime }\right)\le y,x+\frac{{\gamma }_1{T}^{\prime }}{\zeta }+W\left(T^{\prime }\right)\in dz\right)=\hfill \\ \hfill & =\frac{1}{\sqrt{2\pi T^{\prime }}}exp\left\{\frac{{\gamma }_1}{\zeta }(z-x)-\frac{{\gamma }_1^2{T}^{\prime }}{2{\zeta }^2}-\frac{{\left(\right|z-y|+x-y)}^2}{2T^{\prime }}\right\}dz\hfill \\ \hfill & =\frac{1}{\sqrt{2\pi \zeta T}}exp\left\{\frac{{\gamma }_1}{\zeta }(z-x)-\frac{{\gamma }_1^{2}T}{2\zeta }-\frac{{\left(\right|z-y|+x-y)}^2}{2\zeta T}\right\}dz\hfill \\ \hfill & =g(x,y,z,{\gamma }_1/\zeta ,\zeta T)dz.\hfill \end{array}$$

(14)

Therefore, we can rewrite the second probability under the integral in Equation (11) as:

$$\begin{array}{cc}\hfill & \mathbf{P}\left(\underset{0\le s\le T}{inf}{X}_1\left(s\right)>0,X_1\left(T\right)\in dz\right)\approx \hfill \\ & \approx \mathbf{P}\left(\underset{0\le s\le T}{inf}{x}_1+{\gamma }_{1}s-\sqrt{\zeta }W\left(s\right)>0,x_1+{\gamma }_{1}T-\sqrt{\zeta }W\left(T\right)\in dz\right)=\hfill \\ \hfill & =\mathbf{P}\left(x_1+{\gamma }_{1}T-\sqrt{\zeta }W\left(T\right)\in dz\right)-\hfill \\ \hfill & -\mathbf{P}\left(\underset{0\le s\le T}{inf}{x}_1+{\gamma }_{1}s-\sqrt{\zeta }W\left(s\right)\le 0,x_1+{\gamma }_{1}T-\sqrt{\zeta }W\left(T\right)\in dz\right)=\hfill \\ \hfill & =\frac{1}{\sqrt{\zeta T}}{\mathsf{\Phi }}^{\prime }\left(\frac{z-(x_1+{\gamma }_{1}T)}{\sqrt{\zeta T}}\right)-g(x_1,0,z,{\gamma }_1/\zeta ,\zeta T)dz,\hfill \end{array}$$

(15)

where ${\mathsf{\Phi }}^{\prime }(·)$ is the pdf of the standard normal distribution.

Finally, we substitute approximations (12) and (15) into Equation (11) to obtain the formula for the ruin probability in finite time t:

$$\begin{array}{ccc}\hfill & \psi (x_1,x_2,t)\approx \hfill & \hfill \\ \hfill & =1-{\int }_0^{\infty }\left\{\mathsf{\Phi }\left(\frac{z+{\gamma }_2(t-T)}{\sqrt{\zeta (t-T)}}\right)-exp\left\{-\frac{2{\gamma }_{2}z}{\zeta }\right\}\mathsf{\Phi }\left(\frac{-z+{\gamma }_2(t-T)}{\sqrt{{\gamma }_2(t-T)}}\right)\right\}·\hfill & \hfill \\ \hfill ·& \left\{\frac{1}{\sqrt{\zeta T}}{\mathsf{\Phi }}^{\prime }\left(\frac{z-(x_1+{\gamma }_{1}T)}{\sqrt{\zeta T}}\right)-g\left(x_1,0,z,{\gamma }_1/\zeta ,\zeta T\right)\right\}dz.\hfill \end{array}$$

(16)

In the infinite-horizon case, we substitute the classical approximation Formula (8) into Equation (11) and the first probability under the integral (or equivalently let $t\to \infty$ in Equation (16)) and obtain:

$$\begin{array}{ccc}\hfill & \psi (x_1,x_2)\approx \hfill & \\ \hfill & =1-{\int }_0^{\infty }\left(1-exp\left\{-\frac{2{\gamma }_{2}z}{\zeta }\right\}\right)\left\{\frac{1}{\sqrt{\zeta T}}{\mathsf{\Phi }}^{\prime }\left(\frac{z-(x_1+{\gamma }_{1}T)}{\sqrt{\zeta T}}\right)-g\left(x_1,0,z,{\gamma }_1/\zeta ,\zeta T\right)\right\}dz.\hfill \end{array}$$

(17)

□

Remark 1.

For the Poisson process with intensity $\lambda >0$, Theorem 1 holds with ${\gamma }_i={\theta }_i\lambda \mu$ ($i=1,2$) and $\zeta =\lambda ({\sigma }^2+{\mu }^2)$.

5. Results—Simulation Study

In this section, we present a simulation study that gives us insight into the quality of the diffusion approximation for the probability of the ruin in the insurer–reinsurer model that we obtained in previous parts of this work. The main idea is to compare the results of ruin probability approximations given by Theorem 1 with the results of the Monte Carlo simulations for both infinite and finite time horizons. Furthermore, to assess the accuracy of the ruin probability approximation introduced in this work, we compare it with the results of the De Vylder-type approximation proposed in a previous study Burnecki et al. (2019).

Empirical studies prove that insurance loss data are of different nature and can be described using light- or heavy-tailed distributions. Therefore, in our simulation study, we investigate the quality of the proposed approximation by considering several cases of both claim-distribution classes that were fitted to real-world data Burnecki et al. (2019); Burnecki and Teuerle (2011); Ma and Ma (2013). We distinguish between these two classes by their tail behaviour Burnecki and Teuerle (2011); Mikosch (2009). Namely, the distribution X is said to be light-tailed if there exist constants $a>0$ and $b>0$ such that ${\overline{F}}_X\left(x\right)=1-F_X\left(x\right)\le a{e}^{-bx}$, and is heavy-tailed if for all $a>0$ and $b>0$ such that ${\overline{F}}_X\left(x\right)>a{e}^{-bx}$, where $F_X\left(x\right)$ is a cumulative distribution function of X.

In the simulation study, we choose the following light-tailed distributions that are defined using their probability density functions:

• mixture of two exponential distributions:

  $$\begin{array}{c}\hfill f\left(x\right)=a_1{\beta }_1{e}^{-{\beta }_{1}x}+a_2{\beta }_2{e}^{-{\beta }_{2}x},x\ge 0,\end{array}$$

  (18)

  where ${\beta }_1=9.63$, ${\beta }_2=0.77$ and weights $a_1=0.25$ and $a_2=0.75$;
• gamma distribution:

  $$\begin{array}{c}\hfill f\left(x\right)=\frac{{\beta }^{\alpha }}{\mathsf{\Gamma }\left(\alpha \right)}{x}^{\alpha -1}{e}^{-\beta x},x\ge 0,\end{array}$$

  (19)

  where $\alpha =1$ and $\beta =1$;
• and Weibull distribution:

  $$\begin{array}{c}\hfill f\left(x\right)=c\tau x^{\tau -1}{e}^{-c{x}^{\tau }},x\ge 0,\end{array}$$

  (20)

  with $c=0.929$, $\tau =1.201$.

Among the heavy-tailed distributions, we propose the following three distributions:

• log-normal distribution:

  $$\begin{array}{c}\hfill f\left(x\right)=\frac{1}{\sqrt{2\pi }\sigma x}{e}^{-\frac{{(ln\left(x\right)-\mu )}^2}{2{\sigma }^2}},x\ge 0,\end{array}$$

  (21)

  where $\mu =-0.6$ and $\sigma =\sqrt{1.2}$;
• generalized Pareto distribution:

  $$\begin{array}{c}\hfill f\left(x\right)=\frac{1}{\sigma }{\left(1+k\frac{x-\theta }{\sigma }\right)}^{-1-\frac{1}{k}},x\ge 0,\end{array}$$

  (22)

  where $k=0.05$, $\sigma =0.42$ and $\theta =0.56$;
• and Burr type XII distribution:

  $$\begin{array}{c}\hfill f\left(x\right)=\frac{\frac{kc}{\alpha }{\left(\frac{x}{\alpha }\right)}^{c-1}}{{\left(1+{\left(\frac{x}{\alpha }\right)}^c\right)}^{k+1}},x\ge 0,\end{array}$$

  (23)

  with $\alpha =1.65$, $c=1.85$, $k=2.75$.

The parameters of the considered heavy-tailed distributions are chosen in such a way that the variance exists. Furthermore, similarly to Burnecki et al. (2019), we choose the parameters of all distributions to ensure that their means are approximately equal to 1.

In our numerical study, we compare the ruin probabilities in finite- and infinite-time horizons calculated by means of the diffusion approximations from Theorem 1 with the Monte Carlo (MC) results. It is worth clarifying that the MC results for the considered ruin probabilities in the infinite-time horizon are obtained for a long enough time horizon and are used as a reference ’exact’ value for the comparison.

Furthermore, we assume that the initial capital of the insurer $u_1$ takes the values: $u_1=0,1,\dots ,50$ and the capital of the reinsurer is equal to $u_2=u_1/3$. We also assume that the proportion parameter of the quota-share $\delta$ is equal $0.8$, while the safety loadings for the insurer and the reinsurer take the values ${\theta }_1=0.3$ and ${\theta }_2=0.03$, respectively. The accuracy of the diffusion approximations of the ruin probability ${\psi }_{·}(·)$ is measured by the relative error with respect to ${\psi }_{MC}(·)$ using the following formula $({\psi }_{MC}(·)-{\psi }_{·}(·))/{\psi }_{MC}(·)$.

In Figure 1, Figure 2, Figure 3 and Figure 4 we present the results for the simulated renewal process with inter-arrival times being log-normal with mean $1/\alpha =0.1$ and variance ${\sigma }_M=0.1$. We calculate the relative errors of the diffusion approximations for the renewal process and the diffusion approximations for the Poisson process with $\lambda =10$ (corresponding to parameters of the simulated renewal process) with respect to the MC results. In Figure 1 and Figure 2 the results for the finite-time horizon $t=100$ are depicted. We can see that the absolute relative errors for all claim amount distributions are on average not larger than $10\%$, although for larger capitals in the case of gamma, lognormal, generalized Pareto and Burr XII distributions they reach a level of $20\%$. For all cases, except for the Weibull distribution, the relative error gradually increases with the initial capital. It is worth noting that in the case of generalized Pareto distribution the relative error for $u_1\in \{40,\dots ,50\}$ increases when ruin probability values are close to 0. Figure 3 and Figure 4 depict the results for the infinite-time horizon, where for the MC simulations we take $t=1000$. In general, the obtained values of the ruin probabilities are slightly higher compared to the finite-time horizon case. However, the conclusions are analogous. It is worth mentioning that for both finite- and infinite-time horizons, the accuracy of both diffusion approximations is quite similar.

In our last analyses presented in Figure 5 and Figure 6 we simulate the Poisson process with $\lambda =10$) and compare the diffusion approximation for the Poisson process with the so-called De Vylder-type approximation Burnecki et al. (2019). We can observe clearly that in all cases, for the log-normal, mixture of exponential and gamma distributions, the diffusion approximation leads to smaller relative errors. In the case of the Weibull distribution, both approximations give comparable results. For larger initial capital $u_1\in \{20,\dots ,50\}$ and the generalized Pareto and Burr XII distributions, the De Vylder-type approximation outperforms the diffusion approximation.

To conclude, we found that the diffusion approximation works very well for both light- and heavy-tailed distributions that are widely used in non-life insurance mathematics. For many cases, the error was almost negligible from a practical point of view. This is somehow in contrast to performance analyses performed for the classical one-dimensional diffusion approximation which found the approximation relatively poor, see Burnecki and Weron (2005); Grandell (1991). The results were obtained for the insurer–reinsurer model which corresponds to the situation of an insurer and a reinsurer under the quota share contact in which the ceding insurer cedes an agreed-on percentage of every risk it insures that falls within a class of business subject to the treaty. A quota share treaty is utilized when an insurer wants to lower the financial risk to be able to underwrite more policies.

6. Conclusions

Ruin theory can be viewed as the theoretical foundation for modeling insolvency risk. In this paper, we addressed the issue of efficient approximation of the ruin probability within the insurer–reinsurer model for a general counting process, namely the renewal. The considered two-dimensional model classically describes the situation of insurance and reinsurance companies with a proportional quota share reinsurance contract or two branches of the same insurance company (Avram et al. 2008a, 2008b). Since only in the special case of the Poisson process and exponentially distributed claims, a simple ruin probability formula is known in infinite time Avram et al. (2008b); Burnecki et al. (2021), easy-to-use analytical approximations are needed.

In the literature, a De Vylder-type approximation for the insurer–reinsurer model was recently introduced for the infinite-time case Burnecki et al. (2019). The approximation requires only the first three finite moments of the claim amount distribution, and similarly to the above-mentioned results for the exponential claims, the counting process is assumed to be Poisson.

The idea of the diffusion approximation presented here is based on the weak convergence of stochastic processes, which enables the replacement of the initial risk process with an arithmetic Brownian motion. The counting process is more general than Poisson, namely, it is a renewal process. The only assumption for the claim amount distributions is that they should satisfy the central limit theorem, so their variance is finite. We applied this idea to the insurer–reinsurer model to derive simple ruin probability approximations for both finite and infinite time. The special case of the Poisson process was also presented.

Next, we checked the usefulness of the introduced approximations by studying several claim amount distributions (light- and heavy-tailed) and performing a Monte Carlo simulation study for finite and infinite time. We considered the Poisson process and renewal process with log-normally distributed waiting times as counting processes. To make the results comparable, the parameters of the claim amount distributions were chosen to make their means similar. Similarly, for the counting processes, the mean waiting times match. For the simulation study in infinite time, we took sufficiently large time horizons. In the infinite-time case, for the Poisson case, we also compared our results with the De Vylder-type approximation.

We found the introduced approximations to be very accurate for both light- and heavy-tailed distributions. The absolute relative error was quite small and gradually increased with initial capital. The quality of the Poisson process approximation in infinite time was usually better than the De Vylder-type especially in the light-tailed case. We also note that the diffusion approximation works for a more general class of distributions and for a more general counting process.