Learning · Actuarial programming

Actuarial science in R and Python: open source

Programming turns a valuation into something that can be reviewed, repeated and audited. Here you will find a verified catalog of open-source actuarial packages in R, Python and Julia, four recipes we ran to confirm that they work, with their actual output, and the practices that model standards require.

20
verified packages
4
recipes run
52,135
IBNR of the RAA triangle
Why programR, Python or JuliaCatalogChain ladderPoisson GLMAnnuitiesVaR and TVaREnvironmentBooks and coursesGood practiceFAQGo deeper
01

Why program in actuarial work

An actuarial valuation is a chain of decisions: which data are cleaned, which assumptions are adopted, which method is applied and how the results are presented. In a spreadsheet that chain is scattered across hundreds of cells and is hard to review. In an R or Python script every step is written down, in order, with comments. A reviewer can read it from top to bottom, run it again and get the same number.

That is reproducibility: with the same data, the same code and the same random seed, the result repeats. It is also the foundation of auditability. An auditor or supervisor who receives the code and the data can verify the IBNR reserve or the employee benefit liability without relying on the memory of whoever calculated it, and changes from one closing to the next are recorded in the version history.

Standards of practice push in that direction. ASOP No. 56 of the U.S. Actuarial Standards Board, on modeling, has been in effect since October 1, 2020 and requires the actuary to understand the model, test it, apply governance controls to mitigate model risk and document it in enough detail for another qualified actuary to evaluate the work. The IAA has a model standard devoted to the subject, ISAP 1A on the governance of models (dated December 1, 2018 in its catalog), and CONAC (Colegio Nacional de Actuarios) promotes standards substantially consistent with international ones (see actuarial standards). A model written as code, with tests and version control, meets those principles almost by construction.

There is also a practical reason: the open ecosystem has already solved much of the repetitive work. The development triangle, Mack's chain ladder, the GLM with an exposure offset and the parametric mortality table are all implemented, tested by thousands of users and documented. The actuary's value lies in choosing well, validating and explaining.

02

R, Python or Julia: how to choose

R

Statistics and reserving

It was born for statistical analysis and has the oldest and most mature actuarial packages: ChainLadder, actuar, lifecontingencies, StMoMo. Its y ~ x formula and the glm() function are the standard in linear model textbooks. Ideal for reserving, experience studies and reports with R Markdown or Quarto.

Python

Integration and automation

A general-purpose language that dominates data engineering and machine learning. The CAS maintains chainladder in Python; statsmodels and scikit-learn cover GLMs and predictive models; lifelib offers complete life models. The right choice when the calculation must connect to databases, APIs or production processes.

Julia

Speed in large models

It compiles to native code and reaches speeds close to C without losing readability. The JuliaActuary community publishes packages for mortality tables, life contingencies, yield curves and experience studies. Useful for seriatim projections and massive stochastic scenario runs.

Rule of thumb: pick the language your team uses and learn to read the other two. The recipes below show the same calculation in R and in Python so you can compare syntax with identical results.
04

Recipe 1: chain ladder on the RAA triangle

The RAA triangle contains cumulative losses from automatic facultative general liability business, accident years 1981 to 1990, published by the Reinsurance Association of America in 1991. It is a classic example in the stochastic reserving literature, used by Thomas Mack and by England and Verrall. It ships with both R's ChainLadder and Python's chainladder, so it is a good way to check that the two implementations agree. The method is explained in detail in technical reserves.

Python 3.13 · chainladder 0.10.1
import chainladder as cl

# RAA triangle shipped with the package (cumulative losses, accident years 1981 to 1990)
raa = cl.load_sample("raa")
print(raa)

# Deterministic chain ladder: volume-weighted development factors
dev = cl.Development(average="volume").fit(raa)
print(dev.ldf_.round(4))

modelo = cl.Chainladder().fit(raa)
print(modelo.ultimate_.round(0))
print("IBNR total:", round(float(modelo.ibnr_.sum()), 0))

# Mack standard error
mack = cl.MackChainladder().fit(raa)
print(mack.summary_.round(0))
print("Error estándar total:", round(float(mack.total_mack_std_err_.values.ravel()[0]), 0))
Output produced by running the code
         12       24       36       48       60       72       84       96       108      120
1981  5012.0   8269.0  10907.0  11805.0  13539.0  16181.0  18009.0  18608.0  18662.0  18834.0
1982   106.0   4285.0   5396.0  10666.0  13782.0  15599.0  15496.0  16169.0  16704.0      NaN
1983  3410.0   8992.0  13873.0  16141.0  18735.0  22214.0  22863.0  23466.0      NaN      NaN
1984  5655.0  11555.0  15766.0  21266.0  23425.0  26083.0  27067.0      NaN      NaN      NaN
1985  1092.0   9565.0  15836.0  22169.0  25955.0  26180.0      NaN      NaN      NaN      NaN
1986  1513.0   6445.0  11702.0  12935.0  15852.0      NaN      NaN      NaN      NaN      NaN
1987   557.0   4020.0  10946.0  12314.0      NaN      NaN      NaN      NaN      NaN      NaN
1988  1351.0   6947.0  13112.0      NaN      NaN      NaN      NaN      NaN      NaN      NaN
1989  3133.0   5395.0      NaN      NaN      NaN      NaN      NaN      NaN      NaN      NaN
1990  2063.0      NaN      NaN      NaN      NaN      NaN      NaN      NaN      NaN      NaN
        12-24   24-36   36-48   48-60   60-72   72-84   84-96  96-108  108-120
(All)  2.9994  1.6235  1.2709  1.1717  1.1134  1.0419  1.0333  1.0169   1.0092
         2261
1981  18834.0
1982  16858.0
1983  24083.0
1984  28703.0
1985  28927.0
1986  19501.0
1987  17749.0
1988  24019.0
1989  16045.0
1990  18402.0
IBNR total: 52135.0
       Latest     IBNR  Ultimate  Mack Std Err
1981  18834.0      NaN   18834.0           NaN
1982  16704.0    154.0   16858.0         143.0
1983  23466.0    617.0   24083.0         592.0
1984  27067.0   1636.0   28703.0         713.0
1985  26180.0   2747.0   28927.0        1452.0
1986  15852.0   3649.0   19501.0        1995.0
1987  12314.0   5435.0   17749.0        2204.0
1988  13112.0  10907.0   24019.0        5354.0
1989   5395.0  10650.0   16045.0        6332.0
1990   2063.0  16339.0   18402.0       24566.0
Error estándar total: 26881.0
R 4.3.3 · ChainLadder 0.2.21
suppressPackageStartupMessages(library(ChainLadder))

# RAA triangle shipped with the package (Reinsurance Association of America, 1991)
data(RAA)
print(RAA)

# Volume-weighted development factors
f <- sapply(1:9, function(k)
  sum(RAA[1:(10 - k), k + 1]) / sum(RAA[1:(10 - k), k]))
print(round(f, 4))

# Mack stochastic chain ladder: ultimate, IBNR and standard error
mack <- MackChainLadder(RAA)   # final sigma extrapolated log-linearly (default option)
print(mack)
Output produced by running the code
      dev
origin    1     2     3     4     5     6     7     8     9    10
  1981 5012  8269 10907 11805 13539 16181 18009 18608 18662 18834
  1982  106  4285  5396 10666 13782 15599 15496 16169 16704    NA
  1983 3410  8992 13873 16141 18735 22214 22863 23466    NA    NA
  1984 5655 11555 15766 21266 23425 26083 27067    NA    NA    NA
  1985 1092  9565 15836 22169 25955 26180    NA    NA    NA    NA
  1986 1513  6445 11702 12935 15852    NA    NA    NA    NA    NA
  1987  557  4020 10946 12314    NA    NA    NA    NA    NA    NA
  1988 1351  6947 13112    NA    NA    NA    NA    NA    NA    NA
  1989 3133  5395    NA    NA    NA    NA    NA    NA    NA    NA
  1990 2063    NA    NA    NA    NA    NA    NA    NA    NA    NA
[1] 2.9994 1.6235 1.2709 1.1717 1.1134 1.0419 1.0333 1.0169 1.0092
MackChainLadder(Triangle = RAA)

     Latest Dev.To.Date Ultimate   IBNR Mack.S.E CV(IBNR)
1981 18,834       1.000   18,834      0        0      NaN
1982 16,704       0.991   16,858    154      143    0.928
1983 23,466       0.974   24,083    617      592    0.959
1984 27,067       0.943   28,703  1,636      713    0.436
1985 26,180       0.905   28,927  2,747    1,452    0.529
1986 15,852       0.813   19,501  3,649    1,995    0.547
1987 12,314       0.694   17,749  5,435    2,204    0.405
1988 13,112       0.546   24,019 10,907    5,354    0.491
1989  5,395       0.336   16,045 10,650    6,332    0.595
1990  2,063       0.112   18,402 16,339   24,566    1.503

              Totals
Latest:   160,987.00
Dev:            0.76
Ultimate: 213,122.23
IBNR:      52,135.23
Mack.S.E   26,880.74
CV(IBNR):       0.52

How to read the result. The nine development factors agree to four decimal places in both languages (2.9994 from the first to the second year, 1.0092 from the ninth to the tenth). Total IBNR is 52,135 and the estimated ultimate is 213,122. The total Mack standard error is 26,881, equivalent to a coefficient of variation of 52%: the uncertainty is huge because 1990 has only one observed value and its IBNR depends on the product of all the factors. With Mack's original approximation for the volatility of the last period (est.sigma = "Mack" in R) the total rises to 26,909.

05

Recipe 2: Poisson frequency GLM with an exposure offset

In pricing, frequency is modeled with a Poisson GLM with a log link. Because each policy was exposed for a different length of time, the log of exposure enters as an offset: a term with its coefficient fixed at 1 that turns the expected count into a rate per policy year. Here we simulate 20,000 policies with a fixed seed and known frequencies, so we can check whether the model recovers the true parameters.

Python 3.13 · statsmodels 0.15.0
import numpy as np
import pandas as pd
import statsmodels.api as sm
import statsmodels.formula.api as smf

rng = np.random.default_rng(2026)          # fixed seed: reproducible results
n = 20_000

# Simulated auto portfolio: exposure in policy years between 0.1 and 1
df = pd.DataFrame({
    "exposicion": rng.uniform(0.1, 1.0, n),
    "edad": rng.choice(["18-25", "26-45", "46-70"], n, p=[0.2, 0.5, 0.3]),
    "zona": rng.choice(["urbana", "rural"], n, p=[0.7, 0.3]),
})

# "True" frequencies per policy year (what the GLM should recover)
base = 0.08
f_edad = df["edad"].map({"18-25": 1.8, "26-45": 1.0, "46-70": 0.85})
f_zona = df["zona"].map({"urbana": 1.25, "rural": 1.0})
lam = base * f_edad * f_zona * df["exposicion"]
df["siniestros"] = rng.poisson(lam)

# Poisson GLM with log link and offset = log(exposure)
modelo = smf.glm(
    "siniestros ~ C(edad, Treatment('26-45')) + C(zona, Treatment('rural'))",
    data=df,
    family=sm.families.Poisson(),
    offset=np.log(df["exposicion"]),
).fit()

print(modelo.summary().tables[1])
relatividades = np.exp(modelo.params).round(3)
print(relatividades)
print("Pearson chi2 / g.l.:", round(modelo.pearson_chi2 / modelo.df_resid, 3))
print("Siniestros observados:", int(df["siniestros"].sum()),
      " Exposición total:", round(df["exposicion"].sum(), 1))
Output produced by running the code
=========================================================================================================
                                            coef    std err          z      P>|z|      [0.025      0.975]
---------------------------------------------------------------------------------------------------------
Intercept                                -2.4959      0.066    -37.606      0.000      -2.626      -2.366
C(edad, Treatment('26-45'))[T.18-25]      0.5022      0.068      7.424      0.000       0.370       0.635
C(edad, Treatment('26-45'))[T.46-70]     -0.1823      0.075     -2.422      0.015      -0.330      -0.035
C(zona, Treatment('rural'))[T.urbana]     0.2343      0.067      3.479      0.001       0.102       0.366
=========================================================================================================
Intercept                                0.082
C(edad, Treatment('26-45'))[T.18-25]     1.652
C(edad, Treatment('26-45'))[T.46-70]     0.833
C(zona, Treatment('rural'))[T.urbana]    1.264
dtype: float64
Pearson chi2 / g.l.: 0.965
Siniestros observados: 1169  Exposición total: 11023.0

Validation. The estimated relativities are 1.652 for drivers aged 18 to 25 (true value 1.8), 0.833 for ages 46 to 70 (true 0.85) and 1.264 for the urban zone (true 1.25); the estimated base frequency is 0.082 against 0.08. All the true values fall within the 95% confidence intervals reported by the model. The ratio of the Pearson statistic to degrees of freedom, 0.965, is close to 1, a sign of no overdispersion, as expected because the data were generated as Poisson. With real data, a value well above 1 suggests using quasi-Poisson or negative binomial. The reference guide for this workflow is CAS Monograph 5, cited in the books section.

Fitted model
\[\log \mathbb{E}[N_i] = \log(e_i) + \beta_0 + \beta_{\text{age}} + \beta_{\text{zone}}, \qquad \text{relativity} = e^{\beta}\]

Ni is the number of claims for policy i, ei its exposure in policy years and the β are the coefficients for each level relative to the base level (ages 26 to 45, rural zone).

06

Recipe 3: life annuity and insurance under Gompertz-Makeham

The Gompertz-Makeham law, μ(x) = A + Bcx, has a closed-form survival function, which makes it possible to calculate annuities and insurance benefits without a table. We use the parameters A = 0.00022, B = 2.7 × 10-6 and c = 1.124 of the Standard Ultimate Survival Model from the Dickson, Hardy and Waters textbook, used in life contingencies texts and exams, with 5% interest. We did the calculation directly in Python and replicated it with lifecontingencies in R from a life table built with the same law. For Mexican tables, see mortality tables.

Python 3.13 · NumPy
import numpy as np

# Gompertz-Makeham law: mu(x) = A + B * c**x
# Parameters of the Standard Ultimate Survival Model (Dickson, Hardy and Waters)
A, B, c = 0.00022, 2.7e-6, 1.124
i = 0.05
v = 1 / (1 + i)

def tpx(x, t):
    """Probability that (x) survives t years (closed form)."""
    return np.exp(-A * t - B * c**x * (c**t - 1) / np.log(c))

def anualidad_anticipada(x, omega=130):
    """ä_x = sum_{k>=0} v^k * kpx (payments of 1 at the start of each year)."""
    k = np.arange(0, omega - x)
    return np.sum(v**k * tpx(x, k))

def seguro_vida_entera(x, omega=130):
    """A_x = sum_{k>=0} v^(k+1) * kpx * q_{x+k} (paid at the end of the year of death)."""
    k = np.arange(0, omega - x)
    q = 1 - tpx(x + k, 1)
    return np.sum(v**(k + 1) * tpx(x, k) * q)

for x in (40, 60, 65):
    a = anualidad_anticipada(x)
    Ax = seguro_vida_entera(x)
    d = i / (1 + i)
    print(f"x={x}: ä_x={a:.4f}  A_x={Ax:.5f}  control 1 - d*ä_x={1 - d*a:.5f}")

# Annual net premium for a 1,000,000 whole life policy at age 40
P = 1_000_000 * seguro_vida_entera(40) / anualidad_anticipada(40)
print(f"Prima neta nivelada, suma asegurada 1,000,000, edad 40: {P:,.2f}")
Output produced by running the code
x=40: ä_x=18.4578  A_x=0.12106  control 1 - d*ä_x=0.12106
x=60: ä_x=14.9041  A_x=0.29028  control 1 - d*ä_x=0.29028
x=65: ä_x=13.5498  A_x=0.35477  control 1 - d*ä_x=0.35477
Prima neta nivelada, suma asegurada 1,000,000, edad 40: 6,558.72
R 4.3.3 · lifecontingencies 1.6.0
suppressPackageStartupMessages(library(lifecontingencies))

# Gompertz-Makeham: mu(x) = A + B c^x (same parameters as the Python recipe)
A <- 0.00022; B <- 2.7e-6; c <- 1.124
edades <- 0:130
tpx0 <- exp(-A * edades - B * (c^edades - 1) / log(c))   # survival from age 0
tabla <- new("lifetable", x = edades, lx = 100000 * tpx0, name = "GM-SUSM")

i <- 0.05
for (x in c(40, 60, 65)) {
  a <- axn(tabla, x = x, i = i)            # life annuity-due
  Ax <- Axn(tabla, x = x, i = i)           # whole life insurance, paid at end of year
  cat(sprintf("x=%d  a..x=%.4f  Ax=%.5f\n", x, a, Ax))
}
P <- 1e6 * Axn(tabla, x = 40, i = i) / axn(tabla, x = 40, i = i)
cat(sprintf("Prima neta nivelada (SA 1,000,000, edad 40): %.2f\n", P))
Output produced by running the code
x=40  a..x=18.4578  Ax=0.12106
x=60  a..x=14.9041  Ax=0.29028
x=65  a..x=13.5498  Ax=0.35477
Prima neta nivelada (SA 1,000,000, edad 40): 6558.72

Validation. Both languages give exactly the same values: ä60 = 14.9041 and A60 = 0.29028, which match those published for this model. The control column uses the identity Ax = 1 - d·äx, with d = i/(1+i), and confirms that the insurance and the annuity are consistent with each other. The level net premium for a whole life policy of 1,000,000 issued at age 40 comes to 6,558.72 per year. Environment note: the current version on CRAN is 1.6.3, which requires R 4.4; here version 1.6.0 was run on R 4.3.3, with the same axn() and Axn() interface.

Gompertz-Makeham survival
\[{}_tp_x = \exp\!\left(-At - \frac{B\,c^x\,(c^t-1)}{\ln c}\right), \qquad \ddot a_x = \sum_{k\ge 0} v^k\,{}_kp_x, \qquad A_x = 1 - d\,\ddot a_x\]

v = 1/(1+i) is the discount factor and d = iv the discount rate. The sum is truncated at age 130, where the survival probability is already negligible.

07

Recipe 4: VaR and TVaR by Monte Carlo simulation

The annual aggregate loss of a portfolio is modeled as a random sum: a Poisson number of claims, each with lognormal severity. It has no closed form, but it is easy to simulate. The 99.5% VaR is the quantile exceeded only once every 200 years, the measure behind the capital requirement in solvency; the TVaR is the average of the losses that exceed that quantile and answers the question of how much is lost, on average, in the bad years.

Python 3.13 · NumPy
import numpy as np

rng = np.random.default_rng(2026)
n_sim = 200_000

# Annual aggregate loss S = X1 + ... + XN
# N ~ Poisson(lambda = 50) claims; X ~ Lognormal(mu = 11, sigma = 1.2) in pesos
lam, mu, sigma = 50, 11.0, 1.2
N = rng.poisson(lam, n_sim)
X = rng.lognormal(mu, sigma, N.sum())
idx = np.repeat(np.arange(n_sim), N)          # which year each claim belongs to
S = np.bincount(idx, weights=X, minlength=n_sim)

def var_tvar(perdidas, alfa):
    var = np.quantile(perdidas, alfa)
    tvar = perdidas[perdidas >= var].mean()
    return var, tvar

media_teorica = lam * np.exp(mu + sigma**2 / 2)
print(f"Media simulada: {S.mean():,.0f}   media teórica: {media_teorica:,.0f}")
for alfa in (0.95, 0.99, 0.995):
    var, tvar = var_tvar(S, alfa)
    print(f"alfa={alfa:.3f}  VaR={var:,.0f}  TVaR={tvar:,.0f}  VaR - media={var - S.mean():,.0f}")

# Monte Carlo error of the 99.5% VaR: repeat with 20 different seeds
vars_ = []
for s in range(20):
    r = np.random.default_rng(s)
    Nn = r.poisson(lam, n_sim)
    Xx = r.lognormal(mu, sigma, Nn.sum())
    Ss = np.bincount(np.repeat(np.arange(n_sim), Nn), weights=Xx, minlength=n_sim)
    vars_.append(np.quantile(Ss, 0.995))
print(f"VaR 99.5% entre semillas: media {np.mean(vars_):,.0f}, desv. est. {np.std(vars_, ddof=1):,.0f}")
Output produced by running the code
Media simulada: 6,147,626   media teórica: 6,150,371
alfa=0.950  VaR=9,346,492  TVaR=10,746,657  VaR - media=3,198,866
alfa=0.990  VaR=11,490,305  TVaR=13,232,977  VaR - media=5,342,679
alfa=0.995  VaR=12,559,969  TVaR=14,520,485  VaR - media=6,412,342
VaR 99.5% entre semillas: media 12,561,324, desv. est. 70,896

Validation. The simulated mean (6,147,626) differs by 0.04% from the theoretical λ·eμ+σ²/2 = 6,150,371. The 99.5% VaR is 12.56 million and the TVaR 14.52 million; the capital above the expected loss is about 6.4 million. Repeating the simulation with 20 different seeds shows that the 99.5% VaR has a Monte Carlo error of about 71,000, that is, 0.6%: with 200,000 simulated years the quantile is stable, but with 10,000 the error would be about four times larger. Reporting that error is part of good model documentation.

Lab

Poisson-lognormal aggregate loss in your browser

The same logic as the recipe, programmed in JavaScript with a seeded pseudorandom generator. Change the seed and watch how much the VaR moves: that is the simulation error. To avoid freezing the browser, the maximum is 50,000 simulated years.

·
Simulated mean
·
Theoretical mean λeμ+σ²/2
·
VaR α
·
TVaR α
·
VaR minus mean

08

Reproducible environment: versions used here

A result is reproducible only if the environment that produced it is known. These are the versions used to run the recipes on October 16, 2026, on Ubuntu 24.04. Record the same in your own projects: in Python with a requirements.txt or pyproject.toml file with pinned versions and a virtual environment; in R with the renv package, which saves a lockfile with the exact version of each dependency.

ComponentVersionUsed in
Python3.13.16Recipes 1 to 4
NumPy / pandas2.5.3 / 3.0.5Recipes 1 to 4
chainladder0.10.1Recipe 1
statsmodels0.15.0Recipe 2
R4.3.3Recipes 1 and 3
ChainLadder0.2.21Recipe 1
lifecontingencies1.6.0Recipe 3
Terminal · minimal installation
# Python: virtual environment with pinned versions
python3 -m venv .venv
source .venv/bin/activate
pip install chainladder==0.10.1 statsmodels==0.15.0 numpy pandas
pip freeze > requirements.txt

# R: packages and lockfile with renv
Rscript -e 'install.packages(c("renv", "ChainLadder", "lifecontingencies"))'
Rscript -e 'renv::init(); renv::snapshot()'
09

Open books and courses

Loss Data Analytics (Open Actuarial Textbooks)

A free interactive book on modeling loss data: frequency, severity, aggregate models, credibility and reserving, with R code. The second edition came out in October 2024 in HTML and PDF. There is a Spanish version translated mainly by the actuarial group at the Universitat de Barcelona and a short course with 34 videos and 52 exercises in R.

CAS Monograph 5: Generalized Linear Models for Insurance Rating

A practical guide by Goldburd, Khare, Tevet and Guller to building a rating plan with GLMs from start to finish: distributions, offsets, interactions, validation and variable selection. Second edition (2025 revision), free to download and adopted in the syllabus of CAS Exam 8.

Modern Financial Modeling

A book by Alec Loudenback and Yun-Tien Lee for actuaries and financial professionals on why and how to program models: abstractions, performance, software engineering practices and stochastic modeling, with examples in Julia. Available in full in HTML and PDF.

Actuarial Data Science tutorials (Swiss Association of Actuaries)

Notebooks with code and public datasets on neural networks, trees, unsupervised learning and survival analysis applied to insurance, prepared by the Swiss association's data science working group. In September 2025 it added a tutorial on survival models.

CAS AI Tools and Resources

A CAS page that brings together its AI Primer (second edition, July 2026), research on language models in claims analysis and the Almost Nowhere podcast. It complements the technical material with AI governance, risk and ethics.

JuliaActuary: Learn

Tutorials and examples from the JuliaActuary ecosystem: mortality tables, policy valuation, yield curves and performance comparisons with other languages.

Suggested path: Loss Data Analytics for the fundamentals with R, Monograph 5 for GLMs in pricing and Modern Financial Modeling to make the leap to model engineering. If you are interested in applied AI, continue with the site's Actuar con IA section.
10

Good practice and model governance

01 · Version control

Git from day one

Store code, not results: the repository holds scripts, documentation and tests; sensitive data stay out (.gitignore file). Every significant change goes in a commit with a clear message, and changes of method are reviewed by someone else before being merged. The history answers who changed what and when, the first question in any audit.

02 · Testing

Cases with a known answer

An actuarial model has natural tests: reproducing Mack's RAA triangle, the identity Ax = 1 - d·äx, a reserve calculated by hand for a single policy or the theoretical mean of a simulation. Automate them with pytest or testthat so they run on every change. The recipes on this page include checks of that kind.

03 · Documentation

So another actuary can evaluate it

Purpose of the model, data sources and their quality, assumptions and their justification, method, known limitations, tests performed and software version. That is what ASOP 56 asks for in its documentation section and what ASOP 23 requires on data quality. Quarto, R Markdown or Jupyter notebooks let you combine code, results and text in a single document.

04 · Model governance

Controls proportionate to risk

Model inventory, owners, independent validation, access and change control, and periodic review of assumptions. ASOP 56 gives examples in its appendix and the IAA's ISAP 1A reflects them internationally. In Mexico, the actuarial function and the integrated risk management function required by the LISF (Insurance and Surety Institutions Law) are the natural owners of that governance.

Unrecorded seeds

A simulation without a fixed seed cannot be reproduced or audited. Record the seed, the generator and the number of scenarios, and report the Monte Carlo error.

Copying code without understanding it

A package or an AI-generated snippet can run without errors and still be misapplied. Validate against a known case before using it in a closing.

11

Frequently asked questions

Do I need to know how to program to be an actuary in Mexico?
It is not a formal requirement of the CNSF (National Insurance and Surety Commission) or CONAC, but in practice almost every technical area uses R, Python or SQL for reserving, pricing and valuations. Knowing how to program lets you review and reproduce calculations that stay hidden in spreadsheet cells, and makes it easier to document the model as standards of practice require.
R or Python to start?
If your immediate goal is reserving, applied statistics and exams covering linear models, R has mature actuarial packages such as ChainLadder and actuar. If you want to work with data teams or automate processes, Python offers the CAS's chainladder, statsmodels and scikit-learn. Both paths are valid; what matters is learning version control and testing from the start.
Can I use these packages in work regulated by the CNSF?
Yes, as long as the responsible actuary understands the method, validates the results and documents assumptions, data and limitations. An open-source package does not replace professional judgment: it is a tool whose choice and use must also be documented, just as it would be with a commercial model.
Why does the Mack standard error of the RAA triangle appear as 26,881 here and as 26,909 in other texts?
The difference lies in how the volatility of the last development period is extrapolated, since it cannot be estimated from the data. ChainLadder in R and chainladder in Python use a log-linear extrapolation by default and give 26,881; with Mack's original approximation (option est.sigma = "Mack" in R) the total is 26,909.
What should I do if a package is no longer maintained?
Check the date of the latest release and the repository activity, as we did in the catalog. If the package is stable and its method does not change, it can remain useful, but pin its version in your environment, keep a copy of the source code and document the risk. If it is critical to a regulated process, consider a replacement plan.
Are these recipes useful for studying for SOA or CAS exams?
They help you understand the methods, because they let you reproduce classic examples such as the RAA triangle or the values of the standard survival model used in life contingencies texts. Multiple-choice exams are not solved with code, but programming the calculations cements intuition and catches conceptual errors.
12

Go deeper

PreviousOpen data for actuariesNextActuarial video library