There is a genuine syntax error in the printed command: random=list(...) must be a separate argument, following a comma rather than a plus. The intended call is
model2 <- gamm(npres ~ s(age, bs="cr") + s(wbc, bs="cr"), random=list(doctor=~age))
With the default Gaussian family, this generalized additive mixed model is
Here across doctors, is an estimated covariance matrix, and the Gaussian errors have variance and are independent of the random effects. The formula ~age includes a random intercept as well as a random slope; it does not force these two effects to be uncorrelated.
The random intercept allows different prescribing baselines, and the random slope allows doctor-specific linear modifications to the population age curve. Conditional on these effects the observations are independent under the model, whereas two patients of the same doctor have additional covariance . This accounts for clustered data and avoids treating doctor variation as independent patient-level noise. The population smooth and doctor-specific predictions are different targets. With only four doctors, estimates of the random-effect covariance can be imprecise; the Gaussian support issue for counts also remains. The new model accounts for doctor-specific baselines and age slopes.
The random-slope linear mixed model is
with independence between the random effects and errors. The estimates are
Incubators are treated as exchangeable representatives of possible growth environments, so a random slope permits their growth rates to differ while estimating a population mean growth rate. The formula 0+hours deliberately removes a random intercept: larvae are assigned only after hatching, so it is reasonable to assume no incubator-specific baseline at time zero. Random assignment supports this common-baseline assumption in expectation; it does not prove that observed initial sizes or other baseline differences are exactly identical.
For a new incubator, no observations are available to estimate its realized random slope, and its mean effect is zero. Under squared-error loss the best prediction integrates over that effect and the residual error:
This is the plug-in conditional expectation given the population parameters, rather than a prediction using a fitted effect from one of the four existing incubators. Ignoring parameter-estimation uncertainty, the prediction variance for a new individual is ; it is not just the uncertainty in the population mean.
The command fits a Gaussian generalized additive mixed model
The population growth curve is a centered penalized regression spline with a cubic regression basis. The incubator term remains a random slope without a random intercept, because ~0+hours removes the intercept. A population linear component is allowed within the smooth's unpenalized part; curvature is penalized rather than being required.
A gamm fit returns a list with gam and lme components. Plot the population smooth and pointwise uncertainty bands with
plot(fly.model2$gam, se=TRUE, residuals=TRUE)
Inspect whether the fitted curve departs materially from an affine function, particularly in regions supported by observations. The effective degrees of freedom near one suggest that the penalization has selected an approximately linear centered smooth. Bands and partial residuals help distinguish convincing curvature from noisy deviations, although pointwise bands are not a simultaneous test of linearity. A formal comparison would need to account for both smoothing selection and the incubator dependence. The smooth plot tests the visual plausibility of a common straight growth curve while retaining incubator-specific slope variation.
Random slope 2026-10-05
A random slope multiplies a predictor by a group-specific latent coefficient. With , a contribution induces covariance for two observations in group . A random slope does not necessarily include a random intercept.