Skip to content

Instantly share code, notes, and snippets.

@bbolker
Created July 27, 2021 19:47
Show Gist options
  • Save bbolker/00de1e4dfdf0cf7ff56ecc451a8d09e2 to your computer and use it in GitHub Desktop.
Save bbolker/00de1e4dfdf0cf7ff56ecc451a8d09e2 to your computer and use it in GitHub Desktop.
exploring inference for 'separated' trunc distribution
```{r}
dd <- read_excel("Truncated.compois.problem.xlsx")
dd <- transform(dd, pop = factor(pop), morph = factor(morph))
m0 <- glmmTMB(mates ~ pop*morph,
family='truncated_compois', data=dd)
m1 <- update(m0, contrasts = list(pop = my_sum(dd$pop), morph = my_sum(dd$morph)))
tt <- purrr::map_dfr(list(treatment=m0,sum=m1),
tidy, effects = "fixed", .id = "contrasts", conf.int = TRUE) %>%
filter(term != "(Intercept)") %>%
mutate(across(term, ~gsub("(pop\\.?|sum_)","",.)))
```
```{r dwplot1}
ggplot(tt, aes(y = term, x = estimate, colour = contrasts)) +
scale_x_continuous(limits=c(-25,20), oob = scales::squish) +
geom_pointrange(aes(xmin=conf.low, xmax = conf.high), position = position_dodge(width = 0.5)) +
facet_wrap(~contrasts,scales = "free_y")
```
`drop1` gives identical results (overall significance of the interaction):
```{r drop1}
drop1(m0, test="Chisq")
drop1(m1, test="Chisq")
```
These are different even though they're using type II estimates ... ? (Could check and see if they're the same if we drop `ELS`)
```{r anova}
car::Anova(m0)
car::Anova(m1)
```
```{r}
m0B <- update(m0, . ~ morph + pop)
m1B <- update(m1, . ~ morph + pop)
drop1(m0B, test="Chisq")
drop1(m1B, test="Chisq")
```
Again identical (except warning from `m1B`, maybe because of weird custom contrasts?
```{r}
noELS <- droplevels(filter(dd, pop != "ELS"))
m2 <- update(m1, data = noELS,
contrasts = list(pop = my_sum(noELS$pop)))
dwplot(m2, effects = "fixed") + geom_vline(xintercept = 0, lty=2)
```
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment