nma.ab generates a summary file for effect sizes by conducting the arm-based approach, proposed by Zhang et al (2014), in network meta-analysis. The generated summary file contains odds ratio (OR) and patient-centered paramenters, such as risk ratio (RR), risk difference (RD), and absolute risk (AR). Also, it can provide DIC statistics for checking the goodness of fit; give trace plots to check the MCMC convergence; generate posterior density plot for the population-averaged event rates.nma.ab(s.id, t.id, event.n, total.n, o.path = getwd(), f.name = "",
model = "het1", sigma2.a = 0.001, sigma2.b = 0.001, mu.prec = 0.001,
R, param = c("probt","RR","RD","OR","rk","best"), trtname,
n.iter = 200000, n.burnin = floor(n.iter/2), n.chains = 2,
n.thin = max(1, floor((n.iter - n.burnin)/100000)), dic = TRUE,
trace = FALSE, postdens = FALSE)dic = TRUE), trace plot (if trace = TRUE) and posterior density plot (if p"hom" (homogeneous variance, denoted as model HOM), "het1" (model HET1), or "het2" (mosigma2.a and sigma2.b are set as the default, the prior is weatN by tN covariance matrix for Wishart prior in the model HET2 (model = "het2"), where tN is the number of treatments. The default is a matrix with diagonal elements being 1 and off-diagont.id. It is optional, and if specified, the names would be shown in the posterior density plot for event rate (when postn.iter/2.TRUE (the default), the function would generate a file containing the DIC statistics, and a node named "deviance" would be contained in the summary result file; otherwise, the DIC statistics would not be calculated.TRUE, the function would save the trace plots for monitored nodes in "probt", "RR", "RD", and "OR". The default is FALSE.TRUE, the function would save the posterior density plot for population-averaged event rate (the node "probt"). The default is FALSE.nma.ab generates a summary statistics file using the arm-based method. Furthermore, this function would give a DIC statistics file if dic = TRUE, a trace plot file if trace = TRUE, a posterior density file if postdens = TRUE.
In the summary file, each row contains statistics for corresponding OR, RD, RR, population-averaged event rate ("probt"), rank of treatment ("rk"), probability of being the best treatment ("best"), etc. Note that RR[i, j], RD[i, j] or OR[i, j] means that treatment i is compared with treatment j, e.g., RD[i,j] = probt[i] - probt[j]. The columns show the statistics of these nodes, including mean, standard deviance, 2.5% percentile, median, and 97.5% percentile. Also, potential scale reduction factor (PSRF) in Gelman and Rubin's MCMC convergence diagnostic is saved in the column "Rhat" for each node.
The DIC file contains the value of pD and DIC; the trace plot file shows the traces for each node in "probt", "RR", "RD" and "OR"; posterior density plot file contains the combined posterior density for population-averaged event rate.model = "hom") considers a common variance for the random effects, and it assumes that the random effects for different treatments are the same in each study. The JAGS model is given as follows:
model{
for(i in 1:sN){
p[i] <- phi(mu[t[i]] + sigma*vi[s[i]])
r[i] ~ dbin(p[i], totaln[i])
}
for(j in 1:tS){
vi[j] ~ dnorm(0, 1)
}
sigma <- 1/sqrt(tau)
tau ~ dgamma(sigma2.a, sigma2.b)
for(j in 1:tN){
mu[j] ~ dnorm(0, mu.prec)
probt[j] <- phi(mu[j]/sqrt(1 + 1/tau))
}
for(j in 1:tN){
for(k in 1:tN){
RR[j, k] <- probt[j]/probt[k]
RD[j, k] <- probt[j] - probt[k]
OR[j, k] <- probt[j]/(1 - probt[j])/probt[k]*(1 - probt[k])
}
}
rk[1:tN] <- tN + 1 - rank(probt[])
best[1:tN] <- equals(rk[], 1)
}
The first heterogeneous model (model = "het1") accounts for the heterogeneity for the variances of random effects, but it still assumes that the random effects for different treatments are the same in each study. The following shows the corresponding JAGS model:
model{
for(i in 1:sN){
p[i] <- phi(mu[t[i]] + sigma[t[i]]*vi[s[i]])
r[i] ~ dbin(p[i], totaln[i])
}
for(j in 1:tS){
vi[j] ~ dnorm(0, 1)
}
for(j in 1:tN){
mu[j] ~ dnorm(0, mu.prec)
tau[j] ~ dgamma(sigma2.a, sigma2.b)
sigma[j] <- 1/sqrt(tau[j])
probt[j] <- phi(mu[j]/sqrt(1 + 1/tau[j]))
}
for(j in 1:tN){
for(k in 1:tN){
RR[j, k] <- probt[j]/probt[k]
RD[j, k] <- probt[j] - probt[k]
OR[j, k] <- probt[j]/(1 - probt[j])/probt[k]*(1 - probt[k])
}
}
rk[1:tN] <- tN + 1 - rank(probt[])
best[1:tN] <- equals(rk[], 1)
}
The second heterogeneous model (model = "het2") covers the most general cases, and it employs a Wishart prior for the inverse of covariance matrix for random effects. The JAGS model is defined as follows:
model{
for(i in 1:sN){
p[i] <- phi(mu[t[i]] + vi[s[i], t[i]])
r[i] ~ dbin(p[i], totaln[i])
}
for(j in 1:tS){
vi[j, 1:tN] ~ dmnorm(mn[1:tN], T[1:tN, 1:tN])
}
invT[1:tN, 1:tN] <- inverse(T[,])
for(j in 1:tN){
mu[j] ~ dnorm(0, mu.prec)
sigma[j] <- sqrt(invT[j, j])
probt[j] <- phi(mu[j]/sqrt(1 + invT[j, j]))
}
T[1:tN, 1:tN] ~ dwish(R[1:tN, 1:tN], tN)
for(j in 1:tN){
for(k in 1:tN){
RR[j, k] <- probt[j]/probt[k]
RD[j, k] <- probt[j] - probt[k]
OR[j, k] <- probt[j]/(1 - probt[j])/probt[k]*(1 - probt[k])
}
}
rk[1:tN] <- tN + 1 - rank(probt[])
best[1:tN] <- equals(rk[], 1)
}library(lattice)
data(Middleton10)
attach(Middleton10)
set.seed(12345)
nma.ab(s.id = sid, t.id = tid, event.n = r, total.n = n,
model = "hom", f.name = "Middleton10_hom_", n.iter = 500)
nma.ab(s.id = sid, t.id = tid, event.n = r, total.n = n,
model = "het1", f.name = "Middleton10_het1_", n.iter = 500,
trtname = c("FG", "H", "SG", "M"), trace = TRUE, postdens = TRUE)
detach(Middleton10)Run the code above in your browser using DataLab