Fuzzy Markovian Bonus-Malus Systems in Non-Life Insurance

: Markov chains (MCs) are widely used to model a great deal of financial and actuarial problems. Likewise, they are also used in many other fields ranging from economics, management, agricultural sciences, engineering or informatics to medicine. This paper focuses on the use of MCs for the design of non-life bonus-malus systems (BMSs). It proposes quantifying the uncertainty of transition probabilities in BMSs by using fuzzy numbers (FNs). To do so, Fuzzy MCs (FMCs) as defined by Buckley and Eslami in 2002 are used, thus giving rise to the concept of Fuzzy BMSs (FBMSs). More concretely, we describe in detail the common BMS where the number of claims follows a Poisson distribution under the hypothesis that its characteristic parameter is not a real but a triangular FN (TFN). Moreover, we reflect on how to fit that parameter by using several fuzzy data analysis tools and discuss the goodness of triangular approximates to fuzzy transition probabilities, the fuzzy stationary state, and the fuzzy mean asymptotic premium. The use of FMCs in a BMS allows obtaining not only point estimates of all these variables, but also a structured set of their possible values whose reliability is given by means of a possibility measure. Although our analysis is circumscribed to non-life insurance, all of its findings can easily be extended to any of the above-mentioned fields with slight modifications.


Motivation
A bonus-malus system (BMS) is a common method for posteriori ratemaking in nonlife insurance. It is based on partitioning the insurer's portfolio into a finite number of classes: bonus and malus classes. A typical case is automobile third-party liability insurance [1]. In a BMS, policyholders do not have a fixed price for their contracts throughout periods (e.g., the mathematical expectation of claims value per period). Their membership into a concrete BMS class is reviewed each period according to the number of claims in the previous one. Claim-free years are rewarded by discounts or bonuses on a base-premium; at-fault accidents are penalized by surcharges called maluses. Some overviews on how BMS are applied in different countries can be found in [2][3][4].
Following [5], in most commercial BMSs, by knowing the insured's class in the current period and fitting the statistical distribution for the number of claims per period, it is possible to determine the probabilities of the insured's class in the next period. Therefore, these BMSs are Markovian. For that reason, the academic literature on BMSs uses extensively MCs for their modeling ( [1,[5][6][7][8][9][10][11]). Therefore, a key question in a BMS is fitting the value of the one-step transition probability matrix. Following [12,13], if full knowledge of the probabilities of this matrix is not available, they have to be estimated somehow with the uncertainty that any estimation procedure involves. Uncertainty may be due to randomness, hazard, vagueness, incomplete information, etc. In our paper, we consider that the claiming process is probabilistic, but the uncertainty about the parameter that governs this random behavior is captured by means of a fuzzy number (FN) and, as a consequence, fuzzy Markov chains (FMCs) will be used.

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 is a Poisson random variable (RV). So, ~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][17][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][26][27][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. 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, , 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.

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 classes in such a way that the percentage of the base-premium to be paid by the th class, , satisfies < , = 1,2, … , − 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):


The initial class, where new insureds are assigned, .  The premium scale = ( , … , ).  The transition rules, that is to say, the rules that define the conditions for an insured in one class to be transferred to another class in the next period.

Class Premium Level Class after Claims
Source: Own elaboration based on [10].
Let us model the insured's class at time as a discrete stochastic process ( ) ∈ℕ , being its state space the classes = {1,2, … , } ⊂ ℕ. Furthermore, as it is usually done in the literature (e.g., [10]), we consider that the BMS is a finite MC, i.e., ( = | = , = , … , = ) = ( = | = ), ∀ , , … , ∈ . An insurer uses a finite Markovian BMS when the following conditions hold [1]: There exists a finite number of classes such that each insured stays in one class through each period.  The premium for each insured depends only on the class where they stay.  The class for a given period is determined by the class in the preceding period and the number of claims reported in that period.
A finite MC is said to be homogeneous if ( = | = ) does not depend on . In this case, transition probabilities = ( = | = ), i.e., probabilities of moving from class to class in one-step (period), can be collected in a transition matrix = with order . The elements satisfy ≥ 0 and ∑ =1. If are the probabilities of moving from state to state in period. From Equation (1), it follows that: for some continuous functions ( ) , i.e., the elements in are some functions of the elements in .
A homogeneous MC is regular if each state is accessible from any other state, either in one step or more, i.e., there exists ∈ ℕ such that ( ) > 0 ∀ , . 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 where the rows of are identical. So, any regular MC with transition matrix has a stationary distribution, , such that: The vector = , ,…, can be interpreted as the probability that an insured belongs to class , = 1,2, … , after periods, → ∞. That vector does not depend on the insured's initial class, . So, two main outputs in a BMS are:


The stationary distribution of , , as defined in Equation (3).  The mean asymptotic premium, * , i.e., the average premium paid by the insured in that stationary distribution, defined as: The mean asymptotic premium, * 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, , as a discrete RV. In our paper, as it is commonplace in actuarial literature, is supposed to follow a Poisson distribution with parameter , ~Po(l), ( [1,[7][8][9][10][11][12]). Therefore: where (·) 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 is modeled with (5), is a function of the risk parameter , ℎ ( ), in such a way that the BMS probabilities in the one-step transition matrix are: where ( ) = 1 if causes the transition from to and 0 otherwise. Numerical application 1. Let = 0.04 in Equation (5) for the Irish BMS in Table 2. The transition matrix, , that corresponds to this BMS is: From (3), = (0.916232, 0.037394, 0.038921, 0.003861, 0.002523, 0.001069) and, by considering Equation (4), the mean asymptotic premium is * = 51.423. Source: [10].
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.

Fuzzy Numbers
A fuzzy number is a fuzzy set on the referential set ℝ that satisfies (i) is normal, (ii) from where, if needed, the membership of could be obtained. The core of is and can be understood as the most reliable value of this TFN. i.e., the possibility of is 1. The support of is [ , ]. TFNs are used in countless practical applications including actuarial ones [22] because they are easy to handle arithmetically and they are well adapted to the way humans think of uncertain quantities. Moreover, when the information about a variable is vague and imprecise, the parsimony principle leads us to represent that information as simply as possible. The linear shape of TFNs meets that requirement. For instance, the uncertain quantity "approximately 0.04" can be represented in a very natural way as the TFN (0.038/0.04/0.042).
Likewise, let it be a TFN : Let be a continuous real-valued function of -real variables = ( , , … , ), ( ) = ( , , … , ). If are not crisp numbers, but FNs with -cuts ( ) = ( ), ( ) , = 1,2, … , , a FN is induced via such that = , , … , . It is often difficult to obtain a closed expression for the membership function of . However, following [33], the -cuts of , ( ), in the usual case where ( ), … , ( ) are not interactive, i.e., the variables that they quantify have an independent behavior, can be obtained as: where ( ) stands for the rectangular domain: So, the lower (upper) bounds of ( ), ( ) ( ( )), are the global minimum (maximum) of ( ) within the rectangular domain in Equation (13), that is to say: being , = 1,2, … , 2 , a vertex of the domain (13) and ( ), = 1,2, … , , an extreme value of the function within this domain that takes this value at point = . Therefore, if ( , , … , ) is monotonic, the lower and upper bounds of ( ), ( ) and ( ), are in one of the 2 vertexes of (13). Without loss of generality, let us suppose that increases with respect to , = 1,2, … , , ≤ , and decreases in the last − variables, [34] demonstrates that: (15) cannot be used to evaluate ( ). However, according to [35], the general formulation to obtain the lower and upper bounds of ( ) from Equations (12)-(14) is still valid but now the number of vertexes, , is less than 2 . In [35], the authors study the role of interactive fuzzy variables in decision-making problems and analyze some particular cases. Concretely, when ( ) is the mathematical expectation function, is the probability of the th outcome, and it is quantified as a FN, the domain in Equation (13) turns into: and Equation (14) becomes: due to the fact that ( ) is a linear function. It is worth noting that the result of evaluating a non-linear with the TFNs is not necessarily a TFN. However, often admits a good triangular approximation through the secant approach. It builds up the shape of the triangular approximate FN to by means of the secant lines that unite the 0-cut and the 1-cut of . Such that is a TFN as: This approximation, as shown in [36], works pretty well for nonlinear monotonic functions of TFNs such as product, division, power, etc. Likewise, [37,38] show that this approach fits satisfactorily common actuarial and financial calculations with TFN parameters, e.g., the present value of a stream of fuzzy cash-flows. Keeping the triangular shape of the initial data when handling FNs is quite interesting. According to [39], complex shapes of FNs can generate problems with calculations in computer work or interpreting results intuitively. [40] state that a triangular approximate is a kind of defuzzification that is richer than just transforming a FN into a crisp representative value. If defuzzification is carried out too early, a great loss of information occurs, so it is preferable to drag all the fuzzy information in the calculations for as long as possible. The triangular approximation involves a compromise between simplification in computation and interpretation, and not oversimplifying the value of fuzzy parameters. In addition, TFNs have a very intuitive interpretation and, therefore, from the insurance industry point of view, a triangular approximate to actuarial variables and parameters could be very useful in decision-making processes.

Fuzzy Markov Chains
Fuzzy set literature has provided three approaches to MC under fuzziness, namely FMC. The first one, due to [41], supposes fuzzy probabilities and proposes calculating the matrices , ∀ > 1, by applying Zadeh's extension principle [42]. The second approach, in [43], consists of defining the matrix that governs the transition between states by means of a fuzzy relation. In the third one, [2], like [41], suppose that the probabilities of the onestep transition matrix are FNs. However, Buckley and Eslami's framework of FMCs uses restricted matrix multiplication to operate with probabilities in such a way that the constraint of being a well-formed probability distribution always holds. That is to say, they take into account the interdependence between the probabilities of a distribution function, similarly to Equations (16) and (17). This paper follows this last approach.
Let us assume that some probabilities in the one-step transition matrix = are uncertain and are quantified by means of the FNs , with -cuts ( ). Now we have a fuzzy transition matrix = , with 0 ≤ ≤ 1. See Equations (10) and (11). Of course, some elements in may be crisp since crisp numbers are a particular case of a FN. FMCs defined by [2] have uncertainty in the transition probabilities but not in the set of outcomes, that is discrete. So, the following constraint on is added: To compute the -period transition matrix, ( ) , and the fuzzy stationary distribution, , the following process is implemented: Step 1. For a given ∈ [0,1] , obtain the matrix of intervals ( ) = ( ) = ( ), ( ) .
Step 3. Define the domain of the matrix ( ), for the given ∈ [0,1], as: This domain defines a set of matrices that satisfy that each row sums up 1 with a possibility level of at least . So, each matrix = , ∈ ( ), is a crisp MC.
Step 4. Since ( ) in Equation (19) are compact sets, ( ) in Equation (20) is also compact. So, any continuous function applied to its elements has a compact image. Then, if Equation (2) is applied, such an image for that is a compact interval. This interval is set as the -cut of the FN ( ) , i.e., ( ) ( ) = ( ) ( ) , which is surely normal [2].
To determine ( ) ( ) we must find the lower and upper bounds, ( ) ( ) and ( ) ( ) by solving: and Notice that Equations (21) and (22) can be easily solved in low-dimensional problems. However, in more complex problems it is necessary to use an algorithm (see, e.g., [44]) or a heuristic constrained optimization technique [45].

Implementing a Markovian Fuzzy Bonus-Malus System Governed by a Fuzzy Poisson Discrete Random Variable
In this Section, we propose an integral methodology to develop a FBMS under the hypothesis that ~Po . It embeds the fitting of the risk parameter as a TFN, the obtaining of fuzzy transition matrix and the triangular approximate of the stationary distribution calculated by using Equations (23) and (24), and also the determination of the fuzzy mean asymptotic premium in Equation (4).
Step 1. Fit the risk factor as a TFN. We point out three different options to estimate this parameter. Option 1 Given that has an intuitive interpretation since it is the mean number of claims in one period, it may be quantified as a FN based on experts' opinions. For example, an expert may judge that a concrete type of driver generates approximately one claim every 5 years and so the TFN l = 0.2 = (0.18/0.2/0.22) can be considered. Imprecise or subjective quantitative predictions can often come from a pool of experts, leading to a set of fuzzy quantifications. This set of fuzzy opinions can be aggregated simply by their arithmetic mean or other more sophisticated methods (see [46][47][48] for full details).
Option 2 Papers [49,50] consider a standard 1 − a statistical confidence interval as the observed a-cut of the FN, for some increasing values of ≤ a < 1, where is an arbitrary value near 0 (it is often chosen to be 0.001, 0.005 or 0.01). In [50] it is suggested that by placing those confidence intervals one on top of the other, a FN close to triangular-shaped is obtained. So, we point out two alternatives to apply that idea: b) Papers [49,51] propose making fuzzy predictions from statistical linear regression models. In [49] it is stated that a 1 − a statistical confidence interval of coefficients adjusted with a linear regression may be interpreted as the a-cut of a FN for these coefficients. Therefore, let us suppose that a GLM estimate of is determined, as usual, by: being , = 1, … , , the coefficients and the explanatory variables (e.g., age, gender, and driving experience in a car insurance context) that are crisp nonnegative observations (in fact, they are usually modeled as dichotomic variables). For the estimate of each coefficient, it is possible to generate a FN whose acuts, ( ), = 1, … , , are: where * is the GLM point estimate of and the standard deviation of that estimate. So, from the fuzzy function = ∑ , and bearing in mind Equations (15) and (28), the following FN is induced: A similar approach may be developed from the results in [51]. However, in this case, it must be taken into account that their approach to making fuzzy predictions from a statistical regression is built up from the interval predictions of residuals instead of using interval estimates of coefficients. So = ∑ where ̃ is a fuzzy error term induced from the residuals of the conventional regression and, so: Equations (26), (29) and (30) do not give a TFN. However, can be approximated as a TFN simply by using Equation (18).
Option 3 Fuzzy Regression Methods (FRMs) have been applied in several actuarial issues to fit relevant variables [52] for a comprehensive description of application areas). In this way, [53] fits the term structure of interest rates, [54,55] predicts claim provisions, and [25,28] adjusts the Lee-Carter mortality law.
To fit , the fuzzy extension of the log-Poisson regression by [55] may be used. It combines the conventional Poisson GLM and the minimum fuzziness principle by [56]. In this case, the coefficients in Equation (27)  Once the parameters have been estimated, obtaining = ∑ is straightforward (see Equation (29)).
Let us remark again that although is a TFN, is not. Nevertheless, can be approximated as a TFN ′ = ( / / ) with Equation (18): Step 2. Obtain the fuzzy transition matrix. We now suppose that after performing any of the options in Step 1, and the corresponding triangular approximate, the risk parameter is given as the TFN = ( / / ), Fuzzy transitions probabilities come from the fuzzified version of Equation (6): To obtain the -cuts of by using Equation (15), it is necessary to determine the sign of the first derivative of ℎ ( ). Let us show the case of the Irish BMS whose transition matrix is Expression (7), and is either zero, ( = 0) , ( = 1) , ( ≥ 1) or ( ≥ 2). Then: Similarly, any other probabilities for different FBMSs could be calculated. Notice that the FNs whose -cuts are Equations (37)-(40) do not have a triangular shape but they admit a triangular approximation by using the secant approach described in Section 3.1. If this is done, we obtain:  Table 3 shows that when approximating ( ≥ 1) with Equation (43), the errors incurred on the lower and upper bounds of its -cuts are negligible since they are never over 0.001%. Moreover, notice that we measure the performance of the calculations on a scale of eleven grades of possibility. Following [38], this scale provides sufficient discernment without being excessive since we are using imprecise data and, therefore, more precision is not necessary for a FN representation. Table 3. α-cuts of p (N ≥ 1), it's triangular approximate, p (N ≥ 1), and errors. Step 3. Determine the fuzzy stationary distribution function. Once the fuzzy transition matrix associated with the FBMS has been obtained, to determine the fuzzy stationary state, Steps 1 to 4 in Section 3.2 should be applied and, therefore, optimization problems in Equations (23) and (24) must be solved. Notice that although the probabilities , = 1, … , , are obtained by solving complex optimization problems, the results of the numerical applications 4 and 5, that have been obtained with the R package FuzzyStatProb by [13] (see Figures 1 and 2) suggest that its triangular approximate by using Equation (18) provides a satisfactory fitting. Numerical Application 4. Now we compute the fuzzy stationary distribution ( , , ) for the fuzzy transition matrix in numerical application 2. In order to do so, we use the R package FuzzyStatProb described in [13], which is based on the use of Equations (23) and (24). The pseudo-codes and codes used are included in Appendixes 1 and 2, respectively, and the result of ( , , ) in Figure 2.
It should be remarked that the probabilities obtained by [11] (p. 847) are intervals whose values are the 0-cuts of the probabilities in our FMC.  Table 4 shows the supports and cores of the fuzzy stationary state , = 1, … ,6 when considering fuzzy probabilities in Equations (37)- (40). The pseudo-codes and the codes of the R package FuzzyStatProb that have been used to get these results are included in Appendixes 1 and 2, respectively. Figure 3 depicts the graphical representation of , = 1, … ,6.  The shape of the fuzzy stationary distribution suggests that their triangular approximate must work quite well. Table 5 shows that the relative deviations of the lower and upper bounds of the approximate of , ( ), with respect to the respective bounds of ( ) are always under 1%. Therefore, the above intuition is confirmed in the case of . We have also observed that this fact also applies to the other stationary probabilities in the numerical application. So, it can be written: It is worth pointing out that the intervals fitted for the stationary state in [11] (p. 848) are the 0-cuts of , = 1, … ,6, in the fuzzy version of the Irish BMS.  and , their triangular approximates work pretty well. We consider this result interesting for two reasons: a) The calculations can be done easily with less computational effort. For example, in the two first rows of Table 5, we have performed the calculations on a scale with eleven grades of possibility. So, 20 optimization programs have been solved for a single probability (10 minimizing programs for the lower bounds of -cuts = 0, 0.1, … , 0.9, and other 10 maximizing programs for the respective upper bounds). Likewise, obtaining the 1-cut implies nothing but solving a conventional Markov chain. This computational effort is reduced drastically by using the triangular approximate in Equation (18), which leads us to obtain the results in columns 3 and 4 of Table 5. In this case, it is enough to solve 2 optimization programs (1 minimizing program for the lower bound of the 0-cut and 1 maximizing for the upper one) and also, of course, evaluating a conventional BMS in . The interest in this result is amplified by the fact that the Irish BMS is relatively simple (there are 6 classes) and so it embeds only 36 ( ) and 6 probabilities . However, BMSs often have more than 20 classes (e.g., Belgian or German BMSs). b) From the perspective of an actuary, a triangular approximate of the fuzzy probabilities can be very useful. A TFN provides an estimate of the most feasible, minimum, and maximum probability that can be interpreted intuitively without any knowledge of fuzzy set theory (FST). Therefore, the triangular approximates presented in this paper could facilitate the use of FBMSs in the insurance industry.
In fact, the problem to solve in this case is the same as that in [35]. A triangular approximate for * , * = ( * / * / * ), can be obtained by using the 0-cut and the 1-cut obtained from Equations (46) and (47) In this latter way, 20 linear problems that come from Expressions (46) and (47) when * is performed with a scale of eleven grades of possibility are reduced to 2 linear programs. Moreover: * = Numerical Application 6. Let us consider again the Irish BMS in numerical application 3, i.e., ~Po (0.038/0.04/0.042) . By using the premium level of each class (see Table 2), we obtain the -cuts of the fuzzy mean asymptotic premium, * , from Equations (45)- (47). These results are in Table 6, which also show the -cuts of its triangular approximate, * , by Equations (48)- (50). From that table, it can be seen that * is practically triangular since the errors by ′( ) in fitting ( ) are negligible. Notice that the triangular approximate * provides a straightforward generalization of both the point estimate by a crisp BMS, 51.423, as well as the modal interval estimate in [11] (p. 849), [51.344, 51.498].

Sensitivity Analysis
In this Section, we evaluate the sensitivity of the errors of the triangular approximates seen in Section 4 with respect to the parameter . The following assumptions are considered:


The core of may be low (0.04), medium (0.5), or high (0.96).  The uncertainty of , which can be measured by its spreads, is symmetrical, i.e., left and right spreads are equal. This uncertainty can take two possible values: 0.002 or 0.015.  We use the Irish BMS in Table 2.
Only results for ( ≥ 1) and are shown. Furthermore, in order to avoid very long calculations, we have performed them on a scale with five grades of possibility. However, it can be verified that for other transition and stationary probabilities, and for a greater scale of grades of possibility, the conclusions to be drawn are practically the same. Table 7 shows that:


The goodness of triangular approximates is always better than acceptable as can be checked in Table 7. In the worst case, for , = 0.5 and risk parameter = 0.025/0.04/0.055, errors are below 5%.  Errors (as defined for Tables 3 and 5) increase with respect to the uncertainty of risk parameter and decrease with respect to the value of its core.  Table 8 shows the mean asymptotic premiums for the cores of considered in Table  7 and the most uncertain scenario (left and right spread equal to 0.015). It can be checked that triangular approximates * always reach a practically perfect match to * , i.e., errors (as defined for Table 6) are very close to 0.

Summary and Further Research
BMSs are often modelled by means of MCs with crisp probabilities. In this paper, it is considered that transition probabilities of Markovian BMSs are not crisp but uncertain. This uncertainty is captured by using a FN, thus giving rise to the concept of FBMSs. FBMSs modeling is based on the concept of FMC by Buckley and Eslami in [12]. As a result, conventional BMSs can be understood as a particular case of our model where transition probabilities are singletons. The model in [11] represents the uncertainty by means of modal intervals. Since its results can be interpreted as the 0-cuts of ours, that model can also be seen as a particular case of our FBMS.
We assume, as it is often done in actuarial literature, that the number of claims in a period is a Poisson RV. Nonetheless, due to uncertainty, its parameter is not a real number but a TFN. So, to implement the model presented in the paper, it is necessary, firstly, to structure available information of the behavior of that RV. From this information, the Poisson parameter can be fitted by means of a TFN. Three alternatives to do so are proposed. Subsequently, by using -cut arithmetic, transition probabilities, the stationary distribution function, and the mean asymptotic premium of the FBMS are obtained by means of their -cuts. The lower and upper bounds of these -cuts can be understood as the result of a sensitivity analysis of the BMS that evaluates two extreme scenarios with possibility . That output can be very useful in actuarial decision-making processes since it provides a set of sensitivity analyses that is structured on the basis of their grade of reliability.
Although the mean number of claims, , is assumed to be a TFN, the outputs from our FBMS do not maintain that shape. However, in the numerical applications developed within the framework of the Irish BMS, we have verified that all the outputs obtained from a triangular are well approximated by a TFN that maintains the support and core of the original FN. This result is quite interesting. On the one hand, other more complex shapes of FNs can produce drawbacks in information modeling, such as problems with calculations in computer implementation. In this regard, we have observed that the number of optimizing problems to be solved in order to obtain transition probabilities, the stationary distribution, and the mean asymptotic premium is reduced drastically. Likewise, TFNs are very attractive from an insurance decision-making perspective since TFNs admit a very intuitive interpretation even without any knowledge of FST. At least, a TFN provides an estimate of the maximum, minimum, and most feasible values of a variable. Therefore, we feel that the triangular approximations introduced in this document would make it easier to use FMCs in the implementation of BMSs by the insurance industry.