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) are supposed to be TFNs *ai* = (*aiL*/*aiC*/*aiU*), *i* = 1, ... , *m*. These coefficients are fitted in two stages. At the first stage, the centres *aic* are adjusted as in a conventional log-Poisson regression for *i* = 1, ... , *m*. At the second stage, the spreads of *ai*, *aiC* − *aiL* and *aiU* − *aiC* and, consequently, *aiL* and *aiU*, *i* = 1, ... , *m*, are fitted by solving a quadratic programming problem that minimizes the fuzziness of the system allowing that estimates on the dependent variable contain its observed values. *m*

Once the parameters *ai* have been estimated, obtaining - *λ* = *e* ∑ *i*=1 *ai xi* is straightforward (see Equation (29)).

Let us remark again that although *ai* is a TFN, *λ* is not. Nevertheless, *λ* can be approximated as a TFN - *λ* = (*<sup>λ</sup>L*/*λC*/*λU*) with Equation (18):

$$
\lambda\_L = \epsilon^{\stackrel{m}{i=1}a\_{i\:L\:X\_i}} \; , \; \lambda\_{\mathbb{C}} = \epsilon^{\stackrel{m}{i=1}a\_{i\:C\:X\_i}} \text{ and } \lambda\_{\!\!II} = \epsilon^{\stackrel{m}{i=1}a\_{i\:\!IX\_i}} \tag{31}
$$


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 *λ* = (*<sup>λ</sup>L*/*λC*/*λU*), with *α*-cuts, ∀*α* ∈ [0, 1], *<sup>λ</sup>*(*α*) = *λ*(*α*), *<sup>λ</sup>*(*α*) = [*<sup>λ</sup>L* + (*<sup>λ</sup>C* − *<sup>λ</sup>L*)*<sup>α</sup>*, *λU* − (*<sup>λ</sup>U* − *<sup>λ</sup>C*)*<sup>α</sup>* ]. FuzzytransitionsprobabilitiescomefromthefuzzifiedversionofEquation(6):

$$h\_{ij}\left(\overline{\lambda}\right) = \widetilde{p}\_{ij} = \sum\_{k=0}^{\infty} \frac{\widetilde{\lambda}^k}{k!} e^{-\overline{\lambda}} t\_{ij}(k) \tag{32}$$

To obtain the *α*-cuts of *pij* by using Equation (15), it is necessary to determine the sign of the first derivative of *hij*(*λ*). Let us show the case of the Irish BMS whose transition matrix is Expression (7), and *pij* is either zero, *p*(*<sup>N</sup>* = <sup>0</sup>), *p*(*<sup>N</sup>* = <sup>1</sup>), *p*(*<sup>N</sup>* ≥ 1) or *p*(*<sup>N</sup>* ≥ <sup>2</sup>). Then:

$$
v(N=0) = e^{-\lambda}\tag{33}$$

$$p(N=1) = \lambda \sigma^{-\lambda} \tag{34}$$

$$p(N \ge 1) = 1 - e^{-\lambda} \tag{35}$$

$$p(N \ge 2) = 1 - (1 + \lambda)e^{-\lambda} \tag{36}$$

So, in Equations (33)–(36), *∂p*(*<sup>N</sup>*=<sup>0</sup>) *∂λ* < 0, *∂p*(*<sup>N</sup>*=<sup>1</sup>) *∂λ* > 0, *∂p*(*<sup>N</sup>*≥<sup>1</sup>) *∂λ* > 0, *∂p*(*<sup>N</sup>*≥<sup>2</sup>) *∂λ* > 0 and therefore: 

$$p(N=0)(\mathfrak{a}) = \left[\mathfrak{e}^{-\overline{\lambda}(a)}, \mathfrak{e}^{-\underline{\lambda}(a)}\right] \tag{37}$$

$$p(N=1)(\alpha) = \left[\underline{\lambda}(a)e^{-\underline{\lambda}(a)}, \overline{\lambda}(a)e^{-\overline{\lambda}(a)}\right] \tag{38}$$

$$p(N \ge 1)(a) = \left[1 - e^{-\underline{\lambda}(a)}, 1 - e^{-\overline{\lambda}(a)}\right] \tag{39}$$

$$p(N \ge 2)(a) = \left[1 - (1 + \underline{\lambda}(a))e^{-\underline{\lambda}(a)}, 1 - \left(1 + \overline{\lambda}(a)\right)e^{-\overline{\lambda}(a)}\right] \tag{40}$$

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:

$$
\overrightarrow{p}(N=0) \approx \left( e^{-\lambda\_{l\bar{l}}} / e^{-\lambda\_{\bar{c}}} / e^{-\lambda\_{\bar{L}}} \right) \tag{41}
$$

$$\tilde{p}(N=1) \approx \left(\lambda\_L \varepsilon^{-\lambda\_L} / \lambda\_\mathbb{C} \varepsilon^{-\lambda\_\mathbb{C}} / \lambda\_{ll} \varepsilon^{-\lambda\_{ll}}\right) \tag{42}$$

$$\hat{p}(N \ge 1) \approx \left(1 - \varepsilon^{-\lambda\_L}/1 - \varepsilon^{-\lambda\_C}/1 - \varepsilon^{-\lambda\_{lI}}\right) \tag{43}$$

$$\hat{p}(N \ge 2) \approx \left(1 - \left(1 + \lambda\_L\right)e^{-\lambda\_L}/1 - (1 + \lambda\_C)e^{-\lambda\_C}/1 - (1 + \lambda\_{\bar{U}})e^{-\lambda\_{\bar{U}}}\right) \tag{44}$$

Numerical Application 2. Example 3 in [11] (p. 846) considers *λ* = [0.038, 0.042] in Equation (6), and obtains the modal interval version of this crisp transition matrix:

$$P = \begin{pmatrix} p(N \ge 1) & p(N = 0) & 0\\ p(N \ge 1) & 0 & p(N = 0) \\ p(N \ge 1) & 0 & p(N = 0) \end{pmatrix}$$

Let us suppose that this interval is the 0-cut of the fuzzy estimate of a triangular - *λ* in a Poisson FBMS (i.e., *λ*(0) = [0.038, 0.042]) and *λ*(1) = 0.04, that is to say, - *λ* = (0.038/0.04/0.042). By considering Expression (32) and using Equations (41)–(44), the fuzzy transition matrix, *P* - , which corresponds to a FMC, is:

$$
\bar{P} = \begin{pmatrix}
(0.037287/0.039211/0.041130) & (0.958870/0.960789/0.962713) & 0 \\
(0.037287/0.039211/0.041130) & 0 & (0.958870/0.960789/0.962713) \\
(0.037287/0.039211/0.041130) & 0 & (0.958870/0.960789/0.962713)
\end{pmatrix}.
$$

From this matrix, elements different from 0 in the associated matrix *<sup>P</sup>*(*α*) = *pij*(*α*) are, from Equation (9):

$$p(N=0)(a) = p\_{12}(a) = p\_{23}(a) = p\_{33}(a) = [0.958870 + 0.001919a, 0.962713 - 0.001924a] \text{ (from)} = [0.958870 + 0.001919a, 0.962713 - 0.001919a] \text{ (from)} = [0.958870 + 0.001919a, 0.962713 - 0.001919a] \text{ (from)} = [0.958870 + 0.001919a, 0.962713 - 0.001919a] \text{ (from)}$$

$$p(N \ge 1)(\mathfrak{a}) = p\_{11}(\mathfrak{a}) = p\_{21}(\mathfrak{a}) = p\_{31}(\mathfrak{a}) = [0.037287 + 0.001924\mathfrak{a}, 0.041190 - 0.001919\mathfrak{a}].$$

Numerical Application 3. In example 4 of [11] (p. 848), it is considered the risk factor *λ* = [0.038, 0.042] for an Irish BMS. Like in our numerical application above, again, this interval is the 0-cut of the triangular fuzzy estimate for - *λ* and *λ*(1) = 0.04, i.e., - *λ* = (0.038/0.04/0.042). So, the triangular approximates by Equations (41)–(44) to the probabilities of the transition matrix in Expression (7) and induced by Equation (32) are:

$$
\tilde{p}\_{11} = \tilde{p}\_{21} = \tilde{p}\_{32} = \tilde{p}\_{43} = \tilde{p}\_{54} = \tilde{p}\_{65} = \tilde{p}(N=0) \approx (0.958870/0.960789/0.962713)
$$

$$
\tilde{p}\_{13} = \tilde{p}\_{24} = \tilde{p}\_{35} = \tilde{p}(N=1) \approx (0.036583/0.038432/0.040273)
$$

$$
\tilde{p}\_{46} = \tilde{p}\_{56} = \tilde{p}\_{66} = \tilde{p}(N \ge 1) \approx (0.037287/0.039211/0.041130)
$$

$$
\tilde{p}\_{16} = \tilde{p}\_{26} = \tilde{p}\_{36} = \tilde{p}(N \ge 2) \approx (0.000704/0.000779/0.000858)
$$

$$
p\_{ij} = 0 \text{ otherwise}
$$

Let us remark that approximates in Equations (41)–(44) produce small errors of the real values by Equations (37)–(40). Table 3 shows that when approximating *p*-(*N* ≥ 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* ≥ <sup>1</sup>), it's triangular approximate, *p*- (*N* ≥ <sup>1</sup>), and errors.

Source: Own elaboration. *err*(*α*) = *p*(*<sup>N</sup>*≥<sup>1</sup>)(*α*) and *err*(*α*) = *p*(*<sup>N</sup>*≥<sup>1</sup>) (*α*)

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 *πj*, *j* = 1, ... ,*r*, 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) sugges<sup>t</sup> that its triangular approximate by using Equation (18), *π*- *j* = *<sup>π</sup>jL*/*<sup>π</sup>jC*/*<sup>π</sup>jU* = *<sup>π</sup>j*(0)/*<sup>π</sup>j*(1) = *<sup>π</sup>j*(1)/*<sup>π</sup>j*(0) provides a satisfactory fitting.

**Figure 2.** Results of Numerical Application 2. Source: Own elaboration.

Numerical Application 4. Now we compute the fuzzy stationary distribution (*π*-1, *π*-2, *π*-3) 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 Appendices A and B, respectively, and the result of (*π*-1, *π*-2, *π*-3) 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.

Numerical Application 5. Let us consider the Irish BMS in numerical application 3. Table 4 shows the supports and cores of the fuzzy stationary state *πj*, *j* = 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 ge<sup>t</sup> these results are included in Appendices A and B, respectively. Figure 3 depicts the graphical representation of *πj*, *j* = 1, . . . , 6.

**Stationary Probabilities** *α* **= 0** *α* **= 1** *π* - 1 0.912318, 0.920394 0.916232 *π* - 2 0.035705, 0.039075 0.037394 *π* - 3 0.037080, 0.040717 0.038921 *π* - 4 0.003519, 0.004186 0.003861 *π* - 5 0.002275, 0.002758 0.002523 *π* - 6 0.000954, 0.001190 0.001069

**Table 4.** Results of the Irish BMS when - *λ* = (0.038/0.04/0.042)—supports and cores.

**Figure 3.** Results of the Irish BMS when - *λ* = (0.038/0.04/0.042)—graphical representation. Source: Own elaboration.

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 *π*-4, *<sup>π</sup>* <sup>4</sup>(*α*), with respect to the respective bounds of *<sup>π</sup>*4(*α*) are always under 1%. Therefore, the above intuition is confirmed in the case of *π*-4. We have also observed that this fact also applies to the other stationary probabilities in the numerical application. So, it can be written:

$$
\hat{\pi}\_1 \approx (0.912318/0.916232/0.920394).
$$

*π*-2 ≈ (0.035705/0.037394/0.039075). *π*-3 ≈ (0.037080/0.038921/0.040717). *π*-4 ≈ (0.003519/0.003861/0.004186). *π*-5 ≈ (0.002275/0.002523/0.002758). *π*-6 ≈ (0.000954/0.001069/0.001190).

**Table 5.** *α*-cuts of *π*-4, it's triangular approximate, *π*- 4, and errors.


Source: Own elaboration. \* *err*(*α*) = *π*4(*α*)−*<sup>π</sup>* 4(*α*) *<sup>π</sup>*4(*α*) and *err*(*α*) = *π*4(*α*)−*<sup>π</sup>* 4(*α*) *<sup>π</sup>*4(*α*) .

It is worth pointing out that the intervals fitted for the stationary state in [11] (p. 848) are the 0-cuts of *πj*, *j* = 1, . . . , 6, in the fuzzy version of the Irish BMS.

Two considerations are worth highlighting:

	- (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 *λC*. 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 *p*-(*n*) *ij* and 6 probabilities *πj*. 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.

Step 4. Obtain the mean asymptotic premium, *b*<sup>∗</sup>.

In order to obtain the asymptotic mean premium, we have to evaluate the fuzzy version of Equation (4), - *b*∗ = *r* ∑ *j*=1 *bjπj*. Bearing in mind Equations (12)–(14), (16) and (17), wefirstconsiderthedomain:

$$Dom(a) = \left\{ \pi\_{\dot{\jmath}} \in \left[ \underline{\pi\_{\dot{\jmath}}(a)}, \overline{\pi\_{\dot{\jmath}}}(a) \right], \sum\_{j=1}^{r} \pi\_{\dot{\jmath}} = 1 \right\} \tag{45}$$

and then the set \$*b*∗ = *r*∑*j*=1 *bj<sup>π</sup>j πj* ∈ *Dom*(*α*)'. The *α*-cuts of *b*<sup>∗</sup>, *b*∗(*α*) = *<sup>b</sup>*∗(*α*), *<sup>b</sup>*∗(*α*), are obtained by solving:

$$\underline{b^\*} (a) = \min \left\{ b^\* = \sum\_{j=1}^r b\_j \pi\_j \, \middle| \, \pi\_j \in \left[ \underline{\pi\_j} (a), \overline{\pi\_j} (a) \right], \sum\_{j=1}^r \pi\_j = 1 \right\} \tag{46}$$

$$\overline{b^\*} (a) = \max \left\{ b^\* = \sum\_{j=1}^r b\_j \pi\_j \, \middle| \, \pi\_j \in \left[ \underline{\pi\_j} (a), \overline{\pi\_j} (a) \right], \sum\_{j=1}^r \pi\_j = 1 \right\} \tag{47}$$

which are linear programming problems and so solvable, e.g., with the simplex algorithm. In fact, the problem to solve in this case is the same as that in [35].



A triangular approximate for *b*<sup>∗</sup>, *b*∗ = *<sup>b</sup>*<sup>∗</sup>*L*/*b*<sup>∗</sup>*C*/*b*<sup>∗</sup>*U*, can be obtained by using the 0-cut and the 1-cut obtained from Equations (46) and (47) or, alternatively, if these results have not been previously calculated, by considering the TFNs *<sup>π</sup>jL*/*<sup>π</sup>jC*/*<sup>π</sup>jU*, *j* = 1, ... , *r*, and solving the following linear problems:

$$b\_L^\* = \min \left\{ b^\* = \sum\_{j=1}^r b\_j \pi\_j \, \middle| \, \pi\_j \in \left[ \pi\_{jL}, \pi\_{jL} \right], \sum\_{j=1}^r \pi\_j = 1 \right\} \tag{48}$$

$$b\_{l\bar{l}}^{\*} = \max \left\{ b^{\*} = \sum\_{j=1}^{r} b\_{j} \pi\_{j} \, \middle| \, \pi\_{\bar{j}} \in \left[ \pi\_{\bar{j}L,r} \pi\_{\bar{j}L} \right], \sum\_{j=1}^{r} \pi\_{\bar{j}} = 1 \right\} \tag{49}$$

In this latter way, 20 linear problems that come from Expressions (46) and (47) when *b*∗ is performed with a scale of eleven grades of possibility are reduced to 2 linear programs. Moreover:

$$b\_{\mathbb{C}}^{\*} = \sum\_{j=1}^{r} b\_{j} \pi\_{j\mathbb{C}} \tag{50}$$

Numerical Application 6. Let us consider again the Irish BMS in numerical application 3, i.e., *N* ∼ 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, - *b*∗ , from Equations (45)–(47). These results are in Table 6, which also show the *α*-cuts of its triangular approximate, - *b*∗ , by Equations (48)–(50). From that table, it can be seen that *b*∗ is practically triangular since the errors by *<sup>b</sup>*(*α*) in fitting *b*(*α*) are negligible. Notice that the triangular approximate - *b*∗ 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].


**Table 6.** *α*-cuts of *b*<sup>∗</sup>, it's triangular approximate *b*∗ and its errors.

Source: Own elaboration. \* *err*(*α*) = *b*<sup>∗</sup>(*α*)−*b*<sup>∗</sup> (*α*) *b*∗(*α*) and *err*(*α*) = *b*∗(*α*)−*b*<sup>∗</sup> (*α*) *b*∗(*α*).
