**1. Introduction**

The classical ruin problem formulated by Lundberg [1] and Cramér [2] was reconsidered by Taylor [3], who incorporated inflationary conditions into it. He determined that the inflation rate should be taken into consideration and increased. This result is perhaps not surprising since the presence of inflation means that if the effects of inflation are not offset by the premium adjustment, for example, it is the insurer who will suffer prejudice when claims are ultimately initiated. To generalize the problem of ruin, Gerber and Shiu [4] introduced the expected discounted penalty function, which allowed them to analyze the joint distribution of the time of ruin, the surplus immediately before ruin, and the deficit at ruin through the Laplace transform of the defective renewal equation. Cai [5] built upon Gerber and Shiu [4]'s work by including a constant rate of interest and allowing the initial surplus to be negative. Based on these assumptions, he studied absolute ruin using the defective renewal equation of the expected discounted penalty function. Originally studied in the insurance industry, the ruin problem and discounted penalty functions have been found to be applicable to the financial sector as well. We refer readers to Adekambi and Essiomle [6] and Gerber et al. [7] and the references therein for further information on the use of ruin theory in the banking sector.

Sundt and Teugels [8] extended the classical compound Poisson model by assuming that both the premium and the initial reserve are invested in a risk-free asset with a constant interest rate. They derived a differential equation of the ruin probability using the Poisson risk process. In addition, they discussed two special cases where there is no initial reserve and assumed that the claim amounts are exponentially distributed. In the same vein, Kalashnikov and Konstantinides [9] extended Taylor [3]'s work, under a subexponential claims distribution assumption. Similarly, Kalashnikov and Konstantinides [9] and Yuen [10] studied the Gerber–Shiu discounted penalty function when they extended Taylor [3]'s work by incorporating a force of interest and a constant barrier into the classic compound risk model.

Li and Lu [11], in turn, extended Yuen [10]'s work, when they considered a risk process, which incorporates credit and debt interest, and studied the generalized Gerber–Shiu

**Citation:** Essiomle, K.; Adekambi, F. A Note on Gerber–Shiu Function with Delayed Claim Reporting under Constant Force of Interest. *Math. Comput. Appl.* **2022**, *27*, 51. https:// doi.org/10.3390/mca27030051

Received: 7 May 2022 Accepted: 13 June 2022 Published: 20 June 2022

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2022 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

233

penalty function. They derived the exact expressions for the penalty functions when credit and debt interest rates are equal and when claim sizes are exponentially distributed.

Several studies have shown that the assumption of independence between claim sizes and inter-claim times is unrealistic and is created purely for mathematical tractability. Scholars assumed that the risk would have a certain dependence structure to make the model more realistic. Furthermore, Adekambi and Essiomle [12] and Essiomle and Adekambi [13] examined the asymptotic tail probability and the ruin probability of the discounted claim. In their study, they developed a closed analytical form of the ruin probability that assumes a specific distribution of claims.

Cheung et al. [14] examined the structure of various Gerber–Shiu functions using Sparre Andersen's model and allowed a dependence structure between claim sizes and inter-claim times. Cheung et al. [14] analyzed Lundberg's fundamental equation and the generalized adjustment coefficient using the defective renewal equation. Based on Cheung et al. [14]'s research, Willmot and Woo [15] assumed a more general type of dependency; they analyzed the expected discounted penalty function, the surplus prior to ruin, the deficit at ruin, and the surplus after the second claim prior to the moment of ruin.

Furthermore, Gerber argues that when interest yield is present, the surplus can fall below zero and then rise above zero. As a result, the ruin may not occur. As a result of this scenario, Gerber defined what he termed the absolute ruin time, which states that a ruin occurs if the surplus falls below a certain threshold, say <sup>−</sup>*<sup>c</sup> <sup>δ</sup>* , where *c* represents the constant premium rate and *δ* > 0 represents the force of interest. In alignment with this idea, Konstantinides et al. [16] assumed a constant force of interest and studied the probabilities of finite and infinite absolute ruin time. Their study derived explicit asymptotic expressions for the finite and infinite absolute ruin probabilities for the compound Poisson model.

In most of the actuarial literature, it is assumed that claims are not delayed, which means they are reported or settled as soon as they occur. As a result, the model is unrealistic since the insured reported the claim after some days of the event. Furthermore, due to some legislative constraints (for example, claim investigation by experts), the insurers may take some months to settle claims. To address these observations, Cai [17] examined the effects of of the timing of claim payments and interest rate on ruin probability. In his study, the rate of interest was assumed to be dependent on the auto-regressive structure, as well as on the Lundberg inequality. The same perspective was adopted by Zou and Jie [18], who assumed that the claim occurrence may be delayed, and so, the total loss within the time horizon can be divided into two sub-claims: the main claim and the by-claim. In their analyze, they used the mean of Lagrange interpolation and the defective renewal equation of the discounted penalty function. Through a compound geometric distribution, they were able to establish an exact representation of the solution.

This paper incorporates the reporting delay time (when the claim occurs, it takes a while for the insured to declare that claim or for the insurance company to report that claim). Put differently, we assume that the exact amount of claim can only be determined after a reporting delay. We extend Li and Lu [11]'s work by assuming that the loss sizes are only known after a random delay time, during which the insurer continues to receive the premium and investment return. With these assumptions, we establish an integrodifferential equation of the discounted penalty function and find that the solution to this equation can be represented by an infinite series. We derive an explicit solution for the exponential claim distribution.

To our best knowledge, there is no such study in the literature. Thus, our result contributes to the body of knowledge on the effects of delayed reporting times. The rest of the paper is structured as follows: Section 2 provides an overview of the model and a brief explanation of its assumptions. In Section 3, we present the main results of our work and discuss some special cases. The following section provides some numerical illustrations. We conclude the paper in Section 5 with suggested ways in which the work can be extended.

#### **2. Model Setting and Rationality**

Suppose that the premiums are continuously received and invested in the market with a constant force of interest (*δ* > 0) and that **U***δ*(*t*) is the surplus at time *t*. The initial reserve is also subject to investment and generates a constant force of interest. Let us further assume that the *n*-th claim occurs at *Tn* = ∑*<sup>n</sup> <sup>k</sup>*=<sup>0</sup> *τ<sup>k</sup>* and it takes *Vn* time for the exact claim size or for the claim to be reported. Then, the surplus process (**U***δ*(*t*), *t* ≥ 0) is as follows:

$$\mathbf{U}\_{\delta}(t) = \mu \mathbf{e}^{\delta t} + c \int\_{0}^{t} \mathbf{e}^{\delta(t-y)} \mathbf{d}y - \int\_{0}^{t} \mathbf{e}^{\delta(t-y)} \mathbf{d}S\_{y}, \quad \text{t} \ge \mathbf{0}, \tag{1}$$

where *St* = **N**(**t**) ∑ *k*=1 *Xk*{*Tk*+*Vk* <sup>≤</sup> *<sup>t</sup>*} and:


Making things simple, we assume that *Tk* and *Vk* are independent, an assumption that can be justified by the fact that the delay time may depend mostly on the claim size. Thus, the joint distribution of the claim size, the claim arrival time, and the delay time are given by

$$f\_{\mathbf{X}\_k T\_k, V\_k}(\mathbf{x}, t, z) \;=\; f\_{T\_k}(t) f\_{V\_k}(z) f\_{\mathbf{X}\_k \mid T\_k, V\_k}(\mathbf{x}) .$$

Under some circumstances, the claim reporting and settlement process may take more time due to bad economics in the company (insufficient liquidity during wartime or natural disasters). In addition, the delay time is more relevant in a health insurance contract because one has to determine if the illness began within the period that the insurance covers before payments can be made. Additionally, during that holding period, the company will continue to receive premiums from policyholders, along with a compounded yield on the initial reserve. By delaying the claim payment time, companies can partially resolve the liquidity problem. Thus, they will minimize the absolute ruin probability, and therefore reduce the bankruptcy rate.

Gerber and Shiu [4] and Cai [5] argue that when the surplus drops below <sup>−</sup>*<sup>c</sup> <sup>δ</sup>* , the surplus cannot be positive due to the debts of the insurer, because the present value at that time is less than or equal to *c*, which is the present value of all premium income available after that time, and the delay in reporting and payment. A good investment strategy may be required in some types of contracts; for example, when the insurer covers only low risks and invests heavily in high asset returns (we refer to this scenario as a <sup>−</sup>*<sup>c</sup> <sup>δ</sup>* small enough risk pool), the surplus may be positive as the company will continue to receive premiums and returns at the time of the ruin. In light of this, the time to absolute ruin defined by Gerber and Shiu [4] and Cai [5] is considered to be an absolute (relative) time to ruin. Throughout this paper, we assume that *<sup>c</sup> <sup>δ</sup>* is small enough. That is, the insurer covers only low risk and invests in high return assets (*δ* >> *c*).

#### **3. Main Results**

We present our main results in this section. To start with, we use the Gerber–Shiu function and the formula for the absolute (relative) time to ruin. The absolute ruin probability *T<sup>δ</sup>* is calculated as follows:

$$\mathbf{T}\_{\delta} = \begin{cases} \inf \{ \mathbf{t} \, | \, \mathbf{U}\_{\delta}(t) < -\frac{c}{\delta} \} \\ \text{as} \quad \text{if } \mathbf{U}\_{\delta}(t) \ge -\frac{\delta}{\delta} \end{cases} \tag{2}$$

Here, *<sup>c</sup>*<sup>0</sup> <sup>=</sup> <sup>−</sup> *<sup>c</sup> <sup>δ</sup>* and is defined through the Gerber–Shiu function by

$$\Phi\_{\mathbf{a},\delta}(\mathbf{u},\mathbf{c}\_0) = \mathbb{E}\left[e^{-\mathbf{a}\mathbf{T}\_{\delta}}\mathbf{W}\left(\mathbf{U}\_{\delta}(\mathbf{T}\_{\delta}^-),\mathbf{c}\_0 - \mathbf{U}\_{\delta}(\mathbf{T}\_{\delta})\right)\mathbf{1}\_{\{\mathbf{T}\_{\delta} < \infty\}}\right],\tag{3}$$

*α* is interpreted as the discounted rate, and the function W() represents the penalty function and is assumed to be bounded.

#### *3.1. Integro-Differential Equation of the Gerber–Shiu Function*

We obtain the following proposition.

**Proposition 1.** *Assume from the financial surplus given by Equation* (1) *that the delay time is exponentially distributed with parameter β. Then, the Gerber–Shiu function satisfies the integrodifferential equation below:*

$$\begin{split} \left(\delta u + c\right)^{2} \Phi\_{a,\delta}^{\prime\prime}(\mu, c\_{0}) &= \eta\_{0} (\delta u + c) \Phi\_{a,\delta}^{\prime}(\mu, c\_{0}) - \eta\_{1} \Phi\_{a,\delta}(\mu, c\_{0}) \\ &+ \lambda \beta \int\_{0}^{u - \varepsilon\_{0}} \Phi\_{a,\delta}(u - x, c\_{0}) dF\_{\mathbb{X}[\tau, V]}(x) + \lambda \beta \omega(u), \end{split} \tag{4}$$

*where <sup>η</sup>*<sup>0</sup> <sup>=</sup> <sup>2</sup>*<sup>α</sup>* <sup>+</sup> *<sup>β</sup>* <sup>+</sup> *<sup>λ</sup>* <sup>−</sup> *<sup>δ</sup>, <sup>η</sup>*<sup>1</sup> <sup>=</sup> (*<sup>λ</sup>* <sup>+</sup> *<sup>α</sup>*)(*<sup>α</sup>* <sup>+</sup> *<sup>β</sup>*)*, <sup>c</sup>*<sup>0</sup> <sup>=</sup> <sup>−</sup>*<sup>c</sup> <sup>δ</sup> , λ is the claim number process rate in a Poisson process, and*

$$
\omega(\mathfrak{u}) = \int\_{\mathfrak{u} - \mathfrak{c}\_0}^{\infty} \mathfrak{w}(\mathfrak{u}, \mathfrak{c}\_0 + \mathfrak{x} - \mathfrak{u}) \mathrm{d}F\_{\mathfrak{X}|\mathfrak{r}, V}(\mathfrak{x}) .
$$

**Proof.** Let *<sup>s</sup>*(*t*, *<sup>v</sup>*) = *ueδ*(*t*+*v*) <sup>+</sup> *<sup>c</sup> <sup>e</sup>δ*(*t*+*v*)−<sup>1</sup> *<sup>δ</sup>* and *<sup>c</sup>*<sup>0</sup> <sup>=</sup> <sup>−</sup>*<sup>c</sup> δ* .

By conditioning on the first claim occurrence time, the corresponding claim size, and the delay time, we have the equation:

<sup>Φ</sup>*δ*,*α*(*u*, *<sup>c</sup>*0) = <sup>∞</sup> 0 ∞ 0 *s*(*t*,*v*)−*c*<sup>0</sup> 0 *e* −*α*(*t*+*v*) <sup>Φ</sup>*δ*,*α*(*s*(*t*, *<sup>v</sup>*) <sup>−</sup> *<sup>x</sup>*, *<sup>c</sup>*0)d*FX*|*τ*,*V*(*x*)d*Fτ*(*t*)d*FV*(*v*) + ∞ 0 ∞ 0 ∞ *s*(*t*,*v*)−*c*<sup>0</sup> *e* −*α*(*t*+*v*) <sup>W</sup>(*s*(*t*, *<sup>v</sup>*), *<sup>c</sup>*<sup>0</sup> <sup>+</sup> *<sup>x</sup>* <sup>−</sup> *<sup>s</sup>*(*t*, *<sup>v</sup>*))d*FX*|*τ*,*V*(*x*)d*Fτ*(*t*)d*FV*(*v*) = ∞ 0 ∞ 0 *s*(*t*,*v*)−*c*<sup>0</sup> 0 *e* <sup>−</sup>*α*(*t*+*v*)*βe* −*βt λe* −*λt* <sup>Φ</sup>*δ*,*α*(*s*(*t*, *<sup>v</sup>*) <sup>−</sup> *<sup>x</sup>*, *<sup>c</sup>*0)d*FX*|*τ*,*V*(*x*)d*t*d*<sup>v</sup>* + ∞ 0 ∞ 0 ∞ *s*(*t*,*v*)−*c*<sup>0</sup> *e* <sup>−</sup>*α*(*t*+*v*)*βe* −*βt λe* −*λt* <sup>W</sup>(*s*(*t*, *<sup>v</sup>*), *<sup>c</sup>*<sup>0</sup> <sup>+</sup> *<sup>x</sup>* <sup>−</sup> *<sup>s</sup>*(*t*, *<sup>v</sup>*))d*FX*|*τ*,*V*(*x*)d*t*d*v*.

Let *z* = *s*(*t*, *v*). Then,

$$\begin{array}{rcl} \mathfrak{c}^{t} &=& \left(\frac{\delta z + c}{\delta u + c}\right)^{1/\delta} \mathfrak{e}^{-v}; \text{and} \\ \mathbf{d}t &=& \frac{1}{\delta z + c} \mathbf{d}z. \end{array}$$

Plugging these relations in the above equation yields

<sup>Φ</sup>*δ*,*α*(*u*, *<sup>c</sup>*0) = <sup>∞</sup> 0 ∞ *s*(0,*v*) *z*−*c*<sup>0</sup> 0 *λβ <sup>δ</sup><sup>z</sup>* <sup>+</sup> *<sup>c</sup> δu* + *c* −(*α*+*λ*)/*<sup>δ</sup> e* <sup>−</sup>(*β*−*λ*)*<sup>v</sup>* <sup>1</sup> *δz* + *c* Φ*δ*,*α*(*z* − *x*, *c*0) <sup>×</sup> <sup>d</sup>*FX*|*τ*,*V*(*x*)d*z*d*<sup>v</sup>* + ∞ 0 ∞ *s*(0,*v*) ∞ *z*−*c*<sup>0</sup> *λβ <sup>δ</sup><sup>z</sup>* <sup>+</sup> *<sup>c</sup> δu* + *c* −(*α*+*λ*)/*<sup>δ</sup> e* <sup>−</sup>(*β*−*λ*) <sup>1</sup> *δz* + *c* <sup>W</sup>(*z*, *c*<sup>0</sup> + *x* − *z*) <sup>×</sup> <sup>d</sup>*FX*|*τ*,*V*(*x*)d*z*d*v*, = *λβ*(*δu* + *c*) (*α*+*λ*)/*δ* ∞ 0 ∞ *s*(0,*v*) (*δz* + *c*) <sup>−</sup> *<sup>α</sup>* <sup>+</sup> *<sup>λ</sup> <sup>δ</sup>* <sup>−</sup><sup>1</sup> *e* <sup>−</sup>(*β*−*λ*)*<sup>v</sup>* (5) × *z*−*c*<sup>0</sup> 0 <sup>Φ</sup>*δ*,*α*(*<sup>z</sup>* <sup>−</sup> *<sup>x</sup>*, *<sup>c</sup>*0)d*FX*|*τ*,*V*(*x*) + <sup>∞</sup> *z*−*c*<sup>0</sup> <sup>W</sup>(*z*, *<sup>c</sup>*<sup>0</sup> <sup>+</sup> *<sup>x</sup>* <sup>−</sup> *<sup>z</sup>*)d*FX*|*τ*,*V*(*x*) d*z*d*v*.

Now, we define *ω*(*z*) := ∞ *z*−*c*<sup>0</sup> <sup>W</sup>(*z*, *<sup>c</sup>*<sup>0</sup> + *<sup>x</sup>* − *<sup>z</sup>*)d*FX*|*τ*,*V*(*x*) and *<sup>y</sup>* = *<sup>s</sup>*(0, *<sup>v</sup>*) so that

$$\mathfrak{e}^{\upsilon} \quad = \left(\frac{\delta y + c}{\delta u + c}\right)^{1/\delta} \implies \mathrm{d}\upsilon = \frac{1}{\delta y + c} \mathrm{d}y.$$

Then, inverting the above relations into Equation (5) gives

$$\begin{split} \Phi\_{\delta,\mathbb{A}}(\mu,c\_{0}) &= \ &\lambda\beta(\delta u+c)^{(a+\lambda)/\delta} \int\_{\mu}^{\infty} \int\_{y}^{\infty} (\delta z+c)^{-\frac{\alpha+\lambda}{\delta}-1} \frac{1}{\delta y+c} \left(\frac{\delta y+c}{\delta u+c}\right)^{-(\beta-\lambda)/\delta} \\ &\times \ &\left[\int\_{0}^{z-c\_{0}} \Phi\_{\delta,\mathbb{A}}(z-x,c\_{0}) \mathrm{d}F\_{\mathbb{X}[\tau,V}(x)+\omega(z)\right] \mathrm{d}zdy, \\ &= &\lambda\beta(\delta u+c)^{(a+\beta)/\delta} \int\_{\mu}^{\infty} \int\_{y}^{\infty} (\delta z+c)^{-\frac{\alpha+\lambda}{\delta}-1} (\delta y+c)^{-\frac{\beta-\lambda}{\delta}-1} \\ &\times \ &\left[\int\_{0}^{z-c\_{0}} \Phi\_{\delta,\mathbb{A}}(z-x,c\_{0}) \mathrm{d}F\_{\mathbb{X}[\tau,V}(x)+\omega(z)\right] \mathrm{d}zdy, \\ &= &\lambda\beta(\delta u+c)^{(a+\beta)/\delta} \int\_{\mu}^{\infty} (\delta y+c)^{-\frac{\beta-\lambda}{\delta}-1} \int\_{y}^{\infty} (\delta z+c)^{-\frac{\alpha+\lambda}{\delta}-1} \\ &\times \ &\left[\int\_{0}^{z-c\_{0}} \Phi\_{\delta,\mathbb{A}}(z-x,c\_{0}) \mathrm{d}F\_{\mathbb{X}[\tau,V}(x)+\omega(z)\right] \mathrm{d}zdy. \tag{6} \end{split}$$

Taking the first derivative of Equation (6):

$$
\begin{split}
\Phi\_{\delta,\mu}'(\boldsymbol{u},\boldsymbol{c}\_{0}) &= \quad \frac{\boldsymbol{\alpha}+\beta}{\delta\boldsymbol{u}+\boldsymbol{c}}\Phi\_{\delta,\mu}(\boldsymbol{u},\boldsymbol{c}\_{0}) - \lambda\beta \frac{(\delta\boldsymbol{u}+\boldsymbol{c})^{(\boldsymbol{\alpha}+\boldsymbol{\lambda})/\delta}}{\delta\boldsymbol{u}+\boldsymbol{c}} \int\_{\boldsymbol{u}}^{\infty} (\delta\boldsymbol{z}+\boldsymbol{c})^{-\frac{\boldsymbol{\alpha}+\boldsymbol{\lambda}}{\delta}-1} \\ & \quad \times \left[\int\_{0}^{\boldsymbol{z}-\boldsymbol{c}\_{0}} \Phi\_{\delta,\mu}(\boldsymbol{z}-\boldsymbol{x},\boldsymbol{c}\_{0}) \mathrm{d}\mathcal{F}\_{\boldsymbol{X}|\boldsymbol{r},\boldsymbol{V}}(\boldsymbol{x}) + \omega(\boldsymbol{z}) \right] \mathrm{d}\boldsymbol{z}, \\\ (\delta\boldsymbol{u}+\boldsymbol{c})\Phi\_{\delta,\mu}'(\boldsymbol{u},\boldsymbol{c}\_{0}) &= \quad (\boldsymbol{a}+\beta)\Phi\_{\delta,\mu}(\boldsymbol{u},\boldsymbol{c}\_{0}) - \lambda\beta(\delta\boldsymbol{u}+\boldsymbol{c})^{(\boldsymbol{a}+\boldsymbol{\lambda})/\delta} \int\_{\boldsymbol{u}}^{\infty} (\delta\boldsymbol{z}+\boldsymbol{c})^{-\frac{\boldsymbol{a}+\boldsymbol{\lambda}}{\delta}-1} \\ & \quad \times \left[\int\_{0}^{\boldsymbol{z}-\boldsymbol{c}\_{0}} \Phi\_{\delta,\mu}(\boldsymbol{z}-\boldsymbol{x},\boldsymbol{c}\_{0}) \mathrm{d}\mathcal{F}\_{\boldsymbol{X}|\boldsymbol{r},\boldsymbol{V}}(\boldsymbol{x}) + \omega(\boldsymbol{z}) \right] \mathrm{d}\boldsymbol{z}. \end{split}
\tag{7}$$

Upon rearranging Equation (7), we obtain:

$$\begin{split} & -\lambda \beta (\delta u + c)^{(\alpha + \lambda)/\delta} \int\_{\mu}^{\infty} (\delta z + c)^{-\frac{a + \lambda}{\delta} - 1} \left[ \int\_{0}^{z - \varepsilon\_{0}} \Phi\_{\delta, \mu}(z - x, c\_{0}) \mathrm{d}F\_{\mathbf{X}|\mathbf{r}, \mathcal{V}}(x) + \omega(z) \right] d\mathbf{z} \\ & = (\delta u + c) \Phi\_{\delta, \mu}'(u, c\_{0}) - (a + \beta) \Phi\_{\delta, \mu}(u, c\_{0}). \end{split} \tag{8}$$

By taking the derivative of (7), we have

(*δu* + *c*)Φ *<sup>δ</sup>*,*α*(*u*, *<sup>c</sup>*0) + *<sup>δ</sup>*Φ *<sup>δ</sup>*,*α*(*u*, *<sup>c</sup>*0) = (*<sup>α</sup>* + *<sup>β</sup>*)Φ *<sup>δ</sup>*,*α*(*u*, *<sup>c</sup>*0) <sup>−</sup> (*<sup>α</sup>* <sup>+</sup> *<sup>λ</sup>*) *λβ δu* + *c* (*δu* + *c*) (*α*+*λ*)/*δ* × ∞ *u* (*δz* + *c*) <sup>−</sup> *<sup>α</sup>*+*<sup>λ</sup> <sup>δ</sup>* −1 *z*−*c*<sup>0</sup> 0 <sup>Φ</sup>*δ*,*α*(*<sup>z</sup>* <sup>−</sup> *<sup>x</sup>*, *<sup>c</sup>*0)d*FX*|*τ*,*V*(*x*) + *<sup>ω</sup>*(*z*) d*z* <sup>+</sup> *λβ δu* + *c u*−*c*<sup>0</sup> 0 <sup>Φ</sup>*δ*,*α*(*<sup>u</sup>* <sup>−</sup> *<sup>x</sup>*, *<sup>c</sup>*0)d*FX*|*τ*,*V*(*x*) + *<sup>ω</sup>*(*u*) , (*δu* + *c*) 2 Φ *<sup>δ</sup>*,*α*(*u*, *<sup>c</sup>*0) = (*<sup>α</sup>* <sup>+</sup> *<sup>β</sup>* <sup>−</sup> *<sup>δ</sup>*)(*δ<sup>u</sup>* <sup>+</sup> *<sup>c</sup>*)Φ *<sup>δ</sup>*,*α*(*u*, *c*0) − (*α* + *λ*)*λβ*(*δu* + *c*) (*α*+*λ*)/*δ* × ∞ *u* (*δz* + *c*) <sup>−</sup> *<sup>α</sup>*+*<sup>λ</sup> <sup>δ</sup>* −1 *z*−*c*<sup>0</sup> 0 <sup>Φ</sup>*δ*,*α*(*<sup>z</sup>* <sup>−</sup> *<sup>x</sup>*, *<sup>c</sup>*0)d*FX*|*τ*,*V*(*x*) + *<sup>ω</sup>*(*z*) d*z* <sup>+</sup> *λβ <sup>u</sup>*−*c*<sup>0</sup> 0 <sup>Φ</sup>*δ*,*α*(*<sup>u</sup>* <sup>−</sup> *<sup>x</sup>*, *<sup>c</sup>*0)d*FX*|*τ*,*V*(*x*) + *<sup>ω</sup>*(*u*) . (9)

Finally, inserting Equation (8) into (9), the equation:

$$\begin{aligned} \left(\delta\boldsymbol{u} + \boldsymbol{c}\right)^{2} \boldsymbol{\Phi}\_{\delta,\mathfrak{a}}^{\boldsymbol{\mu}}(\boldsymbol{u},\boldsymbol{c}\_{0}) &= \left(2\boldsymbol{a} + \boldsymbol{\lambda} + \boldsymbol{\beta} - \boldsymbol{\delta}\right) (\delta\boldsymbol{u} + \boldsymbol{c}) \boldsymbol{\Phi}\_{\delta,\mathfrak{a}}^{\boldsymbol{\prime}}(\boldsymbol{u},\boldsymbol{c}\_{0}) - (\boldsymbol{\lambda} + \boldsymbol{a})(\boldsymbol{u} + \boldsymbol{\beta}) \boldsymbol{\Phi}\_{\delta,\mathfrak{a}}(\boldsymbol{u},\boldsymbol{c}\_{0}) \\ &+ \boldsymbol{\lambda}\boldsymbol{\beta} \left[\int\_{0}^{\boldsymbol{u} - \boldsymbol{c}\_{0}} \boldsymbol{\Phi}\_{\delta,\mathfrak{a}}(\boldsymbol{u} - \boldsymbol{x},\boldsymbol{c}\_{0}) \mathrm{d}\boldsymbol{F}\_{\boldsymbol{X}|\boldsymbol{\tau},\boldsymbol{V}}(\boldsymbol{x}) + \boldsymbol{\omega}(\boldsymbol{u})\right]\_{\boldsymbol{\prime}} \end{aligned}$$

which proves the statement by setting *η*<sup>0</sup> = (2*α* + *λ* + *β* − *δ*) and *η*<sup>1</sup> = (*λ* + *α*)(*α* + *β*).

**Remark 1.** *In deriving their results, Cai [5] and Li and Lu [11] distinguished between* Φ<sup>+</sup> *<sup>α</sup>*,*δ*(*u*) *for the return interest rate when the initial surplus exceeds zero and* Φ− *<sup>α</sup>*,*δ*(*u*) *for debt interest when the initial surplus drops below zero. However, as pointed out by Cai [5]'s work, when the initial surplus approaches zero, then both quantities are equalized. In this paper, we assume that <sup>c</sup> <sup>δ</sup> approaches* 0*, and thanks to the delayed reporting time, it is reasonable to assume that in the interval* [−*<sup>c</sup> <sup>δ</sup>* , 0]*,* Φ− *<sup>α</sup>*,*δ*(*u*) <sup>≈</sup> <sup>Φ</sup><sup>+</sup> *<sup>α</sup>*,*δ*(*u*)*.*

Equation (4) is difficult to solve. In the following, we use the Laplace transform to derive a more practical equation for the Gerber–Shiu function.

Let us define as in Li and Lu [11] the following transformation:

$$\chi\_{\delta,\mathfrak{a}}(\boldsymbol{u},\boldsymbol{c}\_{0}) = \begin{cases} \frac{\Phi\_{\delta,\mathfrak{a}}(\boldsymbol{c}\_{0},\boldsymbol{c}\_{0}) - \Phi\_{\delta,\mathfrak{a}}(\boldsymbol{u},\boldsymbol{c}\_{0})}{\Phi\_{\delta,\mathfrak{a}}(\boldsymbol{c}\_{0},\boldsymbol{c}\_{0})}, & \text{if } \boldsymbol{u} \ge \boldsymbol{c}\_{0} \\ 0 & \text{if } \quad \boldsymbol{u} < \boldsymbol{c}\_{0}. \end{cases} \tag{10}$$

Clearly, *Yδ*,*α*(*c*0, *c*0) = 0. If in addition, we assume that the penalty function <sup>W</sup>(*x*, *y*) is bounded, then lim *<sup>u</sup>*−→<sup>∞</sup> <sup>Φ</sup>*δ*,*α*(*u*, *<sup>c</sup>*0) = 0 (see Li and Lu [11]), which implies that lim *<sup>u</sup>*−→<sup>∞</sup> *<sup>Y</sup>δ*,*α*(*u*, *<sup>c</sup>*0) = 1.

**Theorem 1.** *In the model given by Equation* (1)*, the Laplace transform of the discounted penalty function satisfies the following non-homogeneous differential equation.*

$$
\tilde{y}\_{\delta,\mathfrak{a}}^{''}(\mathfrak{s}) + f\_1(\mathfrak{s})\tilde{y}\_{\delta,\mathfrak{a}}^{'}(\mathfrak{s}) + f\_0(\mathfrak{s})\tilde{y}\_{\delta,\mathfrak{a}}(\mathfrak{s}) = \mathfrak{g}(\mathfrak{s}),
\tag{11}
$$

*where Yδ*,*α*(*u*, *c*0) *is given by Equation* (10)*,*

$$\begin{split} f\_{1}(s) &= \quad \frac{\eta\_{2} + 2}{\delta} \times \frac{1}{s} + 2c\_{0}; \quad \eta\_{2} = \eta\_{0} + 2\delta; \quad \eta\_{3} = \eta\_{1} + \delta\eta\_{2}, \\\ f\_{0}(s) &= \quad \frac{c^{2}}{\delta^{2}} + \frac{1}{s^{2}} \left(\frac{\eta\_{3}}{\delta^{2}} - \frac{\lambda\beta}{\delta^{2}}\tilde{f}(s)\right) - \frac{1}{s} \left(\frac{2c}{\delta} + \frac{\eta\_{2}c}{\delta^{2}}\right), \\\ g(s) &= \quad -\frac{1}{\left(\delta s\right)^{2}} \left[\frac{c^{-\infty\_{0}}}{s} \left(\lambda\beta\tilde{f}(s) - \eta\_{1}\right) + \frac{\lambda\beta}{\Phi\_{\delta,\mu}(c\_{0})}\tilde{w}(s)\right]. \end{split}$$

*y*˜*δ*,*α*(*s*), ˜ *<sup>f</sup>*(*s*) *and <sup>w</sup>*˜(*s*) *are respectively the Laplace transforms of <sup>Y</sup>δ*,*α*(*u*), d*FX*|*τ*,*V*(*x*)*, and <sup>ω</sup>*(*u*)*.*

**Proof.** Hereafter, let d*F*(*x*) := d*Fτ*,*V*(*x*).

Integrating Equation (5) on both sides from *c*<sup>0</sup> to *u* yields

$$\begin{split} \int\_{\boldsymbol{c}0}^{u} (\delta t + \boldsymbol{c})^2 \boldsymbol{\Phi}\_{\delta,\mathfrak{a}}^{\boldsymbol{\eta}}(t, \boldsymbol{c}\_0) \mathrm{d}t &= \quad \eta\_0 \int\_{\boldsymbol{c}0}^{u} (\delta t + \boldsymbol{c}) \boldsymbol{\Phi}\_{\delta,\mathfrak{a}}^{\boldsymbol{\prime}}(t, \boldsymbol{c}\_0) \mathrm{d}t - \eta\_1 \int\_{\boldsymbol{c}0}^{u} \boldsymbol{\Phi}\_{\delta,\mathfrak{a}}(t, \boldsymbol{c}\_0) \mathrm{d}t \\ &+ \quad \lambda \boldsymbol{\beta} \int\_{\boldsymbol{c}0}^{u} \int\_{0}^{t-\varepsilon\_0} \boldsymbol{\Phi}\_{\delta,\mathfrak{a}}(t-\mathbf{x}, \boldsymbol{c}\_0) \mathrm{d}F(\mathbf{x}) \mathrm{d}t + \lambda \boldsymbol{\beta} \int\_{\boldsymbol{c}0}^{u} \boldsymbol{\omega}(t) \mathrm{d}t. \end{split}$$

That is,

$$\begin{split} \left(\delta\boldsymbol{u} + \boldsymbol{c}\right)^{2} \boldsymbol{\Phi}\_{\delta,\boldsymbol{a}}^{\prime}(\boldsymbol{u},\boldsymbol{c}\_{0}) &= \quad \eta\_{2} \int\_{\boldsymbol{c}\_{0}}^{\boldsymbol{u}} (\delta\boldsymbol{t} + \boldsymbol{c}) \boldsymbol{\Phi}\_{\delta,\boldsymbol{a}}^{\prime}(\boldsymbol{t},\boldsymbol{c}\_{0}) \mathrm{d}t - \eta\_{1} \int\_{\boldsymbol{c}\_{0}}^{\boldsymbol{u}} \boldsymbol{\Phi}\_{\delta,\boldsymbol{a}}(\boldsymbol{t},\boldsymbol{c}\_{0}) \mathrm{d}t \\ &+ \quad \lambda\boldsymbol{\beta} \int\_{0}^{\boldsymbol{u} - \boldsymbol{c}\_{0}} \mathrm{d}F(\boldsymbol{x}) \left(\int\_{\boldsymbol{c}\_{0}}^{\boldsymbol{u} - \boldsymbol{x}} \boldsymbol{\Phi}\_{\delta,\boldsymbol{a}}(\boldsymbol{y}) \mathrm{d}\boldsymbol{y}\right) + \lambda\boldsymbol{\beta} \int\_{\boldsymbol{c}\_{0}}^{\boldsymbol{u}} \boldsymbol{\omega}(\boldsymbol{t}) \mathrm{d}t. \end{split}$$

so that,

$$\begin{split} \left(\delta u + c\right)^{2} \Phi\_{\delta,\mathfrak{a}}^{'} (u, c\_{0}) &= \quad \eta\_{2} (\delta u + c) \Phi\_{\delta,\mathfrak{a}} (u, c\_{0}) - \eta\_{3} \int\_{c\_{0}}^{u} \Phi\_{\delta,\mathfrak{a}} (t, c\_{0}) \mathrm{d}t \\ &+ \quad \lambda \pounds \int\_{0}^{u - \varepsilon\_{0}} \mathrm{d}F(\mathbf{x}) \left( \int\_{c\_{0}}^{u - x} \Phi\_{\delta,\mathfrak{a}} (y, c\_{0}) \mathrm{d}y \right) + \lambda \pounds \int\_{c\_{0}}^{u} \omega(t) \mathrm{d}t. \end{split} \tag{12}$$

Plugging Equation (10) into (12) gives

$$\begin{split} -(\delta u + c)^{2} \Phi\_{\delta,\mathfrak{a}}(c\_{0},c\_{0}) Y\_{\delta,\mathfrak{a}}^{'}(u,c\_{0}) &= \eta\_{2}(\delta u + c) [\Phi\_{\delta,\mathfrak{a}}(c\_{0},c\_{0}) - \Phi\_{\delta,\mathfrak{a}}(c\_{0},c\_{0}) Y\_{\delta,\mathfrak{a}}(u,c\_{0})] \\ &- \eta\_{3} \int\_{c\_{0}}^{u} [\Phi\_{\delta,\mathfrak{a}}(c\_{0},c\_{0}) - \Phi\_{\delta,\mathfrak{a}}(c\_{0},c\_{0}) Y\_{\delta,\mathfrak{a}}(t,c\_{0})] dt \\ &+ \ \lambda \beta \int\_{c\_{0}}^{u} \omega(t) \mathrm{d}t + \lambda \beta \int\_{0}^{u-c\_{0}} \mathrm{d}F(\mathbf{x}) \\ &\times \left( \int\_{c\_{0}}^{u-x} [\Phi\_{\delta,\mathfrak{a}}(c\_{0},c\_{0}) - \Phi\_{\delta,\mathfrak{a}}(c\_{0},c\_{0}) Y\_{\delta,\mathfrak{a}}(y,c\_{0})] d\mathfrak{y} \right). \end{split}$$

Dividing both sides by Φ*δ*,*α*(*c*0, *c*0) yields

−(*δu* + *c*) 2 *Y <sup>δ</sup>*,*α*(*u*, *c*0) = *η*2(*δu* + *c*) − *η*2(*δu* + *c*)*Yδ*,*α*(*u*, *c*0) − *η*3(*u* − *c*0) + *η*<sup>3</sup> *u c*0 *Yδ*,*α*(*t*, *c*0)d*t* <sup>+</sup> *λβ* Φ*δ*,*α*(*c*0, *c*0) *u c*0 *<sup>ω</sup>*(*t*)d*<sup>t</sup>* <sup>+</sup> *λβ <sup>u</sup>*−*c*<sup>0</sup> 0 (*u* − *c*<sup>0</sup> − *x*)d*F*(*x*) <sup>−</sup> *λβ <sup>u</sup>*−*c*<sup>0</sup> 0 d*F*(*x*) *u*−*x c*0 *Yδ*,*α*(*y*, *c*0)d*y* , −(*δu* + *c*) 2 *Y <sup>δ</sup>*,*α*(*u*, *c*0) = −*η*<sup>3</sup> *u c*0 d*t* − *η*2(*δu* + *c*)*Yδ*,*α*(*u*, *c*0) + *η*<sup>3</sup> *u c*0 *Yδ*,*α*(*t*, *c*0)d*t* <sup>+</sup> *λβ* Φ*δ*,*α*(*c*0, *c*0) *u c*0 *<sup>ω</sup>*(*t*)d*<sup>t</sup>* <sup>+</sup> *λβ <sup>u</sup>*−*c*<sup>0</sup> 0 (*u* − *c*<sup>0</sup> − *x*)d*F*(*x*) <sup>−</sup> *λβ <sup>u</sup>*−*c*<sup>0</sup> 0 d*F*(*x*) *u*−*x c*0 *Yδ*,*α*(*y*, *c*0)d*y* .

Finally, we obtain

$$\begin{split} -\delta^{2}\mathbf{Y}\_{\delta,\mu}^{\prime}(\boldsymbol{u},\boldsymbol{c}\_{0}) &= \ 2\delta\boldsymbol{c}\mathbf{Y}\_{\delta,\mu}^{\prime}(\boldsymbol{u},\boldsymbol{c}\_{0}) + \varepsilon^{2}\mathbf{Y}\_{\delta,\mu}^{\prime}(\boldsymbol{u},\boldsymbol{c}\_{0}) - \eta\_{1}\int\_{\boldsymbol{c}\_{0}}^{\boldsymbol{u}}\mathrm{d}t - \eta\_{2}\delta\boldsymbol{u}\boldsymbol{Y}\_{\delta,\mu}(\boldsymbol{u},\boldsymbol{c}\_{0}) \\ &- \quad \eta\_{2}\varepsilon\mathbf{Y}\_{\delta,\mu}(\boldsymbol{u},\boldsymbol{c}\_{0}) + \eta\_{3}\int\_{\boldsymbol{c}\_{0}}^{\boldsymbol{u}}\boldsymbol{Y}\_{\delta,\mu}(\boldsymbol{t},\boldsymbol{c}\_{0})\mathrm{d}t + \frac{\lambda\beta}{\mathsf{\mathsf{B}}\_{\delta,\mu}(\boldsymbol{c}\_{0},\boldsymbol{c}\_{0})}\int\_{\boldsymbol{c}\_{0}}^{\boldsymbol{u}}\boldsymbol{\omega}(\boldsymbol{t})\mathrm{d}t \\ &+ \quad \lambda\beta\int\_{0}^{\boldsymbol{u}-\boldsymbol{c}\_{0}}(\boldsymbol{u}-\boldsymbol{c}\_{0}-\boldsymbol{x})\mathrm{d}\mathcal{F}(\boldsymbol{x}) - \lambda\beta\int\_{0}^{\boldsymbol{u}-\boldsymbol{c}\_{0}}\mathrm{d}\mathcal{F}(\boldsymbol{x})\left(\int\_{\boldsymbol{c}\_{0}}^{\boldsymbol{u}-\boldsymbol{x}}\boldsymbol{Y}\_{\delta,\mu}(\boldsymbol{y},\boldsymbol{c}\_{0})\mathrm{d}\boldsymbol{y}\right). \end{split} \tag{13}$$

Taking the Laplace transform of Equation (13) gives

$$\begin{split}-\delta^{2}\Big(s\overline{\boldsymbol{y}}\_{\delta,a}^{''}(\boldsymbol{s})+2\overline{\boldsymbol{y}}\_{\delta,a}^{'}(\boldsymbol{s})\Big)&=&-2\delta c\Big(s\overline{\boldsymbol{y}}\_{\delta,a}^{'}(\boldsymbol{s})+\overline{\boldsymbol{y}}\_{\delta,a}(\boldsymbol{s})\Big)+c^{2}s\overline{\boldsymbol{y}}\_{\delta,a}(\boldsymbol{s})\\ &+&\eta\_{2}\delta\overline{\boldsymbol{y}}\_{\delta,a}^{'}(\boldsymbol{s})-\eta\_{2}c\overline{\boldsymbol{y}}\_{\delta,a}(\boldsymbol{s})+\frac{\eta\_{3}}{s}\overline{\boldsymbol{y}}\_{\delta,a}(\boldsymbol{s})+\frac{\lambda\beta}{\Phi\_{\delta,a}(\boldsymbol{c}\_{0})}\frac{\overline{\boldsymbol{w}}(\boldsymbol{s})}{s}\\ &-&\eta\_{1}\frac{e^{-s\boldsymbol{c}\_{0}}}{s^{2}}+\frac{\lambda\beta}{s^{2}}e^{-s\boldsymbol{c}\_{0}}f(\boldsymbol{s})-\frac{\lambda\beta}{s}f(\boldsymbol{s})\overline{\boldsymbol{y}}\_{\delta,a}(\boldsymbol{s}),\end{split}$$

which completes the proof.

The differential equation given in (11) is not straightforward to solve. Thus, one may think of transforming this equation to a well-known equation, which can be easily solved analytically or numerically.

For simplicity, let

$$\tilde{y}\_{\delta,\mathfrak{a}}(\mathsf{c}\_{0},\mathsf{c}\_{0}) = b\_{0} = \frac{e^{-\mathfrak{c}\_{0}^{2}}}{\mathfrak{c}\_{0}} - \frac{\tilde{\Phi}\_{\delta,\mathfrak{a}}(\mathsf{c}\_{0},\mathsf{c}\_{0})}{\Phi\_{\delta,\mathfrak{a}}(\mathsf{c}\_{0},\mathsf{c}\_{0})}, \text{ and } \quad \tilde{y}\_{\delta,\mathfrak{a}}'(\mathsf{c}\_{0},\mathsf{c}\_{0}) = b\_{1} = -\left[\frac{\Phi\_{\delta,\mathfrak{a}}'(\mathsf{c}\_{0},\mathsf{c}\_{0})}{\Phi\_{\delta,\mathfrak{a}}(\mathsf{c}\_{0},\mathsf{c}\_{0})} + e^{-\mathfrak{c}\_{0}^{2}} \left(1 + \frac{1}{\mathfrak{c}\_{0}^{2}}\right)\right].$$

In the following corollary, we transform Equation (11) into the Volterra equation of the second kind with a degenerate kernel.

**Corollary 1.** *Under the assumption of Theorem 1, Equation* (11) *becomes*

$$\mathcal{U}(\mathbf{s}) - \int\_{\varepsilon\_0}^{\mathbf{s}} \mathbf{K}(\mathbf{s}, t) \mathcal{U}(t) \mathbf{d}t \,\,=\,\, h(\mathbf{s}), \tag{14}$$

*where* U(*s*)*, K*(*s*, *t*)*, and h*(*s*) *are given by*

$$\begin{array}{rcl} q(\mathbf{s}) &=& \int\_{c\_0}^{\mathbf{s}} (\mathbf{s} - \mathbf{t}) \mathcal{U}(\mathbf{t}) \mathbf{d} \mathbf{t}, \\ \mathfrak{F}\_{\delta, \mathfrak{a}}(\mathbf{s}) &=& q(\mathbf{s}) + b\_1 (\mathbf{s} - \mathbf{c}\_0) + b\_{0\prime} \\ K(\mathbf{s}, t) &=& -[f\_1(\mathbf{s}) + f\_0(\mathbf{s})(\mathbf{s} - \mathbf{t})], \\ h(\mathbf{s}) &=& \mathfrak{g}(\mathbf{s}) - [b\_0 + b\_1(\mathbf{s} - \mathbf{c}\_0) f\_0(\mathbf{s}) + b\_1 f\_1(\mathbf{s})]. \end{array}$$

**Proof.** Consider Equation (11) and introduce the function *q*(*s*) such that *y*˜*δ*,*α*(*s*) = *q*(*s*) + *b*1(*s* − *c*0) + *b*<sup>0</sup> . Clearly, *q*(*c*0) = *q* (*c*0) = 0 and

$$
\tilde{y}\_{\delta,\mathfrak{a}}^{\prime\prime}(\mathfrak{s}) \, = \, q^{\prime\prime}(\mathfrak{s}), \quad \tilde{y}\_{\delta,\mathfrak{a}}^{\prime}(\mathfrak{s}) \, = \, q^{\prime}(\mathfrak{s}) + b\_{1}.
$$

Plugging these expressions into Equation (11) yields

$$q''(s) + f\_1(s)q'(s) + f\_0(s)q(s) \ = \
h(s),\tag{15}$$

where *h*(*s*) = *g*(*s*) − [*b*<sup>0</sup> + *b*1(*s* − *c*0)*f*0(*s*) + *b*<sup>1</sup> *f*1(*s*)]. Let us further introduce another auxiliary function <sup>U</sup>(*t*) such that *<sup>q</sup>*(*s*) = *<sup>s</sup> c*0 (*s* − *t*)U(*t*)d*t*.

U(*t*) satisfies the following Volterra equation:

$$\mathcal{U}(t) + \int\_{c0}^{s} f\_1(s) \mathcal{U}(t) \mathrm{d}t + \int\_{c0}^{s} f\_0(s)(s-t) \mathcal{U}(t) \mathrm{d}t = \mathcal{U}(s),$$

and this completes the proof.

The mathematical literature is replete with several papers that try to solve Equation (14). The book of Polyanin and Manzhirov [19] presents an overview of the analytical solution of (14), which gives a particular non-trivial solution of Equation (15) without the second member. The following remark presents the analytical solution in a general form given that we know a non-trivial solution (*q*1(*s*)).

**Remark 2.** *Let q*1(*s*) *be a non-trivial particular solution of Equation* (15) *without the second member. Then, the other solution of Equation* (15) *without the second member is in the form*

$$q\_2(s) = q\_1(s) \int\_{c\_0}^s \frac{\left| \frac{t}{c\_0} \right|^{\frac{\eta\_2 + 2}{\delta}} \exp\left[\frac{2c}{\delta}(c\_0 - t)\right]}{q\_1(t)^2} dt,$$

*from which the final solution for Equation* (15) *is given by*

$$q(s) = q\_2(s) \int\_{\ell\_0}^{\delta} \frac{q\_1(t)}{\left|\frac{t}{\varepsilon\_0}\right|^{\frac{\eta\_1+2}{\delta}} \exp\left[\frac{2\varepsilon}{\delta}(\mathbf{c}\_0 - t)\right]} h(t) \mathrm{d}t - q\_1(s) \int\_{\ell\_0}^{\delta} \frac{q\_2(t)}{\left|\frac{t}{\varepsilon\_0}\right|^{\frac{\eta\_1+2}{\delta}} \exp\left[\frac{2\varepsilon}{\delta}(\mathbf{c}\_0 - t)\right]} h(t) \mathrm{d}t.$$

Although this solution is interesting, the main problem is the determination of *q*1(*s*), which is not straightforward. An alternative way to solve this problem is to find the solution numerically. In recent years, scholars in applied mathematics have investigated the numerical solution of Equation (14). For an overview of these method, the readers are referred to Maleknejad et al. [20] and Nadir and Rahmoune [21] and the references therein. The following remark presents the numerical way to solve (14) by Simpson's quadrature rule as given in Nadir and Rahmoune [21].

**Remark 3.** *Let us consider the transformed Equation* (14) *given below:*

$$\mathcal{U}(s) + \int\_{\varepsilon\_0}^{s} \mathcal{K}\_0(s, t) \mathcal{U}(t) \mathrm{d}t \, = \,\, h(s), \tag{16}$$

*where K*0(*s*, *t*) = *f*1(*s*) + *f*0(*s*)(*s* − *t*)*.*

*For s*<sup>0</sup> = *c*<sup>0</sup> < *s*<sup>1</sup> < ··· < *s*2*<sup>n</sup> and, moreover, for notation simplicity, we set* U(2*j*) = U(*s*2*j*)*, K*0(2*j*, 2*j*) = *K*0(*s*2*j*,*s*2*j*), *h*(2*j*) = *h*(*s*2*j*)*. Hence, for j* = 1, ··· , *n, it was proven in Nadir and Rahmoune [21] that*

$$\begin{aligned} l\mathcal{U}(2j)\left(1-\frac{h}{3}\left[2K\_0(2j,2j-1)+K\_0(2j,2j)\right]\right) &=& l(2j)+\frac{h}{3}\left(K\_0(2j,0)+2K\_0(2j,1)\right)\mathcal{U}(\varepsilon\_0)+\frac{2}{3}h\\ &\quad \times \quad \sum\_{i=1}^{j-1} \left(K\_0(2j,2i-1)+K\_0(2j,2i)+K\_0(2j,2i+1)\right) \\ &\quad \times \quad l\mathcal{U}(2i),\end{aligned}$$

*where h is the step of subdivision,* <sup>U</sup>(2*<sup>j</sup>* <sup>+</sup> <sup>1</sup>) = <sup>U</sup>(2*j*) + <sup>U</sup>(2*<sup>j</sup>* <sup>+</sup> <sup>2</sup>) <sup>2</sup> *,* <sup>U</sup>(*c*0) = *<sup>h</sup>*(*c*0)*.*

This formula provides us with a powerful tool to approximate the Gerber–Shiu function. Section 4 aims to investigate the numerical solution of *y*˜*δ*,*α*(*s*) and, thus, its inverse Laplace transform in order to find the Gerber–Shiu function in particular cases where a given distribution of conditional claim sizes is assumed with a particular penalty function.

#### *3.2. Exponentially Distributed Conditional Claim*

Equation (4) can be transformed to a non-homogeneous differential equation of order three under an exponentially distributed conditional claim. The following proposition states the form of the differential equation satisfied by the discounted penalty function under this assumption.

**Proposition 2.** *In the model given in Equation* (1)*, under exponentially distributed conditional claims with parameter ρ, Equation* (4) *becomes*

$$\begin{split} \left(\mu - \varepsilon\_{0}\right)^{2} \Phi\_{\delta,\mathfrak{a}}^{\prime\prime\prime}(\mu, \mathfrak{c}\_{0}) &= \left[\mu\_{0}(\mu - \varepsilon\_{0}) - \rho(\mu - \varepsilon\_{0})^{2}\right] \Phi\_{\delta,\mathfrak{a}}^{\prime\prime}(\mu, \mathfrak{c}\_{0}) + \frac{\lambda\beta}{\delta^{2}} \left[\rho\omega(\mu) - \omega^{\prime}(\mathfrak{u})\right] \\ &+ \mu\_{3} \Phi\_{\delta,\mathfrak{a}}(\mu, \mathfrak{c}\_{0}) + \left[\mu\_{1} + \mu\_{2}(\mathfrak{u} - \mathfrak{c}\_{0})\right] \Phi\_{\delta,\mathfrak{a}}^{\prime}(\mu, \mathfrak{c}\_{0}), \end{split} \tag{17}$$

*where*

$$
\mu\_0 = \begin{array}{c c c c c} \frac{\eta\_0}{\delta} - 2; & \mu\_1 & = & \frac{\eta\_0 \delta - \eta\_1}{\delta^2}, \quad \mu\_2 & = & \frac{\rho \eta\_0}{\delta} \text{ and } \quad \mu\_3 & = & \frac{\rho}{\delta^2} (\lambda \beta - \eta\_1).
\end{array}
$$

**Proof.** Differentiating Equation (4) gives:

$$\begin{split} \left(\delta\mu + c\right)^{2} \Phi\_{\delta,\mathfrak{a}}^{\prime\prime\prime}(\mathfrak{u}, c\_{0}) &= \left(\eta\_{0} - 2\delta\right) (\delta\mathfrak{u} + c) \Phi\_{\delta,\mathfrak{a}}^{\prime\prime}(\mathfrak{u}, c\_{0}) + (\eta\_{0}\delta - \eta\_{1}) \Phi\_{\delta,\mathfrak{a}}^{\prime}(\mathfrak{u}, c\_{0}) \\ &+ \quad \lambda\beta \frac{\mathbf{d}}{\mathrm{d}u} \bigg[ e^{-\rho(u)} \int\_{c\_{0}}^{u} \Phi\_{\delta,\mathfrak{a}}(y) \rho e^{\rho y} \mathrm{d}y \bigg] + \lambda\beta \omega^{\prime}(u) . \end{split} \tag{18}$$

That is,

$$\begin{split} \left(\delta\mu + c\right)^{2} \Phi\_{\delta,\mu}^{\prime\prime\prime}(\mu, c\_{0}) &= \left(\eta\_{0} - 2\delta\right) \left(\delta\mu + c\right) \Phi\_{\delta,\mu}^{\prime\prime}(\mu, c\_{0}) + \left(\eta\delta - \eta\_{1}\right) \Phi\_{\delta,\mu}^{\prime}(\mu, c\_{0}) + \lambda\beta\omega^{\prime}(\mu) \\ &+ \quad \lambda\beta\rho \Big(\Phi\_{\delta,\mu}(\mu, c\_{0}) - \int\_{0}^{\mu - \varepsilon\_{0}} \Phi\_{\delta,\mu}(\mu - \mathbf{x}, c\_{0}) \rho e^{-\rho x} \mathbf{dx}\Big). \end{split} \tag{19}$$

From Equation (4), we have

$$\begin{split} -\lambda\beta \int\_{0}^{u-\varepsilon\_{0}} \Phi\_{\delta,\mathfrak{a}}(u-\mathbf{x},\mathfrak{c}\_{0}) \rho \mathbf{e}^{-\rho \mathbf{x}} \mathrm{d}\mathbf{x} &= \begin{split} -(\delta u + \varepsilon)^{2} \Phi\_{\delta,\mathfrak{a}}^{\eta}(\mathfrak{a},\mathfrak{c}\_{0}) + \eta\_{0}(\delta u + \varepsilon) \Phi\_{\delta,\mathfrak{a}}^{\prime}(\mathfrak{a},\mathfrak{c}\_{0}) \\ -\eta\_{1} \Phi\_{\delta,\mathfrak{a}}(\mathfrak{a},\mathfrak{c}\_{0}) + \lambda\beta \omega(u) .\end{split} \end{split} \tag{20}$$

Plugging Equation (20) into (19) yields

$$\begin{split} \left(\delta\mu + c\right)^{2} \Phi\_{\delta,\mathfrak{a}}^{\prime\prime\prime}(\mathfrak{u}, c\_{0}) &= \left[ \left(\eta\_{0} - 2\delta\right)(\delta\mathfrak{u} + c) - \rho(\delta\mathfrak{u} + c)^{2} \right] \Phi\_{\delta,\mathfrak{a}}^{\prime\prime}(\mathfrak{u}, c\_{0}) \\ &+ \quad \lambda\beta \left[ \rho\omega(\mathfrak{u}) + \omega^{\prime}(\mathfrak{u}) \right] + \left[ (\eta\_{0}\delta - \eta\_{1}) + \rho\eta\_{0}(\delta\mathfrak{u} + c) \right] \Phi\_{\delta,\mathfrak{a}}^{\prime}(\mathfrak{u}, c\_{0}) \\ &- \quad \rho\eta\_{1}\Phi\_{\delta,\mathfrak{a}}(\mathfrak{u}, c\_{0}) + \lambda\beta\rho\Phi\_{\delta,\mathfrak{a}}(\mathfrak{u}, c\_{0}). \end{split}$$

Dividing both sides of the above formula by *δ*<sup>2</sup> completes the proof.

Under some assumptions on the penalty function (a penalty function depends only on the deficit at ruin), Equation (17) can be solved via the sequence expansion technique.

**Theorem 2.** *Suppose that the penalty function depends only on the deficit at ruin. Then, Equation* (17)*, with initial condition* Φ*δ*,*α*(*c*0, *c*0)*, has the following solution:*

$$\Phi\_{\delta,\mathfrak{a}}(\mathfrak{u},\mathfrak{c}\_0) = \Phi\_{\delta,\mathfrak{a}}(\mathfrak{c}\_0,\mathfrak{c}\_0) \sum\_{n=0}^{\infty} \left[ \prod\_{i=1}^n (\nu)\_i \right] (\mathfrak{u} - \mathfrak{c}\_0)^{\mathfrak{u}}\_{\ \prime} \tag{21}$$

*where* ∏<sup>0</sup> *<sup>i</sup>*=<sup>1</sup> = 1 *and*

$$(\nu)\_n = \frac{(n-1)[\mu\_2 - \rho(n-2)] + \mu\_3}{n[(n-1)(n-2-\mu\_0) - \mu\_1]}.$$

**Proof.** Under the assumption of Theorem 2, Equation (17) becomes

$$\begin{split} \left(\mu - \mathfrak{c}\_{0}\right)^{2} \Phi\_{\delta,\mathfrak{a}}^{\prime\prime}(\mathsf{u}, \mathsf{c}\_{0}) &= \left[\mu\_{0}(\mathsf{u} - \mathsf{c}\_{0}) - \rho(\mathsf{u} - \mathsf{c}\_{0})^{2}\right] \Phi\_{\delta,\mathfrak{a}}^{\prime\prime}(\mathsf{u}, \mathsf{c}\_{0}) \\ &+ \left[\mu\_{1} + \mu\_{2}(\mathsf{u} - \mathsf{c}\_{0})\right] \Phi\_{\delta,\mathfrak{a}}^{\prime}(\mathsf{u}, \mathsf{c}\_{0}) + \mu\_{3} \Phi\_{\delta,\mathfrak{a}}(\mathsf{u}, \mathsf{c}\_{0}). \end{split} \tag{22}$$

Let Φ*δ*,*α*(*u*, *c*0) be defined by

$$\Phi\_{\delta,\mathfrak{a}}(u,\mathfrak{c}\_0) \quad := \sum\_{n=0}^{\infty} a\_n (u - \mathfrak{c}\_0)^n \text{ , hence}$$

$$\Phi\_{\delta,\mathfrak{a}}^{(i)}(u,\mathfrak{c}\_0) \quad = \sum\_{n=i}^{\infty} n(n-1) \cdot \cdots \cdot (n-i+1) a\_n \left(u - \mathfrak{c}\_0\right)^{n-i} \text{ }\tag{23}$$

where Φ(*i*) *<sup>δ</sup>*,*α*(*u*, *c*0) is the *i*-th derivative of Φ*δ*,*α*(*u*, *c*0), for *i* = 1, 2, 3. Plugging Equation (23) into (22) yields

$$\begin{aligned} \sum\_{n=1}^{\infty} n(n-1)(n-2)a\_{\mathfrak{n}}(\mu - c\_{0})^{n-1} &= \quad \mu\_{0} \sum\_{n=1}^{\infty} n(n-1)a\_{\mathfrak{n}}\left(\mu - c\_{0}\right)^{n-1} \\ &- \quad \rho \sum\_{n=0}^{\infty} n(n-1)a\_{\mathfrak{n}}\left(\mu - c\_{0}\right)^{n} + \mu\_{3} \sum\_{n=0}^{\mathfrak{n}} a\_{\mathfrak{n}}\left(\mu - c\_{0}\right)^{n} \\ &+ \quad \mu\_{1} \sum\_{n=1}^{\infty} n a\_{\mathfrak{n}}\left(\mu - c\_{0}\right)^{n-1} + \mu\_{2} \sum\_{n=0}^{\infty} n a\_{\mathfrak{n}}\left(\mu - c\_{0}\right)^{n}. \end{aligned}$$

Letting *n* be *n* + 1, then

$$\begin{aligned} \left(\sum\_{n=0}^{\infty} (n+1)(n)(n-1)a\_{n+1}\left(u - c\_0\right)^n \right. &= \left. \mu\_0 \sum\_{n=0}^{\infty} (n+1)(n)a\_{n+1}\left(u - c\_0\right)^n \right. \\ &\left. - \left. \rho \sum\_{n=0}^{\infty} n(n-1)a\_n \left(u - c\_0\right)^n + \mu\_3 \sum\_{n=0}^{\mu} a\_n \left(u - c\_0\right)^n \right. \\ &\left. + \left. \mu\_1 \sum\_{n=0}^{\infty} (n+1)a\_{n+1}\left(u - c\_0\right)^n + \mu\_2 \sum\_{n=0}^{\infty} n a\_n \left(u - c\_0\right)^n \right. \end{aligned}$$

Thus,

$$\begin{aligned} &\sum\_{n=0}^{\infty} \left[ (n+1)(n)(n-1) - \mu\_0(n+1)(n) - \mu\_1(n+1) \right] a\_{n+1} \left( \mu - c\_0 \right)^n \\ &= \quad \sum\_{n=0}^{\infty} \left[ \mu\_5 + \mu\_2(n) - \rho(n)(n-1) \right] a\_n \left( \mu - c\_0 \right)^n \end{aligned}$$

from which it follows that

$$a\_{n+1} = \frac{n\left(\mu\_2 - \rho(n-1)\right) + \mu\_3}{(n+1)\left[n\left((n-1) - \mu\_0\right) - \mu\_1\right]} a\_n.$$

By expanding the above relation recursively, we have

$$\begin{aligned} a\_n &= \, \_0\prod\_{i=1}^n \left[ \frac{(i-1)\left(\mu\_2 - \rho(i-2)\right) + \mu\_3}{i\left[ (i-1)\left((i-2) - \mu\_0 \right) - \mu\_1 \right]} \right] \\ &= \, \_0\prod\_{i=1}^n (\nu)\_{i\prime} \end{aligned}$$

which proves the statement by inserting this relation back into Equation (23) and by the fact that Φ*δ*,*α*(*c*0, *c*0) = *a*0.

**Remark 4.** *In Proposition 1, if we let α* = 0*, we obtain the integro-differentiating of the (relative) absolute ruin. Under the exponential claim distribution, one obtains the explicit expression of the (relative) absolute ruin probability by letting α* = 0 *in Theorem 2.*

#### **4. Numerical Illustration**

This section describes how to solve the Volterra equation given in Equation (16) and from there calculate the ruin probability using the inverse Laplace transform of Equation (10).

As outlined by Nadir and Rahmoune [21], we first find the solution to Equation (16) by numerical approximation. Moreover, since we must integrate the latter function to obtain the function *q*(*s*) given by Equation (14) and, hence, the Laplace transform function *y*˜*δ*,*α*(*s*), we approximate the solution to Equation (16) by a function using spline approximation (the spline approximation was chosen over polynomial, cubic spline, or Lagrange polynomial, linear approximations because it yields the best accuracy). Although the rational approximation is also a good approximation, we prefer the late one since the rational approximation produces a singular function that does not allow the integral to converge. *Yδ*,*α*(*u*) is then obtained by applying the inverse Laplace to the function *y*˜*δ*,*α*(*s*) derived above.

In the absence of real data, we chose our parameters arbitrarily. For comparison purposes, we varied the distribution parameter of the reporting delay time to examine its impact on the (relative) absolute ruin. Numerical calculations show that the ruin probability decreases as the reporting delay increases on average, regardless of the initial capital, which confirms that the delay in claim reporting can positively affect the ruin probability.

The numerical setting assumes a Gamma distribution with parameter *a* = 10.2, *b* = 0.008, which means that each claim average is 1275.

Table 1, as well as Table 2 show that the ruin probability is impacted by the initial surplus; however, the ruin probability is normally not influenced by only the initial surplus (increasing the initial reserve will not always decrease the ruin probability), as this may depend on the delay in reporting or settlement times as well. In comparing Tables 1 and 2, it turns out that an increase in the average reporting delay decreases the ruin probability. This finding potentially provides risk managers and actuaries with a powerful instrument to manage risks, in particular the risk of insolvency.


**Table 1.** Ruin probability for *δ* = 0.02, *λ* = 0.01, *c* = 0.1%, *β* = 0.05.

**Table 2.** Ruin probability for *δ* = 0.02, *λ* = 0.01, *c* = 0.1%, *β* = 0.02.


#### **5. Discussion and Conclusions**

This study analyzed the impact of the reporting delay time on the Gerber–Shiu discounted penalty function under an interest rate environment, where the claim effective's cost distribution is determined by the occurrence of the claim and the reporting delays (this somehow incorporated indirectly the inflation on the claim cost). We showed that the Gerber–Shiu function follows a degenerate Volterra equation with an arbitrary distribution of conditional claim sizes. A numerical illustration of our result using an approximation

method (prescribed by Nadir and Rahmoune [21]) confirmed the effect of the delay in reporting time on the ruin probability. In addition, we derived a closed analytical expression of the Gerber–Shiu function under the exponential conditional claim distribution assumption. The numerical results indicated that the ruin probability as a function of the initial surplus is not a monotonic function as one usually expects (increases of the initial surplus decrease the ruin probability).

Although this study contributes to the existing literature, it has some drawbacks in that we assumed that when the surplus belongs to −*c <sup>δ</sup>* , 0 , the insurer could continue to operate without borrowing (in contrast with Gerber and Shiu [4] and Cai [5]), which may be the case when *<sup>c</sup> <sup>δ</sup>* is small enough and there is flexibility in the delayed reporting times. As a result, insurers must seek assets with high risk-free rates and protect only low-risk contracts (*δ* >> *c*). A good distribution of the reporting delayed time must be bounded from above, due to regulation constraints (anti-selection ban and maximum accepting reporting time). In future studies, we will extend our results by analyzing the discounted penalty function with stochastic interest rates, since evidence in finance studies has demonstrated that there is no such rate as a risk-free rate for reporting.

**Author Contributions:** The authors contributed equally to this work from the Introduction to Conclusions. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by the Global Excellence Stature Fellowship 4.0 and NRF INCENTIVE GRANT from the University of Johannesburg.

**Conflicts of Interest:** The authors declare no conflict of interest.

## **References**

