library(pcnetmeta)
## Loading required package: rjags
## Loading required package: coda
## Linked to JAGS 4.2.0
## Loaded modules: basemod,bugs
data(smoke)
# increase n.iter to reach convergence
# increase n.adapt to enhance efficiency
set.seed(1234)
nma.out <- nma.ab.bin(s.id, t.id, r, n, data = smoke,
trtname = c("NC", "SH", "IC", "GC"), param= "AR",
model = "het_cor", n.adapt = 400, n.iter = 100, n.chains = 1)
## Start running MCMC...
## Compiling model graph
## Resolving undeclared variables
## Allocating nodes
## Graph information:
## Observed stochastic nodes: 50
## Unobserved stochastic nodes: 29
## Total graph size: 611
##
## Initializing model
absolute.plot(nma.out)
## The absolute plot is saved in your working directory:
## C:/Users/HP/Documents
absolute.plot(nma.out, alphabetic = FALSE)
## The absolute plot is saved in your working directory:
## C:/Users/HP/Documents
data(smoke)
# increase n.iter to reach convergence
# increase n.adapt to enhance efficiency
set.seed(1234)
nma.out <- nma.ab.bin(s.id, t.id, r, n, data = smoke,
trtname = c("NC", "SH", "IC", "GC"), param= "LOR",
model = "het_cor", n.adapt = 400, n.iter = 100, n.chains = 1)
## Start running MCMC...
## Compiling model graph
## Resolving undeclared variables
## Allocating nodes
## Graph information:
## Observed stochastic nodes: 50
## Unobserved stochastic nodes: 29
## Total graph size: 611
##
## Initializing model
contrast.plot(nma.out)
## The contrast plot is saved in your working directory:
## C:/Users/HP/Documents
data(smoke)
# For the smoke cessation data,
# higher event rate indicates better treatment
# use the model = "het_cor"
#set.seed(1234)
#het.cor.out <- nma.ab.bin(s.id, t.id, r, n, data = smoke,
# trtname = c("NC", "SH", "IC", "GC"), param = c("AR", "OR", "RR", "LOR",
# "LRR", "RD", "rank.prob"), model = "het_cor", higher.better = TRUE,
# n.iter = 200000, n.thin = 1, conv.diag = TRUE, dic = TRUE,
# trace = c("AR", "LOR"), postdens = TRUE)
# use the model = "hom_eqcor"
# increase n.iter to reach convergence
# increase n.adapt to enhance efficiency
set.seed(1234)
hom.eqcor.out <- nma.ab.bin(s.id, t.id, r, n, data = smoke,
param = c("AR", "LRR"), model = "hom_eqcor", prior.type = "unif", c = 10,
higher.better = TRUE, n.adapt = 400, n.iter = 100, n.chains = 1)
## Start running MCMC...
## Compiling model graph
## Resolving undeclared variables
## Allocating nodes
## Graph information:
## Observed stochastic nodes: 50
## Unobserved stochastic nodes: 30
## Total graph size: 760
##
## Initializing model
data(parkinson)
# increase n.iter to reach convergence of MCMC
# increase n.adapt to enhance efficiency
set.seed(1234)
cont.out <- nma.ab.cont(s.id, t.id, mean, sd, n, data = parkinson,
param = c("mu", "diff"), model = "hom_eqcor", prior.type = "unif",
n.adapt = 200, n.iter = 100, n.chains = 1)
## Start running MCMC...
## Compiling model graph
## Resolving undeclared variables
## Allocating nodes
## Graph information:
## Observed stochastic nodes: 15
## Unobserved stochastic nodes: 14
## Total graph size: 362
##
## Initializing model
data(diabetes)
# increase n.iter to reach convergence of MCMC
# increase n.adapt to enhance efficiency
set.seed(1234)
followup.out <- nma.ab.followup(s.id, t.id, r, n, folup, data = diabetes,
model = "het_cor", n.adapt = 500, n.iter = 100, n.chains = 1)
## Start running MCMC...
## Compiling model graph
## Resolving undeclared variables
## Allocating nodes
## Graph information:
## Observed stochastic nodes: 48
## Unobserved stochastic nodes: 29
## Total graph size: 791
##
## Initializing model
data(dietaryfat)
# increase n.iter to reach convergence of MCMC
# increase n.adapt to enhance efficiency
set.seed(1234)
py.out <- nma.ab.py(s.id, t.id, r, py, data = dietaryfat, model = "het_cor",
n.adapt = 300, n.iter = 100, n.chains = 1)
## Start running MCMC...
## Compiling model graph
## Resolving undeclared variables
## Allocating nodes
## Graph information:
## Observed stochastic nodes: 21
## Unobserved stochastic nodes: 13
## Total graph size: 245
##
## Initializing model
data(smoke)
nma.networkplot(s.id, t.id, data = smoke, title = "Smoke Cessation",
trtname = c("NC", "SH", "IC", "GC"))
# NC: No contact; SH: Self-help
# IC: individual counselling; GC: group counselling
data(diabetes)
nma.networkplot(s.id, t.id, data = diabetes, title = "Diabetes",
trtname = c("Diuretic", "Placebo", "b-blocker", "CCB", "ACE inhibitor",
"ARB"))
data(smoke)
# increase n.iter to reach convergence
# increase n.adapt to enhance efficiency
set.seed(1234)
nma.out <- nma.ab.bin(s.id, t.id, r, n, data = smoke,
trtname = c("NC", "SH", "IC", "GC"), param= "rank.prob", model = "het_cor",
higher.better = TRUE, n.adapt = 400, n.iter = 100, n.chains = 1)
## Start running MCMC...
## Compiling model graph
## Resolving undeclared variables
## Allocating nodes
## Graph information:
## Observed stochastic nodes: 50
## Unobserved stochastic nodes: 29
## Total graph size: 646
##
## Initializing model
rank.prob(nma.out)
library("pcnetmeta")
data("smoke", package = "pcnetmeta")
head(smoke)
## s.id t.id r n
## 1 1 1 9 140
## 2 1 3 23 140
## 3 1 4 10 138
## 4 2 2 11 78
## 5 2 3 12 85
## 6 2 4 29 170
data("parkinson", package = "pcnetmeta")
parkinson
## s.id t.id mean sd n
## 1 1 1 -1.22 3.70 54
## 2 1 3 -1.53 4.28 95
## 3 2 1 -0.70 3.70 172
## 4 2 2 -2.40 3.40 173
## 5 3 1 -0.30 4.40 76
## 6 3 2 -2.60 4.30 71
## 7 3 4 -1.20 4.30 81
## 8 4 3 -0.24 3.00 128
## 9 4 4 -0.59 3.00 72
## 10 5 3 -0.73 3.00 80
## 11 5 4 -0.18 3.00 46
## 12 6 4 -2.20 2.31 137
## 13 6 5 -2.50 2.18 131
## 14 7 4 -1.80 2.48 154
## 15 7 5 -2.10 2.99 143
data("smoke", package = "pcnetmeta")
nma.networkplot(s.id, t.id, data = smoke, trtname = c("NC", "SH", "IC", "GC"))
data("parkinson", package = "pcnetmeta")
nma.networkplot(s.id, t.id, data = parkinson)
data("diabetes", package = "pcnetmeta")
nma.networkplot(s.id, t.id, data = diabetes, trtname = c("Diuretic", "Placebo", "BB", "CCB", "ACEI", "ARB"))
data("smoke", package = "pcnetmeta")
set.seed(12345)
smoke.out <- nma.ab.bin(s.id, t.id, r, n, data = smoke, trtname = c("NC", "SH", "IC", "GC"),
param = c("AR", "OR", "RR", "LOR", "LRR", "RD", "rank.prob"),
model = "het_cor", higher.better = TRUE, digits = 3,
n.adapt = 10, n.iter = 20, n.thin = 1, conv.diag = TRUE,
dic = TRUE, trace = "LOR", postdens = TRUE)
## Start running MCMC...
## Compiling model graph
## Resolving undeclared variables
## Allocating nodes
## Graph information:
## Observed stochastic nodes: 50
## Unobserved stochastic nodes: 29
## Total graph size: 646
##
## Initializing model
##
## Adaptation incomplete; users may increase n.adapt.
## NOTE: Stopping adaptation
##
##
## Start calculating MCMC convergence diagnostic statistics...
## Start calculating deviance information criterion statistics...
## Start saving trace plots...
## Start saving posterior density plot for absolute risk...
smoke.out$AbsoluteRisk
## $Mean_SD
## Mean (SD)
## NC 0.102 (0.012)
## SH 0.145 (0.017)
## IC 0.172 (0.014)
## GC 0.213 (0.025)
##
## $Median_CI
## Median (95% CI)
## NC 0.102 (0.083, 0.122)
## SH 0.146 (0.119, 0.173)
## IC 0.173 (0.150, 0.194)
## GC 0.213 (0.166, 0.260)
smoke.out$LogOddsRatio$Median_CI
## NC SH IC
## NC -- -0.402 (-0.729, -0.062) -0.566 (-0.982, -0.335)
## SH 0.402 (0.062, 0.729) -- -0.189 (-0.504, 0.077)
## IC 0.566 (0.335, 0.982) 0.189 (-0.077, 0.504) --
## GC 0.861 (0.524, 1.250) 0.475 (0.088, 0.843) 0.269 (-0.064, 0.494)
## GC
## NC -0.861 (-1.250, -0.524)
## SH -0.475 (-0.843, -0.088)
## IC -0.269 (-0.494, 0.064)
## GC --
smoke.out$TrtRankProb
## rank1 rank2 rank3 rank4
## NC 0.0000 0.0000 0.0000 1.0000
## SH 0.0000 0.1000 0.9000 0.0000
## IC 0.0500 0.8500 0.1000 0.0000
## GC 0.9500 0.0500 0.0000 0.0000
smoke.out$DIC
##
## D.bar 315.4477
## pD 75.4222
## DIC 390.8699
data("parkinson", package = "pcnetmeta")
set.seed(12345)
parkinson.out <- nma.ab.cont(s.id, t.id, mean, sd, n, data = parkinson,
model = "hom_eqcor", prior.type = "unif", digits = 3, n.adapt = 100,
n.iter = 100, n.thin = 1, conv.diag = TRUE, trace = "mu",
postdens = TRUE)
## Start running MCMC...
## Compiling model graph
## Resolving undeclared variables
## Allocating nodes
## Graph information:
## Observed stochastic nodes: 15
## Unobserved stochastic nodes: 14
## Total graph size: 402
##
## Initializing model
##
## Start calculating MCMC convergence diagnostic statistics...
## Start saving trace plots...
## Start saving posterior density plot for treatment effect...
parkinson.out$TrtEffect$Median_CI
## Median (95% CI)
## Trt1 -0.699 (-1.420, 0.135)
## Trt2 -2.580 (-3.520, -1.830)
## Trt3 -1.070 (-1.830, -0.358)
## Trt4 -1.410 (-1.980, -0.857)
## Trt5 -1.910 (-2.760, -1.140)
parkinson.out$EffectDiff$Median_CI
## Trt1 Trt2 Trt3
## Trt1 -- 1.830 (0.846, 3.520) 0.372 (-0.631, 1.370)
## Trt2 -1.830 (-3.520, -0.846) -- -1.560 (-2.750, -0.293)
## Trt3 -0.372 (-1.370, 0.631) 1.560 (0.293, 2.750) --
## Trt4 -0.729 (-1.680, 0.133) 1.200 (0.346, 2.090) -0.330 (-1.260, 0.488)
## Trt5 -1.260 (-2.640, -0.084) 0.646 (-0.257, 1.740) -0.918 (-2.040, 0.448)
## Trt4 Trt5
## Trt1 0.729 (-0.133, 1.680) 1.260 (0.084, 2.640)
## Trt2 -1.200 (-2.090, -0.346) -0.646 (-1.740, 0.257)
## Trt3 0.330 (-0.488, 1.260) 0.918 (-0.448, 2.040)
## Trt4 -- 0.522 (-0.143, 1.260)
## Trt5 -0.522 (-1.260, 0.143) --
data("diabetes", package = "pcnetmeta")
set.seed(12345)
diabetes.out <- nma.ab.followup(s.id, t.id, r, n, folup, data = diabetes,
trtname = c("Diuretic", "Placebo", "BB", "CCB", "ACEI", "ARB"),
model = "het_cor", digits = 3, n.adapt = 100, n.iter = 200,
n.thin = 2, conv.diag = TRUE, trace = "lograte", postdens = TRUE)
## Start running MCMC...
## Compiling model graph
## Resolving undeclared variables
## Allocating nodes
## Graph information:
## Observed stochastic nodes: 48
## Unobserved stochastic nodes: 29
## Total graph size: 791
##
## Initializing model
##
## Adaptation incomplete; users may increase n.adapt.
## NOTE: Stopping adaptation
##
##
## Start calculating MCMC convergence diagnostic statistics...
## Start saving trace plots...
## Start saving posterior density plot for log rate...
diabetes.out$LogRateRatio$Median_CI
## Diuretic Placebo
## Diuretic -- 0.571 (-0.325, 2.600)
## Placebo -0.571 (-2.600, 0.325) --
## BB -0.454 (-2.960, 0.649) -0.012 (-1.230, 1.240)
## CCB -1.890 (-3.260, 0.741) -1.450 (-2.400, 1.680)
## ACEI -1.160 (-3.130, -0.236) -0.476 (-1.840, 0.227)
## ARB -0.443 (-2.850, 5.660) -0.156 (-2.000, 6.630)
## BB CCB
## Diuretic 0.454 (-0.649, 2.960) 1.890 (-0.741, 3.260)
## Placebo 0.012 (-1.240, 1.230) 1.450 (-1.680, 2.400)
## BB -- 1.710 (-2.300, 2.810)
## CCB -1.710 (-2.810, 2.300) --
## ACEI -0.571 (-2.250, 1.010) 0.440 (-1.860, 1.580)
## ARB -0.353 (-2.450, 7.090) 1.510 (-0.132, 6.330)
## ACEI ARB
## Diuretic 1.160 (0.236, 3.130) 0.443 (-5.660, 2.850)
## Placebo 0.476 (-0.227, 1.840) 0.156 (-6.630, 2.000)
## BB 0.571 (-1.010, 2.250) 0.353 (-7.090, 2.450)
## CCB -0.440 (-1.580, 1.860) -1.510 (-6.330, 0.132)
## ACEI -- -0.272 (-6.700, 0.712)
## ARB 0.272 (-0.712, 6.700) --
absolute.plot(smoke.out, width = 5, height = 1.5)
## The absolute plot is saved in your working directory:
## C:/Users/HP/Documents
absolute.plot(parkinson.out, width = 5, height = 1.5)
## The absolute plot is saved in your working directory:
## C:/Users/HP/Documents
absolute.plot(diabetes.out, width = 8, height = 2.5)
## The absolute plot is saved in your working directory:
## C:/Users/HP/Documents
contrast.plot(smoke.out, reference = "NC", width = 5, height = 1.5)
## The contrast plot is saved in your working directory:
## C:/Users/HP/Documents
contrast.plot(parkinson.out, reference = "Trt1", width = 5, height = 1.5)
## The contrast plot is saved in your working directory:
## C:/Users/HP/Documents
contrast.plot(diabetes.out, reference = "Placebo", width = 8, height = 2.5)
## The contrast plot is saved in your working directory:
## C:/Users/HP/Documents
rank.prob(smoke.out, cex.axis = 2, cex.lab = 2)
rank.prob(parkinson.out, cex.axis = 2, cex.lab = 2)
rank.prob(diabetes.out, cex.axis = 1, cex.lab = 2)
