When Chain-Ladder Factors Misbehave: Smoothing Reserves with Polars and Whittaker-Henderson
Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.
TD;DR
- We run a full Chain-Ladder reserving calculation for 233 insurers in a handful of
polarsexpressions, cross-checked against R’sChainLadderpackage. - One insurer’s incurred development pattern turns out non-monotonic — the kind of pattern that trips up standard parametric curve fits.
- Whittaker-Henderson smoothing, newly available as
scipy.signal.whittaker_henderson, handles it gracefully and shifts the estimated reserve by about 3%.
A Reserving Problem, Solved Without Excel
Ask a reserving actuary how they run a Chain-Ladder and you’ll usually hear “Excel” or the name of a pricey specialized tool. It turns out a modern dataframe library handles it just as well — in a few lines, for hundreds of companies at once.
We use the CAS loss reserving data, specifically the other liability line of business (LoB): 233 US insurers (“GRNAME”), 10 accident years (1998–2007), paid and incurred losses at every development lag (1-10). We treat 2007 as our reporting year, i.e. we simulate a year-end closing.
After reading in the data, we prepare it a bit
import polars as pl
# Other Liability Data Set, December 2025
df = pl.read_csv(
"https://www.casact.org/sites/default/files/2026-03/othliab_pos_98-07.csv"
)
df_triangle = (
df
.group_by("GRNAME", "AccidentYear", "DevelopmentLag")
.agg(pl.sum("IncurredLosses").alias("Incurred"), pl.sum("CumPaidLoss").alias("Paid"))
.sort("GRNAME", "AccidentYear", "DevelopmentLag")
)
This is still the “full triangle”, i.e. all 10 development lags for all 10 accident years. Note that df_triangle is not in a triangle format but in a long data format, see Tidy Data by Hadley Wickham.
Chain-Ladder in a Few Lines of Polars
The mechanics are the textbook ones: aggregate claims by accident year i and development lag j (starts at 1, not 0), keep only the upper-left (already observed) triangle, and compute CL development factors
f^{CL}_j = \frac{\sum_{i=1998}^{2008-j} C_{i,j}}{\sum_{i=1998}^{2008-j} C_{i,j-1}}
separately for every company. In polars this is just a group_by plus a couple of window (.over(...)) expressions — no manual loops over triangles required:
def df_2_cl_factors(df, group_by=None):
g = [] if group_by is None else group_by
cl_factors = (
df
.group_by(g + ["AccidentYear", "DevelopmentLag"])
.agg(pl.sum("Paid"), pl.sum("Incurred"))
# Filter upper left triangle
.filter(REPORTING_YEAR >= pl.col("AccidentYear") + pl.col("DevelopmentLag") - 1)
.sort(g + ["AccidentYear", "DevelopmentLag"])
.with_columns(
PreviousPaid=pl.col("Paid").shift(1).over(g + ["AccidentYear"]),
PreviousIncurred=pl.col("Incurred").shift(1).over(g + ["AccidentYear"]),
)
# Calculate volume-weighted factors per lag
.group_by(g + ["DevelopmentLag"])
.agg(
number_of_years=pl.len(),
Paid=pl.col("Paid").sum(),
PreviousPaid=pl.col("PreviousPaid").sum(),
Incurred=pl.col("Incurred").sum(),
PreviousIncurred=pl.col("PreviousIncurred").sum(),
)
.with_columns(
f_CL_paid=pl.when(pl.col("PreviousPaid") == 0).then(1).otherwise(pl.col("Paid") / pl.col("PreviousPaid")),
f_CL_inc=pl.when(pl.col("PreviousIncurred") == 0).then(1).otherwise(pl.col("Incurred") / pl.col("PreviousIncurred")),
)
.sort(g + ["DevelopmentLag"])
)
return cl_factors
cl_factors = df_2_cl_factors(df_triangle, group_by=["GRNAME"])
We validated the result against R’s ChainLadder package for one insurer, Grinnell Mut Grp, and the incurred factors matched.
Next, we calculate the Chain-Ladder development factors, again for both paid and incurred, this time separately for each company. Note that we take care to only account for the upper left triangle, which is what’s usually available in practice.
def df_2_cl_factors(df, group_by=None):
g = [] if group_by is None else group_by
cl_factors = (
df
.group_by(g + ["AccidentYear", "DevelopmentLag"])
.agg(pl.sum("Paid"), pl.sum("Incurred"))
# Filter upper left triangle
.filter(REPORTING_YEAR >= pl.col("AccidentYear") + pl.col("DevelopmentLag") - 1)
.sort(g + ["AccidentYear", "DevelopmentLag"])
.with_columns(
PreviousPaid=pl.col("Paid").shift(1).over(g + ["AccidentYear"]),
PreviousIncurred=pl.col("Incurred").shift(1).over(g + ["AccidentYear"]),
)
# Calculate volume-weighted factors per lag
.group_by(g + ["DevelopmentLag"])
.agg(
number_of_years=pl.len(),
Paid=pl.col("Paid").sum(),
PreviousPaid=pl.col("PreviousPaid").sum(),
Incurred=pl.col("Incurred").sum(),
PreviousIncurred=pl.col("PreviousIncurred").sum(),
)
.with_columns(
f_CL_paid=pl.when(pl.col("PreviousPaid") == 0).then(1).otherwise(pl.col("Paid") / pl.col("PreviousPaid")),
f_CL_inc=pl.when(pl.col("PreviousIncurred") == 0).then(1).otherwise(pl.col("Incurred") / pl.col("PreviousIncurred")),
)
.sort(g + ["DevelopmentLag"])
)
return cl_factors
cl_factors = df_2_cl_factors(df_triangle, group_by=["GRNAME"])
cl_factors.filter(pl.col("GRNAME") == "Grinnell Mut Grp")
┌────────────┬────────────┬────────────┬────────┬───┬──────────┬────────────┬───────────┬──────────┐
│ GRNAME ┆ Developmen ┆ number_of_ ┆ Paid ┆ … ┆ Incurred ┆ PreviousIn ┆ f_CL_paid ┆ f_CL_inc │
│ --- ┆ tLag ┆ years ┆ --- ┆ ┆ --- ┆ curred ┆ --- ┆ --- │
│ str ┆ --- ┆ --- ┆ i64 ┆ ┆ i64 ┆ --- ┆ f64 ┆ f64 │
│ ┆ i64 ┆ u32 ┆ ┆ ┆ ┆ i64 ┆ ┆ │
╞════════════╪════════════╪════════════╪════════╪═══╪══════════╪════════════╪═══════════╪══════════╡
│ Grinnell ┆ 1 ┆ 10 ┆ 59563 ┆ … ┆ 191044 ┆ 0 ┆ 1.0 ┆ 1.0 │
│ Mut Grp ┆ ┆ ┆ ┆ ┆ ┆ ┆ ┆ │
│ Grinnell ┆ 2 ┆ 9 ┆ 89554 ┆ … ┆ 172605 ┆ 166116 ┆ 1.727141 ┆ 1.039063 │
│ Mut Grp ┆ ┆ ┆ ┆ ┆ ┆ ┆ ┆ │
│ Grinnell ┆ 3 ┆ 8 ┆ 107323 ┆ … ┆ 151543 ┆ 151052 ┆ 1.38351 ┆ 1.003251 │
│ Mut Grp ┆ ┆ ┆ ┆ ┆ ┆ ┆ ┆ │
│ Grinnell ┆ 4 ┆ 7 ┆ 106828 ┆ … ┆ 129664 ┆ 131101 ┆ 1.135997 ┆ 0.989039 │
│ Mut Grp ┆ ┆ ┆ ┆ ┆ ┆ ┆ ┆ │
│ Grinnell ┆ 5 ┆ 6 ┆ 96179 ┆ … ┆ 105967 ┆ 106851 ┆ 1.088687 ┆ 0.991727 │
│ Mut Grp ┆ ┆ ┆ ┆ ┆ ┆ ┆ ┆ │
│ Grinnell ┆ 6 ┆ 5 ┆ 82541 ┆ … ┆ 86134 ┆ 86800 ┆ 1.05192 ┆ 0.992327 │
│ Mut Grp ┆ ┆ ┆ ┆ ┆ ┆ ┆ ┆ │
│ Grinnell ┆ 7 ┆ 4 ┆ 64703 ┆ … ┆ 66151 ┆ 66159 ┆ 1.019475 ┆ 0.999879 │
│ Mut Grp ┆ ┆ ┆ ┆ ┆ ┆ ┆ ┆ │
│ Grinnell ┆ 8 ┆ 3 ┆ 48924 ┆ … ┆ 49409 ┆ 49457 ┆ 1.008015 ┆ 0.999029 │
│ Mut Grp ┆ ┆ ┆ ┆ ┆ ┆ ┆ ┆ │
│ Grinnell ┆ 9 ┆ 2 ┆ 32840 ┆ … ┆ 33029 ┆ 32969 ┆ 1.006343 ┆ 1.00182 │
│ Mut Grp ┆ ┆ ┆ ┆ ┆ ┆ ┆ ┆ │
│ Grinnell ┆ 10 ┆ 1 ┆ 15785 ┆ … ┆ 15915 ┆ 15908 ┆ 1.002413 ┆ 1.00044 │
│ Mut Grp ┆ ┆ ┆ ┆ ┆ ┆ ┆ ┆ │
└────────────┴────────────┴────────────┴────────┴───┴──────────┴────────────┴───────────┴──────────┘
We validated the result against R’s ChainLadder package for one insurer, Grinnell Mut Grp, and the incurred factors matched, see the linked notebook.

A Pattern That Doesn’t Play Nice
Zooming into the CL factors of individual companies, two stood out: the paid pattern of Virginia Mut Ins Co and the incurred pattern of Grinnell Mut Grp. Neither is monotonic — a red flag, since most commercial reserving tools only offer parametric curve fits that are strictly monotonic (to be fair, parametric curve fits are more for tail factor estimation).

Grinnell’s incurred CL factors rise above 1 early on, dip below 1 (case reserves getting released, maybe subrogation), then climb back above 1 with some zig-zag in later years. A parametric curve simply can’t represent that saddle-shaped pattern.
Enter Whittaker-Henderson
Whittaker-Henderson (WH) smoothing is a non-parametric, penalized smoother — no functional form assumed, just a trade-off between fitting the data and penalizing roughness. That trade-off is controlled by two knobs: the penalty order (we use the standard order=2, i.e. penalizing curvature) and the penalty strength lamb, which we set by eye.
It happens to be a perfect match here: development factors are a discrete-time signal with equal time steps, exactly what WH smoothing was designed for (it dates back to Georg Bohlmann in 1899 — arguably it should be called Bohlmann-Whittaker-Henderson). As of scipy 1.18, it ships out of the box as scipy.signal.whittaker_henderson — full disclosure, I contributed that implementation, so I might be a little biased towards finding excuses to use it 
from scipy.signal import whittaker_henderson f_smooth = whittaker_henderson( signal=f_cl_inc, weights=previous_incurred, lamb=1e4 ).x

One neat property of WH smoothing: it preserves the weighted in-sample sum of the signal, so the smoothed factors reproduce the observed incurred losses exactly on the fitted range. Out-of-sample — i.e. for the not-yet-observed lower-right triangle — the smoothed and raw factors diverge, which is exactly where it matters for reserving.
Does It Change the Reserve?
Yes, a bit:
| accident year | reserve CL incurred | reserve smoothed |
| 1998 | 0.0 | 0.0 |
| 1999 | 7.5 | 20.6 |
| 2000 | 37.2 | 39.0 |
| 2001 | 21.5 | 38.3 |
| 2002 | 23.3 | 13.7 |
| 2003 | -124.9 | -118.9 |
| 2004 | -336.1 | -358.5 |
| 2005 | -522.0 | -525.9 |
| 2006 | -482.1 | -456.4 |
| 2007 | 394.4 | 398.2 |
| TOTAL | -981.1 | -950.0 |
The total reserve shifts by about 3% (+31.1) — modest for the total, though individual accident years move more (some by double digits in percentage terms), since smoothing lets a factor’s neighbors pull it away from its own noisy ratio.
If we compare against the incurred losses after all 10 years of development — the closest proxy we have to the true ultimate loss — both CL and the WH-smoothed CL turn out to overestimate it. The un-smoothed CL is only marginally closer.
| incurred | ultimate CL | ultimate smoothed |
| 189,901 | 194,067 | 194,098 |
So for this one company, smoothing made the development pattern easier to reason about, but it didn’t make the forecast more accurate — the difference (about 31) is small compared to the 4,582 standard error Mack’s method reports for this triangle, so neither method is clearly better here. A good reminder to check against a holdout whenever one is available.
Takeaways
- polars makes Chain-Ladder wrangling compact and fast, even across hundreds of companies at once.
- Volume-weighted CL factors computed in polars match R’s
ChainLadderpackage exactly — always reassuring when switching tools. - Whittaker-Henderson smoothing is a flexible alternative to parametric curve fitting whenever development factors are noisy or non-monotonic, and it’s now built into scipy. A smoother pattern isn’t automatically a more accurate one, though — always check against a holdout when you can.
Natural next steps could be to add a tail factor, estimate reserve uncertainty (à la Mack’s method), or let REML pick lamb automatically instead of choosing it by eye.
The full notebook — all the polars code, the charts, and the R comparison — is on GitHub. This post as well as the notebook was AI reviewed.
Spotted a bug, or have a favorite way to smooth development factors? Let me know in the comments!
R-bloggers.com offers daily e-mail updates about R news and tutorials about learning R and many other topics. Click here if you're looking to post or find an R/data-science job.
Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.