Some resources:
- https://bookdown.org/content/4857/adventures-in-covariance.html#summary-bonus-multilevel-growth-models-and-the-melsm
- Williams, D. R., Liu, S., Martin, S. R., & Rast, P. (2019). Bayesian Multivariate Mixed-Effects Location Scale Modeling of Longitudinal Relations among Affective Traits, States, and Physical Activity. https://doi.org/10.31234/osf.io/4kfjp
library(ggplot2)
library(brms)
library(rstan)
library(dplyr)
library(stringr)
library(tidyverse)
library(tidybayes)
library(modelr)
# some STAN settings
options(mc.cores = parallel::detectCores())
rstan_options(auto_write = TRUE)## X alias craving signal day_of_week day_index day_index01 stressor
## 1 1 AEMR0604 7 0 Wednesday 1 0 0
## 2 2 AEMR0604 74 1 Wednesday 1 0 0
## 3 3 AEMR0604 66 2 Wednesday 1 0 0
## 4 4 AEMR0604 13 3 Wednesday 1 0 1
## 5 5 AEMR0604 10 4 Wednesday 1 0 1
## 6 6 AEMR0604 70 5 Wednesday 1 0 0
## stressor_prev stress_gmc stress_smc mood_gmc mood_smc mood alias_fct
## 1 NA -7.832119 6.66242 11.479576 17.519108 74.83333 AEMR0604
## 2 0 -9.832119 4.66242 22.812909 28.852442 86.16667 AEMR0604
## 3 0 -3.832119 10.66242 20.479576 26.519108 83.83333 AEMR0604
## 4 0 56.167881 70.66242 -10.520424 -4.480892 52.83333 AEMR0604
## 5 1 31.167881 45.66242 16.979576 23.019108 80.33333 AEMR0604
## 6 1 1.167881 15.66242 -7.520424 -1.480892 55.83333 AEMR0604
## id
## 1 1
## 2 1
## 3 1
## 4 1
## 5 1
## 6 1## [1] 64## [1] 10920data exploration
plotting craving fluctuations for participants
#single participants
ggplot(d2,aes(x=signal,y=craving,group=day_index,colour=day_index))+geom_line()+ scale_colour_gradient2(low = "blue", mid = "red", high = "yellow", midpoint = 14)+theme_minimal()
ggplot(d1,aes(x=signal,y=craving,group=day_index,colour=day_index))+geom_line()+ scale_colour_gradient2(low = "blue", mid = "red", high = "yellow", midpoint = 14)+theme_minimal()
typical Multi-level Model
for comparison
starting with measurement repetition as only factor and craving as outcome
i.e., fixed intercept overall, fixed effect of time; varying intercept, varying effect of time.
MLM1 = brm(
bf(craving ~ 1 + signal + (1 + signal | alias)),
family = gaussian(),
prior = c(prior(normal(50,10), class = Intercept),
prior(normal(0,5), class = b),
prior(cauchy(0,0.5), class = sigma),
prior(cauchy(0,0.5), class = sd),
prior(lkj_corr_cholesky(1), class = cor)),
data=d,
iter=2000,
file="MLM1",file_refit="on_change")
summary(MLM1,prob=0.9)## Family: gaussian
## Links: mu = identity; sigma = identity
## Formula: craving ~ 1 + signal + (1 + signal | alias)
## Data: d (Number of observations: 9882)
## Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
## total post-warmup draws = 4000
##
## Multilevel Hyperparameters:
## ~alias (Number of levels: 64)
## Estimate Est.Error l-90% CI u-90% CI Rhat Bulk_ESS
## sd(Intercept) 17.80 1.69 15.27 20.67 1.00 1026
## sd(signal) 2.64 0.28 2.22 3.13 1.00 1510
## cor(Intercept,signal) -0.46 0.10 -0.62 -0.28 1.00 2110
## Tail_ESS
## sd(Intercept) 1445
## sd(signal) 2384
## cor(Intercept,signal) 2830
##
## Regression Coefficients:
## Estimate Est.Error l-90% CI u-90% CI Rhat Bulk_ESS Tail_ESS
## Intercept 30.19 2.26 26.46 34.05 1.00 979 1840
## signal -0.07 0.37 -0.66 0.54 1.00 2033 2727
##
## Further Distributional Parameters:
## Estimate Est.Error l-90% CI u-90% CI Rhat Bulk_ESS Tail_ESS
## sigma 24.26 0.17 23.97 24.54 1.00 5850 2557
##
## Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).nd <-
d %>%
filter(alias %in% c("SLHR1002","UERL0203")) %>%
select(alias,id, craving, signal)
MLM1_preds = bind_cols(
# fitted
fitted(MLM1, newdata = nd,probs = c(0.05, 0.95)) %>%
data.frame(),
# predict
predict(MLM1, newdata = nd,probs = c(0.05, 0.95)) %>%
data.frame() %>%
select(Q5:Q95) %>%
setNames(c("p_lower", "p_upper"))) %>%
bind_cols(nd)
ggplot(data=MLM1_preds,aes(x = signal)) +
geom_smooth(aes(y = Estimate, ymin = Q5, ymax = Q95),
stat = "identity",
fill = "#115740", color = "#115740", alpha = 1/2, linewidth = 1/2) +
geom_ribbon(aes(ymin = p_lower, ymax = p_upper),
fill = "#115740", alpha = 1/2) +
geom_point(aes(y = craving),
color = "#115740",alpha=0.6,position=position_jitter(width=0.1)) +
ylab("Craving") + xlab("Time of day")+coord_cartesian(ylim = c(-30, 120)) +
scale_y_continuous(breaks = seq(0, 100, 20))+
facet_wrap(~ id)+theme_minimal(base_size=32)
both participants have the same 90% credible interval since the estimated variance for both is the same!
Mixed-Effects Location Scale Model
Note:
- Link_sigma = “log” is default in family, just added here for
readability.
- Using |i| allows for covariance between random intercept for sigma and
random intercept for mean + random slope for day_index for mean.
bMELS1 = brm(
formula=bf(
craving ~ 1 + signal + (1 + signal |i| alias),
sigma ~ 1 + signal + (1 + signal |i| alias)),
prior = c( #for mean intercept, factors and varying effects (expressed as variation from mean)
prior(normal(20, 10), class = Intercept),
prior(normal(0, 5), class = b),
prior(exponential(2), class = sd), # rate=1 might be a bit strict
# for sigma-specific parameters
prior(normal(0, 10), class = Intercept, dpar = sigma),
prior(exponential(2), class = sd, dpar = sigma),
#prior on varying effects correlations
prior(lkj(2), class = cor)),
family = brmsfamily("gaussian",link_sigma = "log"),
data=d, iter=2000,
file="bMELS1",file_refit = "on_change")
bMELS1## Family: gaussian
## Links: mu = identity; sigma = log
## Formula: craving ~ 1 + signal + (1 + signal | i | alias)
## sigma ~ 1 + signal + (1 + signal | i | alias)
## Data: d (Number of observations: 9882)
## Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
## total post-warmup draws = 4000
##
## Multilevel Hyperparameters:
## ~alias (Number of levels: 64)
## Estimate Est.Error l-95% CI u-95% CI Rhat
## sd(Intercept) 14.77 1.06 12.89 17.03 1.00
## sd(signal) 2.29 0.24 1.86 2.78 1.00
## sd(sigma_Intercept) 0.41 0.04 0.34 0.50 1.00
## sd(sigma_signal) 0.07 0.01 0.06 0.09 1.00
## cor(Intercept,signal) -0.38 0.09 -0.56 -0.19 1.00
## cor(Intercept,sigma_Intercept) 0.16 0.10 -0.04 0.36 1.00
## cor(signal,sigma_Intercept) -0.09 0.12 -0.32 0.15 1.00
## cor(Intercept,sigma_signal) -0.09 0.11 -0.31 0.13 1.00
## cor(signal,sigma_signal) 0.55 0.10 0.34 0.73 1.00
## cor(sigma_Intercept,sigma_signal) -0.48 0.11 -0.66 -0.24 1.00
## Bulk_ESS Tail_ESS
## sd(Intercept) 1016 1441
## sd(signal) 2138 3042
## sd(sigma_Intercept) 1538 2229
## sd(sigma_signal) 1832 2829
## cor(Intercept,signal) 1369 2068
## cor(Intercept,sigma_Intercept) 1030 1772
## cor(signal,sigma_Intercept) 879 1637
## cor(Intercept,sigma_signal) 1676 2748
## cor(signal,sigma_signal) 1887 2309
## cor(sigma_Intercept,sigma_signal) 2147 2717
##
## Regression Coefficients:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## Intercept 28.93 1.93 25.24 32.73 1.01 462 893
## sigma_Intercept 3.01 0.05 2.90 3.11 1.00 784 1105
## signal -0.00 0.32 -0.62 0.65 1.00 1224 1878
## sigma_signal 0.03 0.01 0.01 0.05 1.00 1473 2575
##
## Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).nd <-
d %>%
filter(alias %in% c("SLHR1002","UERL0203")) %>%
select(alias,id, craving, signal)
bMELS1_preds = bind_cols(
# fitted
fitted(bMELS1, newdata = nd,probs = c(0.05, 0.95)) %>%
data.frame(),
# predict
predict(bMELS1, newdata = nd,probs = c(0.05, 0.95)) %>%
data.frame() %>%
select(Q5:Q95) %>%
setNames(c("p_lower", "p_upper"))) %>%
bind_cols(nd)
ggplot(data=bMELS1_preds,aes(x = signal)) +
geom_smooth(aes(y = Estimate, ymin = Q5, ymax = Q95),
stat = "identity",
fill = "#115740", color = "#115740", alpha = 1/2, linewidth = 1/2) +
geom_ribbon(aes(ymin = p_lower, ymax = p_upper),
fill = "#115740", alpha = 1/2) +
geom_point(aes(y = craving),
color = "#115740",alpha=0.6,position=position_jitter(width=0.1)) +
ylab("Craving") + xlab("Time of day")+ coord_cartesian(ylim = c(-30, 120)) +
scale_y_continuous(breaks = seq(0, 100, 20))+
facet_wrap(~ id)+theme_minimal(base_size=32)
-> the credible interval for the second participant is more narrow and closer to the actual data, the predictive accuracy is improved.
adding stressor presence to model for both
bMELS2 = brm(
formula=bf(
craving ~ 1 + signal + stressor + (1 + stressor + signal |i| alias),
sigma ~ 1 + signal + stressor + (1 + stressor + signal |i| alias)),
prior = c( #for mean intercept, factors and varying effects (expressed as variation from mean)
prior(normal(20, 10), class = Intercept),
prior(normal(0, 5), class = b),
prior(exponential(2), class = sd), # rate=1 might be a bit strict
# for sigma-specific parameters
prior(normal(0, 10), class = Intercept, dpar = sigma),
prior(exponential(2), class = sd, dpar = sigma),
#prior on varying effects correlations
prior(lkj(2), class = cor)),
family = brmsfamily("gaussian",link_sigma = "log"),
data=d, iter=4000,chains=4,cores=4,
file="bMELS2",file_refit = "on_change")
summary(bMELS2,prob=0.9)## Family: gaussian
## Links: mu = identity; sigma = log
## Formula: craving ~ 1 + signal + stressor + (1 + stressor + signal | i | alias)
## sigma ~ 1 + signal + stressor + (1 + stressor + signal | i | alias)
## Data: d (Number of observations: 9648)
## Draws: 4 chains, each with iter = 4000; warmup = 2000; thin = 1;
## total post-warmup draws = 8000
##
## Multilevel Hyperparameters:
## ~alias (Number of levels: 64)
## Estimate Est.Error l-90% CI u-90% CI Rhat
## sd(Intercept) 14.87 1.05 13.22 16.67 1.00
## sd(stressor) 2.99 0.81 1.66 4.30 1.00
## sd(signal) 2.24 0.24 1.87 2.64 1.00
## sd(sigma_Intercept) 0.41 0.04 0.35 0.48 1.00
## sd(sigma_stressor) 0.15 0.03 0.11 0.20 1.00
## sd(sigma_signal) 0.07 0.01 0.06 0.08 1.00
## cor(Intercept,stressor) -0.14 0.16 -0.40 0.13 1.00
## cor(Intercept,signal) -0.36 0.09 -0.50 -0.20 1.00
## cor(stressor,signal) -0.20 0.20 -0.52 0.13 1.00
## cor(Intercept,sigma_Intercept) 0.16 0.10 -0.01 0.32 1.00
## cor(stressor,sigma_Intercept) -0.10 0.17 -0.36 0.19 1.01
## cor(signal,sigma_Intercept) -0.03 0.13 -0.24 0.18 1.00
## cor(Intercept,sigma_stressor) -0.12 0.14 -0.35 0.12 1.00
## cor(stressor,sigma_stressor) 0.66 0.16 0.38 0.87 1.00
## cor(signal,sigma_stressor) -0.14 0.16 -0.40 0.13 1.00
## cor(sigma_Intercept,sigma_stressor) -0.09 0.16 -0.35 0.19 1.00
## cor(Intercept,sigma_signal) -0.07 0.11 -0.25 0.11 1.00
## cor(stressor,sigma_signal) -0.16 0.19 -0.47 0.14 1.00
## cor(signal,sigma_signal) 0.52 0.10 0.34 0.68 1.00
## cor(sigma_Intercept,sigma_signal) -0.40 0.12 -0.58 -0.20 1.00
## cor(sigma_stressor,sigma_signal) -0.39 0.17 -0.66 -0.09 1.00
## Bulk_ESS Tail_ESS
## sd(Intercept) 1539 2662
## sd(stressor) 1888 1254
## sd(signal) 3961 5576
## sd(sigma_Intercept) 2424 4027
## sd(sigma_stressor) 4820 5354
## sd(sigma_signal) 4106 5880
## cor(Intercept,stressor) 6514 3709
## cor(Intercept,signal) 2525 4196
## cor(stressor,signal) 700 1608
## cor(Intercept,sigma_Intercept) 1965 3282
## cor(stressor,sigma_Intercept) 534 1035
## cor(signal,sigma_Intercept) 1451 2513
## cor(Intercept,sigma_stressor) 6343 6077
## cor(stressor,sigma_stressor) 1884 2095
## cor(signal,sigma_stressor) 6850 5555
## cor(sigma_Intercept,sigma_stressor) 6354 6195
## cor(Intercept,sigma_signal) 3351 4482
## cor(stressor,sigma_signal) 736 1499
## cor(signal,sigma_signal) 3326 4851
## cor(sigma_Intercept,sigma_signal) 4007 5495
## cor(sigma_stressor,sigma_signal) 2797 5178
##
## Regression Coefficients:
## Estimate Est.Error l-90% CI u-90% CI Rhat Bulk_ESS Tail_ESS
## Intercept 29.04 1.92 25.88 32.10 1.01 657 1301
## sigma_Intercept 3.00 0.05 2.91 3.09 1.00 1449 2902
## signal 0.01 0.32 -0.51 0.54 1.00 1932 3322
## stressor -0.44 0.72 -1.64 0.74 1.00 3289 5286
## sigma_signal 0.03 0.01 0.01 0.05 1.00 3103 4480
## sigma_stressor 0.02 0.03 -0.02 0.07 1.00 5408 5708
##
## Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).fixed effects for stressor and time of day
with dpar=TRUE, we receive a column for sigma and for mu
bMELS2_epred <- d %>%
data_grid(signal,stressor)%>% add_epred_rvars(bMELS2,dpar=TRUE,allow_new_levels=TRUE)ggplot(bMELS2_epred, aes(xdist = sigma, y = signal, fill = factor(stressor),colour=factor(stressor))) +
stat_halfeye(alpha=0.5) +
scale_y_reverse(breaks=seq(0,5,1),labels=c("0" = "08:30 AM","1" = "11:00 AM","2" = "01:30 PM","3" = "04:00 PM","4" = "06:30 PM","5" = "09:30 PM"))+
scale_fill_manual(
values = c("0" = "#9C9C9C", "1" = "#115740"),
labels = c("0" = "no stressor", "1" = "stressor"),
name = "Stressor") +
scale_colour_manual(
values = c("0" = "#9C9C9C", "1" = "#115740"),
labels = c("0" = "no stressor", "1" = "stressor"),
name = "Stressor") +
xlim(0,60)+
ylab("Time-of-day")+xlab(expression(sigma))+
theme_minimal()
ggplot(bMELS2_epred, aes(xdist = mu, y = signal, fill = factor(stressor),colour=factor(stressor))) +
#stat_dotsinterval(quantiles = 100, alpha = 0.6) +
stat_halfeye(alpha=0.5) +
scale_y_reverse(breaks=seq(0,5,1),labels=c("0" = "08:30 AM","1" = "11:00 AM","2" = "01:30 PM","3" = "04:00 PM","4" = "06:30 PM","5" = "09:30 PM"))+
scale_fill_manual(
values = c("0" = "#9C9C9C", "1" = "#115740"),
labels = c("0" = "no stressor", "1" = "stressor"),
name = "Stressor") +
scale_colour_manual(
values = c("0" = "#9C9C9C", "1" = "#115740"),
labels = c("0" = "no stressor", "1" = "stressor"),
name = "Stressor") +
#xlim(0,60)+
ylab("Time-of-day")+xlab(expression(mu))+
theme_minimal()
covariance matrix
rho <-
posterior_summary(bMELS2) %>%
data.frame() %>%
rownames_to_column("param") %>%
filter(str_detect(param, "cor_")) %>%
mutate(param = str_remove(param, "cor_alias__")) %>%
separate(param, into = c("left", "right"), sep = "__") %>%
mutate(label = formatC(Estimate, digits = 2, format = "f") %>% str_replace(., "0.", "."))
rho_ord = rho
rho_ord$left=factor(rho_ord$left,
levels=c("Intercept","stressor","signal","sigma_Intercept","sigma_stressor"),ordered=T,
labels=c("alpha","S[mu]","ToD[mu]","alpha[sigma]","S[sigma]"))
rho_ord$right=factor(rho_ord$right,
levels=c("sigma_signal","sigma_stressor","sigma_Intercept","signal","stressor"),ordered=T,
labels=c("ToD[sigma]","S[sigma]","alpha[sigma]","ToD[mu]","S[mu]"))plotting
rho_ord %>%
ggplot(aes(x = left, y = right)) +
geom_tile(aes(fill = Estimate)) +
geom_text(aes(label = label),size=5,
family = "Franklin Gothic Book") +
scale_fill_gradient2(expression(rho),
low = "#59708b", mid = "#FCF9F0", high = "#A65141", midpoint = 0,
labels = c(-1, "", 0, "", 1), limits = c(-1, 1)) +
scale_x_discrete(NULL, drop = F, labels = ggplot2:::parse_safe, position = "top") +
scale_y_discrete(NULL, drop = F, labels = ggplot2:::parse_safe) +
theme(axis.text = element_text(),
axis.ticks = element_blank(),
legend.text = element_text(hjust = 1))+theme_minimal()
Bonus: with interaction between stressor and signal
bMELS2b = brm(
formula=bf(
craving ~ 1 + signal * stressor + (1 + stressor + signal |i| alias),
sigma ~ 1 + signal * stressor + (1 + stressor + signal |i| alias)),
prior = c( #for mean intercept, factors and varying effects (expressed as variation from mean)
prior(normal(20, 10), class = Intercept),
prior(normal(0, 5), class = b),
prior(exponential(2), class = sd), # rate=1 might be a bit strict
# for sigma-specific parameters
prior(normal(0, 10), class = Intercept, dpar = sigma),
prior(exponential(2), class = sd, dpar = sigma),
#prior on varying effects correlations
prior(lkj(2), class = cor)),
family = brmsfamily("gaussian",link_sigma = "log"),
data=d, iter=2000,chains=4,cores=4,
file="bMELS2b",file_refit = "on_change")
summary(bMELS2b,prob=0.9)## Family: gaussian
## Links: mu = identity; sigma = log
## Formula: craving ~ 1 + signal * stressor + (1 + stressor + signal | i | alias)
## sigma ~ 1 + signal * stressor + (1 + stressor + signal | i | alias)
## Data: d (Number of observations: 9648)
## Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
## total post-warmup draws = 4000
##
## Multilevel Hyperparameters:
## ~alias (Number of levels: 64)
## Estimate Est.Error l-90% CI u-90% CI Rhat
## sd(Intercept) 14.92 1.06 13.29 16.76 1.00
## sd(stressor) 2.91 0.89 1.37 4.32 1.01
## sd(signal) 2.17 0.23 1.80 2.58 1.00
## sd(sigma_Intercept) 0.40 0.04 0.34 0.47 1.00
## sd(sigma_stressor) 0.16 0.03 0.11 0.20 1.00
## sd(sigma_signal) 0.07 0.01 0.05 0.08 1.00
## cor(Intercept,stressor) -0.13 0.16 -0.38 0.15 1.00
## cor(Intercept,signal) -0.37 0.09 -0.52 -0.21 1.00
## cor(stressor,signal) -0.19 0.20 -0.52 0.16 1.00
## cor(Intercept,sigma_Intercept) 0.15 0.10 -0.02 0.32 1.00
## cor(stressor,sigma_Intercept) -0.10 0.17 -0.37 0.18 1.02
## cor(signal,sigma_Intercept) -0.01 0.13 -0.21 0.20 1.00
## cor(Intercept,sigma_stressor) -0.11 0.14 -0.34 0.13 1.00
## cor(stressor,sigma_stressor) 0.63 0.18 0.34 0.85 1.01
## cor(signal,sigma_stressor) -0.11 0.17 -0.39 0.17 1.00
## cor(sigma_Intercept,sigma_stressor) -0.09 0.17 -0.36 0.19 1.00
## cor(Intercept,sigma_signal) -0.07 0.11 -0.26 0.11 1.00
## cor(stressor,sigma_signal) -0.12 0.21 -0.47 0.22 1.00
## cor(signal,sigma_signal) 0.50 0.11 0.31 0.67 1.00
## cor(sigma_Intercept,sigma_signal) -0.39 0.12 -0.57 -0.19 1.00
## cor(sigma_stressor,sigma_signal) -0.38 0.18 -0.65 -0.07 1.00
## Bulk_ESS Tail_ESS
## sd(Intercept) 656 1631
## sd(stressor) 302 139
## sd(signal) 1771 2374
## sd(sigma_Intercept) 1457 2401
## sd(sigma_stressor) 2075 2642
## sd(sigma_signal) 2066 2831
## cor(Intercept,stressor) 2742 1136
## cor(Intercept,signal) 1597 2646
## cor(stressor,signal) 405 571
## cor(Intercept,sigma_Intercept) 838 1783
## cor(stressor,sigma_Intercept) 235 451
## cor(signal,sigma_Intercept) 666 1502
## cor(Intercept,sigma_stressor) 3712 2862
## cor(stressor,sigma_stressor) 452 209
## cor(signal,sigma_stressor) 4004 3417
## cor(sigma_Intercept,sigma_stressor) 3469 3096
## cor(Intercept,sigma_signal) 1762 2855
## cor(stressor,sigma_signal) 380 333
## cor(signal,sigma_signal) 2009 2774
## cor(sigma_Intercept,sigma_signal) 1863 2661
## cor(sigma_stressor,sigma_signal) 1514 2073
##
## Regression Coefficients:
## Estimate Est.Error l-90% CI u-90% CI Rhat Bulk_ESS
## Intercept 29.56 1.92 26.39 32.68 1.01 340
## sigma_Intercept 3.01 0.05 2.92 3.10 1.00 749
## signal -0.22 0.32 -0.76 0.32 1.00 1114
## stressor -2.59 0.96 -4.16 -1.02 1.00 2422
## signal:stressor 0.97 0.31 0.45 1.48 1.00 2979
## sigma_signal 0.02 0.01 0.01 0.04 1.00 1718
## sigma_stressor -0.04 0.04 -0.11 0.02 1.00 2804
## sigma_signal:stressor 0.03 0.01 0.01 0.05 1.00 3862
## Tail_ESS
## Intercept 806
## sigma_Intercept 1553
## signal 2044
## stressor 2616
## signal:stressor 2871
## sigma_signal 2781
## sigma_stressor 2905
## sigma_signal:stressor 2812
##
## Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).interesting interaction for both mean and sigma
## R version 4.5.1 (2025-06-13 ucrt)
## Platform: x86_64-w64-mingw32/x64
## Running under: Windows 10 x64 (build 19045)
##
## Matrix products: default
## LAPACK version 3.12.1
##
## locale:
## [1] LC_COLLATE=German_Austria.utf8 LC_CTYPE=German_Austria.utf8
## [3] LC_MONETARY=German_Austria.utf8 LC_NUMERIC=C
## [5] LC_TIME=German_Austria.utf8
##
## time zone: Europe/Vienna
## tzcode source: internal
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] modelr_0.1.11 tidybayes_3.0.7 lubridate_1.9.4
## [4] forcats_1.0.0 purrr_1.1.0 readr_2.1.5
## [7] tidyr_1.3.1 tibble_3.3.0 tidyverse_2.0.0
## [10] stringr_1.5.1 dplyr_1.1.4 rstan_2.32.7
## [13] StanHeaders_2.32.10 brms_2.22.0 Rcpp_1.1.0
## [16] ggplot2_3.5.2
##
## loaded via a namespace (and not attached):
## [1] gtable_0.3.6 tensorA_0.36.2.1 xfun_0.52
## [4] bslib_0.9.0 QuickJSR_1.8.0 inline_0.3.21
## [7] lattice_0.22-7 tzdb_0.5.0 vctrs_0.6.5
## [10] tools_4.5.1 generics_0.1.4 stats4_4.5.1
## [13] parallel_4.5.1 pkgconfig_2.0.3 Matrix_1.7-3
## [16] checkmate_2.3.2 RColorBrewer_1.1-3 distributional_0.5.0
## [19] RcppParallel_5.1.10 lifecycle_1.0.4 compiler_4.5.1
## [22] farver_2.1.2 Brobdingnag_1.2-9 prettydoc_0.4.1
## [25] codetools_0.2-20 htmltools_0.5.8.1 sass_0.4.10
## [28] bayesplot_1.13.0 yaml_2.3.10 pillar_1.11.0
## [31] jquerylib_0.1.4 arrayhelpers_1.1-0 cachem_1.1.0
## [34] bridgesampling_1.1-2 abind_1.4-8 nlme_3.1-168
## [37] posterior_1.6.1 svUnit_1.0.6 tidyselect_1.2.1
## [40] digest_0.6.37 mvtnorm_1.3-3 stringi_1.8.7
## [43] reshape2_1.4.4 labeling_0.4.3 fastmap_1.2.0
## [46] grid_4.5.1 cli_3.6.5 magrittr_2.0.3
## [49] loo_2.8.0 pkgbuild_1.4.8 broom_1.0.9
## [52] withr_3.0.2 scales_1.4.0 backports_1.5.0
## [55] timechange_0.3.0 estimability_1.5.1 rmarkdown_2.29
## [58] matrixStats_1.5.0 emmeans_1.11.2 gridExtra_2.3
## [61] hms_1.1.3 coda_0.19-4.1 evaluate_1.0.4
## [64] knitr_1.50 ggdist_3.3.3 rstantools_2.4.0
## [67] rlang_1.1.6 glue_1.8.0 rstudioapi_0.17.1
## [70] jsonlite_2.0.0 plyr_1.8.9 R6_2.6.1