Chapter 98 · Solutions
Chapter 98 · Solutions
Five challenges, each verified in code.
Solutions to the five challenges from Chapter 98, with statsmodels GLMs.
In [1]:
import numpy as np, pandas as pd, matplotlib.pyplot as plt
from scipy import stats
import seaborn as sns # seaborn = high-level statistical plots (heatmaps, regplots, pairplots)
import statsmodels.api as sm
from statsmodels.formula.api import ols
from statsmodels.stats.outliers_influence import variance_inflation_factor
from statsmodels.nonparametric.smoothers_lowess import lowess
from sklearn.linear_model import LinearRegression, Ridge, Lasso, LogisticRegression
from sklearn.preprocessing import StandardScaler, PolynomialFeatures
from sklearn.pipeline import make_pipeline
from sklearn.model_selection import train_test_split, cross_val_score
from sklearn.metrics import r2_score, mean_squared_error
ORG="#ea580c"; DEEP="#c2410c"; LIGHT="#fdba74"; INK="#1a2138"; GRID="#e6e9f2"; GREEN="#059669"; RED="#ef4444"; AMBER="#d97706"; BLUE="#2563eb"; PUR="#9333ea"; GREY="#94a3b8"; SLATE="#475569"
plt.rcParams.update({"figure.facecolor":"white","axes.facecolor":"white","figure.dpi":110,"font.size":11,
"axes.edgecolor":GRID,"axes.grid":True,"grid.color":GRID,"axes.axisbelow":True,"axes.spines.top":False,
"axes.spines.right":False,"axes.titlesize":12,"axes.titleweight":"bold","legend.frameon":False})
sns.set_style("whitegrid")
BASE="https://raw.githubusercontent.com/johnfisher-ai/Statistics-Data-Science-AI-Visual-Book/main/data/"
In [2]:
import warnings; warnings.filterwarnings('ignore')
import statsmodels.formula.api as smf
try:
claims = pd.read_excel('../../data/generalized-linear-models--claims.xlsx', sheet_name='Policies')
except FileNotFoundError:
claims = pd.read_excel(BASE + 'generalized-linear-models--claims.xlsx', sheet_name='Policies')
claims['off'] = np.log(claims['exposure_years'])
CHALLENGE 1
Poisson GLM with offset
fit claims with the exposure offset.
In [3]:
pois = smf.glm('claims ~ driver_age + car_age + C(region)', claims, family=sm.families.Poisson(), offset=claims['off']).fit()
print(pois.summary().tables[1])
=========================================================================================
coef std err z P>|z| [0.025 0.975]
-----------------------------------------------------------------------------------------
Intercept 1.6069 0.114 14.050 0.000 1.383 1.831
C(region)[T.suburban] 0.3995 0.097 4.104 0.000 0.209 0.590
C(region)[T.urban] 0.8611 0.089 9.666 0.000 0.686 1.036
driver_age -0.0240 0.002 -14.500 0.000 -0.027 -0.021
car_age 0.0275 0.005 5.823 0.000 0.018 0.037
=========================================================================================
CHALLENGE 2
Read a rate ratio
driver-age rate ratio per 10 years.
In [4]:
rr = np.exp(pois.params['driver_age']*10)
print(f'Each +10 years of driver age multiplies the expected claim rate by {rr:.2f} (a {100*(1-rr):.0f}% drop).')
Each +10 years of driver age multiplies the expected claim rate by 0.79 (a 21% drop).
CHALLENGE 3
Detect overdispersion
Pearson chi2/df and variance vs mean.
In [5]:
print(f'variance/mean = {claims.claims.var()/claims.claims.mean():.1f}x')
print(f'Pearson chi2/df = {pois.pearson_chi2/pois.df_resid:.2f} -> overdispersed (should be ~1)')
variance/mean = 7.3x Pearson chi2/df = 4.83 -> overdispersed (should be ~1)
CHALLENGE 4
Fix with negative binomial
compare AIC and SEs.
In [6]:
nb = smf.glm('claims ~ driver_age + car_age + C(region)', claims, family=sm.families.NegativeBinomial(alpha=1.2), offset=claims['off']).fit()
print(pd.DataFrame({'SE_poisson':pois.bse.round(3),'SE_negbin':nb.bse.round(3)}))
print(f'AIC Poisson={pois.aic:.0f} NegBin={nb.aic:.0f}')
SE_poisson SE_negbin Intercept 0.114 0.227 C(region)[T.suburban] 0.097 0.176 C(region)[T.urban] 0.089 0.166 driver_age 0.002 0.003 car_age 0.005 0.010 AIC Poisson=2881 NegBin=2068
CHALLENGE 5
Logistic is a GLM too
Binomial+logit matches smf.logit.
In [7]:
try:
loans = pd.read_excel('../../data/logistic-regression--loans.xlsx', sheet_name='Loans')
except FileNotFoundError:
loans = pd.read_excel(BASE + 'logistic-regression--loans.xlsx', sheet_name='Loans')
glm_logit = smf.glm('default ~ credit_score + dti', loans, family=sm.families.Binomial()).fit()
logit = smf.logit('default ~ credit_score + dti', loans).fit(disp=0)
print('GLM (Binomial+logit) vs logit() coefficients match:', np.allclose(glm_logit.params, logit.params))
GLM (Binomial+logit) vs logit() coefficients match: True
Statistics, Data Science and AI: A Visual Handbook · © 2026 John Fisher