1. INTRODUCTION
This paper measures the risk from the volatility in a triangle of estimated ultimate losses rather than deriving metrics of risk from the volatility in the underlying paid or incurred triangles. Take a hypothetical omniscient actuary using the perfect crystal ball method. The resulting triangle of estimated ultimate losses would be constant along each accident year row. The development factors on the ultimate losses at each age would be 1. However, in the real world, the resulting estimates will change at each evaluation and the ultimate development factors will be different from 1. Intuitively, the degree of variation in the ultimate losses can be used to measure reserve volatility.
The correlation among the agetoage factors could have a significant impact on reserve risk. Although the correlations among the agetoage paid or incurred development factors tend to be positive, correlations among the agetoage ultimate development factors can be negative. In the context of reserve risk, this means that when using paid or incurred loss triangles, reserve risk with correlation included is often higher than reserve risk calculated without correlation terms. On the other hand, when based on ultimate loss triangles, reserve risk with correlation terms in the formula can sometimes be larger, and at other times lower, than when the correlation terms are omitted.
This paper will also highlight and diagnose irregularities that can arise from use of the variancecovariance matrix when it is derived from columns of development factors. Typically, each entry in an empirical variancecovariance matrix is based on taking the variance and covariance of vectors of the same size. Such a matrix will be wellbehaved and have desirable mathematical properties. However, all bets are off when the vectors are of different sizes. In a loss triangle, there are fewer observations in the more mature development age columns. As will be shown, a nominal variancecovariance matrix based on such triangular structure can exhibit surprising pathologies. Another aspect of the “triangle structure problem” is that the later age development factors are more leveraged, yet more impactful. Hence the resulting risk estimates are inherently more unstable. Tentative solutions will be proposed to these problems based on an approach of filling in the triangle of factors and recalibrating the covariance matrix to offset any bias introduced.
1.1. Research Context
There are several ways to estimate reserve risk using historical loss data. Mack (1993) provided the first example, estimating the total runoff reserve risk using paid or incurred loss triangles. Merz and Wuthrich (2008) extended this methodology to evaluate the oneyear reserve risk. Robbin (2012) provided a refinement to the Standard Formula used in Solvency II. Rehman and Klugman (2010) proposed using the triangle of estimated ultimate losses to estimate ultimate reserve risk. They assumed that the development factors, rather than losses, follow a lognormal distribution. Siegenthaler (2019) also used a triangle of estimated ultimate losses to derive both a oneyear and a total runoff formula, but he used different assumptions and obtained different resulting formulas than Rehman and Klugman.
1.2. Objectives
This paper will introduce the FengRobbin method, which makes several specific algorithmic enhancements in estimating reserve risk.^{[1]}
The approach starts with the triangle model promulgated by Rehman and Klugman. In this model, agetoage factors of ultimate losses are assumed to be lognormally distributed. The Rehman and Klugman model will be enhanced and extended in several ways:

Develop a oneyear reserve risk estimate.

Examine the covariance terms and derive a slightly different formula for the total runoff reserve risk.

Propose solutions to address the “triangle structure problem” that will eliminate potentially anomalous behaviors when raw variancecovariance matrices are used.

Propose an estimate of the parameter risk component of reserve risk.
1.3. Outline
The remainder of the paper proceeds as follows. Section 2 starts with a general discussion of reserve risk. Section 3 is devoted to presenting the formulas for oneyear and total runoff reserve risk estimates. Section 4 simplifies the formulas using matrix notation. Section 5 demonstrates how the raw variancecovariance matrix can go awry and yield a negative overall variance. It then presents several alternatives for addressing the “triangle structure problem.” Section 6 proposes an estimate for parameter error. Section 7 provides a comparison of the various formulabased methods as well as insights gained from examining the application of the methods to numerical data.
This paper will include simple examples to demonstrate the methodology described in each section. In the more substantial examples in Section 7 and Appendix C, different sets of ultimate triangles will be analyzed using a variety of reserving methods, and different methods will be applied to these ultimate loss triangles to calculate reserve risk standard errors.
2. Reserve Risk
What is reserve risk? In the context of this paper, it is the potential for adverse and favorable deviations between the eventual ultimate losses and the latest estimates of ultimate. The very mature accident years should have very little deviation, while the newer accident years could have sizeable volatility.
There are two main drivers of reserve risk. First is the volatility in the underlying loss development data. Such volatility of the underlying development process in turn translates into volatility in the difference between the eventual ultimate losses and the latest estimates on the diagonal. Second is the volatility that stems from the method used to project the ultimate estimates. While some methods are mathematically more stable than others, this does not necessarily imply they produce more accurate estimates of ultimate more quickly.
This paper will develop formulae for ultimate and for oneyear reserve risk. The oneyear reserve risk measures the volatility of the estimate of ultimate losses over the next year, whereas the ultimate reserve risk, also referred to as total runoff reserve risk, measures the variability of the estimated ultimate losses until all outstanding contract obligations are fulfilled. In the context of Solvency II, the expected unpaid loss is called the undiscounted best estimate and it is assumed to have no builtin prudential margin.
3. One Year and Ultimate Runoff Standard Errors
The main results will be presented in this section, following a summary of the notation shown in Table 1. The notation is largely consistent with that used by Mack, Merz and Wuthrich, and Siegenthaler,^{[2]} but differs from that used by Rehman and Klugman.^{[3]}
In addition, Table 2 has the notation for various lognormal and normal distribution parameters:
The estimators for these parameters will be denoted by an accent character \((\widehat{\ \ })\) on top of the variable. The formula for the estimators will be given in Sections 3.3 and 3.4. Additionally, let \({\sigma_{i,j}}^{2}\) denote the covariance between \(ln(g_{i})\) and \(ln(g_{j}).\) If \(i = j,\) then \({\sigma_{i,i}}^{2} = {\sigma_{i}}^{2}\) is the variance of \(ln(g_{i}).\) Even though the covariance terms could be negative, we used \({\sigma_{i,j}}^{2}\) instead of the more common \(\sigma_{i,j}\) to simplify the collection of terms within the variancecovariance matrix.
Using the above notation, the main result for the oneyear standard deviation for the entire triangle is as follows.
\[\begin{align} {\text{StDev}}\left( \sum_{y = 1}^{Y}U_{y,D  y + 2} \right) &= \sum_{y = 1}^{Y}U_{y,D  y + 1\ }\ \\ &\quad \times \ \sqrt{\left( e^{2\widehat{\alpha} + 2{\widehat{\theta}}^{2}}  e^{2\widehat{\alpha} + {\widehat{\theta}}^{2}} \right)} \end{align} \tag{3.01}\]
where \(\widehat{\alpha} ≔ \sum_{y = 1}^{Y}r_{y}{\widehat{\mu}}_{D  y + 1}\) and\(\ {\widehat{\theta}}^{2} ≔ \sum_{i = 1}^{Y}{\sum_{j = 1}^{Y}r_{i}r_{j}}{{\widehat{\sigma}}_{D  i + 1,D  j + 1}}^{2}\)
The main result for the total runoff standard deviation for the entire triangle is as follows.
\[\begin{align} {\text{StDev}}\left( \sum_{y = 1}^{Y}U_{y,D} \right) &= \sum_{y = 1}^{Y}U_{y,D  y + 1}\ \\ &\quad \times \ \sqrt{\left( e^{2\widehat{\omega} + 2{\widehat{\lambda}}^{2}}  e^{2\widehat{\omega} + {\widehat{\lambda}}^{2}} \right)} \end{align} \tag{3.02}\]
where \(\widehat{\omega} ≔ \sum_{d = 1}^{D}R_{d}{\widehat{\mu}}_{d}\) and \(\hat{\lambda}^2:=\sum_{i=1}^D \sum_{j=1}^D R_i R_j \hat{\sigma}_{i, j}{ }^2\)
3.1. Development Model Setup
The setup of this paper is largely consistent with that used by Rehman and Klugman.
For a given line of business with \(Y\) accident years and \(D\) development periods, historical estimated ultimate losses are arrayed in the usual triangle format as shown in Table 3.
A frequently used quantity is the sum over the latest diagonal \(\sum_{y = 1}^{Y}U_{y,\ D  y + 1},\) which represents the sum of the best available estimate of ultimate over all accident years at the time of evaluation.^{[4]} Estimates in all future development periods are uncertain.
Throughout this paper, the illustrative triangle of ultimate losses in Table 4 will be used to demonstrate the methodology developed.
3.1.1. Ultimate Loss Development Factors
Consistent with Rehman and Klugman’s, and Siegenthaler’s, definitions, the agetoage ultimate development factor, or UDF, is defined as
\[\begin{align} g_{y,d} ≔ \ \frac{U_{y,\ d + 1}}{U_{y,\ d}}\tag{3.03} \end{align}\]
Similar to the agetoage development factors for paid or incurred losses, a triangle of such ultimate development factors shown in Table 5 could be computed using the observed ultimate estimates, as illustrated using the same structure as Table 3.
Define the random variable \(g_{d}\) as the agetoage ultimate development factor applicable for ultimate losses from age \(d\) to \(d\)+1. This makes each \(g_{y,\ \ d}\) an observed value for \(g_{d}\) for \(y \leq D  d.\) As later sections of this paper would demonstrate, if the distribution of \(g_{d}\) is known, variance of the ultimate losses in future development periods would follow.
The observed agetoage ultimate development factors for the illustrative triangle are shown in Table 6. For example, 1.0502 in accident year 1 age 23 is calculated as 921÷877.
A related concept is the agetoultimate ultimate^{[5]} development factor, \(G_{d},\) defined as the product of all agetoage development from age \(d\) up to maturity age \(D.\) That is,
\[\begin{matrix} G_{d} = \prod_{j = d}^{D}g_{j}\tag{3.04} \end{matrix}\]
The estimation of agetoultimate UDF is computationally similar to the agetoultimate paid or incurred loss development factors and could be used similarly.
3.1.2. Accident Year Weight Factors
To address the different volume of exposure in each accident year, a weight factor \(r_{y}\) is calculated for each accident year \(y\) as a proportion of the latest estimate of ultimate losses:
\[\begin{matrix} r_{y} = \ \frac{U_{y,D  y + 1}}{\sum_{j = 1}^{Y}U_{j,\ D  j + 1\ }}\tag{3.05} \end{matrix}\]
It is easy to check that \(\sum_{y = 1}^{Y}r_{y} = 1.\)
For total runoff risks, a cumulative weight is defined for each age as follows:
\[\begin{matrix} R_{d} ≔ \sum_{j = Y  d + 1}^{Y}{r_{j}\ }\tag{3.06} \end{matrix}\]
For the illustrative triangle, the weight factors are shown in Table 7.
3.2. Model Assumptions
Model assumptions are entirely consistent with those used in Rehman and Klugman (2010). The main assumption is that agetoage UDFs for each development age follow a lognormal distribution for all origin periods. It follows that the logarithm of UDFs follow a normal distribution. No assumption is made regarding the independence of the agetoage ultimate development factors, and no assumption is required regarding the method of determining the ultimate loss estimates at each calendar year. Mathematically,
\[g_{d}\ \sim\ LogNormal\ (\mu_{d},\ {\sigma_{d}}^{2})\]
The mean, variance and covariance of \(g_{d}\) can be initially estimated by taking the logarithm of the triangle of UDFs.^{[6]}
Following the methodology developed in Rehman and Klugman, \(\mu_{d}\) and\(\ {\sigma_{d}}^{2}\) could be estimated using the following sample mean and variance formulae.
\[\begin{matrix} {\widehat{\mu}}_{d} ≔ \ \frac{\sum_{y = 1}^{Y  d}{{\text{ln}(g}_{y,d})\ }}{Y  d}\ \tag{3.07} \end{matrix}\]
\[\begin{matrix} {{\widehat{\sigma}}_{d}}^{2} ≔ \ \frac{\sum_{y = 1}^{Y  d}{{{(\text{ln}(g}_{y,d})  {\widehat{\mu}}_{d})}^{2}\ }}{Y  d  1}\ \tag{3.08} \end{matrix}\]
The observed logarithms of the agetoage ultimate development factors for the illustrative triangle are presented in Table 8.
The estimated unadjusted mean and variance of each log UDF factor is shown in Table 9.
The variancecovariance matrix presented in Table 10 is adjusted according to the fulladjustment procedure described in Section 5.2.3. This guarantees a positive calculated variance.
It follows from the definition of \(g_{d}\) that
\[\begin{matrix} {\widehat{U}}_{y,\ d + 1} = \ U_{y,\ d}\ \times \ g_{d}\ \tag{3.09} \end{matrix}\]
The above equation implicitly assumes that future ultimate losses can be predicted using past ultimate losses. This relationship would fail to hold if there is extraneous information indicating future changes in reserving methodology, and such future changes are not reflected in the past data triangle of \(U_{y,\ d}\)'s. In such situations, computed risk amounts must be manually adjusted. Retroactive adjustments of historical ultimate losses could also be made.
3.2.1. Simplifying Assumptions
This paper follows Rehman and Klugman’s use of a firstdegree Taylor series approximation for the exponential and logarithmic functions. It is assumed that if \(x \approx 0,\) \(e^{x} \approx 1 + x;\) if \(x \approx 1,\) \(\ln(x) \approx x  1.\) For a triangle of ultimate losses, agetoage\(\ g\) factors are approximately 1, thus \(\ln(g)\) is approximately 0. The exponential approximation is used where \(x = \ln(g)\) while the logarithmic approximation is used where \(x = g.\) Consequently, these approximations are appropriate. These approximations are utilized in the proof in Appendix A.
3.3. Standard Error of Ultimate Runoff Reserve Risk
In this section, a slightly modified version of Rehman and Klugman’s formula is presented due to a more detailed analysis of covariance terms. The result is similar to Rehman (2016).
3.3.1. Total Runoff Reserve Risk – Single Accident Year
The ultimate estimate at maturity for a single accident year is \(U_{y,D}.\) It has variance
\[ \begin{align} \text{Var}\left\lbrack U_{y,D} \right\rbrack &= {U_{y,D  y + 1}}^{2} \times \ \text{Var}\left\lbrack \frac{U_{y,D}}{U_{y,D  y + 1}} \right\rbrack \\&= {U_{y,D  y + 1}}^{2}\ \times \text{Var}\ \left\lbrack G_{D  y + 1} \right\rbrack \end{align} \]
Here \(U_{y,D  y + 1}\) is the latest known ultimate estimate for accident year \(y.\) Since \(G_{D  y + 1} = \ \prod_{d = D  y + 1}^{D}g_{d}\) is distributed lognormally, \(\ln{(G_{D  y + 1}}) = \ \sum_{d = D  y + 1}^{D}{ln(}g_{d})\) is distributed normally with associated mean and variance given by
\[\begin{matrix} {\widehat{\omega}}_{y} ≔ E\left\lbrack \sum_{d = D  y + 1}^{D}{ln(}g_{d}) \right\rbrack = \sum_{d = D  y + 1}^{D}{\widehat{\mu}}_{d}\ \tag{3.10} \end{matrix}\]
\[\begin{align} {{\widehat{\lambda}}_{y}}^{2} &≔ Var\left\lbrack \sum_{d = D  y + 1}^{D}{ln(}g_{d}) \right\rbrack \\ &= \ \sum_{i = D  y + 1}^{D}{\sum_{j = D  y + 1}^{D}{{\widehat{\sigma}}_{i,j}}^{2}} \end{align} \tag{3.11}\]
In the above equations, \({\widehat{\mu}}_{d}\) is the mean of the \(\text{ln}(g_{d})\) and \({{\widehat{\sigma}}_{i,j}}^{2}\) is the covariance between \(\text{ln}(g_{i})\) and \(ln(g_{j}).\) If \(i = j,\) \({{\widehat{\sigma}}_{i,i}}^{2} = {{\widehat{\sigma}}_{i}}^{2}\) is the sample variance of \(\text{ln}(g_{i}).\)
Thus, the standard error of total runoff reserve risk for a single accident year is given by
\[\begin{align} \text{StDev}\left( U_{y,D} \right) &= U_{y,D  y + 1}\ \\ &\quad \times \ \sqrt{\left( e^{2{\widehat{\omega}}_{y} + 2{{\widehat{\lambda}}_{y}}^{2}}  e^{2{\widehat{\omega}}_{y} + {{\widehat{\lambda}}_{y}}^{2}} \right)} \tag{3.12} \end{align}\]
Equations 3.10 through 3.12 produce identical results as the RehmanKlugman method for a single accident year.
Using the illustrative triangle, the results for each accident year are given in Table 11.
3.3.2. Total Runoff Reserve Risk – All Accident Years Combined
The ultimate estimate at maturity for the sum of all accident year is \(\sum_{y = 1}^{Y}U_{y,D}.\) This quantity has variance
\[ \begin{align} \operatorname{Var}\left[\sum_{y=1}^Y U_{y, D}\right]&=\left(\sum_{y=1}^Y U_{y, Dy+1}\right)^2 \\ &\quad \times \operatorname{Var}\left[\frac{\sum_{y=1}^Y U_{y, D}}{\sum_{y=1}^Y U_{y, Dy+1}}\right] \end{align} \]
Lemma A.1 in the Appendix demonstrates that
\[\ln\left( \frac{\sum_{y = 1}^{Y}U_{y,\ D}}{\sum_{y = 1}^{Y}U_{y,D  y + 1\ }} \right) \approx \ \sum_{d = 1}^{D}{R_{d}\text{ln}(g_{d})}\]
Since \(g_{d}\) has a lognormal distribution, \(\ln\left( g_{d} \right)\) is normally distributed with mean \({\widehat{\mu}}_{d}\) and variance \({{\widehat{\sigma}}_{d}}^{2}.\) The linear combination of normally distributed random variables with coefficients \(R_{d}\) is a normally distributed random variable with mean and variance given by:
\[\begin{matrix} \widehat{\omega} ≔ E\left\lbrack \ \sum_{d = 1}^{D}R_{d}ln\left( g_{d} \right) \right\rbrack = \sum_{d = 1}^{D}R_{d}{\widehat{\mu}}_{d} \tag{3.13} \end{matrix}\]
\[\begin{matrix} {\widehat{\lambda}}^{2} ≔ Var\left\lbrack \sum_{d = 1}^{D}R_{d}ln\left( g_{d} \right) \right\rbrack = \ \sum_{i = 1}^{D}{\sum_{j = 1}^{D}R_{i}R_{j}}{{\widehat{\sigma}}_{i,j}}^{2} \tag{3.14} \end{matrix}\]
Note that \(\widehat{\omega}\) and \({\widehat{\lambda}}^{2}\) are total runoff estimators for the entire triangle, whereas \({\widehat{\mu}}_{d}\) and \({{\widehat{\sigma}}_{i,j}}^{2}\) are estimators for specific ages. It is also important to distinguish \({\widehat{\omega}}_{y}\) from \(\widehat{\omega}\) and \({\widehat{\lambda}}_{y}\) from \({\widehat{\lambda}}^{2}.\)
In some cases, the formula for \({\widehat{\lambda}}^{2}\) above yields an anomalous negative variance if the raw variancecovariance matrix is used. This issue is further explored in Section 5. Starting from Equation 3.14, the FengRobbin method diverges from the RehmanKlugman method.
Recovering the variance of the lognormal distribution yields
\[\text{Var}\left( \frac{\sum_{y = 1}^{Y}U_{y,\ D}}{\sum_{y = 1}^{Y}U_{y,D  y + 1\ }} \right) = e^{2\widehat{\omega} + 2{\widehat{\lambda}}^{2}}  e^{2\widehat{\omega} + {\widehat{\lambda}}^{2}}\]
The completed formula for the standard error of total runoff reserve risk is
\[\begin{align} \text{StDev}\left( \sum_{y = 1}^{Y}U_{y,D} \right) &= \sum_{y = 1}^{Y}U_{y,D  y + 1}\ \\ &\quad \times \ \sqrt{\left( e^{2\widehat{\omega} + 2{\widehat{\lambda}}^{2}}  e^{2\widehat{\omega} + {\widehat{\lambda}}^{2}} \right)}\tag{3.15} \end{align}\]
Using the illustrative triangle, the total runoff results for the sum of all accident years is given in Table 12.
3.4. Standard Error of OneYear Reserve Risk
This is a new result not covered by either Rehman and Klugman (2010) or Rehman (2016).
3.4.1. Oneyear Reserve Risk – Single Accident Year
The ultimate estimate for a single accident year evaluated 12 months in the future is \(U_{y,D  y + 2}.\) The variance of this quantity is
\[ \begin{align} \text{Var}\left\lbrack U_{y,D  y + 2} \right\rbrack &= {U_{y,D  y