MixedModels.jl cheatsheet

Authors

Phillip Alday

Reinhold Kliegl

Published

2026-08-23

This page is a lookup reference, not a lesson. Use it to:

  • find the extractor function for a particular piece of information about a fitted model;
  • understand what each of the commonly used extractors returns.
NoteBefore you start

Prerequisites: Analysis of the sleepstudy data.

Datasets used: sleepstudy (see the dataset catalog).

Here we introduce most of the commands available in the MixedModels.jl package for the immediate inspection and analysis of results returned in a fitted linear mixed-effect model. All the examples use the random-slopes model for the sleepstudy data, which is small enough that everything here runs in seconds.

Code 
using DataFrames
using MixedModels
using MixedModelsDatasets: dataset
using SMLP2026: fit_or_restore
m1 = fit(MixedModel,
         @formula(reaction ~ 1 + days + (1 + days | subj)),
         dataset(:sleepstudy))

1 Overall summary statistics

m1.optsum            # MixedModels.OptSummary: all info about the optimization
Initialization
Initial parameter vector [1.0, 0.0, 1.0]
Initial objective value 1784.6423001241371
Optimizer settings
Optimizer LN_NEWUOA
Backend nlopt
ftol_rel 1.0e-12
ftol_abs 1.0e-8
xtol_rel 0.0
xtol_abs [1.0e-10, 1.0e-10, 1.0e-10]
initial_step [1.0, 1.0, 1.0]
maxfeval -1
maxtime -1.0
xtol_zero_abs 0.001
ftol_zero_abs 1.0e-5
pirls_maxiter 10
pirls_ftol_rel 1.4901161193847656e-8
pirls_ftol_abs 1.0e-5
pirls_maxhalfstep 10
Result
Function evaluations 82
Final parameter vector [0.9292, 0.0182, 0.2226]
Final objective value 1751.9393
Return code FTOL_REACHED 
loglikelihood(m1) # StatsBase.loglikelihood: log-likelihood of the model
-875.9696734578358
deviance(m1)      # StatsBase.deviance: negative twice the log-likelihood relative to saturated model
1751.9393469156716
objective(m1)    # MixedModels.objective: saturated model not clear: negative twice the log-likelihood
1751.9393469156716
nobs(m1) # n of observations; they are not independent
180
n, p, q, k = size(m1)  # nobs, fixed-effects params, VCs+CPs, grouping factors
(180, 2, 36, 1)
dof(m1)  # n of degrees of freedom is number of model parameters
6
geom_df = sum(leverage(m1)) # trace of hat / rank of model matrix / geom dof
28.611612940625573
resid_df = nobs(m1) - geom_df  # effective residual degrees of freedom
151.38838705937442
aic(m1)  # objective(m1) + 2*dof(m1)
1763.9393469156716
bic(m1)  # objective(m1) + dof(m1)*log(nobs(m1))
1783.097088021013

2 Fixed-effect statistics

coeftable(m1) # StatsBase.coeftable: fixed-effects statistics; default level=0.95
Coef. Std. Error z Pr(>|z|)
(Intercept) 251.405 6.63228 37.91 <1e-99
days 10.4673 1.50224 6.97 <1e-11 
coef(m1)              # StatsBase.coef; parts of the table
2-element Vector{Float64}:
 251.4051060532072
  10.46728550560944
fixef(m1)    # MixedModels.fixef: not the same as coef() for rank-deficient case
2-element Vector{Float64}:
 251.4051060532072
  10.46728550560944
m1.β                  # alternative extractor
2-element Vector{Float64}:
 251.4051060532072
  10.46728550560944
fixefnames(m1)        # works also for coefnames(m1)
2-element Vector{String}:
 "(Intercept)"
 "days"
vcov(m1)   # StatsBase.vcov: var-cov matrix of fixed-effects coefficients
2×2 Matrix{Float64}:
 43.9871  -1.3705
 -1.3705   2.25672
vcov(m1; corr=true) # StatsBase.vcov: correlation matrix of fixed-effects coefficients
2×2 Matrix{Float64}:
  1.0       -0.137555
 -0.137555   1.0
stderror(m1)       # StatsBase.stderror: SE for fixed-effects coefficients
2-element Vector{Float64}:
 6.632279012514935
 1.5022367701697787
propertynames(m1)  # names of available extractors
(:formula, :reterms, :Xymat, :feterm, :sqrtwts, :parmap, :dims, :A, :L, :optsum, :θ, :theta, :β, :beta, :βs, :betas, :λ, :lambda, :stderror, :σ, :sigma, :σs, :sigmas, :σρs, :sigmarhos, :b, :u, :X, :y, :corr, :vcov, :PCA, :rePCA, :objective, :pvalues)

3 Covariance parameter estimates

These commands inform us about the model parameters associated with the random-effects structure (RES).

issingular(m1) # Test if model is singular for parameter vector m1.theta (default)
false
VarCorr(m1) # MixedModels.VarCorr: estimates of random-effect structure (RES)
Column Variance Std.Dev Corr.
subj (Intercept) 565.51589 23.78058
days 32.68220 5.71683 +0.08
Residual 654.94103 25.59182 
m1.σs      # VCs; m1.sigmas
(subj = (var"(Intercept)" = 23.7805780726916, days = 5.716834485317427),)
m1.θ       # Parameter vector for RES (w/o residual); m1.theta
3-element Vector{Float64}:
 0.9292259065867788
 0.018165149307614335
 0.2226454658393477
BlockDescription(m1) #  Description of blocks of A and L in a LinearMixedModel
rows subj fixed
36 BlkDiag
3 Dense Dense

4 Model “predictions”

These commands extract the conditional modes/means and (co-)variances of the random effects, that is the model’s “predictions” for the units (levels) of the grouping (random) factors.

condVar(m1)
1-element Vector{Array{Float64, 3}}:
 [140.9671928188997 -20.604478452614167; -20.604478452614167 5.157784594104255;;; 140.9671928188997 -20.604478452614167; -20.604478452614167 5.157784594104255;;; 140.9671928188997 -20.604478452614167; -20.604478452614167 5.157784594104255;;; … ;;; 140.9671928188997 -20.604478452614167; -20.604478452614167 5.157784594104255;;; 140.9671928188997 -20.604478452614167; -20.604478452614167 5.157784594104255;;; 140.9671928188997 -20.604478452614167; -20.604478452614167 5.157784594104255]

The conditional covariances are hard to look at directly; they are usually visualized instead:

  • Caterpillar plots (caterpillar, qqcaterpillar from MixedModelsMakie) display prediction intervals of the random effects per level; they are introduced in the sleepstudy analysis.
  • Shrinkage plots (shrinkageplot from MixedModelsMakie) visualize how the conditional modes are “shrunk” toward the grand mean relative to per-level estimates; they are the topic of Shrinkage of random effects.

5 See Also

This page was rendered from git revision 29e8d33 using Quarto 1.10.18 and Julia 1.12.7.

Back to top