Adapting inlaNMA to BRMS


#1

Hello everyone,

I am trying to adapt code from nmaINLA to allow for an easy interface for running and expanding these models from stan. Network meta-analysis is an extension of pairwise meta-analysis for more than two treatments. The code for these models is primarily based off code from NICE TSD 2 with some extensions to allow for inconsistency between direct (A vs B) and indirect (A vs B through C vs A and C vs B).

Stepping through the nma_inla function from their example code

library(nmaINLA)
library(INLA)
SmokdatINLA <- create_INLA_dat(dat = Smokdat, armVars = c(‘treatment’ = ‘t’, ‘responders’ = ‘r’
,‘sampleSize’ = ‘n’), nArmsVar = ‘na’)
if(requireNamespace(‘INLA’, quietly = TRUE)){
require(‘INLA’, quietly = TRUE)
fit.Smok.cons.INLA <- nma_inla(SmokdatINLA, likelihood = ‘binomial’, type = ‘consistency’,
tau.prior = ‘uniform’, tau.par = c(0, 5))
}

They create a long format data frame with study indicators for non-baseline arms (stored in het below) and indicators for multi-arm trials (g) and the basic parameters (d12 …).

study treatment responders sampleSize na baseline mu d12 d13 d14 g het inc
1 1 9 140 3 1 1 0 0 0 NA NA NA
1 3 23 140 3 1 1 0 1 0 1 1 1
1 4 10 138 3 1 1 0 0 1 2 1 1
2 2 11 78 3 2 2 0 0 0 NA NA NA
2 3 12 85 3 2 2 -1 1 0 1 2 2
2 4 29 170 3 2 2 -1 0 1 2 2 2
3 1 75 731 2 1 3 0 0 0 NA NA NA
3 3 363 714 2 1 3 0 1 0 1 3 3

The multi-arm bit is used to re-build equation 14 from NICE TSD2 (p.35). Ultimately, the model call to INLA is:

INLA::inla(stats::as.formula(inla.form), data = datINLA,
family = “binomial”, verbose = verbose, control.fixed = list(expand.factor.strategy = “inla”,
mean = fixed.par[1], prec = 1/fixed.par[2]), Ntrials = Ntrials,
control.compute = list(dic = TRUE, waic = TRUE, cpo = TRUE,
mlik = TRUE, config = TRUE), control.inla = list(strategy = inla.strategy,
correct = correct, correct.factor = correct.factor))

Where fixed.par is a vague prior on treatment effects, Ntrials is sample size in the binomial case and inla.form is

“Y ~ -1 + mu + d12+d13+d14 + f(het, model=“iid”, hyper = list(theta = list(prior = unif.prior.table)), group = g, control.group = list(model = “exchangeable”, hyper = list(rho = list(fixed = TRUE, initial = cor.inla.init))))”

Where Y is the number of events in arm i of study j and mu is the outcome in the baseline treatment of that arm. Where I’m stuck is how to translate the f() portion of the call to something for BRMS. I had originally assumed 1 | studyid would give me the single between-study variance parameter I was looking for but the multi-arm piece has me stumped.


#2

I am not used to reading INLA code. Could you write down the statistical model itself? Then I can better tell you whether such a model is possible in brms.


#3

I’m on a phone but it’s a random effects model here the civariance matrix is the kroneker product of an iid random effect and an exchangeable random effect (sigma^2R) where R is a correlation matrix where all the off diagonal elements are the same. The fixed=TRUE part means that the correlation parameter is fixed rather than inferred.

I’m travelling today but let me know if you need more.


#4

Thanks Dan! I don’t yet fully understand all the details, but based on what Dan writes,
you might be able to solve this via argument cov_ranef of brm, which allows you to specify fixed correlation matrices between random effects of the same level. For usage of cov_ranef you may also look at vignette("brms_phylogenetics") although the application presented there is very different from your use case.