Simulating an Experiment from Scratch
Here is a worked example of simulating a partially crossed design from scratch.
First, some setup:
using DataFrames
using MixedModels, MixedModelsSim
using Random
Assemble the Design
We're going to do a 2 x 2 x 2 design. For concreteness, let's think of this as a linguistic design:
- age
old
vs.young
, between subjects - frequency
high
vs.low
, between items - context
matched
vs.unmatched
, within both.
Further, let's assume we want 40 subjects and 80 items. We can specify this design as follows:
n_subj = 40
n_item = 80
subj_btwn = Dict(:age => ["old", "young"])
item_btwn = Dict(:frequency => ["high", "low"])
both_win = Dict(:context => ["matched", "unmatched"])
Dict{Symbol, Vector{String}} with 1 entry:
:context => ["matched", "unmatched"]
and then generate it:
rng = MersenneTwister(42) # specify our random number generator for reproducibility
design = simdat_crossed(rng, n_subj, n_item;
subj_btwn = subj_btwn,
item_btwn = item_btwn,
both_win = both_win)
6400-element Vector{NamedTuple{(:subj, :age, :item, :frequency, :context, :dv), Tuple{String, String, String, String, String, Float64}}}:
(subj = "S01", age = "old", item = "I01", frequency = "high", context = "matched", dv = -0.5560268761463861)
(subj = "S02", age = "young", item = "I01", frequency = "high", context = "matched", dv = -0.444383357109696)
(subj = "S03", age = "old", item = "I01", frequency = "high", context = "matched", dv = 0.027155338009193845)
(subj = "S04", age = "young", item = "I01", frequency = "high", context = "matched", dv = -0.29948409035891055)
(subj = "S05", age = "old", item = "I01", frequency = "high", context = "matched", dv = 1.7778610980573246)
(subj = "S06", age = "young", item = "I01", frequency = "high", context = "matched", dv = -1.14490153172882)
(subj = "S07", age = "old", item = "I01", frequency = "high", context = "matched", dv = -0.46860588216767457)
(subj = "S08", age = "young", item = "I01", frequency = "high", context = "matched", dv = 0.15614346264074028)
(subj = "S09", age = "old", item = "I01", frequency = "high", context = "matched", dv = -2.641991008076796)
(subj = "S10", age = "young", item = "I01", frequency = "high", context = "matched", dv = 1.0033099014594844)
⋮
(subj = "S32", age = "young", item = "I80", frequency = "low", context = "unmatched", dv = -0.0720635904250138)
(subj = "S33", age = "old", item = "I80", frequency = "low", context = "unmatched", dv = -0.5636260084084584)
(subj = "S34", age = "young", item = "I80", frequency = "low", context = "unmatched", dv = -0.48302653966207154)
(subj = "S35", age = "old", item = "I80", frequency = "low", context = "unmatched", dv = 0.5635515836664501)
(subj = "S36", age = "young", item = "I80", frequency = "low", context = "unmatched", dv = -1.5342115330081054)
(subj = "S37", age = "old", item = "I80", frequency = "low", context = "unmatched", dv = -0.7925886663164663)
(subj = "S38", age = "young", item = "I80", frequency = "low", context = "unmatched", dv = 1.1344191438635207)
(subj = "S39", age = "old", item = "I80", frequency = "low", context = "unmatched", dv = 1.3406462534138157)
(subj = "S40", age = "young", item = "I80", frequency = "low", context = "unmatched", dv = 0.7969024201220799)
Note that simdat_crossed
returns a row table, which MixedModels.jl
can process directly. For nicely displaying it, we can again use pretty_table
:
pretty_table(first(design, 5))
┌────────┬────────┬────────┬───────────┬─────────┬───────────┐
│ subj │ age │ item │ frequency │ context │ dv │
│ String │ String │ String │ String │ String │ Float64 │
├────────┼────────┼────────┼───────────┼─────────┼───────────┤
│ S01 │ old │ I01 │ high │ matched │ -0.556027 │
│ S02 │ young │ I01 │ high │ matched │ -0.444383 │
│ S03 │ old │ I01 │ high │ matched │ 0.0271553 │
│ S04 │ young │ I01 │ high │ matched │ -0.299484 │
│ S05 │ old │ I01 │ high │ matched │ 1.77786 │
└────────┴────────┴────────┴───────────┴─────────┴───────────┘
We can also convert it to a DataFrame and change the factors to use pooled arrays to save a bit of memory.
df = pooled!(DataFrame(design))
first(df, 5)
5 rows × 6 columns
subj | age | item | frequency | context | dv | |
---|---|---|---|---|---|---|
String | String | String | String | String | Float64 | |
1 | S01 | old | I01 | high | matched | -0.556027 |
2 | S02 | young | I01 | high | matched | -0.444383 |
3 | S03 | old | I01 | high | matched | 0.0271553 |
4 | S04 | young | I01 | high | matched | -0.299484 |
5 | S05 | old | I01 | high | matched | 1.77786 |
Note that simdat_crossed
generated a column dv
for our dependent variable that has been pre-populated with noise from a standard normal distribution ($N(0,1)$). Typically, we will want to scale that, but we can do that in the simulation step. Also, this dependent variable is pure noise: we haven't yet added in effects. Adding in effects also comes in the simulation step.
But before we get to simulating, let's fit the model to the noise, just to see how things look. We're going to use effects coding for our contrasts.
contrasts = Dict(:age => EffectsCoding(base="young"),
:frequency => EffectsCoding(base="high"),
:context => EffectsCoding(base="matched"))
form = @formula(dv ~ 1 + age * frequency * context +
(1 + frequency + context | subj) +
(1 + age + context | item))
m0 = fit(MixedModel, form, design; contrasts=contrasts)
Est. | SE | z | p | σ_item | σ_subj | |
---|---|---|---|---|---|---|
(Intercept) | -0.0168 | 0.0128 | -1.31 | 0.1895 | 0.0094 | 0.0111 |
age: old | -0.0004 | 0.0128 | -0.03 | 0.9734 | 0.0107 | |
frequency: low | 0.0180 | 0.0131 | 1.37 | 0.1696 | 0.0217 | |
context: unmatched | 0.0099 | 0.0132 | 0.75 | 0.4534 | 0.0062 | 0.0245 |
age: old & frequency: low | 0.0014 | 0.0131 | 0.11 | 0.9147 | ||
age: old & context: unmatched | 0.0068 | 0.0132 | 0.51 | 0.6089 | ||
frequency: low & context: unmatched | -0.0084 | 0.0126 | -0.67 | 0.5043 | ||
age: old & frequency: low & context: unmatched | 0.0000 | 0.0126 | 0.00 | 0.9972 | ||
Residual | 1.0095 |
Assemble the Random Effects
The hard part in simulating right now is specifying the random effects. We're working on making this bit easier, but you need to specify the variance-covariance matrix of the random effects. You can see what this looks like:
vc = VarCorr(m0)
Column | Variance | Std.Dev | Corr. | ||
---|---|---|---|---|---|
item | (Intercept) | 0.00008876 | 0.00942136 | ||
age: old | 0.00011356 | 0.01065636 | +1.00 | ||
context: unmatched | 0.00003836 | 0.00619327 | -1.00 | -1.00 | |
subj | (Intercept) | 0.00012387 | 0.01112973 | ||
frequency: low | 0.00047252 | 0.02173754 | +1.00 | ||
context: unmatched | 0.00060127 | 0.02452087 | +1.00 | +1.00 | |
Residual | 1.01911782 | 1.00951366 |
For each grouping variable (subjects and items), there are two major components: the standard deviations ahd the correlations.
For this example, we'll just assume all the correlations and thus the covariances are 0 in order to make things simple. Then we only have to worry about the standard deviations.
Let's assume that the variability
- between items
- in the intercept is 1.3 times the residual variability
- in age is 0.35 times the residual variability
- in context is 0.75 times the residual variability
- between subjects
- in the intercept is 1.5 times the residual variability
- in frequency is 0.5 times the residual variability
- in context is 0.75 times the residual variability
Note these are always specified relative to the residual standard deviation. In other words, we think about how big the between-subject and between-item differences are relative to the between-observation differences.
We can now create the associated covariance matrices.[cholesky]
re_item = create_re(1.3, 0.35, 0.75)
3×3 LinearAlgebra.LowerTriangular{Float64, LinearAlgebra.Diagonal{Float64, Vector{Float64}}}:
1.3 ⋅ ⋅
0.0 0.35 ⋅
0.0 0.0 0.75
re_subj = create_re(1.5, 0.5, 0.75)
3×3 LinearAlgebra.LowerTriangular{Float64, LinearAlgebra.Diagonal{Float64, Vector{Float64}}}:
1.5 ⋅ ⋅
0.0 0.5 ⋅
0.0 0.0 0.75
We can check that we got these right by installing these parameter values into the model. Note that we have to specify them in the same order as in the output from VarCorr
.
update!(m0, re_item, re_subj)
VarCorr(m0)
Column | Variance | Std.Dev | Corr. | ||
---|---|---|---|---|---|
item | (Intercept) | 1.645411 | 1.282736 | ||
age: old | 0.119268 | 0.345352 | +0.00 | ||
context: unmatched | 0.547659 | 0.740040 | +0.00 | +0.00 | |
subj | (Intercept) | 2.190636 | 1.480080 | ||
frequency: low | 0.243404 | 0.493360 | +0.00 | ||
context: unmatched | 0.547659 | 0.740040 | +0.00 | +0.00 | |
Residual | 0.973616 | 0.986720 |
Looks good. The values don't exactly match the values in our parameter vector because the residual standard deviation isn't exactly 1.0.
For the actual simulation, we'll need the compact form of these covariance matrices that MixedModels.jl stores uses internally. This compact form is the parameter vector θ and we can get it back out of the model where we just installed it:
θ = m0.θ
12-element Vector{Float64}:
1.3
0.0
0.0
0.35
0.0
0.75
1.5
0.0
0.0
0.5
0.0
0.75
Alternatively, we could also create it directly from the covariance matrices we created:
vcat( flatlowertri(re_item), flatlowertri(re_subj) )
12-element Vector{Float64}:
1.3
0.0
0.0
0.35
0.0
0.75
1.5
0.0
0.0
0.5
0.0
0.75
Assemble the Fixed Effects
The last two components we need are the residual variance and the effect sizes for the fixed effects.
σ = 5
β = [1.0, -1.0, 2.0, -1.5, 0.3, -1.3, 1.4, 0]
8-element Vector{Float64}:
1.0
-1.0
2.0
-1.5
0.3
-1.3
1.4
0.0
The entries in the β correspond to the coefficients in the model given by
coefnames(m0)
8-element Vector{String}:
"(Intercept)"
"age: old"
"frequency: low"
"context: unmatched"
"age: old & frequency: low"
"age: old & context: unmatched"
"frequency: low & context: unmatched"
"age: old & frequency: low & context: unmatched"
Simulate
Now we're ready to actually simulate our data. We can use parametricbootstrap
to do this: the parametric bootstrap actually works by simulating new data from an existing model and then looking at how the estimates fit to that new data look. In MixedModels.jl, you can specify different parameter values, such as the ones we made up for our fake data.
# typically we would use a lot more simulations
# but we want to be quick in this example
sim = parametricbootstrap(MersenneTwister(12321), 20, m0; β=β, σ=σ, θ=θ)
MixedModels.MixedModelBootstrap{Float64}(NamedTuple{(:objective, :σ, :β, :se, :θ), Tuple{Float64, Float64, NamedTuple{(Symbol("(Intercept)"), Symbol("age: old"), Symbol("frequency: low"), Symbol("context: unmatched"), Symbol("age: old & frequency: low"), Symbol("age: old & context: unmatched"), Symbol("frequency: low & context: unmatched"), Symbol("age: old & frequency: low & context: unmatched")), NTuple{8, Float64}}, StaticArrays.SVector{8, Float64}, StaticArrays.SVector{12, Float64}}}[(objective = 40018.60227027449, σ = 4.9390780702522346, β = ((Intercept) = 2.4058378521150523, age: old = -1.267088325617489, frequency: low = 2.053371478579786, context: unmatched = -1.7401078608612939, age: old & frequency: low = -0.5190328604292482, age: old & context: unmatched = 0.14511206537539048, frequency: low & context: unmatched = 2.1177660632838906, age: old & frequency: low & context: unmatched = -0.021858881860061054), se = [1.1359055417354553, 0.8670382787646409, 0.8530669025820241, 0.6846553199876846, 0.43496794983017434, 0.5427270835918829, 0.42191451661096074, 0.06173847587815293], θ = [1.3802464337291895, 0.07365149518293726, -0.15134478035345333, 0.36551050750798697, -0.016421822668797173, 0.7403393125738269, 1.0755941968172842, -0.15058476624161662, -0.016288590467011376, 0.4602088962308222, -0.0915926856706881, 0.6841620738033822]), (objective = 40186.70348729695, σ = 5.00126891870514, β = ((Intercept) = -0.17367934549309064, age: old = -1.688473739967146, frequency: low = 2.1030922433750097, context: unmatched = -1.8662651448223795, age: old & frequency: low = 0.44977298527597276, age: old & context: unmatched = -1.2812982019846053, frequency: low & context: unmatched = 1.0553300004935695, age: old & frequency: low & context: unmatched = -0.0565016986833358), se = [1.5490219570688515, 1.4207833710957387, 0.763330985172816, 0.7301427274660097, 0.44925555870110595, 0.626758754474116, 0.3797237141639414, 0.06251586148381426], θ = [1.164709705433686, -0.11327566462780597, 0.05414035161103328, 0.35445172275011816, 0.02704259725655156, 0.6670919956354646, 1.7755782726288738, -0.02704646704556169, 0.01112285191709347, 0.4965378475507712, 0.15190280210620577, 0.7737928552352146]), (objective = 40293.01558970997, σ = 5.041203955981211, β = ((Intercept) = 0.6887229123579159, age: old = -1.80371667666853, frequency: low = 0.65171258828356, context: unmatched = -1.26260749343244, age: old & frequency: low = 0.052042797867456705, age: old & context: unmatched = -1.9319507254949613, frequency: low & context: unmatched = 2.0210983526239397, age: old & frequency: low & context: unmatched = -0.05305680821987915), se = [1.30634249016466, 1.0993894965113127, 0.8260418356774109, 0.7585853515697719, 0.42950166206834384, 0.6168007486040792, 0.44607114745841975, 0.06301504944976515], θ = [1.2952383781067975, 0.05727946538026769, 0.03988855889728099, 0.3272672332335201, 0.01719205317709882, 0.7822925022897953, 1.3568072825830189, 0.05758468623507887, 0.06692080661069301, 0.4749643785733052, -0.0023185530621895947, 0.766854267581006]), (objective = 40086.829209165866, σ = 4.959168242629561, β = ((Intercept) = -0.013630950517976102, age: old = -0.4693672440772701, frequency: low = 2.5381968697550534, context: unmatched = -2.6641251153849588, age: old & frequency: low = 0.8833080201462482, age: old & context: unmatched = -0.77433076863467, frequency: low & context: unmatched = 1.4427059716861463, age: old & frequency: low & context: unmatched = -0.050425034028918604), se = [1.3514927163439765, 1.1374039584873545, 0.8571639694697353, 0.6086790059599934, 0.4493164508589499, 0.39822873009808757, 0.46448543755974936, 0.06198960303286952], θ = [1.3706093453824226, 0.042819324874020855, -0.10367061178738415, 0.37873380520371747, 0.038730660512294554, 0.8228349899366034, 1.4231091924879438, 0.1111641876229951, -0.1468954033928837, 0.4869431516790371, 0.1332212611404095, 0.4608222842474606]), (objective = 40383.62186880702, σ = 5.089727586745207, β = ((Intercept) = 0.2743862976012784, age: old = 0.2544221366676008, frequency: low = 1.426328166433766, context: unmatched = -0.40254682628310534, age: old & frequency: low = -0.2741273897667547, age: old & context: unmatched = -0.7396304540571064, frequency: low & context: unmatched = 1.6351514975027364, age: old & frequency: low & context: unmatched = 0.07347458213196383), se = [1.231942312164965, 0.9875562209614842, 0.8879739733036307, 0.7023354414587631, 0.4960677436826271, 0.5281395206769931, 0.46732368478103814, 0.06362159483431509], θ = [1.323177260133653, 0.0014443863840713912, -0.13108683472873728, 0.27518845895838023, 0.09199618906295373, 0.7976730808740481, 1.2090406116450896, 0.043565322304604207, 0.26118830601308163, 0.5778925588890808, -0.10533925452646775, 0.5874757369094146]), (objective = 40073.141971807585, σ = 4.958400093243766, β = ((Intercept) = 2.17318487448217, age: old = -1.6792985233724493, frequency: low = 2.104991176984805, context: unmatched = -2.637335057515553, age: old & frequency: low = 0.3402084319552244, age: old & context: unmatched = -0.8572717880772645, frequency: low & context: unmatched = 1.257663925265949, age: old & frequency: low & context: unmatched = -0.012957628040293365), se = [1.3402504368915737, 1.0989644234559093, 0.8627652994423082, 0.6905103765852335, 0.3947347620305976, 0.5804552055878615, 0.37910137829114915, 0.061980001165547076], θ = [1.4298369807169993, -0.03195914304585946, -0.06221708077209279, 0.35822197012725354, -0.03652746969098447, 0.6707765364115633, 1.3762247006028288, -0.03877081685872741, -0.037979768451371694, 0.4255352725418778, -0.16409178814645953, 0.7166241536670316]), (objective = 40021.69032262348, σ = 4.931840406782373, β = ((Intercept) = 2.8984620242479515, age: old = 0.3403186082642826, frequency: low = 3.9239223815247883, context: unmatched = -0.8812095444105639, age: old & frequency: low = 0.7198568091003207, age: old & context: unmatched = -1.0617893533396892, frequency: low & context: unmatched = 1.8846416345831252, age: old & frequency: low & context: unmatched = 0.08042828251315388), se = [1.3673501883262156, 1.179909446332868, 0.8210496820301841, 0.7116591753145971, 0.44346628326124204, 0.5934652440821501, 0.3975654190277912, 0.06164800508477967], θ = [1.3066180353313481, -0.005588111442935936, 0.050218144267271475, 0.3698966260267475, -0.12647740360548582, 0.6991744465930981, 1.4882256622405836, 0.10164979528498896, -0.08768693595816449, 0.4882695919980609, -0.0017653901738130054, 0.7518398795516547]), (objective = 40104.20130957483, σ = 4.9630802170834905, β = ((Intercept) = 0.2026519481497094, age: old = -0.2095129064018696, frequency: low = 3.4979679693674974, context: unmatched = -2.418734277831136, age: old & frequency: low = 0.009123554115715242, age: old & context: unmatched = -0.9711464983699467, frequency: low & context: unmatched = 1.8612376860128037, age: old & frequency: low & context: unmatched = 0.04828469357173075), se = [1.3475046769326209, 1.1380312577212153, 0.8475520859160499, 0.8232176358285205, 0.44462436678120765, 0.7257688085873244, 0.39344083194232093, 0.062038502713543636], θ = [1.3466011710339703, 0.0376255364135533, 0.003570577402944517, 0.34777468288447505, 0.09532242923198266, 0.693645547569269, 1.4267785122258858, -0.0277671202334314, -0.2419753622369778, 0.5028183338062486, -0.3152167433396314, 0.831388178133279]), (objective = 40387.830926683644, σ = 5.0669107095427055, β = ((Intercept) = 2.3449323663767223, age: old = -1.669282177401879, frequency: low = 1.5321542492740932, context: unmatched = -1.868705054365737, age: old & frequency: low = 0.4920182881212734, age: old & context: unmatched = -0.5723301717676685, frequency: low & context: unmatched = 0.9305139607120997, age: old & frequency: low & context: unmatched = -0.10765078579583916), se = [1.3995536568114044, 1.2170540138237533, 0.8160811972347296, 0.7668825645047568, 0.43411813455443277, 0.6252873389007002, 0.44848200529485066, 0.06333638386928382], θ = [1.288342426918208, -0.08139129777312587, -0.012002516384396924, 0.40644509153651337, 0.06626862237996227, 0.780841363911027, 1.4884933628839911, 0.1431103990439471, 0.046403869499552566, 0.4254184349549342, 0.03293403656978233, 0.7743861874339877]), (objective = 40129.37192878041, σ = 4.9601436938629755, β = ((Intercept) = 0.35102886386257515, age: old = -3.1126861700706376, frequency: low = 2.201605203321825, context: unmatched = -0.3265068623468301, age: old & frequency: low = 0.0010668170622991225, age: old & context: unmatched = -0.8500950946170217, frequency: low & context: unmatched = 1.4940348325319006, age: old & frequency: low & context: unmatched = 0.029642551622900887), se = [1.608423304325992, 1.406102896583399, 0.9157979550120309, 0.801825471906974, 0.4783155070479197, 0.7019017076989357, 0.39255866167567566, 0.062001796173287195], θ = [1.4500663511027247, -0.02080841495927896, -0.017471626463679696, 0.3450915067665163, -0.06634261106272989, 0.6956130203500185, 1.7743821672008246, 0.15987672746487913, -0.20132372330871764, 0.5295209384620321, 0.2220872152609007, 0.8395715906135471]), (objective = 40068.51334759965, σ = 4.967692014102738, β = ((Intercept) = 3.813398857200134, age: old = -0.18991731215750401, frequency: low = 3.83652857813078, context: unmatched = -0.6362405163702498, age: old & frequency: low = 0.7633862303340991, age: old & context: unmatched = -1.4975619739090549, frequency: low & context: unmatched = 1.4152386154087642, age: old & frequency: low & context: unmatched = 0.01324980056836753), se = [1.419956699197643, 1.240271938597472, 0.7908293013609314, 0.6480369636589242, 0.3839380652744086, 0.5312215225323715, 0.3763130772819809, 0.062096150176284226], θ = [1.286762690968881, 0.045410602016334095, -0.04330005221207892, 0.3226875764023277, -0.058274611751774216, 0.6643038210097459, 1.5601322090952847, -0.020754275607275073, 0.03252592180651187, 0.42326805085843877, -0.026885476375149688, 0.6703546642603438]), (objective = 40521.07801663625, σ = 5.1422256976259515, β = ((Intercept) = 1.7647791150697372, age: old = -1.1438873763291777, frequency: low = 2.4591096960549623, context: unmatched = -2.052479961517255, age: old & frequency: low = 0.6282670328265395, age: old & context: unmatched = -1.8590186710804875, frequency: low & context: unmatched = 2.142001241301375, age: old & frequency: low & context: unmatched = -0.011610630595578653), se = [1.3213444569151551, 1.167020007276592, 0.7319470441160765, 0.7893404027897121, 0.3895266344597332, 0.6667882623562345, 0.4272976982862567, 0.0642778212203244], θ = [1.139699056557817, 0.028345065215859083, -0.019122762645094098, 0.36918413302958436, 0.044611987840668346, 0.733169781289181, 1.409050483564325, 0.05357643805831229, -0.29499831489678674, 0.3896863635722312, -0.07946437631211684, 0.7569513878587282]), (objective = 40327.508949122996, σ = 5.064450041521953, β = ((Intercept) = 0.38441260941249195, age: old = -2.3522764221398123, frequency: low = 2.257616659969222, context: unmatched = -2.494841669098911, age: old & frequency: low = -0.15430399826230215, age: old & context: unmatched = -0.7243257822422107, frequency: low & context: unmatched = 0.991206333069214, age: old & frequency: low & context: unmatched = 0.07881787371977923), se = [1.1834451894932145, 0.96907952454716, 0.7950506213029304, 0.7481385466765316, 0.4131320597528232, 0.6166905689470021, 0.42826583965762577, 0.06330562551902441], θ = [1.2518766541772395, -0.01593451888354515, -0.1529015191803156, 0.3573737219554262, 0.06166953697433065, 0.7296520245812996, 1.18082552334913, 0.012923356599381621, -0.20069305362491796, 0.44246659297441454, -0.1492646861241533, 0.7240823275663995]), (objective = 40195.50050374671, σ = 5.017569748497208, β = ((Intercept) = -0.9419887762199102, age: old = -1.915755937979817, frequency: low = 4.0192759407182015, context: unmatched = -0.035210139113138465, age: old & frequency: low = 0.5068280939465525, age: old & context: unmatched = -1.3640010050398665, frequency: low & context: unmatched = 1.907119943780961, age: old & frequency: low & context: unmatched = 0.08520908798123412), se = [1.2278919025664783, 1.0333587488922342, 0.7914338247303105, 0.7287412921216617, 0.431832465720685, 0.6341588174311766, 0.3644725175965394, 0.0627196218562151], θ = [1.2372347552532037, -0.01703741962596983, -0.17970819124195986, 0.3642226353890782, -0.08145103692529616, 0.6088411770081393, 1.274307621309965, -0.10771716049013025, 0.1796447525716321, 0.46038436852966963, 0.07550453274353051, 0.7711875823973475]), (objective = 40052.15374912469, σ = 4.964586345235049, β = ((Intercept) = 0.5613955381679298, age: old = -0.9105894983625396, frequency: low = 2.436326451692267, context: unmatched = -1.4867157520009648, age: old & frequency: low = 0.29020457179096476, age: old & context: unmatched = -1.4412376617439726, frequency: low & context: unmatched = 1.4912701187493083, age: old & frequency: low & context: unmatched = 0.0713097513561998), se = [1.3121052129574677, 1.1376826387312584, 0.7407155936766693, 0.6360860868257907, 0.3483694695537347, 0.5101106053531039, 0.385024404786732, 0.062057329315438114], θ = [1.2368112608496566, 0.031217841775872945, 0.05269619892329906, 0.37655812346137174, -0.031234088202927175, 0.6818501034379498, 1.4222973209652519, -0.04018250824376225, 0.06596427068524843, 0.3430863768185162, -0.134785628451393, 0.6273220759780028]), (objective = 40363.50619617878, σ = 5.058796667729237, β = ((Intercept) = 0.002264003189656984, age: old = -2.0708710441170344, frequency: low = 1.2263276427679304, context: unmatched = -0.8961769799399376, age: old & frequency: low = 0.6971190204888861, age: old & context: unmatched = -0.31871086161256906, frequency: low & context: unmatched = 2.092878092414033, age: old & frequency: low & context: unmatched = -0.02725510131242307), se = [1.4406821117124793, 1.2460255735146861, 0.8558415995176709, 0.8134476388061065, 0.4576894430829906, 0.6666450103705653, 0.4704042423112601, 0.06323495834661547], θ = [1.3245707005314615, -0.06098075190628558, 0.08591554540748317, 0.3404297046808223, -0.005716486150623569, 0.8196450762399494, 1.5364450158593106, 0.24239961455042983, 0.10750646898263579, 0.45012055580513466, 0.1765134450361844, 0.8035343178263625]), (objective = 40352.9172706416, σ = 5.069466416021452, β = ((Intercept) = 0.14618013717573924, age: old = 0.014007916145419172, frequency: low = 1.5614822145160177, context: unmatched = 0.500097033051294, age: old & frequency: low = 0.07916639522508465, age: old & context: unmatched = -1.307632806509102, frequency: low & context: unmatched = 1.5503568144569446, age: old & frequency: low & context: unmatched = -0.04825743311697761), se = [1.3891119975940656, 1.2166802061505015, 0.8081397178700372, 0.7382023211052826, 0.4514071175517594, 0.6309919014007197, 0.3883393265230018, 0.06336833020026815], θ = [1.2340763711504343, -0.017407531084838647, 0.06532169392507246, 0.3520622181985827, 0.06194049197408412, 0.6699591413162301, 1.4952109864300878, 0.14075589912958164, 0.04465764505175086, 0.4785059022931973, 0.05921689128083644, 0.779712310856859]), (objective = 40239.96502508431, σ = 5.022025504328468, β = ((Intercept) = -0.20636296964775225, age: old = -2.0911493727869037, frequency: low = 1.819330273968871, context: unmatched = -1.2598168770168416, age: old & frequency: low = 0.19363133423597612, age: old & context: unmatched = -0.5923372828396185, frequency: low & context: unmatched = 1.7626785828809246, age: old & frequency: low & context: unmatched = -0.024342188739193048), se = [1.394552087658454, 1.2019464966751663, 0.820935760851337, 0.7897256119158993, 0.416935701420731, 0.6733822004755793, 0.41732924034789376, 0.06277531880410585], θ = [1.3023758945377284, -0.03406158332119959, -0.16726143921143927, 0.32971024464551935, -0.010178578499846165, 0.7154480768301148, 1.4933403017089668, -0.02990312603601277, 0.10890829823489043, 0.46219427588405315, 0.2902506927981767, 0.7853682958029525]), (objective = 40187.567793635666, σ = 4.981173157430227, β = ((Intercept) = -1.1624453910693566, age: old = -1.3422147173533197, frequency: low = 2.083957103408143, context: unmatched = -0.6215757087665388, age: old & frequency: low = 0.17175117822591884, age: old & context: unmatched = -0.8010893303336478, frequency: low & context: unmatched = 1.5582579218646695, age: old & frequency: low & context: unmatched = -0.06538092602603829), se = [1.3889054044447813, 1.1556197828819088, 0.8825074471506774, 0.7507514104513171, 0.430369904112398, 0.5632877248717311, 0.5002114630281158, 0.062264664467877845], θ = [1.4321835295647278, -0.012039062740323979, 0.08265324743453817, 0.37026350165437477, 0.09557083538464473, 0.8821991811840227, 1.441542172892155, -0.03831427699665369, -0.008975895336875683, 0.47143994290298913, -0.18133572655262342, 0.6872412300199247]), (objective = 40425.436556324596, σ = 5.093662683371879, β = ((Intercept) = -1.796637112205905, age: old = -2.4326073181265477, frequency: low = 1.1918804094706645, context: unmatched = -1.837735150592248, age: old & frequency: low = 0.5170195033029257, age: old & context: unmatched = -1.0168219596162196, frequency: low & context: unmatched = 1.462147688758383, age: old & frequency: low & context: unmatched = -0.007559495292853516), se = [1.3413186699656359, 1.1693345750686852, 0.8164871197584805, 0.7381619035769811, 0.4846223176691234, 0.6277618228388212, 0.3935124627702204, 0.06367078354214849], θ = [1.2162122271255167, 0.03933913357934114, -0.023833341888347306, 0.382408937306542, 0.13126248607037777, 0.6687097279217393, 1.424040018871655, -0.09263249027595273, -0.0631130479859689, 0.5228371972415, -0.034241067420944515, 0.7721105175195111])], LinearAlgebra.LowerTriangular{Float64, Matrix{Float64}}[[1.3 0.0 0.0; 0.0 0.35 0.0; 0.0 0.0 0.75], [1.5 0.0 0.0; 0.0 0.5 0.0; 0.0 0.0 0.75]], [[1, 2, 3, 5, 6, 9], [1, 2, 3, 5, 6, 9]], [0.0, -Inf, -Inf, 0.0, -Inf, 0.0, 0.0, -Inf, -Inf, 0.0, -Inf, 0.0], (item = ("(Intercept)", "age: old", "context: unmatched"), subj = ("(Intercept)", "frequency: low", "context: unmatched")))
See your power and profit!
Finally, we can turn this into a power table:
ptbl = power_table(sim)
8-element Vector{NamedTuple{(:coefname, :power), Tuple{String, Float64}}}:
(coefname = "age: old & frequency: low", power = 0.1)
(coefname = "frequency: low & context: unmatched", power = 1.0)
(coefname = "age: old", power = 0.15)
(coefname = "context: unmatched", power = 0.5)
(coefname = "(Intercept)", power = 0.15)
(coefname = "frequency: low", power = 0.7)
(coefname = "age: old & context: unmatched", power = 0.35)
(coefname = "age: old & frequency: low & context: unmatched", power = 0.0)
For nicely displaying it, we can again use pretty_table
:
pretty_table(ptbl)
┌────────────────────────────────────────────────┬─────────┐
│ coefname │ power │
│ String │ Float64 │
├────────────────────────────────────────────────┼─────────┤
│ age: old & frequency: low │ 0.1 │
│ frequency: low & context: unmatched │ 1.0 │
│ age: old │ 0.15 │
│ context: unmatched │ 0.5 │
│ (Intercept) │ 0.15 │
│ frequency: low │ 0.7 │
│ age: old & context: unmatched │ 0.35 │
│ age: old & frequency: low & context: unmatched │ 0.0 │
└────────────────────────────────────────────────┴─────────┘
We can of course convert the row table into a DataFrame:
DataFrame(ptbl)
8 rows × 2 columns
coefname | power | |
---|---|---|
String | Float64 | |
1 | age: old & frequency: low | 0.1 |
2 | frequency: low & context: unmatched | 1.0 |
3 | age: old | 0.15 |
4 | context: unmatched | 0.5 |
5 | (Intercept) | 0.15 |
6 | frequency: low | 0.7 |
7 | age: old & context: unmatched | 0.35 |
8 | age: old & frequency: low & context: unmatched | 0.0 |
- choleskyTechnically, we're creating the lower Cholesky factor of these matrices, which is a bit like the matrix square root. In other words, we're creating the matrix form of standard deviations instead of the matrix form of the variances.