Model Location and Scale Together (GAMLSS)

Fit both the mean and variance as smooth functions of covariates using a Gaussian location-scale model.

A standard GAM models only the mean of the response. Generalized Additive Models for Location, Scale, and Shape (GAMLSS) go further: every distributional parameter gets its own linear predictor. This is essential when the variance itself depends on covariates. Ignoring that structure produces inefficient estimates and misleading confidence intervals.

The climate dataset records monthly temperatures across weather stations that vary widely in altitude and latitude. Mountainous stations and high-latitude stations tend to be not just colder on average but also more variable (a pattern that a standard model cannot capture).

Fit

Load the dataset and fit a GAMLSS with a GaussianLS() family. The mu formula captures the smooth mean trend. The sigma formula captures how spread changes with geography.

import whittaker as wk

data = wk.load_dataset("climate")

# Fit GAMLSS with separate mu and sigma formulas
model = wk.GAMLSS(
    formulas={
        "mu":    "temperature ~ s(altitude) + s(latitude) + month",
        "sigma": "temperature ~ altitude + latitude",
    },
    family=wk.GaussianLS(),
).fit(data, method="GCV")

Summarize

Inspect the fitted model. The summary reports effective degrees of freedom for each smooth in mu, coefficient estimates for sigma, and overall fit statistics.

model.summary()

'GAMLSS fit summary\n========================================\nFamily: GaussianLS(mu=identity, sigma=log)\nN obs: 600\nGlobal deviance: 4614.5107\nAIC: 4646.2817\nBIC: 4716.1292\nLog-likelihood: -2307.2554\nConverged: True (27 iterations)\n\n--- mu ---\n  EDF total: 12.89\n  Smooth 1: edf = 7.90\n  Smooth 2: edf = 2.99\n\n--- sigma ---\n  EDF total: 3.00\n'

Check model.converged and model.n_iter to confirm the alternating estimation algorithm reached a stable solution.

model.converged, model.n_iter

(True, 27)

Predict across an altitude gradient

Build a prediction grid that holds latitude and month fixed while sweeping altitude. model.predict() returns a GAMLSSPrediction object whose .values dict contains one array per parameter.

import numpy as np

# Build prediction grid across altitude gradient
new_data = {
    "altitude": np.linspace(0, 2000, 50),
    "latitude": np.full(50, 50.0),
    "month":    np.full(50, 7),
}

# Predict and inspect mean temperature values
pred = model.predict(new_data)
pred.values["mu"][:5]

array([17.00599998, 17.5139282 , 18.00578257, 18.44388959, 18.77893675])

Inspect predicted scale

The sigma predictions show how temperature variability grows with altitude, independent of the mean.

pred.values["sigma"][:5]

array([5.03063254, 5.13064434, 5.23264444, 5.33667235, 5.4427684 ])

Interpret

Higher altitude drives mu downward (temperatures are colder) while sigma tends to increase, reflecting greater day-to-day variability at elevation. A standard Gaussian GAM would fix the residual variance at a single global value, masking this pattern entirely. By modeling both parameters simultaneously, GAMLSS produces better-calibrated prediction intervals at every altitude level.