Skip to content

Commit bf1471f

Browse files
committed
more compat updates
1 parent 6af2189 commit bf1471f

File tree

3 files changed

+41
-35
lines changed

3 files changed

+41
-35
lines changed

kkl15.qmd

Lines changed: 24 additions & 19 deletions
Original file line numberDiff line numberDiff line change
@@ -31,7 +31,6 @@ Stroup, W. W. (2012, p. 185). _Generalized linear mixed models: Modern concepts,
3131
#| code-fold: true
3232
#| output: false
3333
using AlgebraOfGraphics
34-
using AlgebraOfGraphics: density
3534
using BoxCox
3635
using CairoMakie
3736
using CategoricalArrays
@@ -41,8 +40,12 @@ using DataFrames
4140
using MixedModels
4241
using MixedModelsMakie
4342
using Random
44-
using SMLP2026: dataset
4543
using StatsBase
44+
45+
using AlgebraOfGraphics: density
46+
using SMLP2026: dataset
47+
48+
progress = isinteractive()
4649
```
4750

4851
# Read data, compute and plot means
@@ -101,7 +104,7 @@ contrasts = Dict(
101104

102105
# Maximum LMM
103106

104-
This is the maximum LMM for the design; `size` is a between-subject factor,
107+
This is the maximum LMM for the design; `size` is a between-subject factor,
105108
ignoring other information such as trial number, age and gender of subjects.
106109

107110
```{julia}
@@ -124,7 +127,7 @@ only(MixedModels.PCA(m_max))
124127
VarCorr(m_max)
125128
```
126129

127-
The LMM `m_max` is overparameterized but it is not immediately apparent why.
130+
The LMM `m_max` is overparameterized but it is not immediately apparent why.
128131

129132
# Reduction strategy 1
130133

@@ -152,7 +155,7 @@ only(MixedModels.PCA(m_zcp1))
152155
VarCorr(m_zcp1)
153156
```
154157

155-
The LMM `m_zcp1` is also overparameterized, but now there is clear evidence for absence of evidence for the VC of one of the interactions and the other two interaction-based VCs are also very small.
158+
The LMM `m_zcp1` is also overparameterized, but now there is clear evidence for absence of evidence for the VC of one of the interactions and the other two interaction-based VCs are also very small.
156159

157160
## Reduced zcp LMM
158161

@@ -178,7 +181,7 @@ only(MixedModels.PCA(m_zcp1_rdc))
178181
VarCorr(m_zcp1_rdc)
179182
```
180183

181-
LMM `m_zcp_rdc` is ok . We add in CPs.
184+
LMM `m_zcp_rdc` is ok . We add in CPs.
182185

183186
## Parsimonious LMM (1)
184187

@@ -200,7 +203,7 @@ issingular(m_prm1)
200203
only(MixedModels.PCA(m_prm1))
201204
```
202205

203-
LMM `m_zcp_rdc` is ok . We add in CPs.
206+
LMM `m_zcp_rdc` is ok . We add in CPs.
204207

205208
```{julia}
206209
VarCorr(m_prm1)
@@ -221,21 +224,22 @@ gof_summary = let
221224
deviance=round.(deviance.(mods), digits=0),
222225
AIC=round.(aic.(mods),digits=0),
223226
AICc=round.(aicc.(mods),digits=0),
224-
BIC=round.(bic.(mods),digits=0),
225-
χ²=vcat(:., round.(lrt.tests.deviancediff, digits=0)),
226-
χ²_dof=vcat(:., round.(lrt.tests.dofdiff, digits=0)),
227-
pvalue=vcat(:., round.(lrt.tests.pvalues, digits=3))
227+
BIC=round.(bic.(mods),digits=0),
228+
χ²=vcat(:., round.(Int, diff(collect(lrt.lrt.deviance)))),
229+
χ²_dof=vcat(:., diff(collect(lrt.lrt.dof))),
230+
# StatsBase.PValue includes some pretty-printing methods
231+
pvalue=vcat(:., StatsBase.PValue.(collect(lrt.lrt.pval[2:end])))
228232
)
229233
end
230234
```
231235

232-
AIC prefers LMM `m_prm1` over `m_zcp1_rdc`; BIC LMM `m_zcp1_rdc`.
236+
AIC prefers LMM `m_prm1` over `m_zcp1_rdc`; BIC LMM `m_zcp1_rdc`.
233237
As the CPs were one reason for conducting this experiment, AIC is the criterion of choice.
234238

235239
# Reduction strategy 2
236240

237241
## Complex LMM
238-
Relative to LMM `m_max`, first we take out interaction VCs and associated CPs, because these VCs are very small. This is the same as LMM `m_prm1` above.
242+
Relative to LMM `m_max`, first we take out interaction VCs and associated CPs, because these VCs are very small. This is the same as LMM `m_prm1` above.
239243

240244
```{julia}
241245
m_cpx = let
@@ -290,16 +294,17 @@ gof_summary = let
290294
deviance=round.(deviance.(mods), digits=0),
291295
AIC=round.(aic.(mods),digits=0),
292296
AICc=round.(aicc.(mods),digits=0),
293-
BIC=round.(bic.(mods),digits=0),
294-
χ²=vcat(:., round.(lrt.tests.deviancediff, digits=0)),
295-
χ²_dof=vcat(:., round.(lrt.tests.dofdiff, digits=0)),
296-
pvalue=vcat(:., round.(lrt.tests.pvalues, digits=3))
297+
BIC=round.(bic.(mods),digits=0),
298+
χ²=vcat(:., round.(Int, diff(collect(lrt.lrt.deviance)))),
299+
χ²_dof=vcat(:., diff(collect(lrt.lrt.dof))),
300+
# StatsBase.PValue includes some pretty-printing methods
301+
pvalue=vcat(:., StatsBase.PValue.(collect(lrt.lrt.pval[2:end])))
297302
)
298303
end
299304
```
300305

301-
The cardinal-related CPs could be removed w/o loss of goodness of fit.
302-
However, there is no harm in keeping them in the LMM.
306+
The cardinal-related CPs could be removed w/o loss of goodness of fit.
307+
However, there is no harm in keeping them in the LMM.
303308
The data support both LMM `m_prm2` and `m_cpx` (same as: `m_prm1`).
304309
We keep the slightly more complex LMM `m_cpx` (`m_prm1`).
305310

mrk17.qmd

Lines changed: 11 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -13,10 +13,12 @@ Packages we (might) use.
1313
using DataFrames
1414
using MixedModels
1515
using MixedModelsMakie
16+
using StatsBase
17+
1618
using SMLP2026: dataset
1719
using Statistics: mean, std
1820
19-
const progress=false
21+
const progress = isinteractive()
2022
```
2123

2224
```{julia}
@@ -283,9 +285,10 @@ gof_summary = let
283285
AIC=round.(aic.(mods),digits=0),
284286
AICc=round.(aicc.(mods),digits=0),
285287
BIC=round.(bic.(mods),digits=0),
286-
χ²=vcat(:., round.(lrt.tests.deviancediff, digits=0)),
287-
χ²_dof=vcat(:., round.(lrt.tests.dofdiff, digits=0)),
288-
pvalue=vcat(:., round.(lrt.tests.pvalues, digits=3))
288+
χ²=vcat(:., round.(Int, diff(collect(lrt.lrt.deviance)))),
289+
χ²_dof=vcat(:., diff(collect(lrt.lrt.dof))),
290+
# StatsBase.PValue includes some pretty-printing methods
291+
pvalue=vcat(:., StatsBase.PValue.(collect(lrt.lrt.pval[2:end])))
289292
)
290293
end
291294
```
@@ -378,9 +381,10 @@ gof_summary = let
378381
AIC=aic.(mods),
379382
AICc=aicc.(mods),
380383
BIC=bic.(mods),
381-
χ²=vcat(:.,lrt.tests.deviancediff),
382-
χ²_dof=vcat(:.,lrt.tests.dofdiff),
383-
pvalue=vcat(:., round.(lrt.tests.pvalues, digits=3))
384+
χ²=vcat(:., round.(Int, diff(collect(lrt.lrt.deviance)))),
385+
χ²_dof=vcat(:., diff(collect(lrt.lrt.dof))),
386+
# StatsBase.PValue includes some pretty-printing methods
387+
pvalue=vcat(:., StatsBase.PValue.(collect(lrt.lrt.pval[2:end])))
384388
)
385389
end
386390
```

singularity.qmd

Lines changed: 6 additions & 9 deletions
Original file line numberDiff line numberDiff line change
@@ -19,15 +19,14 @@ using LinearAlgebra
1919
using MixedModels
2020
using MixedModelsMakie
2121
22-
const progress=false
22+
const progress = isinteractive()
2323
```
2424

25-
Fit a model for reaction time in the sleepstudy example, preserving information on the estimation progress (the `thin=1` optional argument)
25+
Fit a model for reaction time in the sleepstudy example.
2626

2727
```{julia}
28-
m01 = let f = @formula reaction ~ 1 + days + (1 + days|subj)
29-
fit(MixedModel, f, MixedModels.dataset(:sleepstudy); thin=1, progress)
30-
end
28+
m01 = lmm(@formula(reaction ~ 1 + days + (1 + days|subj)), MixedModels.dataset(:sleepstudy); progress)
29+
3130
print(m01)
3231
```
3332

@@ -82,16 +81,14 @@ The elements of this parameter vector are subject to constraints.
8281
In particular, two of the three elements have a lower bound of zero.
8382

8483
```{julia}
85-
m01.lowerbd
84+
lowerbd(m01)
8685
```
8786

8887
That is, the first and third elements of $\theta$, corresponding to diagonal elements of $\lambda$, must be non-negative, whereas the second component is unconstrained (has a lower bound of $-\infty$).
8988

9089
# Progress of the iterations
9190

92-
The `optsum.fitlog` property of the model is a vector of tuples where each tuple contains the value of the $\theta$ vector and the value of the objective at that $\theta$.
93-
The `fitlog` always contains the first and the last evaluation.
94-
When the `thin` named argument is set, this property has a row for every thin'th evaluation.
91+
The `optsum.fitlog` property of the model is table that contains the value of the $\theta$ vector and the value of the objective at that $\theta$.
9592

9693
```{julia}
9794
m01.optsum.fitlog

0 commit comments

Comments
 (0)