This repository has been archived by the owner on Sep 20, 2023. It is now read-only.
-
Notifications
You must be signed in to change notification settings - Fork 3
/
Copy pathcheck_emotikon_transform.qmd
141 lines (107 loc) · 3.13 KB
/
check_emotikon_transform.qmd
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
---
title: Transformed and original metrics in Emotikon
---
In @Fuehner2021 the original metric of two tasks (Star, S20) is time, but they were transformed to speed scores in the publication prior to computing z-scores.
The critical result is the absence of evidence for the age x Sex x Test interaction.
Is this interaction significant if we analyse all tasks in their original metric?
Fitting the LMM of the publication takes time, roughly 1 hour.
However, if you save the model parameters (and other relevant information), you can restore the fitted model object very quickly.
The notebook also illustrates this procedure.
## Getting the packages and data
```{julia}
#| code-fold: true
using AlgebraOfGraphics
using Arrow
using CairoMakie
using DataFrames
using DataFrameMacros
using MixedModels
using MixedModelsMakie
using RCall
using Serialization
using StatsBase
CairoMakie.activate!(; type="svg")
datadir = joinpath(@__DIR__, "data");
```
### Data and figure in publication
```{julia}
dat = DataFrame(Arrow.Table(joinpath(datadir, "fggk21.arrow")))
@transform!(dat, :a1 = :age - 8.5);
select!(groupby(dat, :Test), :, :score => zscore => :zScore);
describe(dat)
```
### Data and figure with z-scores based on original metric
```{julia}
# dat_om = rcopy(R"readRDS('./data/fggk21_om.rds')"); #Don't know what the _om is
# @transform!(dat_om, :a1 = :age - 8.5);
# select!(groupby(dat_om, :Test), :, :score => zscore => :zScore);
# describe(dat_om)
```
## LMMs
### Contrasts
```{julia}
contrasts = Dict(
:Test => SeqDiffCoding(),
:Sex => HelmertCoding(),
:School => Grouping(),
:Child => Grouping(),
:Cohort => Grouping(),
);
```
### Formula
```{julia}
f1 = @formula zScore ~
1 +
Test * a1 * Sex +
(1 + Test + a1 + Sex | School) +
(1 + Test | Child) +
zerocorr(1 + Test | Cohort);
```
### Restore LMM m1 from publication
- Command for fitting LMM m1 = fit(MixedModel, f1, dat, contrasts=contr)
- Fit statistics for LMM m1: Minimizing 5179 Time: 0 Time: 1:00:38 ( 0.70 s/it)
```{julia}
m1x = LinearMixedModel(f1, dat; contrasts)
restoreoptsum!(m1x, "./fits/fggk21_m1_optsum.json")
```
```{julia}
VarCorr(m1x)
```
### Restore new LMM m1_om Star and S20 in original metric
- Command for fitting LMM m1_om = fit(MixedModel, f1, dat_om, contrasts=contr)
- Minimizing 10502 Time: 0 Time: 2:09:40 ( 0.74 s/it)
- Store with: julia> saveoptsum("./fits/fggk21_m1_om_optsum.json", m1_om)
- Only for short-term and when desparate: julia> serialize("./fits/m1_om.jls", m1_om);
#### ... restoreoptsum!()
```{julia}
#| eval: false
m1_om = LinearMixedModel(f1, dat; contrasts=contr);
restoreoptsum!(m1_om, "./fits/fggk21_m1_om_optsum.json");
```
#### ... deserialize()
```{julia}
#| eval: false
m1x_om = deserialize("./fits/m1_om.jls")
```
```{julia}
#| eval: false
VarCorr(m1x_om)
```
### Residual diagnostics for LMM m1
Residual plots for published LMM
```{julia}
#scatter(fitted(m1x), residuals(m1x)
```
```{julia}
#qqnorm(m1x)
```
### Residual diagnostics for LMM m1_om
Residual plots for LMM with Star and Speed in original metric.
```{julia}
#scatter(fitted(m1_om_v2), residuals(m1_om_v2)
```
```{julia}
#qqnorm(m1_om_v2)
```
::: {#refs}
:::