nma.ab generates a summary file, containing estimated event rate ("probt"), risk ratio ("RR"), risk difference ("RD"), odds ratio ("OR") etc., using the arm-based approach, proposed by Zhang et al (2013) for network meta-analysis. Also, it can provide DIC statistics for model selection, and give trace plot and posterior density plot for checking if the Markov chain is convergent.nma.ab(s.id, t.id, event.n, total.n, o.path = getwd(), f.name = "",
hom = FALSE, R, param = c("probt","RR","RD","OR","rk","best"),
ity = "estimate", n.iter = 200000, n.burnin = floor(n.iter/2),
n.chains = 3, 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 pFALSE (the default), this function will use a heterogeneous model; otherwise, a homogeneous model would be applied. See "Details" for the JAGS models.tN by tN covariance matrix for Wishart prior in the heterogeneous model (hom = FALSE), where tN is the number of treatments. The default is a matrix with diagonal elements being 1 and off-mu in the JAGS model. It can be set as "estimate" (the default) and "same". See "Details".n.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.FALSE (the default), the function would not give the trace plot; otherwise, trace plots for estimated event rate ("probt") would be given, and they would be saved in the working directory where o.path indicates.FALSE (the default), the function would not give the posterior density plot; otherwise, posterior density plots for estimated event rate ("probt") would be given, and they would be saved in the working directory where o.pahtnma.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, estimated 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.
The DIC file contains the value of pD and DIC; the trace plot and posterior density plot file contain corresponding plots for the node "probt" (estimated event rate).ity is specified as "estimate", the initial value for fixed effect would be estimated (in probit scale) as the proportion of sum of event numbers for a treatment over sum of total numbers in corresponding study; otherwise, the initial value of mu in the JAGS model would be set as 0. Setting ity as "estimate" may be helpful for the Markov chain converging more quickly.
The heterogeneous model (hom = FALSE) defines the following JAGS model:
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, 0.001)
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)
}
The homogeneous model (hom = TRUE) defines the following JAGS model:
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(1, 1)
for(j in 1:tN){
mu[j] ~ dnorm(0, 0.001)
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)
}data(Ara09)
attach(Ara09)
nma.ab(s.id = Study.ID, t.id = Treatment, event.n = r, total.n = n,
f.name = "Ara09_", n.iter = 2000, dic = FALSE)
detach(Ara09)
data(Lam07)
attach(Lam07$data)
nma.ab(s.id = Study.ID, t.id = Treatment, event.n = r, total.n = n, hom = TRUE,
f.name = "Lam07_", ity = "same", n.iter = 2000, n.chains = 5,
param = c("probt", "RR", "RD", "OR", "rk", "best", "mu", "sigma", "vi"),
trace = TRUE, postdens = TRUE)
detach(Lam07$data)Run the code above in your browser using DataLab