Skip to content

Instantly share code, notes, and snippets.

@bbolker
Last active December 21, 2023 20:01
Show Gist options
  • Star 1 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save bbolker/7694c4e88e2b90c7ae4d4d91142aec04 to your computer and use it in GitHub Desktop.
Save bbolker/7694c4e88e2b90c7ae4d4d91142aec04 to your computer and use it in GitHub Desktop.
attempt at matching metafor results with glmmTMB
library(brms)
library(metafor)
library(broom)
library(broom.mixed)
set.seed(101)
y <- rnorm(10)
sei <- rgamma(10, shape = 1)
## metafor: fixed-effect meta-analysis
mf1 <- rma(y, sei = sei, method = "FE")
## glmmTMB:
## fix all dispersions (standard errors) to known values
## note that we have to specify `betad` as log VARIANCES not
## log SDs (see https://github.com/glmmTMB/glmmTMB/issues/955) --
## this may change in the future
library(glmmTMB)
dd <- data.frame(y, sei, id = factor(1:10))
gt1 <- glmmTMB(y ~ 1,
## different SD value for each `id`
dispformula = ~0+id,
data = dd,
## fix SDs, don't estimate them
map = list(betad = factor(rep(NA,10))),
## specify SD values (actually variances)
start = list(betad = 2*log(dd$sei)))
## brms equivalent
br1 <- brm(y | se(sei) ~ 1, data = dd)
(list(glmmTMB = gt1, brms = br1, metafor = mf1)
|> purrr::map_dfr(tidy, .id = "pkg")
|> dplyr::filter(is.na(effect) | effect == "fixed")
|> dplyr::select(pkg, estimate, std.error)
)
## pkg estimate std.error
## 1 glmmTMB 0.898 0.0389
## 2 brms 0.898 0.0399
## 3 metafor 0.898 0.0389
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment