Data Preparation

data(thrombolytic)

age <- rnorm(nrow(thrombolytic), 60, 10)

dich.slr <- data.prep(arm.data = cbind(thrombolytic, age),
                     varname.t = "treatment",
                     varname.s = "study")

Feasibility Assessment

 net.plot(dich.slr, node.scale = 1.5, edge.scale=0.5)

network.char <- net.tab(data = dich.slr,
                        outcome = "events",
                        N = "sampleSize", 
                        type.outcome = "binomial",
                        time = NULL)

Main analysis

nma.model() creates BUGS code and that will be put into nma.run() and analysed through JAGS (Plummer 2017). The reference parameter indicates the name of the treatment that will be seen as the ‘referent’ comparator, this is often a placebo of some sort. In our case, it is streptokinase (“SK”). Since our outcome is dichotomous, and we are not interested in event rates, we are using the “binomial” family. In our case, we want to compare relative risks, so we are using the \(log\) link. If you are interested in using an odds ratio, set link="logit".

fixed_effects_model <- nma.model(data=dich.slr,
                     outcome="events",
                     N="sampleSize",
                     reference="SK",
                     family="binomial",
                     link="log",
                     effects="fixed")

random_effects_model <- nma.model(data=dich.slr,
                     outcome="events",
                     N="sampleSize",
                     reference="SK",
                     family="binomial",
                     link="log",
                     effects="random")

If you want to review or modify the BUGS code, you can review it by outputting cat(fixed_effects_model$bugs) and cat(random_effects_model$bugs).

The next step is to run the NMA model using nma.run(). Since we are working in a Bayesian framework, we need to specify the number of adaptations, burn-ins, and iterations. A description of Bayesian MCMC is omitted here, we direct the reader to any introductory text on Bayesian Modelling (Lunn et al. 2012).

set.seed(20190829)
fixed_effects_results <- nma.run(fixed_effects_model,
                           n.adapt=1000,
                           n.burnin=1000,
                           n.iter=10000)

random_effects_results <- nma.run(random_effects_model,
                           n.adapt=1000,
                           n.burnin=1000,
                           n.iter=10000)

Assess model fit

par(mfrow = c(1,2))
nma.fit(fixed_effects_results, main = "Fixed Effects Model" )
nma.fit(random_effects_results, main= "Random Effects Model")

sucra.out <- nma.rank(fixed_effects_results, largerbetter=FALSE, sucra.palette= "Set1")
sucra.out$sucraplot

league.out <- nma.league(fixed_effects_results,  
                         central.tdcy="median",
                         order = sucra.out$order)
league.out$heatplot

nma.forest(fixed_effects_results,
           central.tdcy="median",
           comparator = "SK")

fe_inconsistency_model <- nma.model(data=dich.slr,
                     outcome="events",
                     N="sampleSize",
                     reference="SK",
                     family="binomial",
                     link="log",
                     type = "inconsistency",
                     effects="fixed")

fe_inconsistency_results <- nma.run(fe_inconsistency_model,
                           n.adapt=1000,
                           n.burnin=1000,
                           n.iter=10000)

par(mfrow = c(1,2))
fe_model_fit <- nma.fit(fixed_effects_results)
inconsist_model_fit <- nma.fit(fe_inconsistency_results)

nma.compare(fe_model_fit, inconsist_model_fit)