*1.2. Novelties*

Although other hypotheses can be taken, such as, for example, considering that the number of claims in a period, *N*, follows a negative binomial distribution ([14]), academic literature on BMSs usually assumes that *N* is a Poisson random variable (RV). So, *N* ∼ Po(*λ*) and the parameter *λ*, the claim frequency, is perfectly known and can be interpreted as a risk measure of the policy. However, in a more realistic approach, some authors like [15] use intervals to quantify the uncertainty about the parameter of a distribution function that governs a risk variable. This is also the case within a BMS framework of [11], who model the uncertainty about *λ* by means of a modal interval. One extended way to combine randomness and uncertainty of parameters of distribution functions consists in modeling these parameters as FNs. It has been done both for continuous RVs ([16–18]), and in the discrete case ([19]). Following this approach, [20,21] and also [15] model risk financial parameters with FNs. In the actuarial field, FNs have been used to capture the uncertainty of insurance pricing variables ([22]) but also to model parameters that quantify risks. In this regard, we can point out [23] in a non-life insurance context, [24] to interpret the parameter that quantifies the dependence in a Farlie-Gumbel-Morgestein copula, and [25–28] in life insurance pricing. Since any interval can be seen as the *α*-cut of a FN, even in the case of improper intervals ([29]), in this work we consider that *λ* is fitted by means of a FN and, more particularly, by a triangular fuzzy number (TFN). So, this paper builds up a framework to model Markovian BMSs that embed the standard case, where the risk parameter *λ* is crisp, but also the method developed in [11] that quantifies this parameter as a modal interval. Standard BMSs provide point values for the stationary state and the mean asymptotic premium. Modal BMSs as introduced in [11] allow obtaining these variables as modal intervals whose lower and upper bounds may be understood as pessimistic/optimistic scenarios. Our method generalizes both types of BMSs since it quantifies variables related to BMS as FNs. On the one hand, these FNs can be understood as a set of crisp outcomes with an associated possibility measure. On the other, these FNs can be interpreted as a set of intervals that come from pessimistic/optimistic scenarios and are structured by means of possibility levels. Figure 1 shows a graphical synthesis of the methodological framework developed in this paper.

**Figure 1.** Graphical representation of our fuzzy bonus-malus systems (FBMS) model.

Other more complex forms of FNs, such as generalized FNs (GFNs) or intuitionistic FNs (IFNs), could be considered to quantify uncertain probabilities. Tools like GFNs or IFNs provide a more complete capture of uncertainty than FNs. However, their adjustment has a greater cost than in the case of triangular FNs since they incorporate more parameters and their computational handle may be more expensive as well. Therefore, using TFNs supposes a balance between the simplicity of crisp or modal interval probabilities and more complex representations of uncertain quantities such as GFNs or IFNs.

It should be noted that there are several scientific fields in which MCs are in use. In the field of economics and finance, we can observe applications within Leontief's inputoutput model, credit risk measurement, asset price volatility modeling, life insurance, etc. In addition, MCs have shown their usefulness in many other areas: industrial engineering (e.g., queuing theory), computer science (e.g., computer performance evaluation and web search engines), healthcare (e.g., pandemics transmission or evolution of ICU patients), etc. Hence, although our developments are carried out within a non-life insurance context, most of the results can be applied to any problem modeled by means of MCs when the transition probabilities (or the parameters that define them) are not precisely known.

The paper is organized as follows. Section 2 describes briefly how BMSs work. Section 3 shows the basic concepts of FNs and FMCs used throughout the paper. In Section 4, a methodologic approach is proposed to fit a fuzzy BMS (FBMS) when the number of claims within a period, *N*, follows a Poisson distribution with fuzzy parameter *λ*. This methodology is applied to the Irish BMS. A sensitivity analysis is conducted in Section 5. Finally, in Section 6, the work ends with a summary of its main contributions and potential extensions.

#### **2. Markovian Bonus-Malus Systems in Non-Life Insurance**

A BMS is a usual way to deal with risk aversion and moral hazard in some types of insurance, e.g., automobile third-party liability insurance [1]. BMSs classify insureds in *r* classes in such a way that the percentage of the base-premium to be paid by the *i*th class, *bi*, satisfies *bi* < *bi*+1, *i* = 1, 2, ... ,*r* − 1. In a BMS, the transition between classes is governed by a set of rules defined over the insured's number of claims in the current period. To summarize, it can be said that every BMS is determined by three elements (see Table 1):



**Table 1.** Elements of a bonus-malus system (BMS).

Source: Own elaboration based on [10].

Let us model the insured's class at time *t* as a discrete stochastic process (*Xt*)*<sup>t</sup>*∈N, being its state space the classes *S* = {1, 2, . . . ,*<sup>r</sup>*} ⊂ N. Furthermore, as it is usually done in the literature (e.g., [10]), we consider that the BMS is a finite MC, i.e., *p*(*Xn*+<sup>1</sup> = *in*+<sup>1</sup>| *X*0 = *i*0, *X*1 = *i*1,..., *Xn* = *in*) = *p*(*Xn*+<sup>1</sup> = *in*+<sup>1</sup>| *Xn* = *in*), ∀*i*0, *i*1, ... , *in*+<sup>1</sup> ∈ *S*. An insurer uses a finite Markovian BMS when the following conditions hold [1]:


A finite MC is said to be homogeneous if *p*(*Xm*+*<sup>h</sup>* = *j*| *Xm* = *i*) does not depend on *m*. In this case, transition probabilities *pij* = *p*(*Xt*+<sup>1</sup> = *j*| *Xt* = *i*), i.e., probabilities of moving from class *i* to class *j* in one-step (period), can be collected in a transition matrix *P* = *pij r*

with order *r*. The elements *pij* satisfy *pij* ≥ 0 and ∑ *j*=1 *pij* = 1.

If *p*(0) *i i*∈*S* denote the probabilities of initially being in state *i* ∈ *S*, the probabilities of being in state *i* ∈ *S* after *n* periods, *p*(*n*) *i i*∈*S* are:

$$p\_i^{(n)} = p\_i^{(n-1)} \cdot P = p\_i^{(n-2)} \cdot P \cdot P = \dots = p\_i^{(0)} \cdot P^n \tag{1}$$

where *Pn* is represented by *P*(*n*) = *p*(*n*) *ij* and *p*(*n*) *ij* are the probabilities of moving from state *i* to state *j* in *n* period. From Equation (1), it follows that:

$$p\_{ij}^{(n)} = f\_{ij}^{(n)}(p\_{11}, \dots, p\_{1r}, p\_{21}, \dots, p\_{2r}, \dots, p\_{1r}, \dots, p\_{rr}) \tag{2}$$

for some continuous functions *f* (*n*) *ij* , i.e., the elements in *Pn* are some functions of the elements in *P*.

A homogeneous MC is regular if each state is accessible from any other state, either in one step or more, i.e., there exists *n* ∈ N such that *p*(*n*) *ij* > 0 ∀ *i*, *j*. One of the features that characterize regular MCs is its stationary distribution, which represents the probability of the chain being at each state after a large number of periods, namely, lim*n*→∞*Pn* = *π*, where the rows of *π* are identical. So, any regular MC with transition matrix *P* has a stationary distribution, *π*, such that:

$$
\boldsymbol{\pi} = \boldsymbol{\pi} \cdot \boldsymbol{P} \tag{3}
$$

The vector *π* = *<sup>π</sup>j<sup>j</sup>*=1,2,...,*<sup>r</sup>* can be interpreted as the probability that an insured belongs to class *j*, *j* = 1, 2, ... ,*r* after *n* periods, *n* → ∞. That vector does not depend on the insured's initial class, *i*0. So, two main outputs in a BMS are:


$$b^\* = \sum\_{j=1}^r b\_j \pi\_j \tag{4}$$

The mean asymptotic premium, *b*∗ is a concept of the utmost importance because it has been intensively used to assess the efficiency of a BMS (e.g., [1,6,30,31]).

BMSs consider the number of claims, *N*, as a discrete RV. In our paper, as it is commonplace in actuarial literature, *N* is supposed to follow a Poisson distribution with parameter *λ*, *N* ∼ Po(*λ*), ([1,7–12]). Therefore:

$$p(N=k) = \frac{\lambda^k}{k!}e^{-\lambda} \tag{5}$$

where *p*(·) stands for a probability measure.

Poisson RVs are often used in actuarial modeling due to their interesting arithmetical properties. Furthermore, the risk parameter *λ* can be fitted specifically to each insured taking into account relevant rating factors (e.g., gender, age, and social status) by using a generalised linear model (GLM) ([5,8,32]).

If *N* is modeled with (5), *pij* is a function of the risk parameter *λ*, *hij*(*λ*), in such a way that the BMS probabilities in the one-step transition matrix are:

$$h\_{ij}(\lambda) = p\_{ij} = \sum\_{k=0}^{\infty} \frac{\lambda^k}{k!} e^{-\lambda} t\_{ij}(k) \tag{6}$$

where *tij*(*k*) = 1 if *k* causes the transition from *i* to *j* and 0 otherwise.

Numerical application 1. Let *λ* = 0.04 in Equation (5) for the Irish BMS in Table 2. The transition matrix, *P*, that corresponds to this BMS is:


**Table 2.** Irish bonus-malus system.


Source: [10].

From (3), *π* = (0.916232, 0.037394, 0.038921, 0.003861, 0.002523, 0.001069) and, by considering Equation (4), the mean asymptotic premium is *b*∗ = 51.423.

In this paper, we will consider that the risk parameter *λ* cannot be determined precisely. Uncertainty may be the result of different causes: stochastic variability, inaccuracy, incomplete information, etc. Stochastic variability can be described by using RVs or stochastic processes, but inaccuracy and incomplete information can be captured by means of intervals or FNs. Given that any interval can be interpreted as the *α*-cut of a FN, in this work it is assumed that *λ* is a FN.

#### **3. Fuzzy Numbers and Fuzzy Markov Chains**
