Partially-within subjects designs

Authors

Phillip Alday

Reinhold Kliegl

Published

2026-06-21

After working through this page you will be able to:

distinguish within-unit from between-unit factors and identify a partially-within design;
simulate data for such a design with MixedModelsSim;
specify the random-effects structure appropriate to crossed subjects and items.

Before you start

Prerequisites: Analysis of the sleepstudy data.

Datasets used: simulated data (no external dataset).

Begin by loading the packages to be used.

Code

using AlgebraOfGraphics
using CairoMakie
using DataFrames
using MixedModels
using MixedModelsMakie
using MixedModelsSim
using ProgressMeter
using Random
using StableRNGs

const progress = isinteractive()

Code

n_subj = 40
n_item = 3
# things are expressed as "between", so "within subjects" is "between items"
item_btwn = Dict(:frequency => ["high", "medium", "low"])
design = simdat_crossed(StableRNG(42), n_subj, n_item;
                        item_btwn = item_btwn)
design = DataFrame(design)

120×4 DataFrame

95 rows omitted

Row	subj	item	frequency	dv
	String	String	String	Float64
1	S01	I1	high	-0.670252
2	S02	I1	high	0.447122
3	S03	I1	high	1.37363
4	S04	I1	high	1.30954
5	S05	I1	high	0.12607
6	S06	I1	high	0.683948
7	S07	I1	high	-1.0192
8	S08	I1	high	-0.793513
9	S09	I1	high	1.77472
10	S10	I1	high	1.29735
11	S11	I1	high	-1.64385
12	S12	I1	high	0.794439
13	S13	I1	high	-1.30967
⋮	⋮	⋮	⋮	⋮
109	S29	I3	low	1.23469
110	S30	I3	low	-1.12075
111	S31	I3	low	-1.39684
112	S32	I3	low	-0.658154
113	S33	I3	low	-1.51584
114	S34	I3	low	-0.711622
115	S35	I3	low	-0.411443
116	S36	I3	low	-1.25361
117	S37	I3	low	2.08165
118	S38	I3	low	-0.530435
119	S39	I3	low	-1.63993
120	S40	I3	low	-0.768712

Code

unique!(select(design, :item, :frequency))

3×2 DataFrame

Row	item	frequency
	String	String
1	I1	high
2	I2	medium
3	I3	low

Code

m0 = let contrasts, form
    contrasts = Dict(:frequency => HelmertCoding(base="high"))
    form = @formula(dv ~ 1 + frequency +
                    (1 + frequency | subj))
    fit(MixedModel, form, design; contrasts, progress)
end

Code

corrmat = [ 1    0.1 -0.2
            0.1  1    0.1
           -0.2  0.1  1 ]
re_subj = create_re(1.2, 1.5, 1.5; corrmat)

3×3 LinearAlgebra.LowerTriangular{Float64, Matrix{Float64}}:
  1.2    ⋅         ⋅ 
  0.15  1.49248    ⋅ 
 -0.3   0.180907  1.45852

Code

θ = createθ(m0; subj=re_subj)

6-element Vector{Float64}:
  1.2
  0.15000000000000002
 -0.30000000000000004
  1.49248115565993
  0.18090680674665818
  1.4585173044131932

Code

σ = 1;
β = [1.0, -3, -2];

Code

fit!(simulate!(m0; θ, β, σ))

Minimizing 2    Time: 0:00:00 (87.90 ms/it)
   objective: 586.9949278182586


Minimizing 65    Time: 0:00:00 ( 5.61 ms/it)

	Est.	SE	z	p	σ_subj
(Intercept)	1.0744	0.2100	5.12	<1e-06	1.1711
frequency: low	-2.7846	0.2478	-11.24	<1e-28	1.3669
frequency: medium	-2.0779	0.2134	-9.74	<1e-21	1.2748
Residual	1.0846

Code

shrinkageplot(m0)

Code

caterpillar(m0; orderby=nothing, vline_at_zero=true)

Code

design[!, :dv] .= response(m0)

120-element Vector{Float64}:
  7.867093401160148
  5.538829250915448
  5.414289298190056
  7.310706578089037
  6.145051801242777
  4.89132774951212
  7.415483715570554
  1.4665005247908267
  8.664907781396643
  0.8623860676977291
  ⋮
 -1.177447105372094
 -1.5997586226969664
  4.054132094174037
 -3.532926167685095
  4.508469953415926
  0.6623510767242606
 -0.4694143252005385
  4.268023336514854
  2.8708789937181893

Code

design_partial = filter(design) do row
    subj = parse(Int, row.subj[2:end])
    item = parse(Int, row.item[2:end])
    # for even-numbered subjects, we keep all conditions
    # for odd-numbered subjects, we keep only the two "odd" items,
    # i.e. the first and last conditions
    return iseven(subj) || isodd(item)
end
sort!(unique!(select(design_partial, :subj, :frequency)), :subj)

100×2 DataFrame

75 rows omitted

Row	subj	frequency
	String	String
1	S01	high
2	S01	low
3	S02	high
4	S02	medium
5	S02	low
6	S03	high
7	S03	low
8	S04	high
9	S04	medium
10	S04	low
11	S05	high
12	S05	low
13	S06	high
⋮	⋮	⋮
89	S36	medium
90	S36	low
91	S37	high
92	S37	low
93	S38	high
94	S38	medium
95	S38	low
96	S39	high
97	S39	low
98	S40	high
99	S40	medium
100	S40	low

Code

m1 = let contrasts, form
    contrasts = Dict(:frequency => HelmertCoding(base="high"))
    form = @formula(dv ~ 1 + frequency +
                    (1 + frequency | subj))
    fit(MixedModel, form, design_partial; contrasts, progress)
end

Code

shrinkageplot(m1)

Code

caterpillar(m1; orderby=nothing, vline_at_zero=true)

1 Exercises

Classify the factors. For the design simulated here, which experimental factors are manipulated within subjects, and which are within items? Why does that distinction determine where each factor can appear in the random-effects structure?

Solution

A factor can have a by-subject random slope only if it varies within subjects (each subject sees more than one level); likewise a by-item random slope requires the factor to vary within items. A factor that is between subjects (each subject sees only one level) cannot have a by-subject random slope — there is no within-subject variation for the slope to capture. In a partially-within design some factors are within one grouping factor but between the other, which is exactly what shapes the maximal random-effects structure.

Change the design size. Re-run the simulation with more items per subject. How does increasing n_item affect your ability to estimate the by-item random effects?

Solution

With very few items the by-item variance components and any item-level correlations are poorly estimated and the fit is prone to singularity; increasing n_item gives the model more information per item-level effect and generally stabilizes those estimates. Try n_item = 3, 10, 30 and compare issingular and the estimated item variance.

This page was rendered from git revision e563d55 using Quarto 1.9.38 and Julia 1.12.6.