Longitudinal Analyses: 42 months, 48 months, 54 months, Longitudinal
library(brms)
library(readxl)
library(writexl)
library(rmarkdown)
library(tidybayes)
library(DHARMa)
library(ggplot2)
library(tidyr)
library(dplyr)
library(GGally)
library(ggpubr)
library(tidybayes)
library(patchwork)
library(sjPlot)
library(lubridate)
library(zoo)
library(yarrr)
DOD_54mo <- read_xlsx("DOD_54mo_blind.xlsx")
N_total <- nrow(DOD_54mo) # sanity check
DOD_42mo <- read_xlsx("DOD_42mo_prepped.xlsx") %>%
group_by(Blind_ID) %>%
slice(1) %>%
dplyr::select("Blind_ID", "Prod_DOD") %>%
rename("DOD_42" = "Prod_DOD") %>%
ungroup()
DOD_48mo <- read_xlsx("DOD_48mo_prepped.xlsx") %>%
group_by(Blind_ID) %>%
slice(1) %>%
dplyr::select("Blind_ID", "Prod_DOD") %>%
rename("DOD_48" = "Prod_DOD") %>%
ungroup()
Create mismatch variables. I’m also creating some character versions of many factors to make plotting a bit easier later.
DOD_54mo <- DOD_54mo %>%
mutate(
Prime = ifelse(Prime.Type==.5, yes="DOD", no="POD"), # Prime type, character version
Overlap = ifelse(Verb.match==.5, yes="Overlap", no="No Overlap") # overlap
)
DOD_54mo$bias <-"DOD" # set as DOD bias
DOD_54mo$bias[DOD_54mo$Prime.verb=="passed"] = "POD" # Change POD-biased prime verbs.
DOD_54mo$bias[DOD_54mo$Prime.verb=="sent"] = "POD"
DOD_54mo$bias[DOD_54mo$Prime.verb=="threw"] = "POD"
DOD_54mo$biasmatch <-"match" # set to match
DOD_54mo$biasmatch[DOD_54mo$bias=="POD" & DOD_54mo$Prime == "DOD"] = "mismatch" # Change to mismatch.
DOD_54mo$biasmatch[DOD_54mo$bias=="DOD" & DOD_54mo$Prime == "POD"] = "mismatch"
table(DOD_54mo$biasmatch)
##
## match mismatch
## 529 527
DOD_54mo$biasmismatchc <- ifelse(DOD_54mo$biasmatch == "mismatch", yes=.5, no=-.5)
There is one type of trial coded as an admin error that I don’t want to remove at the time being – trials 12r on lists 6 and 8.
DOD_54mo <- DOD_54mo %>%
mutate(
Trial_Drop = ifelse(is.na(RESPONSE.CODE), yes=1, no=0),
Trial_Problem = ifelse(List==6 & Item == "12r", yes = 1, no =ifelse(
List==8 & Item == "12r", yes = 1, no =0
))
)
nrow(filter(DOD_54mo, Trial_Problem==1)) == (nrow(filter(DOD_54mo, List==6))/12) + (nrow(filter(DOD_54mo, List==8))/12)
## [1] TRUE
N_miss <- DOD_54mo %>%
filter(is.na(RESPONSE.CODE)) %>%
nrow()
DOD_54mo %>%
summarise(
Verb_Error = sum(Verb.error, na.rm=TRUE),
Admin_Error = sum(Admin.error, na.rm=TRUE),
Prompts.disallowed = sum(Prompts.disallowed, na.rm=TRUE),
Prep_Error = sum(Prep.error.bin, na.rm=TRUE),
Reversal_Error = sum(Reversal.error, na.rm=TRUE),
NonTargetResponse = sum(NonTarget.response, na.rm=TRUE),
No_Response = sum(No.response, na.rm=TRUE),
GavePOD = sum(Gave.POD, na.rm=TRUE)
)
## # A tibble: 1 x 8
## Verb_Error Admin_Error Prompts.disal~1 Prep_~2 Rever~3 NonTa~4 No_Re~5 GavePOD
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <int> <dbl>
## 1 15 22 1 29 1 33 0 5
## # ... with abbreviated variable names 1: Prompts.disallowed, 2: Prep_Error,
## # 3: Reversal_Error, 4: NonTargetResponse, 5: No_Response
33 non-target responses (this includes lots of stuff, non-dative responses, preposition errors, etc).
DOD_54mo %>%
filter(is.na(NonTarget.response) & is.na(RESPONSE.CODE)) %>%
summarise(
Verb_Error = sum(Verb.error, na.rm=TRUE),
Admin_Error = sum(Admin.error, na.rm=TRUE),
Prompts.disallowed = sum(Prompts.disallowed, na.rm=TRUE),
Prep_Error = sum(Prep.error.bin, na.rm=TRUE),
Reversal_Error = sum(Reversal.error, na.rm=TRUE),
NonTargetResponse = sum(NonTarget.response, na.rm=TRUE),
No_Response = sum(No.response, na.rm=TRUE),
GavePOD = sum(Gave.POD, na.rm=TRUE)
)
## # A tibble: 1 x 8
## Verb_Error Admin_Error Prompts.disal~1 Prep_~2 Rever~3 NonTa~4 No_Re~5 GavePOD
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <int> <dbl>
## 1 10 5 0 0 0 0 0 0
## # ... with abbreviated variable names 1: Prompts.disallowed, 2: Prep_Error,
## # 3: Reversal_Error, 4: NonTargetResponse, 5: No_Response
DOD_54mo %>%
filter(is.na(RESPONSE.CODE) & is.na(NonTarget.response) & Admin.error==1)
## # A tibble: 5 x 41
## Blind_ID Zoom Gender Session.Age List Item Prime.~1 Targe~2 Prime~3 Verb.~4
## <chr> <dbl> <dbl> <chr> <dbl> <chr> <chr> <chr> <dbl> <dbl>
## 1 HcODjkHR 1 1 54 : 27 1 2 showed gave -0.5 -0.5
## 2 HcODjkHR 1 1 54 : 27 1 9 sent gave 0.5 -0.5
## 3 9hKLS9my 1 1 54 : 22 3 1 passed gave 0.5 -0.5
## 4 9hKLS9my 1 1 54 : 22 3 10 brought gave -0.5 -0.5
## 5 LkVZQCTn 1 1 54 : 15 7 6 gave gave -0.5 0.5
## # ... with 31 more variables: Verb.error <dbl>, Admin.error <dbl>,
## # COVID.Zoom.session <dbl>, Repeat.prime <dbl>, Repeat.target.prompt <dbl>,
## # Prime.interference <dbl>, Missing.article <dbl>, Pronoun <dbl>,
## # Prep.error <chr>, Prep.error.bin <dbl>, Reversal.error <dbl>,
## # Prompts.allowable <dbl>, Prompts.disallowed <dbl>,
## # NonTarget.response <dbl>, No.response <lgl>, Gave.POD <dbl>,
## # RESPONSE.CODE <dbl>, Prime.code <dbl>, Transcription <chr>, ...
10 verb errors and 5 Admin Errors. Some of those verb errors are also classified as admin errors. This is because the audio cut out on a number a few Zoom sessions and the child didn’t hear the target verb and
See how many missing values we have by participant
DOD_54mo %>%
group_by(Blind_ID) %>%
mutate(
N_Missing = sum(is.na(RESPONSE.CODE))
) %>%
slice(1) %>%
select("Blind_ID", "N_Missing") %>%
group_by(N_Missing) %>%
count()
## # A tibble: 5 x 2
## # Groups: N_Missing [5]
## N_Missing n
## <int> <int>
## 1 0 67
## 2 1 12
## 3 2 3
## 4 3 5
## 5 10 1
DOD_54mo <- DOD_54mo %>%
group_by(Blind_ID) %>%
mutate(
N_Missing = sum(is.na(RESPONSE.CODE))
) %>%
slice(1) %>%
select("Blind_ID", "N_Missing") %>%
mutate(
Drop = ifelse(N_Missing >= 10, 1, 0),
MightDrop = ifelse(N_Missing > 6, 1, 0)
) %>%
full_join(DOD_54mo, .)
## Joining, by = "Blind_ID"
N_drop = DOD_54mo %>%
filter(Drop == 1 & !is.na(RESPONSE.CODE)) %>%
count()
One participant is missing 10 trials.
And by item.
DOD_54mo %>%
group_by(Target.verb) %>%
mutate(
N_Missing = sum(is.na(RESPONSE.CODE))
) %>%
slice(1) %>%
select("Target.verb", "N_Missing") %>%
group_by(Target.verb)
## # A tibble: 7 x 2
## # Groups: Target.verb [7]
## Target.verb N_Missing
## <chr> <int>
## 1 bought 2
## 2 brought 9
## 3 gave 9
## 4 passed 1
## 5 sent 3
## 6 showed 6
## 7 threw 13
All non-target verbs have been removed. It doesn’t make sense to estimate random effects for a level of a random factor with only a single observation.
And let’s just check the proportion of missing trials across the various experimental conditions.
DOD_54mo %>%
mutate(
Missing = ifelse(is.na(RESPONSE.CODE), yes=1, no = 0)
) %>%
group_by(Prime, Overlap) %>%
summarise(
mean = mean(Missing),
sd = sd(Missing)
)
## `summarise()` has grouped output by 'Prime'. You can override using the
## `.groups` argument.
## # A tibble: 4 x 4
## # Groups: Prime [2]
## Prime Overlap mean sd
## <chr> <chr> <dbl> <dbl>
## 1 DOD No Overlap 0.0478 0.214
## 2 DOD Overlap 0.0288 0.167
## 3 POD No Overlap 0.0586 0.235
## 4 POD Overlap 0.0295 0.170
DOD_54mo %>%
mutate(
Missing = ifelse(is.na(RESPONSE.CODE), yes=1, no = 0)
) %>%
group_by(Prime, biasmatch) %>%
summarise(
mean = mean(Missing),
sd = sd(Missing)
)
## `summarise()` has grouped output by 'Prime'. You can override using the
## `.groups` argument.
## # A tibble: 4 x 4
## # Groups: Prime [2]
## Prime biasmatch mean sd
## <chr> <chr> <dbl> <dbl>
## 1 DOD match 0.0377 0.191
## 2 DOD mismatch 0.0379 0.191
## 3 POD match 0.0379 0.191
## 4 POD mismatch 0.0494 0.217
The differences in missingness between condition don’t look drastic. To sum up things, we have a fair amount of missing data, including 1 participants who are missing 10 or more trials. I’m goig to remove the participant who is missing 10 trials.
DOD_54mo <- filter(DOD_54mo, !is.na(RESPONSE.CODE) & Drop==0) # 1011
nrow(DOD_54mo) == N_total - N_miss - N_drop
## n
## [1,] TRUE
DOD_54mo %>%
summarise(
min = min(Session.Age),
max = max(Session.Age)
)
## # A tibble: 1 x 2
## min max
## <chr> <chr>
## 1 54 : 10 55 : 5
(p1 <- DOD_54mo %>%
group_by(Prime, Overlap) %>%
summarise(mean = mean(RESPONSE.CODE, na.rm=TRUE),
sd = sd(RESPONSE.CODE, na.rm=TRUE)) %>%
ggplot(aes(y=mean, x=Overlap, fill=Prime)) +
geom_bar(stat="identity",position="dodge") +
geom_errorbar(aes(ymin=mean, ymax=mean+sd), width=.2,
position=position_dodge(.9)) +
ylim(0, .6) +
theme_minimal() +
labs(y="Raw Mean") +
theme(legend.position = "none")
)
## `summarise()` has grouped output by 'Prime'. You can override using the
## `.groups` argument.
DOD_54mo %>%
group_by(Prime, Overlap, Target.verb) %>%
summarise(mean = mean(RESPONSE.CODE, na.rm=TRUE),
sd = sd(RESPONSE.CODE, na.rm=TRUE)) %>%
ggplot(aes(y=mean, x=Overlap, fill=Prime)) +
geom_bar(stat="identity",position="dodge") +
geom_errorbar(aes(ymin=mean, ymax=mean+sd), width=.2,
position=position_dodge(.9)) +
facet_wrap(~Target.verb)
## `summarise()` has grouped output by 'Prime', 'Overlap'. You can override using
## the `.groups` argument.
This is totally different than the other time points – it looks like there is a lexical boost for some verbs (showed, brought and threw, but not others (gave, passed or sent).
And let’s look at the pattern of responses when we only consider participants who produced at least 1 DOD.
DOD_54mo <- DOD_54mo %>%
group_by(Blind_ID) %>%
mutate(
Prod_DOD = sum(RESPONSE.CODE)
) %>%
dplyr::select(Blind_ID, Prod_DOD) %>%
slice(1) %>%
ungroup() %>%
full_join(DOD_54mo, ., by="Blind_ID") %>%
left_join(., DOD_42mo, by="Blind_ID") %>%
left_join(., DOD_48mo, by="Blind_ID")
write_xlsx(DOD_54mo, "DOD_54mo_prepped.xlsx")
DOD_54mo2 <- filter(DOD_54mo, !(Prod_DOD==0 & DOD_48==0 & DOD_42==0))
DOD_54mo2 %>%
group_by(Blind_ID) %>%
slice(1) %>%
ungroup() %>%
count()
## # A tibble: 1 x 1
## n
## <int>
## 1 75
(p1b <- DOD_54mo2 %>%
group_by(Prime, Overlap) %>%
summarise(mean = mean(RESPONSE.CODE, na.rm=TRUE),
sd = sd(RESPONSE.CODE, na.rm=TRUE)) %>%
ggplot(aes(y=mean, x=Overlap, fill=Prime)) +
geom_bar(stat="identity",position="dodge") +
geom_errorbar(aes(ymin=mean, ymax=mean+sd), width=.2,
position=position_dodge(.9)) +
ylim(0, .9) +
theme_minimal() +
labs(y="Raw Mean") +
theme(legend.position = "none")
)
## `summarise()` has grouped output by 'Prime'. You can override using the
## `.groups` argument.
DOD_54mo2 %>%
group_by(Prime, Overlap, Target.verb) %>%
summarise(mean = mean(RESPONSE.CODE, na.rm=TRUE),
sd = sd(RESPONSE.CODE, na.rm=TRUE)) %>%
ggplot(aes(y=mean, x=Overlap, fill=Prime)) +
geom_bar(stat="identity",position="dodge") +
geom_errorbar(aes(ymin=mean, ymax=mean+sd), width=.2,
position=position_dodge(.9)) +
facet_wrap(~Target.verb)
## `summarise()` has grouped output by 'Prime', 'Overlap'. You can override using
## the `.groups` argument.
Plot proprortion of DOD responses by each cell in two-way interaction between biasmatch and prime type.
(p1c <- DOD_54mo %>%
filter(Overlap == "No Overlap") %>%
group_by(Prime, biasmatch) %>%
summarise(mean = mean(RESPONSE.CODE, na.rm=TRUE),
sd = sd(RESPONSE.CODE, na.rm=TRUE)) %>%
ggplot(aes(y=mean, x=biasmatch, fill=Prime)) +
geom_bar(stat="identity",position="dodge") +
geom_errorbar(aes(ymin=mean, ymax=mean+sd), width=.2,
position=position_dodge(.75)) +
ylim(0, .9) +
theme_minimal() +
labs(y="Raw Mean")
)
## `summarise()` has grouped output by 'Prime'. You can override using the
## `.groups` argument.
library(ggrepel)
DOD_54mo %>%
filter(Overlap == "No Overlap") %>%
group_by(Prime, biasmatch, Target.verb) %>%
summarise(mean = mean(RESPONSE.CODE, na.rm=TRUE),
sd = sd(RESPONSE.CODE, na.rm=TRUE),
total=n()) %>%
ggplot(aes(y=mean, x=biasmatch, fill=Prime, label=total)) +
geom_bar(stat="identity",position="dodge") +
geom_errorbar(aes(ymin=mean, ymax=mean+sd), width=.2,
position=position_dodge(.9)) +
geom_text_repel() +
facet_wrap(~Target.verb) + theme_minimal()
## `summarise()` has grouped output by 'Prime', 'biasmatch'. You can override using
## the `.groups` argument.
DOD_54mo %>%
filter(Overlap == "No Overlap") %>%
group_by(Target.verb, Prime, biasmatch) %>%
count()
## # A tibble: 22 x 4
## # Groups: Target.verb, Prime, biasmatch [22]
## Target.verb Prime biasmatch n
## <chr> <chr> <chr> <int>
## 1 brought DOD match 19
## 2 brought DOD mismatch 18
## 3 brought POD match 23
## 4 brought POD mismatch 18
## 5 gave DOD match 2
## 6 gave DOD mismatch 46
## 7 gave POD match 41
## 8 gave POD mismatch 1
## 9 passed DOD match 11
## 10 passed DOD mismatch 20
## # ... with 12 more rows
Again, it looks like the pattern is different on the levels of aggregate data and indiviaul verbs, which was the case with the 42 month verbs.
DOD_54mo %>%
filter(Overlap == "No Overlap") %>%
filter(Target.verb=="gave" & Prime=="DOD" & biasmatch == "match")
## # A tibble: 2 x 47
## Blind_ID Zoom Gender Session.Age List Item Prime.~1 Targe~2 Prime~3 Verb.~4
## <chr> <dbl> <dbl> <chr> <dbl> <chr> <chr> <chr> <dbl> <dbl>
## 1 9hKLS9my 1 1 54 : 22 3 2 showed gave 0.5 -0.5
## 2 LkVZQCTn 1 1 54 : 15 7 2 showed gave 0.5 -0.5
## # ... with 37 more variables: Verb.error <dbl>, Admin.error <dbl>,
## # COVID.Zoom.session <dbl>, Repeat.prime <dbl>, Repeat.target.prompt <dbl>,
## # Prime.interference <dbl>, Missing.article <dbl>, Pronoun <dbl>,
## # Prep.error <chr>, Prep.error.bin <dbl>, Reversal.error <dbl>,
## # Prompts.allowable <dbl>, Prompts.disallowed <dbl>,
## # NonTarget.response <dbl>, No.response <lgl>, Gave.POD <dbl>,
## # RESPONSE.CODE <dbl>, Prime.code <dbl>, Transcription <chr>, ...
DOD_54mo %>%
filter(Overlap == "No Overlap") %>%
filter(Target.verb=="gave" & Prime=="POD" & biasmatch == "mismatch")
## # A tibble: 1 x 47
## Blind_ID Zoom Gender Session.Age List Item Prime.~1 Targe~2 Prime~3 Verb.~4
## <chr> <dbl> <dbl> <chr> <dbl> <chr> <chr> <chr> <dbl> <dbl>
## 1 XncBwsN2 NA 1 54 : 22 8 12r gave gave -0.5 -0.5
## # ... with 37 more variables: Verb.error <dbl>, Admin.error <dbl>,
## # COVID.Zoom.session <dbl>, Repeat.prime <dbl>, Repeat.target.prompt <dbl>,
## # Prime.interference <dbl>, Missing.article <dbl>, Pronoun <dbl>,
## # Prep.error <chr>, Prep.error.bin <dbl>, Reversal.error <dbl>,
## # Prompts.allowable <dbl>, Prompts.disallowed <dbl>,
## # NonTarget.response <dbl>, No.response <lgl>, Gave.POD <dbl>,
## # RESPONSE.CODE <dbl>, Prime.code <dbl>, Transcription <chr>, ...
logit_prior = prior(normal(0,2), class=sd)
First lets consider whether we should include random effects by the interaction between ID and Item. In principle this is possible, but I think it will greatly complicated the model. So I’ll fit two variance component models with the two random effects specifications and compare them to one another.
m54mo0<- brm(RESPONSE.CODE ~ 1 + (1 | Blind_ID) + (1 | Target.verb), family="bernoulli", prior=logit_prior, seed=12345, iter=4000, cores=4, data=DOD_54mo, file="54mo/output/m54mo0.rds")
m54mo0b <- brm(RESPONSE.CODE ~ 1 + (1 | Blind_ID) + (1 | Target.verb) + (1 | Blind_ID:Target.verb), seed=12345, family="bernoulli", prior=logit_prior, iter=4000, cores=4, data=DOD_54mo, file="54mo/output/m54mo0b.rds")
Add loo for loo_compare
m54mo0<- add_criterion(m54mo0, "loo")
m54mo0b<- add_criterion(m54mo0b, "loo")
loo(m54mo0, m54mo0b, reloo=TRUE)
## Output of model 'm54mo0':
##
## Computed from 8000 by 1011 log-likelihood matrix
##
## Estimate SE
## elpd_loo -282.7 19.3
## p_loo 51.7 4.7
## looic 565.3 38.6
## ------
## Monte Carlo SE of elpd_loo is 0.1.
##
## Pareto k diagnostic values:
## Count Pct. Min. n_eff
## (-Inf, 0.5] (good) 1003 99.2% 2544
## (0.5, 0.7] (ok) 8 0.8% 1377
## (0.7, 1] (bad) 0 0.0% <NA>
## (1, Inf) (very bad) 0 0.0% <NA>
##
## All Pareto k estimates are ok (k < 0.7).
## See help('pareto-k-diagnostic') for details.
##
## Output of model 'm54mo0b':
##
## Computed from 8000 by 1011 log-likelihood matrix
##
## Estimate SE
## elpd_loo -283.7 19.3
## p_loo 64.6 5.7
## looic 567.4 38.7
## ------
## Monte Carlo SE of elpd_loo is NA.
##
## Pareto k diagnostic values:
## Count Pct. Min. n_eff
## (-Inf, 0.5] (good) 993 98.2% 1153
## (0.5, 0.7] (ok) 17 1.7% 1276
## (0.7, 1] (bad) 1 0.1% 1638
## (1, Inf) (very bad) 0 0.0% <NA>
## See help('pareto-k-diagnostic') for details.
##
## Model comparisons:
## elpd_diff se_diff
## m54mo0 0.0 0.0
## m54mo0b -1.0 0.7
model_weights(m54mo0, m54mo0b)
## Warning: Some Pareto k diagnostic values are slightly high. See help('pareto-k-diagnostic') for details.
## Warning: Some Pareto k diagnostic values are too high. See help('pareto-k-diagnostic') for details.
## m54mo0 m54mo0b
## 9.999983e-01 1.667583e-06
#loo1 <- loo(m54mo0b)
#loo2 <- reloo(m54mo0b, loo1)
#loo3 <- loo(m54mo0)
#loo_compare(loo2, loo3)
Here’s the model of the lexical boost.
m54mo1 <- brm(RESPONSE.CODE ~ 1 + Prime.Type*Overlap + (1 + Prime.Type*Overlap | Blind_ID) + (1 + Prime.Type*Overlap | Target.verb), logit_prior,family="bernoulli", cores=4, seed=12345, control=list(adapt_delta=.9), data=DOD_54mo, file="54mo/output/m54mo1.rds")
m54mo1 <- add_criterion(m54mo1, "loo")
Chains look okay.
#plot(m54mo1, N=4)
summary(m54mo1, prob=.9)
## Family: bernoulli
## Links: mu = logit
## Formula: RESPONSE.CODE ~ 1 + Prime.Type * Overlap + (1 + Prime.Type * Overlap | Blind_ID) + (1 + Prime.Type * Overlap | Target.verb)
## Data: DOD_54mo (Number of observations: 1011)
## Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
## total post-warmup draws = 4000
##
## Group-Level Effects:
## ~Blind_ID (Number of levels: 87)
## Estimate Est.Error l-90% CI
## sd(Intercept) 1.69 0.32 1.20
## sd(Prime.Type) 0.71 0.48 0.06
## sd(OverlapOverlap) 0.62 0.43 0.06
## sd(Prime.Type:OverlapOverlap) 1.04 0.69 0.12
## cor(Intercept,Prime.Type) 0.15 0.40 -0.56
## cor(Intercept,OverlapOverlap) 0.28 0.41 -0.48
## cor(Prime.Type,OverlapOverlap) 0.10 0.44 -0.65
## cor(Intercept,Prime.Type:OverlapOverlap) 0.14 0.41 -0.58
## cor(Prime.Type,Prime.Type:OverlapOverlap) -0.09 0.44 -0.77
## cor(OverlapOverlap,Prime.Type:OverlapOverlap) 0.04 0.45 -0.71
## u-90% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 2.23 1.00 1723 2802
## sd(Prime.Type) 1.57 1.00 1366 1999
## sd(OverlapOverlap) 1.42 1.00 1111 2123
## sd(Prime.Type:OverlapOverlap) 2.29 1.00 1362 1950
## cor(Intercept,Prime.Type) 0.76 1.00 5057 3135
## cor(Intercept,OverlapOverlap) 0.85 1.00 4252 3166
## cor(Prime.Type,OverlapOverlap) 0.78 1.00 2112 2659
## cor(Intercept,Prime.Type:OverlapOverlap) 0.76 1.00 5092 2712
## cor(Prime.Type,Prime.Type:OverlapOverlap) 0.66 1.00 2978 3276
## cor(OverlapOverlap,Prime.Type:OverlapOverlap) 0.76 1.00 2561 3161
##
## ~Target.verb (Number of levels: 6)
## Estimate Est.Error l-90% CI
## sd(Intercept) 1.81 0.62 1.02
## sd(Prime.Type) 0.58 0.51 0.04
## sd(OverlapOverlap) 0.88 0.55 0.13
## sd(Prime.Type:OverlapOverlap) 1.47 1.02 0.14
## cor(Intercept,Prime.Type) 0.00 0.44 -0.71
## cor(Intercept,OverlapOverlap) -0.36 0.38 -0.88
## cor(Prime.Type,OverlapOverlap) 0.01 0.44 -0.71
## cor(Intercept,Prime.Type:OverlapOverlap) -0.09 0.39 -0.69
## cor(Prime.Type,Prime.Type:OverlapOverlap) -0.06 0.45 -0.76
## cor(OverlapOverlap,Prime.Type:OverlapOverlap) 0.06 0.43 -0.66
## u-90% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 2.99 1.00 2038 2561
## sd(Prime.Type) 1.54 1.00 2662 2029
## sd(OverlapOverlap) 1.88 1.00 1927 1998
## sd(Prime.Type:OverlapOverlap) 3.35 1.00 1871 2606
## cor(Intercept,Prime.Type) 0.72 1.00 6343 2947
## cor(Intercept,OverlapOverlap) 0.35 1.00 4546 3046
## cor(Prime.Type,OverlapOverlap) 0.74 1.00 3684 3281
## cor(Intercept,Prime.Type:OverlapOverlap) 0.60 1.00 5935 3263
## cor(Prime.Type,Prime.Type:OverlapOverlap) 0.69 1.00 2913 2794
## cor(OverlapOverlap,Prime.Type:OverlapOverlap) 0.74 1.00 3473 3580
##
## Population-Level Effects:
## Estimate Est.Error l-90% CI u-90% CI Rhat Bulk_ESS
## Intercept -3.19 0.84 -4.51 -1.78 1.00 1739
## Prime.Type 0.84 0.57 -0.06 1.76 1.00 4098
## OverlapOverlap -0.56 0.66 -1.62 0.50 1.00 2148
## Prime.Type:OverlapOverlap 0.72 1.07 -0.97 2.54 1.00 3476
## Tail_ESS
## Intercept 2536
## Prime.Type 3202
## OverlapOverlap 2331
## Prime.Type:OverlapOverlap 2751
##
## 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).
hypothesis(m54mo1, "Prime.Type > 0")
## Hypothesis Tests for class b:
## Hypothesis Estimate Est.Error CI.Lower CI.Upper Evid.Ratio Post.Prob
## 1 (Prime.Type) > 0 0.84 0.57 -0.06 1.76 14.62 0.94
## Star
## 1
## ---
## 'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
## '*': For one-sided hypotheses, the posterior probability exceeds 95%;
## for two-sided hypotheses, the value tested against lies outside the 95%-CI.
## Posterior probabilities of point hypotheses assume equal prior probabilities.
hypothesis(m54mo1, "Prime.Type:OverlapOverlap > 0")
## Hypothesis Tests for class b:
## Hypothesis Estimate Est.Error CI.Lower CI.Upper Evid.Ratio
## 1 (Prime.Type:Overl... > 0 0.72 1.07 -0.97 2.54 3.28
## Post.Prob Star
## 1 0.77
## ---
## 'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
## '*': For one-sided hypotheses, the posterior probability exceeds 95%;
## for two-sided hypotheses, the value tested against lies outside the 95%-CI.
## Posterior probabilities of point hypotheses assume equal prior probabilities.
hypothesis(m54mo1, "Intercept + OverlapOverlap + .5*Prime.Type + .5*Prime.Type:OverlapOverlap >
Intercept + OverlapOverlap - .5*Prime.Type - .5*Prime.Type:OverlapOverlap ")
## Hypothesis Tests for class b:
## Hypothesis Estimate Est.Error CI.Lower CI.Upper Evid.Ratio
## 1 (Intercept+Overla... > 0 1.56 1.05 -0.08 3.32 15.39
## Post.Prob Star
## 1 0.94
## ---
## 'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
## '*': For one-sided hypotheses, the posterior probability exceeds 95%;
## for two-sided hypotheses, the value tested against lies outside the 95%-CI.
## Posterior probabilities of point hypotheses assume equal prior probabilities.
m54mo1.model.check <- createDHARMa(
simulatedResponse = t(posterior_predict(m54mo1)),
observedResponse = DOD_54mo$RESPONSE.CODE,
fittedPredictedResponse = apply(t(posterior_epred(m54mo1)), 1, mean),
integerResponse=TRUE
)
plot(m54mo1.model.check) # standard plots
testDispersion(m54mo1.model.check)
##
## DHARMa nonparametric dispersion test via sd of residuals fitted vs.
## simulated
##
## data: simulationOutput
## dispersion = 0.73983, p-value = 0.0045
## alternative hypothesis: two.sided
Sensitivity analyses (without zoom sessions, without problem trials and without trials with high Pareto’s k) in R Markdown file.
m54mo1b <- DOD_54mo %>%
filter(Trial_Problem == 0) %>%
brm(RESPONSE.CODE ~ 1 + Prime.Type*Overlap + (1 + Prime.Type*Overlap | Blind_ID) + (1 + Prime.Type*Overlap | Target.verb), seed=12345, logit_prior, control=list(adapt_delta=.9), family="bernoulli", cores=4, data=., file="54mo/output/m54mo1b.rds")
summary(m54mo1b, prob=.9)
## Family: bernoulli
## Links: mu = logit
## Formula: RESPONSE.CODE ~ 1 + Prime.Type * Overlap + (1 + Prime.Type * Overlap | Blind_ID) + (1 + Prime.Type * Overlap | Target.verb)
## Data: . (Number of observations: 990)
## Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
## total post-warmup draws = 4000
##
## Group-Level Effects:
## ~Blind_ID (Number of levels: 87)
## Estimate Est.Error l-90% CI
## sd(Intercept) 1.70 0.33 1.17
## sd(Prime.Type) 0.69 0.48 0.06
## sd(OverlapOverlap) 0.57 0.42 0.05
## sd(Prime.Type:OverlapOverlap) 0.86 0.63 0.08
## cor(Intercept,Prime.Type) 0.09 0.41 -0.61
## cor(Intercept,OverlapOverlap) 0.25 0.42 -0.53
## cor(Prime.Type,OverlapOverlap) 0.08 0.45 -0.69
## cor(Intercept,Prime.Type:OverlapOverlap) 0.02 0.42 -0.67
## cor(Prime.Type,Prime.Type:OverlapOverlap) -0.10 0.45 -0.77
## cor(OverlapOverlap,Prime.Type:OverlapOverlap) 0.00 0.45 -0.73
## u-90% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 2.26 1.00 1199 2176
## sd(Prime.Type) 1.58 1.00 1222 1711
## sd(OverlapOverlap) 1.35 1.01 807 1754
## sd(Prime.Type:OverlapOverlap) 2.05 1.00 1280 2266
## cor(Intercept,Prime.Type) 0.74 1.00 3563 2817
## cor(Intercept,OverlapOverlap) 0.84 1.00 2650 2717
## cor(Prime.Type,OverlapOverlap) 0.77 1.00 1883 2707
## cor(Intercept,Prime.Type:OverlapOverlap) 0.70 1.00 4060 2828
## cor(Prime.Type,Prime.Type:OverlapOverlap) 0.67 1.00 2775 2413
## cor(OverlapOverlap,Prime.Type:OverlapOverlap) 0.72 1.00 2715 3163
##
## ~Target.verb (Number of levels: 6)
## Estimate Est.Error l-90% CI
## sd(Intercept) 1.82 0.63 0.99
## sd(Prime.Type) 0.62 0.54 0.05
## sd(OverlapOverlap) 0.98 0.55 0.19
## sd(Prime.Type:OverlapOverlap) 1.81 1.15 0.21
## cor(Intercept,Prime.Type) -0.01 0.44 -0.72
## cor(Intercept,OverlapOverlap) -0.40 0.37 -0.89
## cor(Prime.Type,OverlapOverlap) -0.01 0.44 -0.73
## cor(Intercept,Prime.Type:OverlapOverlap) -0.13 0.37 -0.72
## cor(Prime.Type,Prime.Type:OverlapOverlap) -0.11 0.45 -0.81
## cor(OverlapOverlap,Prime.Type:OverlapOverlap) 0.09 0.41 -0.60
## u-90% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 2.97 1.00 1672 2127
## sd(Prime.Type) 1.69 1.00 1904 2216
## sd(OverlapOverlap) 1.99 1.00 1308 1044
## sd(Prime.Type:OverlapOverlap) 3.87 1.00 951 1472
## cor(Intercept,Prime.Type) 0.71 1.00 4495 2814
## cor(Intercept,OverlapOverlap) 0.28 1.00 2877 2718
## cor(Prime.Type,OverlapOverlap) 0.71 1.00 2568 2771
## cor(Intercept,Prime.Type:OverlapOverlap) 0.51 1.00 3252 2862
## cor(Prime.Type,Prime.Type:OverlapOverlap) 0.65 1.00 1621 2356
## cor(OverlapOverlap,Prime.Type:OverlapOverlap) 0.74 1.00 2558 3393
##
## Population-Level Effects:
## Estimate Est.Error l-90% CI u-90% CI Rhat Bulk_ESS
## Intercept -3.19 0.85 -4.57 -1.79 1.00 1240
## Prime.Type 0.86 0.60 -0.09 1.83 1.00 2126
## OverlapOverlap -0.62 0.70 -1.74 0.49 1.00 1423
## Prime.Type:OverlapOverlap 0.90 1.18 -0.88 2.95 1.00 2031
## Tail_ESS
## Intercept 1867
## Prime.Type 2008
## OverlapOverlap 1967
## Prime.Type:OverlapOverlap 2065
##
## 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).
hypothesis(m54mo1b, "Prime.Type > 0")
## Hypothesis Tests for class b:
## Hypothesis Estimate Est.Error CI.Lower CI.Upper Evid.Ratio Post.Prob
## 1 (Prime.Type) > 0 0.86 0.6 -0.09 1.83 14.75 0.94
## Star
## 1
## ---
## 'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
## '*': For one-sided hypotheses, the posterior probability exceeds 95%;
## for two-sided hypotheses, the value tested against lies outside the 95%-CI.
## Posterior probabilities of point hypotheses assume equal prior probabilities.
hypothesis(m54mo1b, "Prime.Type:OverlapOverlap > 0")
## Hypothesis Tests for class b:
## Hypothesis Estimate Est.Error CI.Lower CI.Upper Evid.Ratio
## 1 (Prime.Type:Overl... > 0 0.9 1.18 -0.88 2.95 3.76
## Post.Prob Star
## 1 0.79
## ---
## 'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
## '*': For one-sided hypotheses, the posterior probability exceeds 95%;
## for two-sided hypotheses, the value tested against lies outside the 95%-CI.
## Posterior probabilities of point hypotheses assume equal prior probabilities.
m54mo1c <- DOD_54mo %>%
filter(is.na(Zoom)) %>%
brm(RESPONSE.CODE ~ 1 + Prime.Type*Overlap + (1 + Prime.Type*Overlap | Blind_ID) + (1 + Prime.Type*Overlap | Target.verb), seed=12345, logit_prior,family="bernoulli", control=list(adapt_delta=.9), cores=4, data=., file="54mo/output/m54mo1c.rds")
summary(m54mo1c, prob=.9)
## Family: bernoulli
## Links: mu = logit
## Formula: RESPONSE.CODE ~ 1 + Prime.Type * Overlap + (1 + Prime.Type * Overlap | Blind_ID) + (1 + Prime.Type * Overlap | Target.verb)
## Data: . (Number of observations: 892)
## Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
## total post-warmup draws = 4000
##
## Group-Level Effects:
## ~Blind_ID (Number of levels: 76)
## Estimate Est.Error l-90% CI
## sd(Intercept) 1.70 0.33 1.19
## sd(Prime.Type) 0.88 0.55 0.08
## sd(OverlapOverlap) 0.65 0.45 0.06
## sd(Prime.Type:OverlapOverlap) 1.16 0.77 0.11
## cor(Intercept,Prime.Type) 0.22 0.39 -0.47
## cor(Intercept,OverlapOverlap) 0.29 0.41 -0.50
## cor(Prime.Type,OverlapOverlap) 0.13 0.43 -0.64
## cor(Intercept,Prime.Type:OverlapOverlap) 0.14 0.40 -0.55
## cor(Prime.Type,Prime.Type:OverlapOverlap) -0.10 0.43 -0.75
## cor(OverlapOverlap,Prime.Type:OverlapOverlap) 0.06 0.44 -0.67
## u-90% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 2.27 1.00 1825 2591
## sd(Prime.Type) 1.84 1.00 1314 1539
## sd(OverlapOverlap) 1.49 1.00 1158 2048
## sd(Prime.Type:OverlapOverlap) 2.57 1.00 1148 2095
## cor(Intercept,Prime.Type) 0.80 1.00 4165 3314
## cor(Intercept,OverlapOverlap) 0.85 1.00 4341 2884
## cor(Prime.Type,OverlapOverlap) 0.78 1.00 2404 3145
## cor(Intercept,Prime.Type:OverlapOverlap) 0.75 1.00 4971 3316
## cor(Prime.Type,Prime.Type:OverlapOverlap) 0.66 1.00 2977 3183
## cor(OverlapOverlap,Prime.Type:OverlapOverlap) 0.76 1.00 2739 3646
##
## ~Target.verb (Number of levels: 6)
## Estimate Est.Error l-90% CI
## sd(Intercept) 1.71 0.59 0.93
## sd(Prime.Type) 0.55 0.49 0.04
## sd(OverlapOverlap) 0.79 0.52 0.11
## sd(Prime.Type:OverlapOverlap) 1.15 0.90 0.08
## cor(Intercept,Prime.Type) -0.01 0.44 -0.72
## cor(Intercept,OverlapOverlap) -0.35 0.39 -0.88
## cor(Prime.Type,OverlapOverlap) 0.02 0.45 -0.72
## cor(Intercept,Prime.Type:OverlapOverlap) -0.04 0.41 -0.69
## cor(Prime.Type,Prime.Type:OverlapOverlap) -0.04 0.46 -0.76
## cor(OverlapOverlap,Prime.Type:OverlapOverlap) -0.02 0.43 -0.71
## u-90% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 2.83 1.00 2760 2942
## sd(Prime.Type) 1.52 1.00 3057 2870
## sd(OverlapOverlap) 1.74 1.00 2089 2150
## sd(Prime.Type:OverlapOverlap) 2.81 1.00 2130 2601
## cor(Intercept,Prime.Type) 0.71 1.00 6629 2998
## cor(Intercept,OverlapOverlap) 0.39 1.00 5038 2943
## cor(Prime.Type,OverlapOverlap) 0.74 1.00 3778 3266
## cor(Intercept,Prime.Type:OverlapOverlap) 0.65 1.00 6355 2582
## cor(Prime.Type,Prime.Type:OverlapOverlap) 0.71 1.00 3616 3304
## cor(OverlapOverlap,Prime.Type:OverlapOverlap) 0.69 1.00 3748 3642
##
## Population-Level Effects:
## Estimate Est.Error l-90% CI u-90% CI Rhat Bulk_ESS
## Intercept -3.26 0.83 -4.57 -1.90 1.00 1632
## Prime.Type 0.74 0.61 -0.28 1.71 1.00 3791
## OverlapOverlap -0.49 0.65 -1.56 0.53 1.00 2167
## Prime.Type:OverlapOverlap 0.55 1.00 -1.07 2.19 1.00 3483
## Tail_ESS
## Intercept 2028
## Prime.Type 2514
## OverlapOverlap 2430
## Prime.Type:OverlapOverlap 2868
##
## 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).
hypothesis(m54mo1c, "Prime.Type > 0")
## Hypothesis Tests for class b:
## Hypothesis Estimate Est.Error CI.Lower CI.Upper Evid.Ratio Post.Prob
## 1 (Prime.Type) > 0 0.74 0.61 -0.28 1.71 8.83 0.9
## Star
## 1
## ---
## 'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
## '*': For one-sided hypotheses, the posterior probability exceeds 95%;
## for two-sided hypotheses, the value tested against lies outside the 95%-CI.
## Posterior probabilities of point hypotheses assume equal prior probabilities.
hypothesis(m54mo1c, "Prime.Type:OverlapOverlap > 0")
## Hypothesis Tests for class b:
## Hypothesis Estimate Est.Error CI.Lower CI.Upper Evid.Ratio
## 1 (Prime.Type:Overl... > 0 0.55 1 -1.07 2.19 2.61
## Post.Prob Star
## 1 0.72
## ---
## 'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
## '*': For one-sided hypotheses, the posterior probability exceeds 95%;
## for two-sided hypotheses, the value tested against lies outside the 95%-CI.
## Posterior probabilities of point hypotheses assume equal prior probabilities.
m54mo1d <- DOD_54mo %>%
mutate(
k = m54mo1$criteria$loo$diagnostics$pareto_k
) %>%
filter(k < .7) %>%
brm(RESPONSE.CODE ~ 1 + Prime.Type*Overlap + (1 + Prime.Type*Overlap | Blind_ID) + (1 + Prime.Type*Overlap | Target.verb), seed=12345, logit_prior,family="bernoulli", control=list(adapt_delta=.9), cores=4, data=., file="54mo/output/m54mo1d.rds")
summary(m54mo1d, prob=.9)
## Family: bernoulli
## Links: mu = logit
## Formula: RESPONSE.CODE ~ 1 + Prime.Type * Overlap + (1 + Prime.Type * Overlap | Blind_ID) + (1 + Prime.Type * Overlap | Target.verb)
## Data: . (Number of observations: 1007)
## Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
## total post-warmup draws = 4000
##
## Group-Level Effects:
## ~Blind_ID (Number of levels: 87)
## Estimate Est.Error l-90% CI
## sd(Intercept) 1.76 0.34 1.25
## sd(Prime.Type) 0.74 0.49 0.07
## sd(OverlapOverlap) 0.74 0.47 0.08
## sd(Prime.Type:OverlapOverlap) 1.19 0.76 0.13
## cor(Intercept,Prime.Type) 0.17 0.40 -0.54
## cor(Intercept,OverlapOverlap) 0.34 0.38 -0.38
## cor(Prime.Type,OverlapOverlap) 0.13 0.43 -0.61
## cor(Intercept,Prime.Type:OverlapOverlap) 0.21 0.40 -0.50
## cor(Prime.Type,Prime.Type:OverlapOverlap) -0.08 0.44 -0.74
## cor(OverlapOverlap,Prime.Type:OverlapOverlap) 0.10 0.44 -0.66
## u-90% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 2.36 1.00 1710 2513
## sd(Prime.Type) 1.63 1.00 1366 1606
## sd(OverlapOverlap) 1.60 1.00 1031 1520
## sd(Prime.Type:OverlapOverlap) 2.54 1.00 1348 2167
## cor(Intercept,Prime.Type) 0.77 1.00 3387 2864
## cor(Intercept,OverlapOverlap) 0.87 1.00 3233 2943
## cor(Prime.Type,OverlapOverlap) 0.78 1.00 1810 2746
## cor(Intercept,Prime.Type:OverlapOverlap) 0.80 1.00 4022 2829
## cor(Prime.Type,Prime.Type:OverlapOverlap) 0.68 1.00 2726 3252
## cor(OverlapOverlap,Prime.Type:OverlapOverlap) 0.78 1.00 2618 3055
##
## ~Target.verb (Number of levels: 6)
## Estimate Est.Error l-90% CI
## sd(Intercept) 1.82 0.64 1.00
## sd(Prime.Type) 0.67 0.57 0.05
## sd(OverlapOverlap) 0.86 0.55 0.12
## sd(Prime.Type:OverlapOverlap) 1.54 1.08 0.12
## cor(Intercept,Prime.Type) 0.02 0.44 -0.68
## cor(Intercept,OverlapOverlap) -0.35 0.40 -0.88
## cor(Prime.Type,OverlapOverlap) 0.03 0.45 -0.72
## cor(Intercept,Prime.Type:OverlapOverlap) -0.08 0.40 -0.70
## cor(Prime.Type,Prime.Type:OverlapOverlap) -0.05 0.45 -0.75
## cor(OverlapOverlap,Prime.Type:OverlapOverlap) 0.06 0.43 -0.64
## u-90% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 2.98 1.00 2103 2663
## sd(Prime.Type) 1.76 1.00 2243 1913
## sd(OverlapOverlap) 1.89 1.00 1614 1424
## sd(Prime.Type:OverlapOverlap) 3.52 1.00 1495 1807
## cor(Intercept,Prime.Type) 0.74 1.00 5030 2935
## cor(Intercept,OverlapOverlap) 0.39 1.00 4192 2644
## cor(Prime.Type,OverlapOverlap) 0.74 1.00 3136 3457
## cor(Intercept,Prime.Type:OverlapOverlap) 0.60 1.00 4229 2937
## cor(Prime.Type,Prime.Type:OverlapOverlap) 0.71 1.00 2947 2988
## cor(OverlapOverlap,Prime.Type:OverlapOverlap) 0.75 1.00 3247 3129
##
## Population-Level Effects:
## Estimate Est.Error l-90% CI u-90% CI Rhat Bulk_ESS
## Intercept -3.24 0.87 -4.64 -1.83 1.00 1507
## Prime.Type 0.77 0.61 -0.21 1.76 1.00 3108
## OverlapOverlap -0.76 0.71 -1.95 0.38 1.00 1858
## Prime.Type:OverlapOverlap 0.58 1.15 -1.16 2.53 1.00 2415
## Tail_ESS
## Intercept 2151
## Prime.Type 2795
## OverlapOverlap 2314
## Prime.Type:OverlapOverlap 2468
##
## 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).
hypothesis(m54mo1d, "Prime.Type > 0")
## Hypothesis Tests for class b:
## Hypothesis Estimate Est.Error CI.Lower CI.Upper Evid.Ratio Post.Prob
## 1 (Prime.Type) > 0 0.77 0.61 -0.21 1.76 9.61 0.91
## Star
## 1
## ---
## 'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
## '*': For one-sided hypotheses, the posterior probability exceeds 95%;
## for two-sided hypotheses, the value tested against lies outside the 95%-CI.
## Posterior probabilities of point hypotheses assume equal prior probabilities.
hypothesis(m54mo1d, "Prime.Type:OverlapOverlap > 0")
## Hypothesis Tests for class b:
## Hypothesis Estimate Est.Error CI.Lower CI.Upper Evid.Ratio
## 1 (Prime.Type:Overl... > 0 0.58 1.15 -1.16 2.53 2.36
## Post.Prob Star
## 1 0.7
## ---
## 'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
## '*': For one-sided hypotheses, the posterior probability exceeds 95%;
## for two-sided hypotheses, the value tested against lies outside the 95%-CI.
## Posterior probabilities of point hypotheses assume equal prior probabilities.
m54mo2 <- brm(RESPONSE.CODE ~ 1 + Prime.Type*Overlap + (1 + Prime.Type*Overlap | Blind_ID) + (1 + Prime.Type*Overlap | Target.verb), logit_prior,family="bernoulli", cores=4, seed=12345, control=list(adapt_delta=.9), data=DOD_54mo2, file="54mo/output/m54mo2.rds")
m54mo2 <- add_criterion(m54mo2, "loo")
summary(m54mo2, prob=.9)
## Family: bernoulli
## Links: mu = logit
## Formula: RESPONSE.CODE ~ 1 + Prime.Type * Overlap + (1 + Prime.Type * Overlap | Blind_ID) + (1 + Prime.Type * Overlap | Target.verb)
## Data: DOD_54mo2 (Number of observations: 872)
## Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
## total post-warmup draws = 4000
##
## Group-Level Effects:
## ~Blind_ID (Number of levels: 75)
## Estimate Est.Error l-90% CI
## sd(Intercept) 1.51 0.31 1.03
## sd(Prime.Type) 0.72 0.49 0.07
## sd(OverlapOverlap) 0.61 0.42 0.06
## sd(Prime.Type:OverlapOverlap) 1.00 0.68 0.09
## cor(Intercept,Prime.Type) 0.14 0.41 -0.57
## cor(Intercept,OverlapOverlap) 0.29 0.40 -0.44
## cor(Prime.Type,OverlapOverlap) 0.11 0.44 -0.66
## cor(Intercept,Prime.Type:OverlapOverlap) 0.12 0.41 -0.60
## cor(Prime.Type,Prime.Type:OverlapOverlap) -0.10 0.46 -0.79
## cor(OverlapOverlap,Prime.Type:OverlapOverlap) 0.03 0.44 -0.69
## u-90% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 2.06 1.00 1735 2613
## sd(Prime.Type) 1.62 1.00 1366 1845
## sd(OverlapOverlap) 1.36 1.00 1050 2128
## sd(Prime.Type:OverlapOverlap) 2.26 1.00 1388 2091
## cor(Intercept,Prime.Type) 0.77 1.00 4482 2596
## cor(Intercept,OverlapOverlap) 0.85 1.00 3763 3080
## cor(Prime.Type,OverlapOverlap) 0.78 1.00 1632 2191
## cor(Intercept,Prime.Type:OverlapOverlap) 0.76 1.00 5009 2700
## cor(Prime.Type,Prime.Type:OverlapOverlap) 0.70 1.00 2872 3136
## cor(OverlapOverlap,Prime.Type:OverlapOverlap) 0.75 1.00 2851 3106
##
## ~Target.verb (Number of levels: 6)
## Estimate Est.Error l-90% CI
## sd(Intercept) 1.82 0.62 1.02
## sd(Prime.Type) 0.57 0.49 0.04
## sd(OverlapOverlap) 0.89 0.56 0.14
## sd(Prime.Type:OverlapOverlap) 1.41 1.03 0.12
## cor(Intercept,Prime.Type) 0.00 0.44 -0.72
## cor(Intercept,OverlapOverlap) -0.36 0.39 -0.88
## cor(Prime.Type,OverlapOverlap) 0.02 0.45 -0.71
## cor(Intercept,Prime.Type:OverlapOverlap) -0.09 0.41 -0.73
## cor(Prime.Type,Prime.Type:OverlapOverlap) -0.05 0.45 -0.77
## cor(OverlapOverlap,Prime.Type:OverlapOverlap) 0.05 0.43 -0.67
## u-90% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 3.01 1.00 2678 3006
## sd(Prime.Type) 1.50 1.00 2628 2497
## sd(OverlapOverlap) 1.91 1.00 1576 1410
## sd(Prime.Type:OverlapOverlap) 3.41 1.00 1473 2318
## cor(Intercept,Prime.Type) 0.72 1.00 6932 3133
## cor(Intercept,OverlapOverlap) 0.36 1.00 4358 3070
## cor(Prime.Type,OverlapOverlap) 0.74 1.00 3174 3233
## cor(Intercept,Prime.Type:OverlapOverlap) 0.61 1.00 5257 2815
## cor(Prime.Type,Prime.Type:OverlapOverlap) 0.69 1.00 2683 3104
## cor(OverlapOverlap,Prime.Type:OverlapOverlap) 0.74 1.00 3221 3459
##
## Population-Level Effects:
## Estimate Est.Error l-90% CI u-90% CI Rhat Bulk_ESS
## Intercept -2.80 0.83 -4.16 -1.43 1.00 1456
## Prime.Type 0.84 0.55 -0.01 1.72 1.00 3608
## OverlapOverlap -0.53 0.63 -1.54 0.52 1.00 2501
## Prime.Type:OverlapOverlap 0.80 1.05 -0.76 2.57 1.00 2816
## Tail_ESS
## Intercept 1959
## Prime.Type 2843
## OverlapOverlap 2744
## Prime.Type:OverlapOverlap 2150
##
## 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).
hypothesis(m54mo2, "Prime.Type > 0")
## Hypothesis Tests for class b:
## Hypothesis Estimate Est.Error CI.Lower CI.Upper Evid.Ratio Post.Prob
## 1 (Prime.Type) > 0 0.84 0.55 -0.01 1.72 17.6 0.95
## Star
## 1
## ---
## 'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
## '*': For one-sided hypotheses, the posterior probability exceeds 95%;
## for two-sided hypotheses, the value tested against lies outside the 95%-CI.
## Posterior probabilities of point hypotheses assume equal prior probabilities.
hypothesis(m54mo2, "Prime.Type:OverlapOverlap > 0")
## Hypothesis Tests for class b:
## Hypothesis Estimate Est.Error CI.Lower CI.Upper Evid.Ratio
## 1 (Prime.Type:Overl... > 0 0.8 1.05 -0.76 2.57 4.12
## Post.Prob Star
## 1 0.8
## ---
## 'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
## '*': For one-sided hypotheses, the posterior probability exceeds 95%;
## for two-sided hypotheses, the value tested against lies outside the 95%-CI.
## Posterior probabilities of point hypotheses assume equal prior probabilities.
hypothesis(m54mo2, "Intercept + OverlapOverlap + .5*Prime.Type + .5*Prime.Type:OverlapOverlap >
Intercept + OverlapOverlap - .5*Prime.Type - .5*Prime.Type:OverlapOverlap ")
## Hypothesis Tests for class b:
## Hypothesis Estimate Est.Error CI.Lower CI.Upper Evid.Ratio
## 1 (Intercept+Overla... > 0 1.64 1 0.11 3.4 23.24
## Post.Prob Star
## 1 0.96 *
## ---
## 'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
## '*': For one-sided hypotheses, the posterior probability exceeds 95%;
## for two-sided hypotheses, the value tested against lies outside the 95%-CI.
## Posterior probabilities of point hypotheses assume equal prior probabilities.
Chains look good.
#plot(m54mo2, N=4)
m54mo2.model.check <- createDHARMa(
simulatedResponse = t(posterior_predict(m54mo2)),
observedResponse = DOD_54mo2$RESPONSE.CODE,
fittedPredictedResponse = apply(t(posterior_epred(m54mo2)), 1, mean),
integerResponse=TRUE
)
plot(m54mo2.model.check) # standard plots
testDispersion(m54mo2.model.check)
##
## DHARMa nonparametric dispersion test via sd of residuals fitted vs.
## simulated
##
## data: simulationOutput
## dispersion = 0.74686, p-value = 0.006
## alternative hypothesis: two.sided
tab_model(m54mo1, m54mo2, show.re.var=FALSE, show.r2=FALSE, show.icc=FALSE, transform=NULL, collapse.ci=TRUE, p.style="scientific", ci.hyphen = " : ", dv.labels=c("Model 1", "Model 2"), string.stat="",
pred.labels=c(
"Intercept",
"Prime (DOD)",
"Overlap",
"Prime * Overlap"
))
| Model 1 | Model 2 | |
|---|---|---|
| Predictors | Log-Odds | Log-Odds |
| Intercept |
-3.21 (-4.80 : -1.42) |
-2.81 (-4.45 : -1.09) |
| Prime (DOD) |
0.83 (-0.27 : 1.99) |
0.83 (-0.23 : 1.94) |
| Overlap |
-0.57 (-1.91 : 0.73) |
-0.55 (-1.80 : 0.73) |
| Prime * Overlap |
0.70 (-1.35 : 3.04) |
0.75 (-1.17 : 3.14) |
| N | 87 Blind_ID | 75 Blind_ID |
| 6 Target.verb | 6 Target.verb | |
| Observations | 1011 | 872 |
DOD_54mo3 <- filter(DOD_54mo, Overlap=="No Overlap")
m54mo3 <- brm(RESPONSE.CODE ~ 1 +Prime.Type*biasmatch + (1 + Prime.Type*biasmatch| Blind_ID) + (1 + Prime.Type*biasmatch| Target.verb), seed=12345, control=list(adapt_delta=.9), logit_prior, family="bernoulli", iter=4000, cores=4, file="54mo/output/m54mo3.rds", data=DOD_54mo3)
summary(m54mo3, prob=.9)
## Family: bernoulli
## Links: mu = logit
## Formula: RESPONSE.CODE ~ 1 + Prime.Type * biasmatch + (1 + Prime.Type * biasmatch | Blind_ID) + (1 + Prime.Type * biasmatch | Target.verb)
## Data: DOD_54mo3 (Number of observations: 479)
## Draws: 4 chains, each with iter = 4000; warmup = 2000; thin = 1;
## total post-warmup draws = 8000
##
## Group-Level Effects:
## ~Blind_ID (Number of levels: 87)
## Estimate Est.Error l-90% CI
## sd(Intercept) 1.32 0.49 0.53
## sd(Prime.Type) 1.31 0.86 0.12
## sd(biasmatchmismatch) 1.05 0.69 0.11
## sd(Prime.Type:biasmatchmismatch) 1.21 0.90 0.10
## cor(Intercept,Prime.Type) 0.00 0.42 -0.69
## cor(Intercept,biasmatchmismatch) 0.13 0.40 -0.56
## cor(Prime.Type,biasmatchmismatch) 0.05 0.43 -0.67
## cor(Intercept,Prime.Type:biasmatchmismatch) 0.12 0.44 -0.64
## cor(Prime.Type,Prime.Type:biasmatchmismatch) -0.13 0.45 -0.80
## cor(biasmatchmismatch,Prime.Type:biasmatchmismatch) 0.05 0.45 -0.70
## u-90% CI Rhat Bulk_ESS
## sd(Intercept) 2.16 1.00 1693
## sd(Prime.Type) 2.88 1.00 1834
## sd(biasmatchmismatch) 2.31 1.00 1635
## sd(Prime.Type:biasmatchmismatch) 2.95 1.00 2880
## cor(Intercept,Prime.Type) 0.68 1.00 4830
## cor(Intercept,biasmatchmismatch) 0.76 1.00 4750
## cor(Prime.Type,biasmatchmismatch) 0.73 1.00 3122
## cor(Intercept,Prime.Type:biasmatchmismatch) 0.77 1.00 7231
## cor(Prime.Type,Prime.Type:biasmatchmismatch) 0.64 1.00 5336
## cor(biasmatchmismatch,Prime.Type:biasmatchmismatch) 0.75 1.00 6242
## Tail_ESS
## sd(Intercept) 1537
## sd(Prime.Type) 3691
## sd(biasmatchmismatch) 2889
## sd(Prime.Type:biasmatchmismatch) 3572
## cor(Intercept,Prime.Type) 4828
## cor(Intercept,biasmatchmismatch) 5632
## cor(Prime.Type,biasmatchmismatch) 4850
## cor(Intercept,Prime.Type:biasmatchmismatch) 5769
## cor(Prime.Type,Prime.Type:biasmatchmismatch) 5873
## cor(biasmatchmismatch,Prime.Type:biasmatchmismatch) 6400
##
## ~Target.verb (Number of levels: 6)
## Estimate Est.Error l-90% CI
## sd(Intercept) 1.99 0.75 1.00
## sd(Prime.Type) 1.50 0.94 0.18
## sd(biasmatchmismatch) 0.86 0.72 0.06
## sd(Prime.Type:biasmatchmismatch) 1.18 0.92 0.08
## cor(Intercept,Prime.Type) 0.11 0.41 -0.60
## cor(Intercept,biasmatchmismatch) 0.04 0.45 -0.69
## cor(Prime.Type,biasmatchmismatch) -0.05 0.44 -0.76
## cor(Intercept,Prime.Type:biasmatchmismatch) -0.02 0.44 -0.73
## cor(Prime.Type,Prime.Type:biasmatchmismatch) 0.02 0.44 -0.71
## cor(biasmatchmismatch,Prime.Type:biasmatchmismatch) -0.06 0.44 -0.76
## u-90% CI Rhat Bulk_ESS
## sd(Intercept) 3.39 1.00 4709
## sd(Prime.Type) 3.24 1.00 3254
## sd(biasmatchmismatch) 2.24 1.00 4062
## sd(Prime.Type:biasmatchmismatch) 3.01 1.00 5048
## cor(Intercept,Prime.Type) 0.73 1.00 8203
## cor(Intercept,biasmatchmismatch) 0.75 1.00 7804
## cor(Prime.Type,biasmatchmismatch) 0.69 1.00 7655
## cor(Intercept,Prime.Type:biasmatchmismatch) 0.71 1.00 10772
## cor(Prime.Type,Prime.Type:biasmatchmismatch) 0.74 1.00 8030
## cor(biasmatchmismatch,Prime.Type:biasmatchmismatch) 0.69 1.00 6905
## Tail_ESS
## sd(Intercept) 5032
## sd(Prime.Type) 3340
## sd(biasmatchmismatch) 4017
## sd(Prime.Type:biasmatchmismatch) 4229
## cor(Intercept,Prime.Type) 6101
## cor(Intercept,biasmatchmismatch) 4789
## cor(Prime.Type,biasmatchmismatch) 6239
## cor(Intercept,Prime.Type:biasmatchmismatch) 5760
## cor(Prime.Type,Prime.Type:biasmatchmismatch) 6046
## cor(biasmatchmismatch,Prime.Type:biasmatchmismatch) 6942
##
## Population-Level Effects:
## Estimate Est.Error l-90% CI u-90% CI Rhat Bulk_ESS
## Intercept -3.03 0.98 -4.66 -1.48 1.00 2717
## Prime.Type 0.87 1.11 -0.92 2.73 1.00 4092
## biasmatchmismatch -0.99 0.92 -2.60 0.33 1.00 3753
## Prime.Type:biasmatchmismatch 0.55 1.51 -1.79 3.16 1.00 4918
## Tail_ESS
## Intercept 4108
## Prime.Type 4534
## biasmatchmismatch 3870
## Prime.Type:biasmatchmismatch 4869
##
## 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).
m54mo3 <- add_criterion(m54mo3, "loo")
Chains look fine.
#plot(m54mo3)
hypothesis(m54mo3, "Prime.Type > 0")
## Hypothesis Tests for class b:
## Hypothesis Estimate Est.Error CI.Lower CI.Upper Evid.Ratio Post.Prob
## 1 (Prime.Type) > 0 0.87 1.11 -0.92 2.73 4.17 0.81
## Star
## 1
## ---
## 'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
## '*': For one-sided hypotheses, the posterior probability exceeds 95%;
## for two-sided hypotheses, the value tested against lies outside the 95%-CI.
## Posterior probabilities of point hypotheses assume equal prior probabilities.
hypothesis(m54mo3, "Prime.Type:biasmatchmismatch > 0") # Bias mismatch
## Hypothesis Tests for class b:
## Hypothesis Estimate Est.Error CI.Lower CI.Upper Evid.Ratio
## 1 (Prime.Type:biasm... > 0 0.55 1.51 -1.79 3.16 1.76
## Post.Prob Star
## 1 0.64
## ---
## 'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
## '*': For one-sided hypotheses, the posterior probability exceeds 95%;
## for two-sided hypotheses, the value tested against lies outside the 95%-CI.
## Posterior probabilities of point hypotheses assume equal prior probabilities.
m54mo3.model.check <- createDHARMa(
simulatedResponse = t(posterior_predict(m54mo3)),
observedResponse = DOD_54mo3$RESPONSE.CODE,
fittedPredictedResponse = apply(t(posterior_epred(m54mo3)), 1, mean),
integerResponse=TRUE
)
plot(m54mo3.model.check) # standard plots
testDispersion(m54mo3.model.check)
##
## DHARMa nonparametric dispersion test via sd of residuals fitted vs.
## simulated
##
## data: simulationOutput
## dispersion = 0.61694, p-value = 0.002
## alternative hypothesis: two.sided
When looking at the data aggregated across verb, it looks like there is the predicted interaction between prime structure and bias mismatch (for which these is no evidence in the model above). However, I think the apparent bias effect reflects aggregation bias (see prime surprisal effects by verb above). I’ll fit models (a) without random effects by verbs, (b) with just random intercepts by verbs, and (c) without random effects but without gave and threw.
m54mo_aggreg1<- brm(RESPONSE.CODE ~ 1 +Prime.Type*biasmatch + (1 + Prime.Type*biasmatch| Blind_ID), seed=12345, logit_prior, family="bernoulli", iter=4000, cores=4, control=list(adapt_delta=.9), file="54mo/output/m54mo_aggreg1.rds", data=DOD_54mo3)
summary(m54mo_aggreg1, prob=.9)
## Family: bernoulli
## Links: mu = logit
## Formula: RESPONSE.CODE ~ 1 + Prime.Type * biasmatch + (1 + Prime.Type * biasmatch | Blind_ID)
## Data: DOD_54mo3 (Number of observations: 479)
## Draws: 4 chains, each with iter = 4000; warmup = 2000; thin = 1;
## total post-warmup draws = 8000
##
## Group-Level Effects:
## ~Blind_ID (Number of levels: 87)
## Estimate Est.Error l-90% CI
## sd(Intercept) 0.85 0.40 0.15
## sd(Prime.Type) 0.69 0.51 0.05
## sd(biasmatchmismatch) 0.78 0.53 0.07
## sd(Prime.Type:biasmatchmismatch) 0.82 0.64 0.06
## cor(Intercept,Prime.Type) 0.02 0.44 -0.71
## cor(Intercept,biasmatchmismatch) 0.04 0.43 -0.65
## cor(Prime.Type,biasmatchmismatch) -0.01 0.44 -0.73
## cor(Intercept,Prime.Type:biasmatchmismatch) 0.02 0.44 -0.71
## cor(Prime.Type,Prime.Type:biasmatchmismatch) -0.15 0.46 -0.82
## cor(biasmatchmismatch,Prime.Type:biasmatchmismatch) -0.07 0.45 -0.76
## u-90% CI Rhat Bulk_ESS
## sd(Intercept) 1.50 1.00 1295
## sd(Prime.Type) 1.65 1.00 2049
## sd(biasmatchmismatch) 1.74 1.00 1154
## sd(Prime.Type:biasmatchmismatch) 2.05 1.00 2317
## cor(Intercept,Prime.Type) 0.73 1.00 4871
## cor(Intercept,biasmatchmismatch) 0.74 1.00 2997
## cor(Prime.Type,biasmatchmismatch) 0.71 1.00 2530
## cor(Intercept,Prime.Type:biasmatchmismatch) 0.74 1.00 6038
## cor(Prime.Type,Prime.Type:biasmatchmismatch) 0.66 1.00 4491
## cor(biasmatchmismatch,Prime.Type:biasmatchmismatch) 0.69 1.00 4919
## Tail_ESS
## sd(Intercept) 1546
## sd(Prime.Type) 3078
## sd(biasmatchmismatch) 2540
## sd(Prime.Type:biasmatchmismatch) 3038
## cor(Intercept,Prime.Type) 4365
## cor(Intercept,biasmatchmismatch) 3711
## cor(Prime.Type,biasmatchmismatch) 4227
## cor(Intercept,Prime.Type:biasmatchmismatch) 5143
## cor(Prime.Type,Prime.Type:biasmatchmismatch) 6185
## cor(biasmatchmismatch,Prime.Type:biasmatchmismatch) 5949
##
## Population-Level Effects:
## Estimate Est.Error l-90% CI u-90% CI Rhat Bulk_ESS
## Intercept -2.53 0.32 -3.10 -2.05 1.00 4183
## Prime.Type -0.35 0.52 -1.24 0.49 1.00 5423
## biasmatchmismatch -0.26 0.50 -1.13 0.50 1.00 2561
## Prime.Type:biasmatchmismatch 2.43 0.84 1.11 3.84 1.00 3418
## Tail_ESS
## Intercept 4181
## Prime.Type 5229
## biasmatchmismatch 3301
## Prime.Type:biasmatchmismatch 3441
##
## 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).
hypothesis(m54mo_aggreg1, "Prime.Type > 0")
## Hypothesis Tests for class b:
## Hypothesis Estimate Est.Error CI.Lower CI.Upper Evid.Ratio Post.Prob
## 1 (Prime.Type) > 0 -0.35 0.52 -1.24 0.49 0.32 0.24
## Star
## 1
## ---
## 'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
## '*': For one-sided hypotheses, the posterior probability exceeds 95%;
## for two-sided hypotheses, the value tested against lies outside the 95%-CI.
## Posterior probabilities of point hypotheses assume equal prior probabilities.
hypothesis(m54mo_aggreg1, "Prime.Type:biasmatchmismatch > 0")
## Hypothesis Tests for class b:
## Hypothesis Estimate Est.Error CI.Lower CI.Upper Evid.Ratio
## 1 (Prime.Type:biasm... > 0 2.43 0.84 1.11 3.84 887.89
## Post.Prob Star
## 1 1 *
## ---
## 'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
## '*': For one-sided hypotheses, the posterior probability exceeds 95%;
## for two-sided hypotheses, the value tested against lies outside the 95%-CI.
## Posterior probabilities of point hypotheses assume equal prior probabilities.
m54mo_aggreg2<- brm(RESPONSE.CODE ~ 1 +Prime.Type*biasmatch + (1 + Prime.Type*biasmatch| Blind_ID) + (1 | Target.verb), seed=12345, logit_prior, family="bernoulli", iter=4000, control=list(adapt_delta=.9), cores=4, file="54mo/output/m54mo_aggreg2.rds", data=DOD_54mo3)
summary(m54mo_aggreg2, prob=.9)
## Family: bernoulli
## Links: mu = logit
## Formula: RESPONSE.CODE ~ 1 + Prime.Type * biasmatch + (1 + Prime.Type * biasmatch | Blind_ID) + (1 | Target.verb)
## Data: DOD_54mo3 (Number of observations: 479)
## Draws: 4 chains, each with iter = 4000; warmup = 2000; thin = 1;
## total post-warmup draws = 8000
##
## Group-Level Effects:
## ~Blind_ID (Number of levels: 87)
## Estimate Est.Error l-90% CI
## sd(Intercept) 1.40 0.51 0.57
## sd(Prime.Type) 1.24 0.84 0.11
## sd(biasmatchmismatch) 0.92 0.65 0.09
## sd(Prime.Type:biasmatchmismatch) 1.10 0.84 0.09
## cor(Intercept,Prime.Type) -0.06 0.42 -0.73
## cor(Intercept,biasmatchmismatch) 0.08 0.42 -0.61
## cor(Prime.Type,biasmatchmismatch) 0.01 0.44 -0.73
## cor(Intercept,Prime.Type:biasmatchmismatch) 0.09 0.44 -0.65
## cor(Prime.Type,Prime.Type:biasmatchmismatch) -0.15 0.45 -0.82
## cor(biasmatchmismatch,Prime.Type:biasmatchmismatch) -0.00 0.44 -0.71
## u-90% CI Rhat Bulk_ESS
## sd(Intercept) 2.26 1.00 1837
## sd(Prime.Type) 2.79 1.00 1889
## sd(biasmatchmismatch) 2.14 1.00 1355
## sd(Prime.Type:biasmatchmismatch) 2.74 1.00 3047
## cor(Intercept,Prime.Type) 0.67 1.00 5462
## cor(Intercept,biasmatchmismatch) 0.74 1.00 5452
## cor(Prime.Type,biasmatchmismatch) 0.73 1.00 3381
## cor(Intercept,Prime.Type:biasmatchmismatch) 0.77 1.00 7549
## cor(Prime.Type,Prime.Type:biasmatchmismatch) 0.65 1.00 6580
## cor(biasmatchmismatch,Prime.Type:biasmatchmismatch) 0.71 1.00 6826
## Tail_ESS
## sd(Intercept) 1829
## sd(Prime.Type) 3548
## sd(biasmatchmismatch) 2663
## sd(Prime.Type:biasmatchmismatch) 4183
## cor(Intercept,Prime.Type) 5380
## cor(Intercept,biasmatchmismatch) 5823
## cor(Prime.Type,biasmatchmismatch) 4802
## cor(Intercept,Prime.Type:biasmatchmismatch) 5844
## cor(Prime.Type,Prime.Type:biasmatchmismatch) 5836
## cor(biasmatchmismatch,Prime.Type:biasmatchmismatch) 6623
##
## ~Target.verb (Number of levels: 6)
## Estimate Est.Error l-90% CI u-90% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 1.89 0.69 1.01 3.19 1.00 4358 6141
##
## Population-Level Effects:
## Estimate Est.Error l-90% CI u-90% CI Rhat Bulk_ESS
## Intercept -3.15 0.95 -4.71 -1.62 1.00 3159
## Prime.Type 1.00 0.74 -0.15 2.23 1.00 4843
## biasmatchmismatch -0.37 0.60 -1.42 0.53 1.00 4488
## Prime.Type:biasmatchmismatch 0.17 1.20 -1.75 2.13 1.00 4258
## Tail_ESS
## Intercept 3904
## Prime.Type 3699
## biasmatchmismatch 3769
## Prime.Type:biasmatchmismatch 4260
##
## 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).
hypothesis(m54mo_aggreg2, "Prime.Type > 0")
## Hypothesis Tests for class b:
## Hypothesis Estimate Est.Error CI.Lower CI.Upper Evid.Ratio Post.Prob
## 1 (Prime.Type) > 0 1 0.74 -0.15 2.23 12.36 0.93
## Star
## 1
## ---
## 'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
## '*': For one-sided hypotheses, the posterior probability exceeds 95%;
## for two-sided hypotheses, the value tested against lies outside the 95%-CI.
## Posterior probabilities of point hypotheses assume equal prior probabilities.
hypothesis(m54mo_aggreg2, "Prime.Type:biasmatchmismatch > 0")
## Hypothesis Tests for class b:
## Hypothesis Estimate Est.Error CI.Lower CI.Upper Evid.Ratio
## 1 (Prime.Type:biasm... > 0 0.17 1.2 -1.75 2.13 1.24
## Post.Prob Star
## 1 0.55
## ---
## 'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
## '*': For one-sided hypotheses, the posterior probability exceeds 95%;
## for two-sided hypotheses, the value tested against lies outside the 95%-CI.
## Posterior probabilities of point hypotheses assume equal prior probabilities.
m54mo_aggreg3<- DOD_54mo3 %>%
filter(Target.verb != "gave" & Target.verb != "threw") %>%
brm(RESPONSE.CODE ~ 1 +Prime.Type*biasmatch + (1 + Prime.Type*biasmatch| Blind_ID), seed=12345, logit_prior, family="bernoulli", iter=4000, cores=4, file="54mo/output/m54mo_aggreg3.rds", data=.)
summary(m54mo_aggreg3, prob=.9)
## Family: bernoulli
## Links: mu = logit
## Formula: RESPONSE.CODE ~ 1 + Prime.Type * biasmatch + (1 + Prime.Type * biasmatch | Blind_ID)
## Data: . (Number of observations: 307)
## Draws: 4 chains, each with iter = 4000; warmup = 2000; thin = 1;
## total post-warmup draws = 8000
##
## Group-Level Effects:
## ~Blind_ID (Number of levels: 87)
## Estimate Est.Error l-90% CI
## sd(Intercept) 1.14 0.63 0.17
## sd(Prime.Type) 1.29 0.91 0.10
## sd(biasmatchmismatch) 1.23 0.87 0.10
## sd(Prime.Type:biasmatchmismatch) 1.80 1.27 0.14
## cor(Intercept,Prime.Type) 0.10 0.44 -0.65
## cor(Intercept,biasmatchmismatch) 0.07 0.43 -0.65
## cor(Prime.Type,biasmatchmismatch) 0.08 0.44 -0.68
## cor(Intercept,Prime.Type:biasmatchmismatch) 0.02 0.43 -0.70
## cor(Prime.Type,Prime.Type:biasmatchmismatch) -0.07 0.45 -0.76
## cor(biasmatchmismatch,Prime.Type:biasmatchmismatch) -0.01 0.45 -0.74
## u-90% CI Rhat Bulk_ESS
## sd(Intercept) 2.24 1.00 1095
## sd(Prime.Type) 2.98 1.00 1871
## sd(biasmatchmismatch) 2.86 1.00 1296
## sd(Prime.Type:biasmatchmismatch) 4.17 1.00 1627
## cor(Intercept,Prime.Type) 0.77 1.00 3083
## cor(Intercept,biasmatchmismatch) 0.76 1.00 3445
## cor(Prime.Type,biasmatchmismatch) 0.77 1.00 2867
## cor(Intercept,Prime.Type:biasmatchmismatch) 0.73 1.00 4115
## cor(Prime.Type,Prime.Type:biasmatchmismatch) 0.69 1.00 3689
## cor(biasmatchmismatch,Prime.Type:biasmatchmismatch) 0.71 1.00 4110
## Tail_ESS
## sd(Intercept) 1426
## sd(Prime.Type) 3127
## sd(biasmatchmismatch) 2377
## sd(Prime.Type:biasmatchmismatch) 2363
## cor(Intercept,Prime.Type) 4144
## cor(Intercept,biasmatchmismatch) 4342
## cor(Prime.Type,biasmatchmismatch) 4562
## cor(Intercept,Prime.Type:biasmatchmismatch) 5042
## cor(Prime.Type,Prime.Type:biasmatchmismatch) 5116
## cor(biasmatchmismatch,Prime.Type:biasmatchmismatch) 5741
##
## Population-Level Effects:
## Estimate Est.Error l-90% CI u-90% CI Rhat Bulk_ESS
## Intercept -3.08 0.60 -4.20 -2.26 1.00 2408
## Prime.Type 0.34 0.82 -1.01 1.66 1.00 3383
## biasmatchmismatch -0.89 0.99 -2.76 0.43 1.00 2062
## Prime.Type:biasmatchmismatch 0.03 1.37 -2.17 2.30 1.00 3116
## Tail_ESS
## Intercept 3419
## Prime.Type 3942
## biasmatchmismatch 2100
## Prime.Type:biasmatchmismatch 3534
##
## 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).
hypothesis(m54mo_aggreg3, "Prime.Type > 0")
## Hypothesis Tests for class b:
## Hypothesis Estimate Est.Error CI.Lower CI.Upper Evid.Ratio Post.Prob
## 1 (Prime.Type) > 0 0.34 0.82 -1.01 1.66 2.15 0.68
## Star
## 1
## ---
## 'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
## '*': For one-sided hypotheses, the posterior probability exceeds 95%;
## for two-sided hypotheses, the value tested against lies outside the 95%-CI.
## Posterior probabilities of point hypotheses assume equal prior probabilities.
hypothesis(m54mo_aggreg3, "Prime.Type:biasmatchmismatch > 0")
## Hypothesis Tests for class b:
## Hypothesis Estimate Est.Error CI.Lower CI.Upper Evid.Ratio
## 1 (Prime.Type:biasm... > 0 0.03 1.37 -2.17 2.3 1.01
## Post.Prob Star
## 1 0.5
## ---
## 'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
## '*': For one-sided hypotheses, the posterior probability exceeds 95%;
## for two-sided hypotheses, the value tested against lies outside the 95%-CI.
## Posterior probabilities of point hypotheses assume equal prior probabilities.
DOD_54mo4 <-
DOD_54mo2 %>%
filter(Overlap=="No Overlap")
m54mo4 <- brm(RESPONSE.CODE ~ 1 +Prime.Type*biasmatch + (1 + Prime.Type*biasmatch| Blind_ID) + (1 + Prime.Type*biasmatch| Target.verb), seed=12345, logit_prior, family="bernoulli", iter=4000, control=list(adapt_delta=.9),cores=4, file="54mo/output/m54mo4.rds", data=DOD_54mo4)
summary(m54mo4, prob=.9)
## Family: bernoulli
## Links: mu = logit
## Formula: RESPONSE.CODE ~ 1 + Prime.Type * biasmatch + (1 + Prime.Type * biasmatch | Blind_ID) + (1 + Prime.Type * biasmatch | Target.verb)
## Data: DOD_54mo4 (Number of observations: 412)
## Draws: 4 chains, each with iter = 4000; warmup = 2000; thin = 1;
## total post-warmup draws = 8000
##
## Group-Level Effects:
## ~Blind_ID (Number of levels: 75)
## Estimate Est.Error l-90% CI
## sd(Intercept) 1.13 0.52 0.25
## sd(Prime.Type) 1.36 0.87 0.15
## sd(biasmatchmismatch) 1.12 0.71 0.12
## sd(Prime.Type:biasmatchmismatch) 1.24 0.92 0.10
## cor(Intercept,Prime.Type) -0.01 0.42 -0.70
## cor(Intercept,biasmatchmismatch) 0.11 0.41 -0.58
## cor(Prime.Type,biasmatchmismatch) 0.06 0.43 -0.67
## cor(Intercept,Prime.Type:biasmatchmismatch) 0.09 0.44 -0.66
## cor(Prime.Type,Prime.Type:biasmatchmismatch) -0.14 0.45 -0.81
## cor(biasmatchmismatch,Prime.Type:biasmatchmismatch) 0.03 0.45 -0.71
## u-90% CI Rhat Bulk_ESS
## sd(Intercept) 2.01 1.00 1249
## sd(Prime.Type) 2.93 1.00 2075
## sd(biasmatchmismatch) 2.42 1.00 1559
## sd(Prime.Type:biasmatchmismatch) 2.98 1.00 3267
## cor(Intercept,Prime.Type) 0.69 1.00 5718
## cor(Intercept,biasmatchmismatch) 0.76 1.00 4226
## cor(Prime.Type,biasmatchmismatch) 0.75 1.00 3402
## cor(Intercept,Prime.Type:biasmatchmismatch) 0.77 1.00 8717
## cor(Prime.Type,Prime.Type:biasmatchmismatch) 0.64 1.00 6976
## cor(biasmatchmismatch,Prime.Type:biasmatchmismatch) 0.75 1.00 6522
## Tail_ESS
## sd(Intercept) 1275
## sd(Prime.Type) 4145
## sd(biasmatchmismatch) 3757
## sd(Prime.Type:biasmatchmismatch) 4841
## cor(Intercept,Prime.Type) 5451
## cor(Intercept,biasmatchmismatch) 4763
## cor(Prime.Type,biasmatchmismatch) 5016
## cor(Intercept,Prime.Type:biasmatchmismatch) 5686
## cor(Prime.Type,Prime.Type:biasmatchmismatch) 6468
## cor(biasmatchmismatch,Prime.Type:biasmatchmismatch) 5967
##
## ~Target.verb (Number of levels: 6)
## Estimate Est.Error l-90% CI
## sd(Intercept) 1.98 0.74 1.00
## sd(Prime.Type) 1.31 0.91 0.13
## sd(biasmatchmismatch) 0.84 0.71 0.06
## sd(Prime.Type:biasmatchmismatch) 1.17 0.91 0.09
## cor(Intercept,Prime.Type) 0.10 0.41 -0.60
## cor(Intercept,biasmatchmismatch) 0.04 0.44 -0.69
## cor(Prime.Type,biasmatchmismatch) -0.05 0.45 -0.75
## cor(Intercept,Prime.Type:biasmatchmismatch) -0.02 0.44 -0.73
## cor(Prime.Type,Prime.Type:biasmatchmismatch) 0.01 0.45 -0.72
## cor(biasmatchmismatch,Prime.Type:biasmatchmismatch) -0.06 0.45 -0.77
## u-90% CI Rhat Bulk_ESS
## sd(Intercept) 3.37 1.00 5486
## sd(Prime.Type) 3.03 1.00 3755
## sd(biasmatchmismatch) 2.24 1.00 5548
## sd(Prime.Type:biasmatchmismatch) 2.93 1.00 6732
## cor(Intercept,Prime.Type) 0.75 1.00 11966
## cor(Intercept,biasmatchmismatch) 0.75 1.00 15059
## cor(Prime.Type,biasmatchmismatch) 0.70 1.00 8959
## cor(Intercept,Prime.Type:biasmatchmismatch) 0.71 1.00 14491
## cor(Prime.Type,Prime.Type:biasmatchmismatch) 0.74 1.00 10468
## cor(biasmatchmismatch,Prime.Type:biasmatchmismatch) 0.70 1.00 7465
## Tail_ESS
## sd(Intercept) 6007
## sd(Prime.Type) 4250
## sd(biasmatchmismatch) 4627
## sd(Prime.Type:biasmatchmismatch) 5155
## cor(Intercept,Prime.Type) 6006
## cor(Intercept,biasmatchmismatch) 5696
## cor(Prime.Type,biasmatchmismatch) 6420
## cor(Intercept,Prime.Type:biasmatchmismatch) 5512
## cor(Prime.Type,Prime.Type:biasmatchmismatch) 6561
## cor(biasmatchmismatch,Prime.Type:biasmatchmismatch) 6717
##
## Population-Level Effects:
## Estimate Est.Error l-90% CI u-90% CI Rhat Bulk_ESS
## Intercept -2.72 0.98 -4.33 -1.16 1.00 3574
## Prime.Type 0.87 1.05 -0.81 2.61 1.00 6421
## biasmatchmismatch -0.99 0.90 -2.56 0.36 1.00 4954
## Prime.Type:biasmatchmismatch 0.63 1.50 -1.65 3.16 1.00 6980
## Tail_ESS
## Intercept 5102
## Prime.Type 5307
## biasmatchmismatch 4619
## Prime.Type:biasmatchmismatch 4703
##
## 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).
m4 <- add_criterion(m54mo4, "loo")
hypothesis(m54mo4, "Prime.Type > 0")
## Hypothesis Tests for class b:
## Hypothesis Estimate Est.Error CI.Lower CI.Upper Evid.Ratio Post.Prob
## 1 (Prime.Type) > 0 0.87 1.05 -0.81 2.61 4.38 0.81
## Star
## 1
## ---
## 'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
## '*': For one-sided hypotheses, the posterior probability exceeds 95%;
## for two-sided hypotheses, the value tested against lies outside the 95%-CI.
## Posterior probabilities of point hypotheses assume equal prior probabilities.
hypothesis(m54mo4, "Prime.Type:biasmatchmismatch > 0") # Bias mismatch
## Hypothesis Tests for class b:
## Hypothesis Estimate Est.Error CI.Lower CI.Upper Evid.Ratio
## 1 (Prime.Type:biasm... > 0 0.63 1.5 -1.65 3.16 1.86
## Post.Prob Star
## 1 0.65
## ---
## 'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
## '*': For one-sided hypotheses, the posterior probability exceeds 95%;
## for two-sided hypotheses, the value tested against lies outside the 95%-CI.
## Posterior probabilities of point hypotheses assume equal prior probabilities.
m54mo4.model.check <- createDHARMa(
simulatedResponse = t(posterior_predict(m54mo4)),
observedResponse = DOD_54mo4$RESPONSE.CODE,
fittedPredictedResponse = apply(t(posterior_epred(m54mo4)), 1, mean),
integerResponse=TRUE
)
plot(m54mo4.model.check) # standard plots
testDispersion(m54mo4.model.check)
##
## DHARMa nonparametric dispersion test via sd of residuals fitted vs.
## simulated
##
## data: simulationOutput
## dispersion = 0.62811, p-value = 0.00175
## alternative hypothesis: two.sided
tab_model(m54mo1, m54mo2, m54mo3, m54mo4, show.ci=.90, show.re.var=FALSE, show.r2=FALSE, show.icc=FALSE, transform=NULL, collapse.ci=TRUE, p.style="scientific", ci.hyphen = " : ", dv.labels=c("Model 1", "Model 2", "Model 3", "Model 4"), string.stat="", file="Table54", bpe="mean",
pred.labels=c(
"Intercept",
"Prime (DOD)",
"Overlap",
"Prime * Overlap",
"Bias Mismatch (POD)",
"Prime * Bias Mismatch "
))
| Model 1 | Model 2 | Model 3 | Model 4 | |
|---|---|---|---|---|
| Predictors | Log-Odds | Log-Odds | Log-Odds | Log-Odds |
| Intercept |
-3.19 (-4.51 : -1.78) |
-2.80 (-4.16 : -1.43) |
-3.03 (-4.66 : -1.48) |
-2.72 (-4.33 : -1.16) |
| Prime (DOD) |
0.84 (-0.06 : 1.76) |
0.84 (-0.01 : 1.72) |
0.87 (-0.92 : 2.73) |
0.87 (-0.81 : 2.61) |
| Overlap |
-0.56 (-1.62 : 0.50) |
-0.53 (-1.54 : 0.52) |
||
| Prime * Overlap |
0.72 (-0.97 : 2.54) |
0.80 (-0.76 : 2.57) |
||
| Bias Mismatch (POD) |
-0.99 (-2.60 : 0.33) |
-0.99 (-2.56 : 0.36) |
||
| Prime * Bias Mismatch |
0.55 (-1.79 : 3.16) |
0.63 (-1.65 : 3.16) |
||
| N | 87 Blind_ID | 75 Blind_ID | 87 Blind_ID | 75 Blind_ID |
| 6 Target.verb | 6 Target.verb | 6 Target.verb | 6 Target.verb | |
| Observations | 1011 | 872 | 479 | 412 |
sessionInfo()
## R version 4.1.0 (2021-05-18)
## Platform: x86_64-w64-mingw32/x64 (64-bit)
## Running under: Windows 10 x64 (build 19045)
##
## Matrix products: default
##
## locale:
## [1] LC_COLLATE=Dutch_Netherlands.1252 LC_CTYPE=Dutch_Netherlands.1252
## [3] LC_MONETARY=Dutch_Netherlands.1252 LC_NUMERIC=C
## [5] LC_TIME=Dutch_Netherlands.1252
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] ggrepel_0.9.1 yarrr_0.1.5 circlize_0.4.15
## [4] BayesFactor_0.9.12-4.3 Matrix_1.4-0 coda_0.19-4
## [7] jpeg_0.1-9 zoo_1.8-9 lubridate_1.8.0
## [10] sjPlot_2.8.10 patchwork_1.1.1 ggpubr_0.4.0
## [13] GGally_2.1.2 dplyr_1.0.7 tidyr_1.2.0
## [16] ggplot2_3.3.5 DHARMa_0.4.6 tidybayes_3.0.2
## [19] rmarkdown_2.13 writexl_1.4.0 readxl_1.3.1
## [22] brms_2.17.0 Rcpp_1.0.8.3
##
## loaded via a namespace (and not attached):
## [1] backports_1.4.1 plyr_1.8.6 igraph_1.2.11
## [4] splines_4.1.0 svUnit_1.0.6 crosstalk_1.2.0
## [7] gap.datasets_0.0.5 TH.data_1.1-0 rstantools_2.1.1
## [10] inline_0.3.19 digest_0.6.29 foreach_1.5.2
## [13] htmltools_0.5.2 fansi_1.0.2 magrittr_2.0.1
## [16] checkmate_2.0.0 doParallel_1.0.17 modelr_0.1.8
## [19] RcppParallel_5.1.5 matrixStats_0.61.0 xts_0.12.1
## [22] sandwich_3.0-1 prettyunits_1.1.1 colorspace_2.0-3
## [25] ggdist_3.1.1 xfun_0.30 callr_3.7.0
## [28] crayon_1.5.0 jsonlite_1.8.0 lme4_1.1-28
## [31] iterators_1.0.14 survival_3.3-1 glue_1.6.2
## [34] gtable_0.3.1 emmeans_1.7.2 MatrixModels_0.5-0
## [37] sjstats_0.18.1 sjmisc_2.8.9 distributional_0.3.0
## [40] car_3.0-12 pkgbuild_1.3.1 rstan_2.21.3
## [43] shape_1.4.6 abind_1.4-5 scales_1.2.1
## [46] mvtnorm_1.1-3 DBI_1.1.2 ggeffects_1.1.1
## [49] rstatix_0.7.0 miniUI_0.1.1.1 performance_0.8.0
## [52] xtable_1.8-4 diffobj_0.3.5 stats4_4.1.0
## [55] StanHeaders_2.21.0-7 DT_0.21 datawizard_0.3.0
## [58] htmlwidgets_1.5.4 threejs_0.3.3 arrayhelpers_1.1-0
## [61] RColorBrewer_1.1-3 posterior_1.2.1 ellipsis_0.3.2
## [64] pkgconfig_2.0.3 reshape_0.8.8 loo_2.5.0
## [67] farver_2.1.0 qgam_1.3.4 sass_0.4.0
## [70] utf8_1.2.2 labeling_0.4.2 effectsize_0.6.0.1
## [73] tidyselect_1.1.2 rlang_1.0.6 reshape2_1.4.4
## [76] later_1.3.0 munsell_0.5.0 cellranger_1.1.0
## [79] tools_4.1.0 cli_3.3.0 generics_0.1.3
## [82] sjlabelled_1.1.8 broom_0.7.12 ggridges_0.5.3
## [85] evaluate_0.15 stringr_1.4.0 fastmap_1.1.0
## [88] yaml_2.3.5 processx_3.5.2 knitr_1.37
## [91] purrr_0.3.4 pbapply_1.5-0 nlme_3.1-155
## [94] mime_0.12 projpred_2.0.2 gap_1.3-1
## [97] compiler_4.1.0 bayesplot_1.9.0 shinythemes_1.2.0
## [100] rstudioapi_0.13 gamm4_0.2-6 ggsignif_0.6.3
## [103] tibble_3.1.6 bslib_0.3.1 stringi_1.7.6
## [106] highr_0.9 parameters_0.17.0 ps_1.6.0
## [109] Brobdingnag_1.2-8 lattice_0.20-45 nloptr_2.0.0
## [112] markdown_1.1 shinyjs_2.1.0 tensorA_0.36.2
## [115] vctrs_0.3.8 pillar_1.8.1 lifecycle_1.0.3
## [118] GlobalOptions_0.1.2 jquerylib_0.1.4 bridgesampling_1.1-2
## [121] estimability_1.3 insight_0.16.0 httpuv_1.6.5
## [124] R6_2.5.1 promises_1.2.0.1 gridExtra_2.3
## [127] codetools_0.2-18 boot_1.3-28 colourpicker_1.1.1
## [130] MASS_7.3-55 gtools_3.9.2 assertthat_0.2.1
## [133] withr_2.5.0 shinystan_2.6.0 multcomp_1.4-18
## [136] bayestestR_0.11.5 mgcv_1.8-39 parallel_4.1.0
## [139] grid_4.1.0 minqa_1.2.4 carData_3.0-5
## [142] shiny_1.7.1 base64enc_0.1-3 dygraphs_1.1.1.6