Chain Ladder in Python

Recently I was trying to implement a simple chain ladder model in Python. I had tried this a few years back when my skills were not as... skilled. Back then I knew there was a chain ladder library for R, but when I tried searching for a similar package in Python, I found nothing notable.

I hadn't approached the problem for a while until a few months back when out of interest I searched again. I found an excellent module developed by a GitHub user by the name of JBogaardt. It seems to be a passion project of his because he is consistently working on and developing it by himself.

https://chainladder-python.readthedocs.io/en/latest/index.html You can view the documentation here

What's great is that it aims to be compatible with the Pandas and Scikit-Learn APIs. This means that if you are used to using either of these modules you'll quickly grasp chainladder. Further, if you are a fan of sklearn pipelines, then you can easily incorporate chainladder into your pipelines, as I will show below.

There is some good documentation available in the docs, but I figured I'd write a bit of an introductory guide here too to demonstrate its power.

The Use Case -- Calculating IBNR for an Insurer

To demonstrate some of the abilities of chainladder we will be calculating the IBNR for an insurer across six lines of business, using a blend of the loss ratio and Cape Cod methods.

(As a note, this post assumes prior knowledge of short-term actuarial reserving. Also, I won't be worrying about the quality of the reserving; the methods are purely for demonstration.)

Loading the data

We will use a dataset made up of claims data from 376 US insurance groups. This dataset is commonly used in research and is built into chainladder. You may recognize it as Schedule P data from the Casualty Actuarial society.

The data is already triangulated for us and in the format required by chainladder. chainladder has a Triangle class that is the required format for its methods. I will explain how the Triangle class works after I load the data.

import chainladder as cl
import pandas as pd
 
data = cl.load_sample('clrd')
data

This is a summary of the triangle data. It gives us a high level snapshot of all the pertinent details.

Valuation tells us the valuation date of the data, in this case December 1997.

The Grain is the granularity of our triangle. This data is loaded by default on an Origin Year Development Year basis (OYDY), but it can be loaded with another grain. We can also easily change grains between monthly, quarterly and annually.

The Shape is the same as the shape we are familiar with from pandas. You will note the shape of the triangle is four dimensional. This might sound strange, since the triangles we are familiar with are two dimensional. But part of the power of the chainladder Triangle data type, is that it allows us to work with multiple triangle at once. The first two items of the shape tuple tell us the number of rows and columns our set has, and the last two tell us the number of origin periods and development periods we have.

If we were to visualise the triangle, it would look something like this:

So each cell of our dataframe is itself a 2-dimensional triangle with 10 origin periods and 10 development periods.

Now, let's look at Index and Columns to complete our understanding of the Triangle data type.

We have a multi-level index. The first level for this dataset is GRNAME which is the name of the insurance group. The second level is LOB; the line of business. For each insurance group, we have claims data for each line of business the group writes. Finally, we have the Columns. This is telling us for each GRNAME and LOB, we have a triangulation of each column. Let's revisit the visualisation, but this time allowing for Index and Columns:

Now, let's zoom into one of the triangles, say the cumulative paid losses (CumPaidLoss) for the Workmen's Compensation line of business (wkcomp) line of business for Allstate Ins Co Grp:

data.loc[('Allstate Ins Co Grp', 'wkcomp'), 'CumPaidLoss']

And there is the triangle we are familiar with! We have 775 * 6 = 4 650 of these in our data, all easily accessible using *pandas* style indexing.

Okay, so that is actually a lot of data, and it is a bit overwhelming. Let's just focus on one insurance group, National American Ins Co. And the only columns we really care about for now are CumPaidLoss and EarnedPremDIR.

naic = data.loc['National American Ins Co', ['CumPaidLoss', 'EarnedPremDIR']]
naic

This is the data we will work with. It has only one index now, LOB, because we selected only one group. And only two columns, the cumulative paid loss and earned premium.

Useful methods of the Triangle class

Triangle's have a number of useful built in methods to make working with them easier. Let's focus on Cumulative Claims Paid for Workmens Compensation for these examples just to simplify some of the indexing.

wkcomp = naic.loc['wkcomp', 'CumPaidLoss']
wkcomp

The Triangle class gives us quick access to many useful methods. We can easily convert from cumulative to incremental, and back:

wkcomp.cum_to_incr()

We can quickly get development factors for each cell of the triangle

wkcomp.link_ratio

And because the API is similar to pandas, we can also plot very easily

wkcomp.T.plot()

wkcomp.link_ratio.plot()

There are numerous other methods and I urge you to explore them in the documentation, or just while playing around with the data in a jupyter workbook.

Using basic chain ladder to calculate IBNR

Now, since we have our triangle, the next step is to apply good old chain ladder. In a standard actuarial workflow, this normally involves a bit of Excel work where we calculate development factors and multiply it out to get our ultimate claims figures. We then subtract the latest developments to get our IBNR portion. It is not terribly difficult if we have standardised workbooks, but there is always some issue that pops up. Not to mention, without really stringent workbook controls, Excel is a very uncontrolled environment prone to human error and operational error. I am sure we have all experienced the darker side of Excel (usually late at night before a deadline).

With chainladder things are a lot easier:

bcl = cl.Chainladder().fit(naic)

And that's pretty much it.

Okay, I am oversimplifying things. But, technically, that one line applied a basic chain ladder to every line of business. We can get an IBNR for our insurance group with the following line:

bcl.ibnr_['CumPaidLoss'].sum()

Let's focus on Workmens Compensation again. We can see the full triangle expectation, as well as back fit the chain ladder to get a full expected triangle:

Full triangle

bcl.full_triangle_.loc['wkcomp', 'CumPaidLoss']

Backfit

bcl.full_expectation_.loc['wkcomp', 'CumPaidLoss']

We can also see development factors:

Cumulative development factors

bcl.cdf_.loc['wkcomp', 'CumPaidLoss']

Incremental development factors

bcl.ldf_.loc['wkcomp', 'CumPaidLoss']

As you can see, chainladder adopts the sklearn API. If you can use sklearn, you can use chainladder.

Creating a simple reserving algorithm

Now, simply fitting a basic chain ladder is not always appropriate. We need to consider the sparsity of the triangle. Where there is too little data to use chain ladder based methods we will use the loss ratio approach. For other cases we will use the Cape Cod approach.

So how do we go about this? Well first, we need an indicator to tell us if the triangle is too sparse to apply chain ladder based methods.

If we convert non zero cells in the triangle to 1, we can compare the sum to the number of elements. If the sum is above a certain threshold we can say that it is too sparse and apply out simple loss ratio based approach.

Indicator for non zero elements

# n_cells is hardcoded for example purposes but should be derived in
# real world workflow
def get_sparsity(triangle, n_cells=55):
    indicator = triangle.unstack() / triangle.unstack()
    density = indicator.sum() / n_cells
    return 1 - density
 
get_sparsity(naic.loc['prodliab', 'CumPaidLoss'])

0.8363636363636364

So what threshold should be chosen? That's quite subjective. So let's just pick 25%, and go with that. Now, if a triangle has more than 25% sparsity, it means that 25% of the cells have no data. In a cumulative triangle, that is a lot (in my subjective opinion).

Now if a triangle is sparse, we want to apply the loss ratio. chainladder does not have a loss ratio method built in (and why would it?). It is not too difficult to create one, but for brevity we will just make a function that returns the IBNR based on a given loss ratio applied for a given number of years.

def loss_ratio_ibnr(triangle, earned_premium, loss_ratio, from_year):
 
    earned_premium = earned_premium.to_frame()
    earned_premium.columns = ['IBNR']
    earned_premium['IBNR'] = loss_ratio * earned_premium
 
    latest_paid = triangle.latest_diagonal.to_frame()
    latest_paid.columns = ['IBNR']
 
    ibnr = earned_premium[['IBNR']] - latest_paid
 
    ibnr.loc[ibnr.index.year < from_year, 'IBNR'] = 0
 
    return ibnr
 
loss_ratio_ibnr(
    naic.loc['prodliab', 'CumPaidLoss'],
    naic.loc['prodliab', 'EarnedPremDIR'].latest_diagonal,
    0.5,
    1995
)

Great, so we have a way to get an IBNR for sparse triangles. Now we can create a simple reserving algorithm that uses loss ratio when the sparsity is above 25%, and Cape Cod otherwise.

def reserving_algorithm(triangle, sparsity_threshold):
 
    lobs = triangle['LOB'].to_numpy()
    ibnr_by_lob = {}
 
    for lob in lobs:
 
        claims = triangle.loc[lob, 'CumPaidLoss']
        premium = triangle.loc[lob, 'EarnedPremDIR'].latest_diagonal
 
        if get_sparsity(claims) >= 0.25:
            ibnr_by_lob[lob] = loss_ratio_ibnr(claims, premium, 1.0, 1995)
            ibnr_by_lob[lob].columns = ['Cape Cod']
        else:
            ibnr_by_lob[lob] = cl.CapeCod(decay=0.75).fit(
                claims,
                sample_weight=premium
            ).ibnr_.to_frame()
            ibnr_by_lob[lob].columns = ['Cape Cod']
 
    ibnr = pd.concat(ibnr_by_lob, axis=1)
    ibnr['Total'] = ibnr.sum(axis=1)
 
    return ibnr.fillna(0)

reserving_algorithm(naic, 0.25)