Detection and treatment of outliers for multivariate robust loss reserving

Abstract

Traditional techniques for calculating outstanding claim liabilities such as the chain-ladder are notoriously at risk of being distorted by outliers in past claims data. Unfortunately, the literature in robust methods of reserving is scant, with notable exceptions such as Verdonck & Debruyne (2011, Insurance: Mathematics and Economics, 48, 85–98) and Verdonck & Van Wouwe (2011, Insurance: Mathematics and Economics, 49, 188–193). In this paper, we put forward two alternative robust bivariate chain-ladder techniques to extend the approach of Verdonck & Van Wouwe (2011, Insurance: Mathematics and Economics, 49, 188–193). The first technique is based on Adjusted Outlyingness (Hubert & Van der Veeken, 2008. Journal of Chemometrics, 22, 235–246) and explicitly incorporates skewness into the analysis while providing a unique measure of outlyingness for each observation. The second technique is based on bagdistance (Hubert et al., 2016. Statistics: Methodology, 1–23) which is derived from the bagplot; however; it is able to provide a unique measure of outlyingness and a means to adjust outlying observations based on this measure.

Furthermore, we extend our robust bivariate chain-ladder approach to an N-dimensional framework. The implementation of the methods, especially beyond bivariate, is not trivial. This is illustrated on a trivariate data set from Australian general insurers and results under the different outlier detection and treatment mechanisms are compared.

1. Introduction

1.1 Motivation

At any moment, the number, timing, and severity of future claims payments for a general insurer is shrouded in uncertainty. Reserves are set up to ensure, with a given degree of confidence, that necessary claims payments are met as they arise. The reserving problem is one of solvency; however, it is also one of capital efficiency; if an insurer holds reserves over and above what is needed, they forfeit the opportunity to utilize this capital elsewhere, and insurance becomes more expensive than necessary. It is critical that the results from models and techniques applied to the loss-reserving problem are as accurate as possible when tasked to a range of different data sets.

Some of these data sets may include abnormal observations; outliers; or deviations from model assumptions. Full inclusion of these data points in an analysis may prove detrimental to the accuracy of reserve estimates and the resulting inference. This is an issue that needs to be addressed if these models and techniques are going to reflect reality and be used to inform decisions. Robustness refers to the ability of a model or estimation procedure to not be overtly influenced by outliers in the data set under investigation and/or deviations from the underlying assumptions of the model. Understanding how outliers may impact the results of a model and having statistically sound procedures to detect and treat abnormal observations will improve the robustness of reserving techniques and ultimately lead to more informed and reliable decisions.

While some authors have explored the issue of robustness in reserving, the body of literature in this area is relatively scant. Of particular importance for this paper is the robust GLM chain-ladder of Verdonck & Debruyne (2011) and the robust bivariate chain-ladder of Verdonck & Van Wouwe (2011); however, there has been notable work that moves away from the chain-ladder technique. For robust non-chain-ladder methods, please also see for example Chan & Choy (2003), Chan et al. (2008), Pitselis et al. (2015).

Verdonck et al. (2009) provide a two-stage deterministic robust chain-ladder technique which fundamentally relies on the analysis of residuals given after fitting an over-dispersed Poisson (ODP) GLM to the cumulative and then incremental claims data as described in England & Verrall (1999). A boxplot is employed on the Pearson residuals after fitting the ODP GLM at each stage to detect outlying observations. Under their robust technique, development factors are calculated as medians which are known to be much more robust than means. Unfortunately, in its original formulation, results were poor as reserves were still being heavily influenced by outliers.

The approach of Verdonck et al. (2009) was refined by Verdonck & Debruyne (2011). In particular, it was found that the standard threshold value of 1.345 to identify outliers was often too low for triangular loss data. To combat this, an additional stage was added to the methodology whereby the threshold point was taken to be the 75%-quantile of the residuals after an initial fit using the threshold value of 1.345 which led to better results in terms of robustness. Verdonck & Debruyne (2011) also showed that the influence function for reserves with respect to incremental claims is unbounded when assuming a Poisson GLM specification. This provides a formal basis for the non-robustness of the chain-ladder technique and related techniques. A comprehensive study of impact functions of central estimates, mean squared errors, and quantiles of reserves can be found in Avanzi et al. (2023).

1.2 Statement of contributions

The refined robust GLM chain-ladder technique is an integral component of the robust bivariate chain-ladder as developed by Verdonck & Van Wouwe (2011). In particular, the residuals given after fitting the robust GLM chain-ladder are used to generate a bagplot (Rousseeuw et al., 1999) which may be considered as a bivariate boxplot. The second approach employed to detect and treat outliers is the minimum covariance determinant (MCD) (Rousseeuw, 1984) Mahalanobis distance, also applied to the residuals.

Each of these two techniques has shortcomings that this paper aims to address. In particular, they fail to take skewness of the data sufficiently into account, which may lead to the misclassification of data points between regular and outliers.

We introduce two alternative methodologies that offer a consistent and structured approach to the detection, measurement and treatment of outlying observations with statistical backing. Our methods are based on Adjusted Outlyingness (AO) (Hubert & Van der Veeken, 2008) and bagdistance (Hubert et al., 2016). AO provides a unique measure of outlyingness for each observation and explicitly incorporates a robust measure of skewness into the detection process. Bagdistance is derived from the bagplot and provides a measure of outlyingness for each observation. Through calculation of the bagdistance, a greater variety of alternative treatments for outliers becomes available than when simply using the bagplot. These methodologies are applied and compared on real data. Improvements over the techniques introduced in Verdonck & Van Wouwe (2011) are shown to be material and significant.

Additionally, we extend the methodologies beyond the bivariate case. Multivariate reserving techniques primarily benefit from the consideration of the dependence structure between loss triangles. Such consideration should provide a more accurate estimation of reserves; however, the impact of outliers can still be significant in this setting. Furthermore, information about other lines of business may inform the identification and adjustment of outliers. Hence, appropriate detection techniques and adjustment procedures help ensure that the benefits of multivariate reserving analysis are not lost due to a lack of robustness. This is also a focus of this paper, which develops a multivariate robust reserving technique. We extend the bivariate chain-ladder to an N-dimensional framework and illustrate the implementation of all four outlier detection and treatment techniques on a trivariate example.

Implementation of the methodology is not trivial, and R codes are freely available on GitHub for any interested user; see section 5.

2. Outlier Detection Techniques

In section 2.1, we start by reviewing some basic concepts related to robust statistics, which are required to introduce four detection techniques in section 2.2.

2.1 Preliminary: general robustness concepts

2.1.1 Influence functions and breakdown points

When considering robust methodologies, two common concepts are influence functions and breakdown points, which we define here for later reference:

Influence functions consider an estimation technique and calculate the potential `influence’ a single data point could have on that estimator. For robustness, it is desirable to have a bounded influence function such that no single data point can have an unlimited impact on estimators. An example of an estimator with an unbounded influence function is the mean, whereas the median has a bounded influence function.

Breakdown points describe the proportion of the data set required to be contaminated (or outlying, however this may be defined) before an estimator provides inaccurate results. For example, when calculating the mean only one data point has to be contaminated to render the technique invalid, however, when calculating the median, 50% of data points need to be contaminated to invalidate the technique.

2.1.2 Halfspace depth and Tukey median

Both the bagplot (Rousseeuw et al., 1999) and the bagdistance (bd) (Hubert et al., 2016) are based on the concept of halfspace depth (Tukey, 1975).

The halfspace depth (Tukey, 1975) of a point is defined as the minimum number of points (from the data sample) in a closed halfplane through the point of interest. We refer to Fig. 1 which illustrates the concept of halfspace depth in the bivariate case (note that this is scalable to higher dimensions). In Fig. 1, to calculate the halfspace depth of the point marked with a green asterisk, we would consider numerous lines through this point (such as the lines denoted as $L_1$ and $L_2$ in the Figure) and look for the minimum number of points on either side of each of these lines. Romanazzi (2001) showed that halfspace depth has a bounded influence function. This means that the impact an outlier may have on the halfspace depth of a given observation is limited and highlights the robustness of this statistic.

Figure 1 Halfspace depth illustration in 2 dimensions.

The depth median (Tukey median) T ^* is the point with the greatest halfspace depth. It represents a central point of the data and is marked by the red asterisks in Fig. 2. If the point is not unique, the center of gravity of the deepest region is used.

Figure 2 Bagplot with fence drawn.

2.1.3 MCD estimation procedure

A popular technique to estimate location and dispersion parameters in a robust fashion is the MCD (Rousseeuw, 1984) procedure, whereby one finds the

(1) \begin{equation} \left \lfloor \frac{n+p+1}{2}\right \rfloor \leq h\leq n \end{equation}

observations whose classical estimator of the covariance matrix has the smallest determinant where $n$ represents the number of observations and $p$ represents the dimension of the data. The variable $h$ represents the number of observations that are ultimately used, and the procedure puts conditions around the minimum number of observations to avoid selection of a small set that does not truly represent the data. The location vector is then the arithmetic mean of these points, and the scale matrix is taken as a multiple of this covariance matrix. Heuristically, this procedure identifies a `core’ set of data points which are mostly related to one another.

2.1.4 Swamping and masking

Inappropriate outlier detection techniques may lead to swamping or masking. Swamping occurs when one classifies too many observations as outliers, that is, the impact outliers have on the technique being applied results in regular observations incorecctly appearing to be outliers. Masking, on the other hand, refers to the failure to classify outliers as such, that is, outliers have impacted the technique being applied in a way that `masks’ them from being detected as outliers.

2.2 Outlier detection techniques

A robust bivariate chain-ladder technique hinges on the detection and adjustment of bivariate outliers. Verdonck & Van Wouwe (2011) put forward two techniques to be used for this task: the MCD (Rousseeuw, 1984) Mahalanobis Distance, and the bagplot (Rousseeuw et al., 1999). In this section, we outline these approaches and two different techniques; the bagdistance and AO – that can be applied to the multivariate reserving problem and which address some of the shortcomings of the former two.

2.2.1 MCD mahalanobis distance

A standard method used to detect outliers in multivariate analysis is to calculate some measure of distance of each data point from the center of the data. A popular measure is the Mahalanobis Distance (MD), given by,

(2) \begin{equation} \text{MD}(\textbf{x}_i)=\sqrt{(\textbf{x}_i-\widehat{\boldsymbol{\mu }})'\widehat{\boldsymbol{\Sigma }}^{(-1)}(\textbf{x}_i-\widehat{\boldsymbol{\mu }})}, \end{equation}

where $\widehat{\boldsymbol{\mu }}$ represents the sample location vector and $\widehat{\boldsymbol{\Sigma }}$ represents the sample scale matrix. If these are not estimated in a robust manner, then outliers may fail to be detected due to masking and swamping effects. Essentially this means that outlying observations will influence the estimate of the central point and dispersion of the data leading to outliers themselves having small MD (masking) and/or non-outlying observations having high MD (swamping). To combat this, the MCD estimates of location and dispersion are used, and outliers are then flagged as observations that have a MD greater than a certain threshold. Note that $\text{MD}^2\sim \chi _p^2$ if the underlying data is normal.

2.2.2 Bagplot

We refer to Fig. 2 to outline the components of a bagplot.

The bag (B) is shown in Fig. 2 by the darker inner area surrounding T ^* and is constructed as follows. First, denote $D_k$ as the region of all data points that have halfspace depth greater than $k$ and $\#D_k$ as the number of data points contained in this region. The bag is given by the linear interpolation with respect to T $^*$ of the two regions that satisfy $\#D_k\leq \lfloor \frac{n}{2}\rfloor \lt \#D_{k-1}$ . Here, linear interpolation refers to the method of constructing new data points that are in the range of both of $\#D_k$ and $\#D_{k-1}$ .

The fence is then given by multiplying the bag by some factor relative to T $^*$ . Typically this factor is three and Rousseeuw et al. (1999) point out that this value was chosen based on simulations. However, in their recent R package, mrfDepth, Segaert, Hubert, Rousseeuw, Raymaekers, and Vakili (2017) have instead opted to use the square root of the 99th percentile of a chi-squared distribution where the degrees of freedom is equal to the number of dimensions being considered. This new fence factor was chosen after derivation of the distribution of the bagdistance (Hubert et al., 2015); see also the next section. In the two-dimensional case, it is very close to three. We will use this newly developed fence factor throughout the paper.

When displaying a bagplot, the fence is generally not drawn (Rousseeuw et al., 1999). Fig. 2, however, shows the same bagplot with the fence drawn in green, and we can see that the four outlying points marked in red are outside this fence. Note that the choice of this fence factor will directly influence both the number of outliers detected and how they are adjusted. Data points outside the fence are considered outliers and are adjusted to facilitate the application of loss-reserving techniques. The adjustment may be done in a purely graphical manner, that is, bringing outliers back to the fence or loop or a weighting function based on bagdistance may be employed.

The perimeter of the lighter blue area is known as the loop and is given by the convex hull of all non-outlying points. This is analogous to the whiskers in a univariate boxplot (Tukey, 1977).

2.2.3 Bagdistance

Now we present the bagdistance (bd) (Hubert et al., 2016). This statistic provides a measure of outlyingness for each observation and hence does not rely solely on graphical representations. It can also handle skewness often found in loss data. However, similarly to the bagplot, the bd utilizes the bag to capture the shape of the data and subsequently detect outliers. As the bag is formulated based on circa 50% of the data points, there is potential that it does not fully encapsulate the skewness in the set. When this is the case, the ensuant outlier detection results may be flawed. To calculate the bd, firstly define c $_x$ as the intersections of the boundary of the bag (B) and the ray from the Tukey median T $^*$ through the points $\textbf{x}$ (red lines in Fig. 2). The bd is defined as follows

(3) \begin{equation} bd(\textbf{x};\,P_n)= \begin{cases}{0}\ \text{if}\ \textbf{x}=\textbf{T}^*;\\ \frac{||\textbf{x}-\textbf{T}^*||}{||\textbf{c}_x-\textbf{T}^*||}\ \text{elsewhere}, \end{cases} \end{equation}

where $P_n$ represents the distribution of the data set and $||.||$ is the Euclidean norm such that $||\textbf{x}||=\sqrt{x_1^2+\ldots +x_n^2}$ . The denominator scales the distance of the data point (x) to T $^*$ , relative to the dispersion of the bag (in that direction of the projected ray). A cutoff point is then set, such that data points with a bd beyond this threshold are considered outliers and are adjusted back to an appropriate point on the ray emanating from T $^*$ passing through x. Note that if the cutoff value is chosen to be the same as the fence factor under the bagplot technique, then the same outliers will be detected.

2.2.4 Adjusted Outlyingness

AO (Hubert & Van der Veeken, 2008) is based on an adjustment of the Stahel-Donoho estimate of outlyingness (Donoho, 1982) to explicitly incorporate skewness. Furthermore, it is based on robust estimates of location, scale, and skewness such that it achieves a theoretical breakdown point of 25%. Further highlighting its robustness, the Adjusted Outlyingness (AO) technique has a bounded influence function (Hubert & Van der Veeken, 2008). The technique is applied by firstly considering a $p$ -dimensional sample $\textbf{X}_n=(\textbf{x}_1, \ldots,\textbf{x}_n)'$ where $\textbf{x}_i=(x_{i1},\ldots,x_{ip})'$ and $\textbf{a}\in \mathbb{R}^p$ . The measure of AO for $\textbf{x}_i$ is given by

(4) \begin{equation} AO_i= \sup _{\textbf{a}\in \mathbb{R}^{p}}AO^{(1)}(\textbf{a}'\textbf{x}_i,\textbf{X}_n\textbf{a}), \end{equation}

where

(5) \begin{equation} AO^{(1)}(x_i,X_n)=\begin{cases} \frac{x_i-\text{med}(X_n)}{w_2-\text{med}(X_n)},& \text{if } x_i\gt \text{med}(X_n), \\ \\ \frac{\text{med}(X_n)-x_i}{\text{med}(X_n)-w_1}, &\text{if } x_i\lt \text{med}(X_n), \end{cases} \end{equation}

where $w_1$ and $w_2$ are the lower and upper whiskers of the skew-adjusted boxplot (Hubert & Vandervieren, 2008) applied to the data set $X_n$ and $\text{med}(X_n)$ denotes the median of that data set. The skew-adjusted boxplot for all the AO values is then constructed and those that are beyond the cutoff value

(6) \begin{equation} \text{cut-off}=Q_3+1.5e^{3MC}IQR \quad \end{equation}

are declared as outliers, where $Q_3$ represents the third quartile of the data, MC represents the medcouple (a robust measure of skewness; see Brys et al., 2004), and where IQR is the interquartile range (i.e., $Q_3-Q_1$ , where $Q_1$ represents the first quartile of the data). We are only concerned with the upper cutoff value as we are performing the skew-adjusted boxplot technique on the measures of outlyingness, and hence, small results in this context are not of interest. An alternative AO cutoff value of $\sqrt{\chi ^2_{\{99,p\}}}\cdot \text{median}(AO)$ is given in mrfDepth (Segaert et al., 2017), where $p$ represents the dimension of the data, $\chi ^2_{\{99,p\}}$ is the 99th percentile of the chi-squared distribution with $p$ degrees of freedom and AO represents the set of all AO values for the data set. This cutoff value aligns closely with the fence factor used for the bagplot approach. In this paper, we consider both cutoff values.

Under the AO methodology, not all univariate vectors a can be considered; however, Hubert & Van der Veeken (2008) note that taking $m=250p$ directions provides a good balance between `efficiency’ and computation time. From here, if $p=2$ , to visualize the bivariate data a version of the bagplot based on the AO values (rather than halfspace depth) may be constructed (Hubert & Van der Veeken, 2008).

When constructing an AO-based bagplot, the bag is given by the convex hull of the 50% of data points with the smallest AO (note this is different from Fig. 2, which was based on halfspace depth). There are three possible approaches, which differ in the mechanism used to detect and hence treat outliers:

1. 1. A fence is drawn by multiplying the AO-based bag by 3 (or some other factor) relative to the point with the lowest AO. Outliers are then flagged as those observations outside the fence, and the loop is the convex hull of all points within the fence.

   Since the fence is generated from the bag which only considers 50% of the data points, it may fail to fully capture the shape of the data and in particular the skewness in the set.

2. 2. Utilize the alternative cutoff value given in mrfDepth where no fence is drawn.

   Here, the cutoff does not incorporate a robust measure of skewness and instead relies on the median value of the AOs and a quantile of the chi-squared distribution.

3. 3. Outliers are flagged using the traditional cutoff value given by equation (6). In this case, the loop will be generated by the convex hull of all points with AO less than this value and no fence will be generated.

   The traditional AO cutoff value incorporates a robust measure of skewness known as the medcouple which considers the whole data set. Hence, it is more equipped to capture the total skewness, and the loop that is generated by this approach more fully captures the skewness in the data in comparison to the fence methodology.

For the reasons explained above and unless stated otherwise, we will use the third approach.

3. An Extension of the Robust Bivariate Chain-ladder

In this section, we implement the AO and bagdistance (see section 2) outlier detection and treatment techniques in a bivariate reserving setting. We show that the four techniques often lead to different results (although the bagplot and bagdistance will be the same if you use the fence factor as the cutoff distance). The differences are not only in the number of outliers flagged but also in which observations are detected and hence treated. In each case, the different techniques should be implemented and results compared.

This section is motivated by some of the shortcomings of the current robust bivariate chain-ladder methodology. Firstly, the MCD Mahalanobis distance approach assumes elliptical symmetry of the multivariate data. If this assumption is not met, we may fail to detect outlying observations as well as falsely declare regular observations as outliers due to masking and swamping effects. The bagplot is better able to effectively visualize bivariate data, highlighting any correlation, skewness, and tail behavior.

We now present comparisons between these four outlier detection techniques when applied to real non-life insurance data.

3.1 The robust bivariate chain-ladder

The robust bivariate chain-ladder technique was put forward by Verdonck & Van Wouwe (2011), and the general steps involved are as follows:

1. 1. Apply the robust Poisson GLM chain-ladder technique (Verdonck & Debruyne, 2011) to each triangle separately. Obtain residuals from each triangle given by

   (7) \begin{equation} r_{ij}=\frac{X_{ij}-\widehat{\mu }_{ij}}{V^{\frac{1}{2}}(\mu _{i,j})}, \end{equation}
   where $X_{ij}$ represents the incremental claim for accident year $i$ and development year $j$ .

2. 2. Store residuals from each triangle as bivariate observations in a $n\times 2$ matrix $\textbf{X}=(\textbf{x}_1,\ldots,\textbf{x}_n)'$ where $\textbf{x}_i=(r_{kj1},r_{kj2})$ . We have that $n=I(I+1)/2$ .

3. 3. Apply either the bagplot or MCD Mahalanobis Distance (see section 2) to these bivariate residuals.

4. 4. Adjust outliers. For the bagplot, outlying observations are brought back to the fence or loop. For the MCD Mahalanobis distance technique, observations are brought back to the tolerance ellipse representing the 95% quantile of the $\chi _2^2$ distribution.

5. 5. Backtransform adjusted residuals to obtain robust incremental claims $X_{i,j}^{Rob}$ .

6. 6. Apply the multivariate time series chain-ladder (Merz & Wüthrich, 2008) to the robust (adjusted) observations.

We now provide an example of the above framework on real bivariate data. We however include two alternative outlier detection and treatment techniques: the AO and bagdistance.

3.2 Bivariate example using belgian data

The data for this example are from Belgian non-life insurers and are given in Appendix A.1 (from Shi et al., 2012).

3.2.1 Bagplot

We begin by illustrating the major shortcoming of the bagplot approach in that it does not provide a measure of outlyingness that is unique for each observation. This is because the bagplot approach is based on halfspace depth whereby multiple different observations can have the same halfspace depth even if they are at different levels of outlyingness. Fig. 3 shows the bagplot for this data, where 6 outliers have been detected. Four of these observations have the lowest halfspace depth of 1 which may be expected for outliers. However, two of these outliers, $X_{6,5}$ and $X_{7,3}$ , have a halfspace depth of 2 and 3, respectively. These points are marked with red crosses. Additionally, one non-outlying observation, $X_{9,2}$ , has a halfspace depth of 1 and is marked with a purple `+.’ Hence without the corresponding bagplot, little inference can be made about how outlying the outliers are and importantly about whether an observation is outlying or not. This will become an even bigger concern when extending this technique to higher dimensions where graphical representations become less available. Indeed, the bagplot is purely based on a form of ranking (halfspace depth) (see section 2.2.2). As a result, without having the graph available (such as for higher dimensions as considered in the next section) it is difficult to communicate how outlying an observation is.

Figure 3 Bagplot.

3.2.2 MCD mahalanobis distance

When employing the MCD Mahalanobis Distance technique, a graphical representation is also available in the bivariate case. In particular, we can plot each data point as well as tolerance ellipses of which have a squared distance to the central estimate of the data equal to a quantile of the $\chi ^2_2$ distribution. The classical tolerance ellipse is constructed when the Mahalanobis Distance is calculated using the classical estimators of the location vector and scale matrix rather than their robust counterparts from the MCD procedure outlined in section 2.1.3. Outliers are adjusted by a technique known as bivariate Winsorization such that an outlying observation x is adjusted according to

(8) \begin{equation} \text{min}\left (\sqrt{\frac{c}{MD^2(\textbf{x})}},1\right )\cdot \textbf{x}, \end{equation}

where $c$ is equal to the 95% quantile of a $\chi ^2_2$ distribution ( $\chi ^2_{0.95,2}$ ). Under this methodology, 7 observations are detected as outliers. These are the same observations as was flagged under the bagplot approach as well as $X_{9,2}$ , the only other observation with a halfspace depth of 1. Fig. 4 gives these plots before and after outliers have been adjusted.

Figure 4 Tolerance ellipses.

3.2.3 AO

We now explore AO as an outlier detection and adjustment technique. AO explicitly incorporates a robust measure of skewness and provides a measure of outlyingness for each multivariate observation, alleviating the aforementioned issue with the bagplot based on halfspace depth. Additionally, it may be better equipped to handle more extreme levels of skewness. A bagplot based on AO is also available in the two-dimensional case as shown in Figs. 5 and 6, which illustrate all three methods outlined in section 2.2.4.

Figure 5 Adjusted outlyingness bagplot using traditional cutoff value (Approach 3).

Figure 6 Adjusted outlyingness bagplots (Approaches 1–2).

For these plots, the red asterisk is the point with the lowest AO and represents a central point of the data analogous to the Tukey median, the dark blue area represents the bag which contains the 50% of points with the lowest AO, and the light blue area represents the loop whose perimeter is constructed by the convex hull of all points not declared as outliers when using the traditional cutoff method (Hubert & Van der Veeken, 2008).

In Fig. 5, we have used approach 3 with the traditional cutoff value given in Hubert & Van der Veeken (2008) ( $\text{cut-off}=Q_3+1.5e^{3MC}IQR$ ) which incorporates a robust measure of skewness calculated from the whole data set and hence more fully considers the shape of the data to declare outliers, which are shown in red. This cutoff value was calculated to be 5.5813. The light blue area represents the convex hull of all points not declared as outliers under this methodology and may be considered as an AO-based loop.

The fence approach (approach 1) only captures the shape of the data from the 50% of observations determined to be least outlying. Fig. 6(a) shows the AO bagplot with the fence drawn, and we see that under this approach an additional 6 observations are detected as outliers. A further issue that presents itself in this situation is the adjustment to the fence or the loop.

Finally, we may use the alternative cutoff value as given in the mrfDepth package (Segaert et al., 2017) (Approach 2). The cutoff value under this approach is calculated to be 7.6227 compared with 5.5813 under the traditional calculation. An AO bagplot for this approach is shown in Fig. 6(b) where no fence is drawn. In this case, we detect 1 less outlier in comparison to the traditional case. Based on our previous arguments regarding the potential lack of total consideration of skewness under the two alternative approaches, we recommend use of the traditional cutoff value (Approach 3).

3.2.4 Bagdistance

The second outlier detection approach we propose is the bagdistance (bd) (see section 2.2.3) which utilizes the bag to capture the shape of the data and provides a distance measure for each observation to represent its outlyingness. A cutoff point is set such that observations with a bd beyond this threshold are classified as outliers and are adjusted back to an appropriate point on the ray emanating from T $^*$ passing through x. An illustration of the bd is given in Fig. 7 where the bd for the outliers is given by taking the ratio of the orange line to the white line (noting that the orange line continues through to T $^*$ ). If we choose the bd cutoff distance to be uniformly the same as the fence factor for a corresponding bagplot, we will detect the same outliers under the two approaches.

Figure 7 Bagdistance illustration.

3.2.5 Summary of detection results

Table 1 summarizes the observations that were detected as outliers by the four techniques discussed here where a tick indicates that the observation was flagged under the relevant technique and a cross indicates that it was not. All observations not listed here were not flagged by any of the techniques. Further, the results in brackets correspond to the halfspace depth, MCD Mahalanobis distance, AO, and bd values for each observation, respectively.

Table 1. Outlier detection results.

^* As detected by Approach 3 (traditional AO cutoff value).

^** Using cutoff distance of $\sqrt{\chi ^2_{0.99,2}}$ .

3.3 Adjustment of outliers

In the previous section, we outlined how outliers may be detected. Of course, an equally crucial task is to determine how to treat those outliers. In this section, we discuss different methods for such treatment.

In some robust analyses, if a data point is extremely outlying, it is simply removed. Unfortunately, this option is not available in most triangular reserving techniques (which require a fully populated triangle), and hence, we must formulate and choose appropriate adjustment mechanisms.

3.3.1 Bagplot based adjustments

Under the bagplot approach, Verdonck & Van Wouwe (2011) suggest that outlying observations should be brought back to the fence however upon inspecting their plots they appear to have brought observations back to the loop.

Fig. 8(a) and (b) show the bagplot with the fence drawn in green and the bagplot after adjusting residuals back to the fence, respectively. In Fig. 8(b), all of the adjusted outliers are within the loop; however, this is not always the case as after adjustment of the outlying observations, the bagplot methodology is performed again and with some data sets, some of the observations will be outlying because the overall shape of the data has changed which will lead to a different bag and hence fence.

Figure 8 Bagplot before and after adjusting outliers to fence.

Fig. 9 shows the bagplot after adjusting outliers to the loop. Note that the adjustment in this case is more drastic in comparison to the fence adjustment.

Figure 9 Bagplot after adjusting outliers to loop.

3.3.2 MCD mahalanobis distance based adjustments

The adjustment methodology under the MCD Mahalanobis Distance approach is outlined in section 3.2.2.

3.3.3 Bagdistance based adjustments

Through calculation of the bd for each observation, a similar adjustment technique as was used in the MCD Mahalanobis example can be employed in that we adjust an outlying observation x according to

(9) \begin{equation} \text{min}\left (\frac{f}{\textit{bd}},1\right )(\textbf{x}-\textbf{T}^*)+\textbf{T}^*, \end{equation}

where $f$ represents the factor of the bag that we wish to adjust outliers back to. For adjustment back to the bagplot-fence, we choose $f=\sqrt{\chi ^2_{0.99,2}}$ . Note that other adjustment functions may be employed.

For illustration purposes, we have considered two other adjustment functions such that, dependant on the level of outlyingness as measured by bd, observations are adjusted to differing degrees.

In the first case, all observations (x) are adjusted according to

(10) \begin{equation} \textbf{x}^{rob}=\begin{cases} \textbf{x}, & \text{if } \textit{bd}_{\textbf{x}}\leq f; \\[6pt] \text{min}\left (\frac{f}{\textit{bd}_{\textbf{x}}},1\right )(\textbf{x}-\textbf{T}^*)+\textbf{T}^*, & \text{if } f\lt \textit{bd}_{\textbf{x}}\leq \frac{1}{2}\left (2f+\sqrt{4f+1}+1\right );\\[6pt] \frac{f+\sqrt{\textit{bd}_{\textbf{x}}}}{\textit{bd}_{\textbf{x}}}(\textbf{x}-\textbf{T}^*)+\textbf{T}^*, & \text{if } \textit{bd}_{\textbf{x}}\gt \frac{1}{2}\left (2f+\sqrt{4f+1}+1\right ). \end{cases} \end{equation}

This means that moderate outliers will be brought back to the fence, whereas more extreme outliers will be adjusted to points beyond the fence according to their levels of outlyingness. Specifically, the $\frac{1}{2}(2f+\sqrt{4f+1}+1)$ constraint is chosen as it is the value of $\textit{bd}_{\textbf{x}}$ where the adjustment function in case 3 of equations 10 and 11 (i.e., $\frac{f+\sqrt{\textit{bd}_{\textbf{x}}}}{\textit{bd}_{\textbf{x}}}\big)$ is equal to 1. For $\textit{bd}_{\textbf{x}}\gt \frac{1}{2}\left (2f+\sqrt{4f+1}+1\right )$ , $\frac{f}{\textit{bd}_{\textbf{x}}}\lt \frac{f+\sqrt{\textit{bd}_{\textbf{x}}}}{\textit{bd}_{\textbf{x}}}\lt 1$ which means that the points that fall into case 3 will be adjusted back towards the center of the data, however, to a point beyond the fence. For completeness, $\frac{f+\sqrt{\textit{bd}_{\textbf{x}}}}{\textit{bd}_{\textbf{x}}}\gt 1$ for $\textit{bd}_{\textbf{x}}\lt \frac{1}{2}\left (2f+\sqrt{4f+1}+1\right )$ , which means that if the constraint was set at a lower value, then case 3 could lead to outliers being adjusted further away from the central point of the data rather than closer to it. This takes one view point on outliers whereby if an observation is very far outlying then it should still represent moderately rare events in the data. Another approach would be to consider outliers beyond a certain limit as misleading or errors and adjust them further towards the center of the data, limiting their influence. Further investigation will be required based on the data set being analyzed, and the adjustment functions can then be designed depending on the conclusions drawn.

The second approach we consider is given by

(11) \begin{equation} \textbf{x}^{rob}=\begin{cases} \textbf{x}, & \text{if } \textit{bd}_{\textbf{x}}\leq f; \\[4pt] \text{min}\left (\frac{f}{\textit{bd}_{\textbf{x}}},1\right )(\textbf{x}-\textbf{T}^*)+\textbf{T}^*, & \text{if } f\lt \textit{bd}_{\textbf{x}}\leq \frac{1}{2}\left (2f+\sqrt{4f+1}+1\right )\text{ OR } \textit{bd}_{\textbf{x}}\gt u\;{\gt f};\\[5pt] \frac{f+\sqrt{\textit{bd}_{\textbf{x}}}}{\textit{bd}_{\textbf{x}}}(\textbf{x}-\textbf{T}^*)+\textbf{T}^*, & \text{if }\frac{1}{2}\left (2f+\sqrt{4f+1}+1\right )\lt \textit{bd}_{\textbf{x}}\;{\leq u}. \end{cases} \end{equation}

Under this adjustment function, we set an additional limit $u$ such that if outliers are beyond this point they are also brought back to the fence. The rationale for this approach is that if an observation is this far outlying, it should be given full adjustment back to the fence as it either represents a data entry error or the information contained in that observation is too irrelevant for the current task to include more fully. Note that $u$ is an arbitrary selection that will be dependent on the data set under review. The key takeaway is that with unique measures of outlyingness for each observation, there is significantly greater flexibility in how observations can be adjusted.

Fig. 10(a) shows the outliers after adjustment according to equation (10), and Fig. 10(b) shows the outliers after adjustment according to equation (11) where we have set $f=\sqrt{\chi ^2_{0.99,2}}$ and $u=15$ . Under the latter approach with $u=15$ , observation $X_{9,1}$ which has bagdistance of 15.8 as shown in Table 1 is brought back to the fence whereas it remains beyond the fence in Fig. 10(a).

Figure 10 Bagdistance adjusted bagplots.

Remark 3.1. Additionally, we may use the bd values as a residual term when fitting a model with a robust loss function (e.g., robust M-estimation). For example, we could utilize the bd in Huber’s loss function such that we have

(12) \begin{equation} L_c^{\textit{bd}}(\textbf{x})=\begin{cases} \frac{1}{2}\textit{bd}^2, & \text{if } \textit{bd}\leq c; \\ c(|\textit{bd}|-\frac{1}{2}c) & \text{otherwise}. \end{cases} \end{equation}

3.3.4 AO-based adjustments

The AO approach provides a measure of outlyingness for each observation. This measure is based on numerous one-dimensional projections of the individual multivariate observations and the total multivariate sample. For each of these projections, the univariate measure of AO is calculated. The maximum of these univariate AO values after all the projections is taken to be the AO measure in the multivariate case. Without knowing which direction has led to the final AO, we are seemingly unable to backtransform the AO measure to lead to adjusted residuals and ultimately adjusted claim amounts. A solution is to use the AO-based bagplot such that we adjust outliers back to the fence or loop. Consideration must be given to which cutoff value or whether a fence has been drawn to detect outliers under the AO approach. In each case, the loop may be different. We recommend using the traditional cutoff value (Approach 3) as it explicitly incorporates skewness.

3.3.5 Summary of adjustment mechanisms

Tables 2 and 3 summarize the claim values for the outliers detected and their adjusted values for Triangle 1 and 2, respectively. Note that in the case of the AO-Fence and AO-Loop adjustments, we are only considering those observations that were flagged as outliers based on the traditional cutoff.

Table 2. Triangle 1 outlier adjustment results.

Table 3. Triangle 2 outlier adjustment results

For the AO-mrfDepth approach, we are considering outliers flagged under this methodology and outliers are adjusted back to the loop in this case. As we adjust residuals back in the direction of the central estimate of the data, this may lead to an upward or downward adjustment for each residual and similar adjustments will be seen for the corresponding claim observations. Note that `-’ values are provided where outliers were not detected under that methodology.

3.4 Final reserves and discussion

Once the outliers have been appropriately identified and treated, the data can be used to compute reserves estimates. Here, we used the multivariate time series chain-ladder technique as described in Merz & Wüthrich (2008).

Table 4 summarizes the final reserve estimates and their associated RMSE for each individual triangle and total reserves under each outlier detection and adjustment technique as well as when we simply apply the multivariate chain-ladder without adjusting any observations. However, we consider the last three development periods as separate univariate triangles. This is because there are few data points for these development periods and applying the multivariate chain-ladder to such periods often leads to highly volatile results or potentially failure in that elements of the estimated correlation matrices may have absolute values greater than one. This in turn may lead to a lack of convergence. Note that some authors (including Merz & Wüthrich, 2008) suggest extrapolation of these correlation variables from the previous periods however as our main focus is on outlier detection and adjustment we have not pursued this option.

Table 4. Bivariate example reserves.

The robust reserves in this example are always modestly less than the original reserves. This phenomenon is well known in the robust literature and is obviously related to all adjustments being reductions; losses are thus decreased and so are the associated calculated reserves. The greatest adjustment is for the bagdistance technique where a limit is set for the adjustment function.

In each case, we saw a much greater reduction in the RMSE than reserves. For example, the bagplot-loop methodology leads to a 0.91% reduction in reserves and a 12.17% reduction in RMSE. This suggests that even if the change in reserve estimates is only minor, the accuracy of such estimates may be enhanced considerably as a result of these robust techniques.

4. Detection and Treatment of Outliers in N-Dimensional Reserving

4.1 An N-dimensional framework

We now put forward a framework for N-dimensional robust chain-ladder reserving. The methodology is as follows. Firstly, perform the robust Poisson GLM chain-ladder technique on each triangle separately to obtain residuals. These residuals are then stored in an N-dimensional matrix. It is these N-dimensional residuals that outlier detection and adjustment techniques are applied to. The techniques are based on AO, halfspace depth, MCD Mahalanobis distance, and bagdistance.

Under the AO approach outliers are detected using the standard cutoff value (see section 2.2.4). In the three-dimensional case, a graphical representation is also constructed. Firstly, declare the point with the lowest AO as the AO-median. Next, form the convex polyhedron that contains the 50% of points with the lowest AO which can be considered as a N-dimensional AO-based bag. We then construct an N-dimensional AO-based loop by forming the convex polyhedron that contains all non-outlying points.

Outliers are adjusted by bringing them back to the intersection point of the ray from the AO-median to the outlying point and the AO-based loop. This intersection is found using the parametric line clipping algorithm (Cyrus & Beck, 1978; Liang & Barsky, 1984). Adjusted residuals are then backtransformed to give robust incremental claims. The multivariate time series chain-ladder (Merz & Wüthrich, 2008) is applied on these robust claims.

For the halfspace depth approach, we calculate residuals in the same fashion; however, the outlier detection and treatment methodology is altered. The halfspace depth methodology outlined in section 2.1.2 can be applied in N dimensions with Fig. 1 providing an illustration in 2 dimensions. Firstly, calculate the halfspace depth of each observation. The Tukey Median is then the observation with the greatest halfspace depth or the center of gravity of the deepest region if this is not a unique point. Next, formulate the N-dimensional bag. Firstly, let $D_s$ represent the region that contains all points with halfspace depth greater than or equal to $s$ . Let $\#D_s$ be the number of data points in $D_s$ . The bag is constructed by linearly interpolating the two convex regions that satisfy $\#D_s\leq \lfloor \frac{n}{2} \rfloor \lt \#D_{s-1}$ with respect to the Tukey Median. This interpolation is done in the same way as for the bivariate bagplot as described in Miller et al. (2003) however now applicable for any polytype rather than only two-dimensional polygons. Firstly, calculate

(13) \begin{equation} \lambda =\frac{\frac{n}{2}-\#D_s}{\#D_{s-1}-\#D_s}. \end{equation}

The bag is given by $\lambda \cdot \text{outer polytype}+(1-\lambda )\cdot \text{inner polytype}$ . Next, the N-dimensional fence is constructed by multiplying the N-dimensional bag by a fence factor with respect to the Tukey Median. Outliers are declared as points beyond the fence. The fence factor we use is $\sqrt{\chi ^2_N}$ . For all $N\gt 1$ , this is greater than the previously used arbitrary value of 3 in the two-dimensional case, and hence, less outliers will be detected and adjusted. Then, construct the N-dimensional loop by forming the convex hull of all non-outlying points. Outliers are treated by bringing them back to either the fence or loop. In three dimensions, the bag, loop, and fence are convex polyhedrons and graphical representation is available. Adjusted residuals are backtransformed and the multivariate time series chain-ladder applied.

The MCD Mahalanobis distance and bagdistance approaches are implemented in the same fashion as they were for the bivariate case. For the trivariate case, a graphical representation of the tolerance ellipsoid for the MCD Mahalanobis technique may be displayed. We now show the implementation of these four techniques on real data and compare the results.

4.2 Illustration with the AUSI dataset

For the purposes of illustration, we use data from three different lines of business of two Australian insurers. A short discussion of this data is provided in Appendix A.2.

4.2.1 Data adjustments

The fitting of the robust Poisson GLM stage of the methodology is only possible if there are no negative incremental claims in each triangle. However, of the 820 observations for each line, the CTP line from insurer 1, the CTP line from insurer 2, and the Home line from insurer 1 have 1, 21, and 57 negative claims, respectively. We have set all these claims to zero. To ensure that this did not significantly impact results, we calculated development factors under the classical chain-ladder technique before and after the adjustment. The largest change in a development factor was 0.231%, and hence, the effect of the adjustment is small and is deemed acceptable.

4.2.2 Detection of outliers

Firstly, we perform the AO technique on these residuals. 34 outliers were flagged under this approach. We then construct a 3D version of the AO-based bagplot such that the 50% of data points with the lowest AO are contained in a convex polyhedron analogous to the bag, and an additional convex polyhedron is formed for all non-outlying observations analogous to the loop. Fig. 11(a) shows the residuals with the 3D AO bag drawn in yellow, and Fig. 11(b) shows the residuals with the 3D AO-loop in orange and the corresponding bag contained inside. Outliers are shown in red in each case.

Figure 11 AO based detection.

We note that this loop extends in the various directions of the data’s dispersion and seemingly captures its general shape. For completeness, we have also investigated the use of the cutoff value for AO of $\sqrt{\chi ^2_{\{99,N\}}}\cdot \text{median}(AO)$ as suggested in Segaert et al. (2017). Under this approach, 48 outliers are detected.

We now illustrate the halfspace depth-based approach. Fig. 12(a) shows the residuals with the bag drawn in yellow, and Fig. 12(b) shows the residuals with the bag loop and fence all drawn. The loop is again in orange and is within the light blue fence. Outliers are declared as the points outside the fence, and under this methodology, we detect 48 outliers. Interestingly, 29 of the 34 observations flagged as outliers under the traditional AO methodology are also detected under the halfspace depth approach. This highlights some agreement between the two methodologies. Note that the bags under both the AO and halfspace depth approaches are relatively similar; however, the fence which is used to declare outliers under the halfspace depth methodology appears unable to capture the three-dimensional skewness in the data. Rather, it remains relatively elliptical in shape whereas the AO-loop stretches in each direction that the data are further reflected in the number of outliers flagged. We conject that the halfspace depth approach may be misclassifying regular observations as outliers in this instance.

Figure 12 Halfspace depth-based detection.

To highlight the differences between a fence-based approach and AO approach, we may also draw an AO-based fence by multiplying the AO bag by three and declare as outliers observations outside this fence. Fig. 13(a) shows the residuals with the AO bag drawn in yellow and AO-fence in light blue. In Fig. 13(b), we have drawn the AO bag, loop, and fence. In this case, the convex polyhedrons representing the loop and fence are overlapping, and in particular, the loop seemingly extends further in each direction that the data are dispersed. Notably, the AO-fence appears reasonably elliptical in shape, similar to its halfspace depth counterpart.

Figure 13 AO based detection (Approach 1 versus Approach 3).

Under this fence-based approach, 66 observations are flagged as outliers. We believe that using an AO-based fence may again fail to fully capture the shape of the data as we suspect has occurred for the halfspace depth approach. In particular, we believe this will be of the greatest concern when there is significant skewness in the data set.

Now, the MCD Mahalanobis Distance technique which does not consider skewness in the data set at all. Analogous to tolerance ellipses as described in section 3.2.2, tolerance ellipsoids have a squared distance to the central estimate of the data equal to a quantile of the $\chi ^2_3$ distribution. This concept is extendable to higher dimensions where the squared distance is equal to a quantile of the $\chi ^2_p$ distribution where $p$ represents the dimension of the data. Fig. 14(a) shows the residuals with the 97.5% robust tolerance ellipsoid drawn in orange, and Fig. 14(b) shows the non-robust tolerance ellipsoid drawn in green with its robust counterpart encapsulated within it.

Figure 14 MCD mahalanobis distance detection.

This is analogous to what we saw in the bivariate case (see Fig. 6). Under the MCD Mahalanobis distance approach, 120 observations were flagged as outliers. This approach saw a substantially greater number of observations detected as outliers than the other available techniques as it does not consider skewness, likely leading to misclassifications, for example, swamping as the tolerance ellipsoid does not extend towards the skewed data.

4.2.3 Treatment of outliers

Now we turn to the treatment of outliers. Under the AO, halfspace depth, and MCD Mahalanobis distance approaches, outlying observations are brought back to a relevant convex polyhedron (i.e., loop or fence) or tolerance ellipsoid. Under the bagdistance methodology, an alternative approach may be used which weights the degree of adjustment according to the level of outlyingness of each observation. In particular, we have again considered the weighting functions given by equations (10) and (11) and set $f=\sqrt{\chi ^2_3}$ and $u=10$ . As noted in section 3.3.3, $u$ is an arbitrary selection that will be dependent on the data set under review. Adjusted residuals under each methodology are backtransformed to give robust incremental observations.

4.2.4 Final reserves and discussion

The time series multivariate chain-ladder technique (Merz & Wüthrich, 2008) is performed on these robust observations. The results for final reserves and estimated RMSE are given in Table 5. Note that we have scaled the presented results within each LOB to ensure commercial data confidentiality.

Table 5. Trivariate reserves (Scaled values with different radix).

Under each methodology, reserves have been reduced. The greatest reduction was seen when outliers were detected and adjusted using the MCD Mahalanobis distance approach. This is likely because there was a much greater proportion of observations detected as outliers under this approach. However as this methodology does not consider the skewness in the data set, it may be inappropriate in this situation (likely excessive). Rather, the AO technique based on the traditional cutoff value is likely more equipped to capture the shape of the data in this regard and we believe it to be the most suitable approach for this data set. The smallest reduction in reserves is given for the bd approach where there was no upper limit on the weighting function. This is likely because under this approach outliers have been adjusted to a smaller degree than in comparison to the alternative mechanisms. Additionally, we have seen the RMSE of reserves decrease in each case. The MCD Mahalanobis distance approach saw the greatest reduction in this metric of 21.71%, whereas the traditional AO approach (loop) saw a 8.17% reduction. Notably, the shift in the RMSE for each approach was significantly greater than the corresponding change in reserves. This highlights that even for small adjustments in terms of the point estimate of reserves when using these robust techniques, we may see significant improvements in accuracy.

5. R Codes

All relevant R codes used for the bivariate and trivariate illustrations can be found at https://github.com/agi-lab/reserving-robust , along with a synthetic dataset (with features inspired by the AUSI data used in section 4) to help replication of results.

6. Conclusion

We put forward two alternative robust bivariate chain-ladder techniques. The first technique is based on AO and explicitly incorporates skewness into the analysis while providing a unique measure of outlyingness for each observation. The second technique is based on bagdistance which is derived from the bagplot; however, it is able to provide a unique measure of outlyingness and a means to adjust outlying observations based on this measure. Using a real bivariate data set, we illustrated how those techniques compared to existing ones. While the (inevitable) reduction in central estimates is (desirably) modest, the reduction in RMSE is much more significant. This highlights that even for small adjustments in terms of the point estimate of reserves when using these robust techniques, we may see significant improvements in accuracy. We showed that our methodology offers material and significant improvements, compared to the techniques introduced in Verdonck & Van Wouwe (2011).

We then extended our framework to N dimensions and implemented all four outlier detection and treatment techniques on a real trivariate data set. The mathematical generalization of the bivariate methodologies was not trivial. The illustration demonstrated good performance again, and results were qualitatively similar to those discussed in the bivariate example.

Through our expositions and extension to higher dimensions, we have added to the toolbox of techniques available to detect and treat outliers in multivariate reserving. We believe that the new techniques applied in this paper address some of the shortcomings of the previous approaches and should be explored as common practice when implementing robust multivariate reserving techniques in practice.

While we focused on the multivariate time series chain-ladder model of Merz & Wüthrich (2008) for our illustrations, it should be noted that the detection and treatment methodologies developed in this paper can be applied to other reserving techniques, including in most recent machine learning actuarial procedures (see, e.g., Richman, 2018). Notably, our work is concerned with the detection and treatment of outliers rather than finding the best reserving model for the data and we have used the chain-ladder model for illustrative purposes. In the case that one proceeds to model the robustified data, back-testing of results compared to those from unadjusted data may provide useful insights into the effectiveness of the techniques.