Research

A permutation-based random forest diagnostic for linear mixed models

Independent research · Sole author · Department of Statistics, University of Washington · 2025

A linear mixed model can be misspecified in several ways at once. Standard diagnostics are often aimed at one particular failure or depend heavily on what an analyst thinks to check. RandomForestSpecCheck takes a broader approach. It asks a random forest to predict the fitted model's residuals from the covariates, then uses a permutation test that respects the clustered design to decide whether the predictive signal is stronger than noise.

The result is a single screening test for systematic mean structure left in the residuals. The analyst specifies the variables to examine, but not the functional form of the possible misspecification.

5,400
simulated datasets across 54 parameter configurations, in the first study alone
1–3%
empirical false-positive rate at a nominal 5% level, under correct specification
80–100%
power against large mean-structure departures in favourable designs
R
packaged as a single function over ranger and dplyr

The problem with checking a mixed model

Linear mixed models are used for correlated data such as repeated measurements on a subject, students within classrooms, or patients within clinics. They model that correlation directly. In doing so, they introduce several assumptions about the fixed effects, random effects, and residual covariance structure.

A model can be wrong in any of these places. Some failures bias the coefficient estimates. Others leave the estimates relatively stable but distort their uncertainty.

The usual diagnostic tools each see only part of the problem. Residual plots rely on visual judgment and can hide structure in a dense cloud of points. Information criteria compare only the models an analyst chose to fit. Formal tests usually target one alternative at a time. In each case, the analyst has to make some prior guess about what may be wrong.

If a mixed model is correctly specified in its mean structure, its conditional residuals should not contain predictable information about the covariates. RandomForestSpecCheck uses that fact directly. It gives the residuals and selected covariates to a random forest and asks whether the forest can recover signal that the fitted model should already have explained.

The method

The procedure begins by fitting the mixed model and extracting its conditional residuals. Those residuals are then regressed on the predictors the analyst wants to examine using a random forest.

The forest is scored with its out-of-bag coefficient of determination,

oob = 1 − MSEoob / var(êij)

where each tree is evaluated only on observations excluded from its own bootstrap sample. This matters because an in-sample forest can fit residual noise extremely well. Out-of-bag scoring gives a much more useful measure of whether the residual structure generalizes.

Under correct mean specification,

H0: E[e | X] = 0   versus   H1: E[e | X] ≠ 0

The analyst supplies the variables in X. The random forest is left to search for the functional form.

The sampling distribution of the statistic is not available in closed form, so the null is generated by permutation. Residuals are permuted within clusters rather than across the full dataset. This breaks the relationship between the residuals and predictors while preserving the clustered structure of the data.

The forest is refit B times to produce the permutation distribution.

Large observed values indicate residual structure. The procedure flags a model only when two conditions are met,

oob > Q0.975   and   R²oob > τ

with τ = 0.05 by default.

The percentile criterion asks whether the observed value is unusual under the permutation null. The practical threshold asks whether the remaining structure is large enough to matter. This second condition is important in large samples, where a very small effect can be statistically detectable without being substantively useful.

The default value of τ is a judgment rather than a theoretical result. It is exposed as an argument for that reason.

The diagnostic, running

The panel below is a working browser implementation of the procedure. It simulates clustered longitudinal data with a quadratic time trend, fits a random-intercepts model that omits the quadratic term, grows a forest on the residuals, constructs the permutation null within clusters, and applies the decision rule.

Nothing is precomputed.

Figure 1 · Computed live

RandomForestSpecCheck applied to an omitted quadratic effect

permutation null distribution, B = 60 Q0.975 practical threshold τ = 0.05 observed R²oob
Running.

Data are generated as yij = 2 + 0.5 tij + γ tij² + bi + εij. The fitted random-intercepts model omits the quadratic term.

The browser version uses 50 trees and B = 60 permutations so it can run interactively. The full procedure used 500 permutations.

The permutation null falls below zero, as it should. Once the residuals are shuffled, there is no relationship left for the forest to predict. A forest scored out of bag can then perform worse than the sample mean. A null centered near zero would be more concerning than a negative one. In checks against the paper at twenty replicates per cell, the reduced browser version reproduces the reported behavior.

What the method detects, and what it does not

Four simulation studies were run. Each generated data under a known model violation and then fit a model that omitted it.

Three studies concern the conditional mean.

  • an omitted quadratic fixed effect
  • an omitted random slope
  • unmodeled AR(1) correlation

The fourth concerns the variance. Residual variance increases with time while the fitted model assumes constant variance.

Each cell was replicated one hundred times.

For the three mean-structure violations, the diagnostic behaves as intended. False positives stay between one and three percent under correct specification. Power rises with the size of the omitted effect and with a stronger signal-to-noise ratio.

Large quadratic departures are detected in roughly eighty to one hundred percent of replicates. Weak departures are often missed. Omitted random slopes show a similar pattern, with additional dependence on cluster size. Unmodeled autocorrelation is the most difficult of the three and reaches only about one quarter of replicates at ρ = 0.6.

Heteroscedasticity is different. It is essentially never flagged.

That result follows from the statistic itself. The forest predicts the residual, not its magnitude. If the conditional variance changes while the conditional mean remains zero, there is no mean signal for the forest to recover.

Scope

RandomForestSpecCheck is a test for unmodeled mean structure. It is not a general test of model adequacy.

Pure variance misspecification, heavy tails, and other distributional departures that leave E[e | X] = 0 are outside its target. Variance functions and residual-spread diagnostics remain necessary when those are the concern.

The method is intended to complement those tools rather than replace them.

Application to the Framingham Heart Study

The diagnostic was applied to a subset of the Framingham Heart Study containing more than 2,600 cholesterol measurements taken at two-year intervals over ten years. Baseline age, sex, and body mass index were included.

Five models were fit. Four were plausible specifications an analyst might reasonably consider. The fifth deliberately omitted time from the fixed effects even though time is known to matter biologically.

Table 1 · Five linear mixed models fitted to the Framingham subset, with diagnostic results
Model AIC BIC Observed R² Permutation 2.5th Permutation 97.5th Verdict
Random intercepts 130708130760 −0.008−0.162−0.127 No strong evidence
Random slopes 130415130482 −0.046−0.164−0.126 No strong evidence
AR(1) residual correlation 130710130770 −0.007−0.162−0.126 No strong evidence
Quadratic in time 130709130770 0.031−0.075−0.050 No strong evidence
Time omitted 133091133136 0.191−0.156−0.123 Likely misspecified
Figure 2 · Framingham

Observed out-of-bag R² against its permutation interval, by model

permutation 95% interval observed R², not flagged observed R², flagged τ = 0.05

The four plausible models sit close to zero at this scale. The model with time omitted reaches an observed R² of 0.191, nearly four times the practical threshold.

The intended separation is therefore clear. The four reasonable models are not flagged. The deliberately misspecified model is.

Residual-versus-fitted plots do not make that distinction nearly as clearly. All five produce dense clouds centered around zero. The misspecification is in the mean structure, and the forest is more sensitive to that structure than the visual residual check.

There is also a less comfortable result. All five observed R² values exceed their own permutation 97.5th percentiles. The percentile criterion alone separates nothing in this application.

The practical threshold is what distinguishes the sabotaged model. This is important because the permutation null is reliably negative when the shuffled residuals contain no predictive signal. On a dataset of this size, almost any real structure can clear the percentile criterion. The operating behavior is therefore strongly influenced by τ.

An analyst should choose that threshold with the scale of a meaningful effect in mind rather than accept the default automatically.

Limitations

The procedure detects conditional mean structure only. Power against subtle effects is modest and depends on the number of clusters, the observations within each cluster, and the noise level.

A null result in a small or noisy study is therefore weak evidence that the model is adequate.

The forest can only detect structure in variables it is given. A relevant omitted variable that never enters the random forest is invisible to the diagnostic.

The within-cluster permutation also creates a specific blind spot. A wholly omitted cluster-level effect will not be flagged because the permutation preserves cluster membership. That failure is often visible through simpler diagnostics, but it still has to be checked elsewhere.

RandomForestSpecCheck is best understood as a screening instrument. It is a useful first question about remaining mean structure, not a final verdict on the model.

The function

The procedure is packaged as a single R function over ranger and dplyr. It takes a fitted model, the data, a formula specifying the predictors to search, the number of permutations, and the practical threshold.

RandomForestSpecCheck(
  fitted_model = my_lmm,
  data         = my_data,
  formula_rf   = resid_model ~ time + cluster,
  B            = 500,
  threshold    = 0.05
)
RandomForestSpecCheck: Random Forest Model Specification Diagnostic for LMMs

Practical threshold:                        0.05
Observed R-squared:                         0.191
Permutation mean R-squared and 95% CI:     -0.009 (-0.011, -0.006)

Result: Model likely misspecified

Interpretation: The observed R-squared exceeds both the
upper bound of the null distribution and the practical
threshold, indicating that the model leaves systematic
structure in the residuals. This suggests that the model
may be misspecified.

The output is intentionally descriptive. A bare logical result would encourage the test to be treated as a verdict.

A flag means that the named predictors still contain information about the residuals. The next step is to inspect the forest's variable-importance measures, fit a plausible parametric extension, and see whether the signal disappears.

Contribution

Permutation inference and tree ensembles are established methods. Machine-learning approaches to residual analysis have also been explored for ordinary and generalized linear models.

The contribution here is their combination in a mixed-model setting. The procedure uses a within-cluster permutation scheme, an out-of-bag random-forest statistic with an interpretable scale, and a two-part decision rule that separates statistical evidence from practical magnitude.

The simulation studies also map where the method has power and where it does not. The negative result for pure variance misspecification is part of that scope rather than something to hide.

The work was written for submission to Computational Statistics & Data Analysis.

← All statistical work