diff --git a/tests/testthat/test-example_smoking.R b/tests/testthat/test-example_smoking.R index 36e95472..74ce8109 100644 --- a/tests/testthat/test-example_smoking.R +++ b/tests/testthat/test-example_smoking.R @@ -8,10 +8,10 @@ skip_on_cran() params <- list(run_tests = FALSE) -## ----code=readLines("children/knitr_setup.R"), include=FALSE-------------------------------------- +## ---- code=readLines("children/knitr_setup.R"), include=FALSE------------------------------------- -## ----eval = FALSE--------------------------------------------------------------------------------- +## ---- eval = FALSE-------------------------------------------------------------------------------- ## library(multinma) ## options(mc.cores = parallel::detectCores()) @@ -37,7 +37,7 @@ smknet <- set_agd_arm(smoking, smknet -## ----eval=FALSE----------------------------------------------------------------------------------- +## ---- eval=FALSE---------------------------------------------------------------------------------- ## plot(smknet, weight_edges = TRUE, weight_nodes = TRUE) ## ----smoking_network_plot, echo=FALSE, fig.width=8, fig.height=6, out.width="100%"---------------- @@ -63,7 +63,7 @@ smkfit <- nma(smknet, smkfit -## ----eval=FALSE----------------------------------------------------------------------------------- +## ---- eval=FALSE---------------------------------------------------------------------------------- ## # Not run ## print(smkfit, pars = c("d", "tau", "mu", "delta")) @@ -122,7 +122,7 @@ summary(smk_nodesplit) ## ----smk_nodesplit, fig.width = 7----------------------------------------------------------------- plot(smk_nodesplit) + - ggplot2::theme(legend.position = "bottom", legend.direct = "horizontal") + ggplot2::theme(legend.position = "bottom", legend.direction = "horizontal") ## ----smoking_releff, fig.height=4.5--------------------------------------------------------------- diff --git a/vignettes/example_smoking.html b/vignettes/example_smoking.html index 8db770ed..9ccfa5fd 100644 --- a/vignettes/example_smoking.html +++ b/vignettes/example_smoking.html @@ -362,8 +362,8 @@

Example: Smoking cessation

-
library(multinma)
-options(mc.cores = parallel::detectCores())
+
library(multinma)
+options(mc.cores = parallel::detectCores())
#> For execution on a local, multicore CPU with excess RAM we recommend calling
 #> options(mc.cores = parallel::detectCores())
 #> 
@@ -371,17 +371,19 @@ 

Example: Smoking cessation

#> The following objects are masked from 'package:stats': #> #> dgamma, pgamma, qgamma
-

This vignette describes the analysis of smoking cessation data (Hasselblad 1998), replicating the -analysis in NICE Technical Support Document 4 (Dias et al. 2011). The -data are available in this package as smoking:

-
head(smoking)
-#>   studyn trtn                   trtc  r   n
-#> 1      1    1        No intervention  9 140
-#> 2      1    3 Individual counselling 23 140
-#> 3      1    4      Group counselling 10 138
-#> 4      2    2              Self-help 11  78
-#> 5      2    3 Individual counselling 12  85
-#> 6      2    4      Group counselling 29 170
+

This vignette describes the analysis of smoking cessation data (Hasselblad +1998), replicating the analysis in NICE Technical Support +Document 4 (Dias et al. +2011). The data are available in this package as +smoking:

+
head(smoking)
+#>   studyn trtn                   trtc  r   n
+#> 1      1    1        No intervention  9 140
+#> 2      1    3 Individual counselling 23 140
+#> 3      1    4      Group counselling 10 138
+#> 4      2    2              Self-help 11  78
+#> 5      2    3 Individual counselling 12  85
+#> 6      2    4      Group counselling 29 170

Setting up the network

We begin by setting up the network. We have arm-level count data @@ -389,37 +391,37 @@

Setting up the network

(n) in each arm, so we use the function set_agd_arm(). Treatment “No intervention” is set as the network reference treatment.

-
smknet <- set_agd_arm(smoking, 
-                      study = studyn,
-                      trt = trtc,
-                      r = r, 
-                      n = n,
-                      trt_ref = "No intervention")
-smknet
-#> A network with 24 AgD studies (arm-based).
-#> 
-#> ------------------------------------------------------- AgD studies (arm-based) ---- 
-#>  Study Treatment arms                                                 
-#>  1     3: No intervention | Group counselling | Individual counselling
-#>  2     3: Group counselling | Individual counselling | Self-help      
-#>  3     2: No intervention | Individual counselling                    
-#>  4     2: No intervention | Individual counselling                    
-#>  5     2: No intervention | Individual counselling                    
-#>  6     2: No intervention | Individual counselling                    
-#>  7     2: No intervention | Individual counselling                    
-#>  8     2: No intervention | Individual counselling                    
-#>  9     2: No intervention | Individual counselling                    
-#>  10    2: No intervention | Self-help                                 
-#>  ... plus 14 more studies
-#> 
-#>  Outcome type: count
-#> ------------------------------------------------------------------------------------
-#> Total number of treatments: 4
-#> Total number of studies: 24
-#> Reference treatment is: No intervention
-#> Network is connected
+
smknet <- set_agd_arm(smoking, 
+                      study = studyn,
+                      trt = trtc,
+                      r = r, 
+                      n = n,
+                      trt_ref = "No intervention")
+smknet
+#> A network with 24 AgD studies (arm-based).
+#> 
+#> ------------------------------------------------------- AgD studies (arm-based) ---- 
+#>  Study Treatment arms                                                 
+#>  1     3: No intervention | Group counselling | Individual counselling
+#>  2     3: Group counselling | Individual counselling | Self-help      
+#>  3     2: No intervention | Individual counselling                    
+#>  4     2: No intervention | Individual counselling                    
+#>  5     2: No intervention | Individual counselling                    
+#>  6     2: No intervention | Individual counselling                    
+#>  7     2: No intervention | Individual counselling                    
+#>  8     2: No intervention | Individual counselling                    
+#>  9     2: No intervention | Individual counselling                    
+#>  10    2: No intervention | Self-help                                 
+#>  ... plus 14 more studies
+#> 
+#>  Outcome type: count
+#> ------------------------------------------------------------------------------------
+#> Total number of treatments: 4
+#> Total number of studies: 24
+#> Reference treatment is: No intervention
+#> Network is connected

Plot the network structure.

-
plot(smknet, weight_edges = TRUE, weight_nodes = TRUE)
+
plot(smknet, weight_edges = TRUE, weight_nodes = TRUE)

@@ -431,180 +433,178 @@

Random effects NMA

for the between-study heterogeneity standard deviation \(\tau\). We can examine the range of parameter values implied by these prior distributions with the summary() method:

-
summary(normal(scale = 100))
-#> A Normal prior distribution: location = 0, scale = 100.
-#> 50% of the prior density lies between -67.45 and 67.45.
-#> 95% of the prior density lies between -196 and 196.
-
summary(half_normal(scale = 5))
-#> A half-Normal prior distribution: location = 0, scale = 5.
-#> 50% of the prior density lies between 0 and 3.37.
-#> 95% of the prior density lies between 0 and 9.8.
+
summary(normal(scale = 100))
+#> A Normal prior distribution: location = 0, scale = 100.
+#> 50% of the prior density lies between -67.45 and 67.45.
+#> 95% of the prior density lies between -196 and 196.
+
summary(half_normal(scale = 5))
+#> A half-Normal prior distribution: location = 0, scale = 5.
+#> 50% of the prior density lies between 0 and 3.37.
+#> 95% of the prior density lies between 0 and 9.8.

The model is fitted using the nma() function. By default, this will use a Binomial likelihood and a logit link function, auto-detected from the data.

-
smkfit <- nma(smknet, 
-              trt_effects = "random",
-              prior_intercept = normal(scale = 100),
-              prior_trt = normal(scale = 100),
-              prior_het = normal(scale = 5))
+
smkfit <- nma(smknet, 
+              trt_effects = "random",
+              prior_intercept = normal(scale = 100),
+              prior_trt = normal(scale = 100),
+              prior_het = normal(scale = 5))

Basic parameter summaries are given by the print() method:

-
smkfit
-#> A random effects NMA with a binomial likelihood (logit link).
-#> Inference for Stan model: binomial_1par.
-#> 4 chains, each with iter=2000; warmup=1000; thin=1; 
-#> post-warmup draws per chain=1000, total post-warmup draws=4000.
-#> 
-#>                               mean se_mean   sd     2.5%      25%      50%      75%    97.5% n_eff
-#> d[Group counselling]          1.11    0.01 0.44     0.27     0.81     1.09     1.40     2.01  1998
-#> d[Individual counselling]     0.84    0.01 0.25     0.39     0.67     0.83     0.99     1.36  1206
-#> d[Self-help]                  0.50    0.01 0.40    -0.27     0.24     0.48     0.75     1.31  1916
-#> lp__                      -5919.72    0.20 6.50 -5933.29 -5923.82 -5919.41 -5915.14 -5907.86  1012
-#> tau                           0.84    0.01 0.19     0.54     0.70     0.82     0.95     1.28   918
-#>                           Rhat
-#> d[Group counselling]         1
-#> d[Individual counselling]    1
-#> d[Self-help]                 1
-#> lp__                         1
-#> tau                          1
-#> 
-#> Samples were drawn using NUTS(diag_e) at Tue Jan  9 17:49:33 2024.
-#> For each parameter, n_eff is a crude measure of effective sample size,
-#> and Rhat is the potential scale reduction factor on split chains (at 
-#> convergence, Rhat=1).
+
smkfit
+#> A random effects NMA with a binomial likelihood (logit link).
+#> Inference for Stan model: binomial_1par.
+#> 4 chains, each with iter=2000; warmup=1000; thin=1; 
+#> post-warmup draws per chain=1000, total post-warmup draws=4000.
+#> 
+#>                               mean se_mean   sd     2.5%      25%      50%      75%    97.5% n_eff
+#> d[Group counselling]          1.11    0.01 0.43     0.32     0.81     1.09     1.38     2.00  1681
+#> d[Individual counselling]     0.84    0.01 0.23     0.41     0.68     0.83     0.99     1.33  1238
+#> d[Self-help]                  0.49    0.01 0.40    -0.29     0.24     0.49     0.74     1.28  1572
+#> lp__                      -5767.91    0.21 6.43 -5781.12 -5772.07 -5767.63 -5763.34 -5756.32   902
+#> tau                           0.84    0.01 0.19     0.54     0.70     0.81     0.95     1.30  1204
+#>                           Rhat
+#> d[Group counselling]         1
+#> d[Individual counselling]    1
+#> d[Self-help]                 1
+#> lp__                         1
+#> tau                          1
+#> 
+#> Samples were drawn using NUTS(diag_e) at Wed Mar  6 12:04:00 2024.
+#> For each parameter, n_eff is a crude measure of effective sample size,
+#> and Rhat is the potential scale reduction factor on split chains (at 
+#> convergence, Rhat=1).

By default, summaries of the study-specific intercepts \(\mu_j\) and study-specific relative effects \(\delta_{jk}\) are hidden, but could be examined by changing the pars argument:

-
# Not run
-print(smkfit, pars = c("d", "tau", "mu", "delta"))
+
# Not run
+print(smkfit, pars = c("d", "tau", "mu", "delta"))

The prior and posterior distributions can be compared visually using the plot_prior_posterior() function:

-
plot_prior_posterior(smkfit)
-

+
plot_prior_posterior(smkfit)
+

By default, this displays all model parameters given prior distributions (in this case \(d_k\), \(\mu_j\), and \(\tau\)), but this may be changed using the prior argument:

-
plot_prior_posterior(smkfit, prior = "het")
-

+
plot_prior_posterior(smkfit, prior = "het")
+

Model fit can be checked using the dic() function

-
(dic_consistency <- dic(smkfit))
-#> Residual deviance: 54 (on 50 data points)
-#>                pD: 43.7
-#>               DIC: 97.7
+
(dic_consistency <- dic(smkfit))
+#> Residual deviance: 53.8 (on 50 data points)
+#>                pD: 43.6
+#>               DIC: 97.5

and the residual deviance contributions examined with the corresponding plot() method

-
plot(dic_consistency)
-

+
plot(dic_consistency)
+

Overall model fit seems to be adequate, with almost all points showing good fit (mean residual deviance contribution of 1). The only two points with higher residual deviance (i.e. worse fit) correspond to the two zero counts in the data:

-
smoking[smoking$r == 0, ]
-#>    studyn trtn            trtc r  n
-#> 13      6    1 No intervention 0 33
-#> 31     15    1 No intervention 0 20
+
smoking[smoking$r == 0, ]
+#>    studyn trtn            trtc r  n
+#> 13      6    1 No intervention 0 33
+#> 31     15    1 No intervention 0 20

Checking for inconsistency

Note: The results of the inconsistency models here -are slightly different to those of Dias et al. (2010, -2011), although the -overall conclusions are the same. This is due to the presence of -multi-arm trials and a different ordering of treatments, meaning that -inconsistency is parameterised differently within the multi-arm trials. -The same results as Dias et al. are obtained if the network is instead -set up with trtn as the treatment variable.

+are slightly different to those of Dias et al. (2010, 2011), although the overall conclusions are +the same. This is due to the presence of multi-arm trials and a +different ordering of treatments, meaning that inconsistency is +parameterised differently within the multi-arm trials. The same results +as Dias et al. are obtained if the network is instead set up with +trtn as the treatment variable.

Unrelated mean effects

-

We first fit an unrelated mean effects (UME) model (Dias et al. -2011) to assess the consistency assumption. Again, we use the -function nma(), but now with the argument +

We first fit an unrelated mean effects (UME) model (Dias et al. 2011) to +assess the consistency assumption. Again, we use the function +nma(), but now with the argument consistency = "ume".

-
smkfit_ume <- nma(smknet, 
-                  consistency = "ume",
-                  trt_effects = "random",
-                  prior_intercept = normal(scale = 100),
-                  prior_trt = normal(scale = 100),
-                  prior_het = normal(scale = 5))
-smkfit_ume
-#> A random effects NMA with a binomial likelihood (logit link).
-#> An inconsistency model ('ume') was fitted.
-#> Inference for Stan model: binomial_1par.
-#> 4 chains, each with iter=2000; warmup=1000; thin=1; 
-#> post-warmup draws per chain=1000, total post-warmup draws=4000.
-#> 
-#>                                                     mean se_mean   sd     2.5%      25%      50%
-#> d[Group counselling vs. No intervention]            1.15    0.02 0.81    -0.35     0.62     1.11
-#> d[Individual counselling vs. No intervention]       0.92    0.01 0.27     0.41     0.75     0.91
-#> d[Self-help vs. No intervention]                    0.33    0.01 0.59    -0.83    -0.04     0.33
-#> d[Individual counselling vs. Group counselling]    -0.26    0.01 0.60    -1.46    -0.65    -0.26
-#> d[Self-help vs. Group counselling]                 -0.61    0.01 0.71    -2.02    -1.05    -0.61
-#> d[Self-help vs. Individual counselling]             0.15    0.02 1.07    -2.02    -0.52     0.15
-#> lp__                                            -5933.38    0.20 6.27 -5946.38 -5937.54 -5933.05
-#> tau                                                 0.94    0.01 0.23     0.59     0.78     0.91
-#>                                                      75%    97.5% n_eff Rhat
-#> d[Group counselling vs. No intervention]            1.64     2.93  2470 1.00
-#> d[Individual counselling vs. No intervention]       1.09     1.51  1094 1.00
-#> d[Self-help vs. No intervention]                    0.70     1.55  2000 1.00
-#> d[Individual counselling vs. Group counselling]     0.13     0.93  2775 1.00
-#> d[Self-help vs. Group counselling]                 -0.16     0.81  2530 1.00
-#> d[Self-help vs. Individual counselling]             0.82     2.27  3983 1.00
-#> lp__                                            -5928.99 -5921.69  1025 1.01
-#> tau                                                 1.07     1.47  1203 1.00
-#> 
-#> Samples were drawn using NUTS(diag_e) at Tue Jan  9 17:49:53 2024.
-#> For each parameter, n_eff is a crude measure of effective sample size,
-#> and Rhat is the potential scale reduction factor on split chains (at 
-#> convergence, Rhat=1).
+
smkfit_ume <- nma(smknet, 
+                  consistency = "ume",
+                  trt_effects = "random",
+                  prior_intercept = normal(scale = 100),
+                  prior_trt = normal(scale = 100),
+                  prior_het = normal(scale = 5))
+smkfit_ume
+#> A random effects NMA with a binomial likelihood (logit link).
+#> An inconsistency model ('ume') was fitted.
+#> Inference for Stan model: binomial_1par.
+#> 4 chains, each with iter=2000; warmup=1000; thin=1; 
+#> post-warmup draws per chain=1000, total post-warmup draws=4000.
+#> 
+#>                                                     mean se_mean   sd     2.5%      25%      50%
+#> d[Group counselling vs. No intervention]            1.11    0.02 0.78    -0.33     0.57     1.08
+#> d[Individual counselling vs. No intervention]       0.89    0.01 0.27     0.37     0.71     0.89
+#> d[Self-help vs. No intervention]                    0.35    0.01 0.58    -0.80     0.00     0.34
+#> d[Individual counselling vs. Group counselling]    -0.28    0.01 0.61    -1.45    -0.68    -0.29
+#> d[Self-help vs. Group counselling]                 -0.60    0.02 0.69    -1.94    -1.05    -0.62
+#> d[Self-help vs. Individual counselling]             0.17    0.02 1.06    -1.97    -0.49     0.16
+#> lp__                                            -5765.20    0.20 6.39 -5778.70 -5769.44 -5764.84
+#> tau                                                 0.93    0.01 0.22     0.58     0.77     0.90
+#>                                                      75%    97.5% n_eff Rhat
+#> d[Group counselling vs. No intervention]            1.60     2.74  2282    1
+#> d[Individual counselling vs. No intervention]       1.07     1.43   961    1
+#> d[Self-help vs. No intervention]                    0.72     1.50  1765    1
+#> d[Individual counselling vs. Group counselling]     0.12     0.95  2376    1
+#> d[Self-help vs. Group counselling]                 -0.16     0.77  1801    1
+#> d[Self-help vs. Individual counselling]             0.83     2.24  3429    1
+#> lp__                                            -5760.67 -5753.72   982    1
+#> tau                                                 1.04     1.44   984    1
+#> 
+#> Samples were drawn using NUTS(diag_e) at Wed Mar  6 12:04:49 2024.
+#> For each parameter, n_eff is a crude measure of effective sample size,
+#> and Rhat is the potential scale reduction factor on split chains (at 
+#> convergence, Rhat=1).

Comparing the model fit statistics

-
dic_consistency
-#> Residual deviance: 54 (on 50 data points)
-#>                pD: 43.7
-#>               DIC: 97.7
-(dic_ume <- dic(smkfit_ume))
-#> Residual deviance: 53.5 (on 50 data points)
-#>                pD: 45
-#>               DIC: 98.5
+
dic_consistency
+#> Residual deviance: 53.8 (on 50 data points)
+#>                pD: 43.6
+#>               DIC: 97.5
+(dic_ume <- dic(smkfit_ume))
+#> Residual deviance: 53.5 (on 50 data points)
+#>                pD: 44.7
+#>               DIC: 98.2

We see that there is little to choose between the two models. However, it is also important to examine the individual contributions to model fit of each data point under the two models (a so-called “dev-dev” plot). Passing two nma_dic objects produced by the dic() function to the plot() method produces this dev-dev plot:

-
plot(dic_consistency, dic_ume, point_alpha = 0.5, interval_alpha = 0.2)
-

+
plot(dic_consistency, dic_ume, point_alpha = 0.5, interval_alpha = 0.2)
+

All points lie roughly on the line of equality, so there is no evidence for inconsistency here.

Node-splitting

-

Another method for assessing inconsistency is node-splitting (Dias et al. -2011, 2010). -Whereas the UME model assesses inconsistency globally, node-splitting -assesses inconsistency locally for each potentially inconsistent -comparison (those with both direct and indirect evidence) in turn.

+

Another method for assessing inconsistency is node-splitting (Dias et al. 2011, 2010). Whereas the UME model assesses +inconsistency globally, node-splitting assesses inconsistency locally +for each potentially inconsistent comparison (those with both direct and +indirect evidence) in turn.

Node-splitting can be performed using the nma() function with the argument consistency = "nodesplit". By default, all possible comparisons will be split (as determined by the get_nodesplits() function). Alternatively, a specific comparison or comparisons to split can be provided to the nodesplit argument.

-
smk_nodesplit <- nma(smknet, 
-                     consistency = "nodesplit",
-                     trt_effects = "random",
-                     prior_intercept = normal(scale = 100),
-                     prior_trt = normal(scale = 100),
-                     prior_het = normal(scale = 5))
-#> Fitting model 1 of 7, node-split: Group counselling vs. No intervention
-#> Fitting model 2 of 7, node-split: Individual counselling vs. No intervention
-#> Fitting model 3 of 7, node-split: Self-help vs. No intervention
-#> Fitting model 4 of 7, node-split: Individual counselling vs. Group counselling
-#> Fitting model 5 of 7, node-split: Self-help vs. Group counselling
-#> Fitting model 6 of 7, node-split: Self-help vs. Individual counselling
-#> Fitting model 7 of 7, consistency model
+
smk_nodesplit <- nma(smknet, 
+                     consistency = "nodesplit",
+                     trt_effects = "random",
+                     prior_intercept = normal(scale = 100),
+                     prior_trt = normal(scale = 100),
+                     prior_het = normal(scale = 5))
+#> Fitting model 1 of 7, node-split: Group counselling vs. No intervention
+#> Fitting model 2 of 7, node-split: Individual counselling vs. No intervention
+#> Fitting model 3 of 7, node-split: Self-help vs. No intervention
+#> Fitting model 4 of 7, node-split: Individual counselling vs. Group counselling
+#> Fitting model 5 of 7, node-split: Self-help vs. Group counselling
+#> Fitting model 6 of 7, node-split: Self-help vs. Individual counselling
+#> Fitting model 7 of 7, consistency model

The summary() method summarises the node-splitting results, displaying the direct and indirect estimates \(d_\mathrm{dir}\) and \(d_\mathrm{ind}\) from each node-split model, the network estimate \(d_\mathrm{net}\) from the consistency @@ -614,104 +614,104 @@

Node-splitting

standard deviation \(\tau\) under each node-split model and under the consistency model is also displayed. The DIC model fit statistics are also provided.

-
summary(smk_nodesplit)
-#> Node-splitting models fitted for 6 comparisons.
-#> 
-#> ------------------------------ Node-split Group counselling vs. No intervention ---- 
-#> 
-#>                  mean   sd  2.5%   25%   50%  75% 97.5% Bulk_ESS Tail_ESS Rhat
-#> d_net            1.09 0.44  0.28  0.79  1.07 1.37  1.97     1793     2275    1
-#> d_dir            1.06 0.75 -0.36  0.55  1.03 1.54  2.63     3096     2875    1
-#> d_ind            1.12 0.55 -0.01  0.76  1.12 1.48  2.19     1517     1709    1
-#> omega           -0.05 0.91 -1.81 -0.66 -0.08 0.54  1.82     2295     2580    1
-#> tau              0.87 0.20  0.56  0.73  0.84 0.98  1.32     1342     1858    1
-#> tau_consistency  0.84 0.19  0.55  0.71  0.82 0.94  1.29     1092     1786    1
-#> 
-#> Residual deviance: 53.9 (on 50 data points)
-#>                pD: 44.1
-#>               DIC: 98
-#> 
-#> Bayesian p-value: 0.93
-#> 
-#> ------------------------- Node-split Individual counselling vs. No intervention ---- 
-#> 
-#>                 mean   sd  2.5%   25%  50%  75% 97.5% Bulk_ESS Tail_ESS Rhat
-#> d_net           0.84 0.24  0.38  0.68 0.83 0.99  1.36     1176     1770 1.00
-#> d_dir           0.89 0.26  0.40  0.72 0.88 1.05  1.42     1218     1860 1.00
-#> d_ind           0.56 0.68 -0.74  0.12 0.55 0.98  1.99     1165     1674 1.00
-#> omega           0.33 0.70 -1.08 -0.10 0.35 0.79  1.65     1163     1632 1.00
-#> tau             0.86 0.19  0.55  0.72 0.83 0.97  1.28     1118     2046 1.01
-#> tau_consistency 0.84 0.19  0.55  0.71 0.82 0.94  1.29     1092     1786 1.00
-#> 
-#> Residual deviance: 54.2 (on 50 data points)
-#>                pD: 44.1
-#>               DIC: 98.3
-#> 
-#> Bayesian p-value: 0.6
-#> 
-#> -------------------------------------- Node-split Self-help vs. No intervention ---- 
-#> 
-#>                  mean   sd  2.5%   25%   50%  75% 97.5% Bulk_ESS Tail_ESS Rhat
-#> d_net            0.48 0.40 -0.28  0.21  0.48 0.73  1.28     1698     1876    1
-#> d_dir            0.33 0.54 -0.71 -0.01  0.33 0.68  1.41     3339     2529    1
-#> d_ind            0.70 0.65 -0.58  0.28  0.69 1.11  2.03     2165     2405    1
-#> omega           -0.36 0.83 -2.09 -0.88 -0.35 0.17  1.26     2351     2589    1
-#> tau              0.87 0.19  0.56  0.73  0.84 0.98  1.31     1311     2303    1
-#> tau_consistency  0.84 0.19  0.55  0.71  0.82 0.94  1.29     1092     1786    1
-#> 
-#> Residual deviance: 53.7 (on 50 data points)
-#>                pD: 44.1
-#>               DIC: 97.8
-#> 
-#> Bayesian p-value: 0.64
-#> 
-#> ----------------------- Node-split Individual counselling vs. Group counselling ---- 
-#> 
-#>                  mean   sd  2.5%   25%   50%   75% 97.5% Bulk_ESS Tail_ESS Rhat
-#> d_net           -0.25 0.41 -1.06 -0.52 -0.24  0.03  0.55     2718     2818    1
-#> d_dir           -0.10 0.51 -1.09 -0.43 -0.10  0.23  0.90     4213     3127    1
-#> d_ind           -0.53 0.61 -1.76 -0.92 -0.51 -0.14  0.67     1866     2406    1
-#> omega            0.43 0.68 -0.88 -0.02  0.42  0.86  1.78     1934     1802    1
-#> tau              0.87 0.20  0.56  0.73  0.84  0.98  1.35     1341     1985    1
-#> tau_consistency  0.84 0.19  0.55  0.71  0.82  0.94  1.29     1092     1786    1
-#> 
-#> Residual deviance: 54 (on 50 data points)
-#>                pD: 44.4
-#>               DIC: 98.5
-#> 
-#> Bayesian p-value: 0.52
-#> 
-#> ------------------------------------ Node-split Self-help vs. Group counselling ---- 
-#> 
-#>                  mean   sd  2.5%   25%   50%   75% 97.5% Bulk_ESS Tail_ESS Rhat
-#> d_net           -0.61 0.49 -1.58 -0.93 -0.61 -0.30  0.36     2737     2551    1
-#> d_dir           -0.63 0.66 -1.98 -1.05 -0.63 -0.21  0.68     3657     3014    1
-#> d_ind           -0.62 0.67 -2.05 -1.04 -0.61 -0.18  0.70     2132     2344    1
-#> omega           -0.02 0.88 -1.75 -0.58 -0.02  0.56  1.72     2421     2701    1
-#> tau              0.86 0.19  0.57  0.72  0.84  0.97  1.29     1013     1808    1
-#> tau_consistency  0.84 0.19  0.55  0.71  0.82  0.94  1.29     1092     1786    1
-#> 
-#> Residual deviance: 54 (on 50 data points)
-#>                pD: 44.2
-#>               DIC: 98.2
-#> 
-#> Bayesian p-value: 0.99
-#> 
-#> ------------------------------- Node-split Self-help vs. Individual counselling ---- 
-#> 
-#>                  mean   sd  2.5%   25%   50%   75% 97.5% Bulk_ESS Tail_ESS Rhat
-#> d_net           -0.36 0.40 -1.15 -0.62 -0.37 -0.11  0.45     2188     2043 1.00
-#> d_dir            0.06 0.65 -1.24 -0.36  0.06  0.50  1.35     3068     2587 1.00
-#> d_ind           -0.61 0.51 -1.62 -0.94 -0.61 -0.28  0.38     1883     2496 1.00
-#> omega            0.68 0.80 -0.91  0.15  0.67  1.20  2.27     2164     2423 1.00
-#> tau              0.86 0.19  0.56  0.72  0.83  0.96  1.30     1066     1156 1.01
-#> tau_consistency  0.84 0.19  0.55  0.71  0.82  0.94  1.29     1092     1786 1.00
-#> 
-#> Residual deviance: 53.8 (on 50 data points)
-#>                pD: 44.1
-#>               DIC: 97.9
-#> 
-#> Bayesian p-value: 0.38
+
summary(smk_nodesplit)
+#> Node-splitting models fitted for 6 comparisons.
+#> 
+#> ------------------------------ Node-split Group counselling vs. No intervention ---- 
+#> 
+#>                  mean   sd  2.5%   25%   50%  75% 97.5% Bulk_ESS Tail_ESS Rhat
+#> d_net            1.10 0.43  0.27  0.82  1.09 1.38  1.98     2360     2515 1.00
+#> d_dir            1.05 0.74 -0.33  0.55  1.01 1.51  2.61     3575     3129 1.00
+#> d_ind            1.13 0.54  0.10  0.78  1.12 1.48  2.22     1793     2010 1.00
+#> omega           -0.08 0.91 -1.81 -0.68 -0.11 0.49  1.77     2659     2417 1.00
+#> tau              0.86 0.20  0.56  0.72  0.83 0.97  1.34     1129     1958 1.01
+#> tau_consistency  0.84 0.18  0.55  0.71  0.82 0.94  1.26     1380     2233 1.00
+#> 
+#> Residual deviance: 54.3 (on 50 data points)
+#>                pD: 44.2
+#>               DIC: 98.5
+#> 
+#> Bayesian p-value: 0.9
+#> 
+#> ------------------------- Node-split Individual counselling vs. No intervention ---- 
+#> 
+#>                 mean   sd  2.5%   25%  50%  75% 97.5% Bulk_ESS Tail_ESS Rhat
+#> d_net           0.85 0.24  0.38  0.68 0.84 1.00  1.35     1430     2196    1
+#> d_dir           0.88 0.25  0.41  0.71 0.87 1.04  1.37     2196     2641    1
+#> d_ind           0.56 0.65 -0.68  0.13 0.54 0.98  1.89     1590     1992    1
+#> omega           0.32 0.67 -1.00 -0.13 0.32 0.78  1.60     1654     2031    1
+#> tau             0.85 0.19  0.55  0.72 0.83 0.96  1.29     1322     1893    1
+#> tau_consistency 0.84 0.18  0.55  0.71 0.82 0.94  1.26     1380     2233    1
+#> 
+#> Residual deviance: 54.4 (on 50 data points)
+#>                pD: 44.3
+#>               DIC: 98.7
+#> 
+#> Bayesian p-value: 0.63
+#> 
+#> -------------------------------------- Node-split Self-help vs. No intervention ---- 
+#> 
+#>                  mean   sd  2.5%   25%   50%  75% 97.5% Bulk_ESS Tail_ESS Rhat
+#> d_net            0.51 0.40 -0.27  0.24  0.50 0.77  1.33     2229     2107    1
+#> d_dir            0.34 0.55 -0.73 -0.01  0.34 0.68  1.46     3389     2939    1
+#> d_ind            0.70 0.63 -0.52  0.30  0.69 1.09  1.98     2487     2488    1
+#> omega           -0.36 0.82 -2.00 -0.88 -0.37 0.17  1.27     2531     2682    1
+#> tau              0.88 0.20  0.56  0.73  0.85 0.99  1.33     1412     1750    1
+#> tau_consistency  0.84 0.18  0.55  0.71  0.82 0.94  1.26     1380     2233    1
+#> 
+#> Residual deviance: 53.5 (on 50 data points)
+#>                pD: 44
+#>               DIC: 97.5
+#> 
+#> Bayesian p-value: 0.65
+#> 
+#> ----------------------- Node-split Individual counselling vs. Group counselling ---- 
+#> 
+#>                  mean   sd  2.5%   25%   50%   75% 97.5% Bulk_ESS Tail_ESS Rhat
+#> d_net           -0.25 0.40 -1.05 -0.51 -0.26  0.01  0.53     3147     2866    1
+#> d_dir           -0.10 0.48 -1.05 -0.40 -0.10  0.22  0.82     3613     3179    1
+#> d_ind           -0.55 0.61 -1.78 -0.94 -0.54 -0.16  0.62     1543     1646    1
+#> omega            0.45 0.66 -0.85  0.02  0.45  0.87  1.79     1692     1832    1
+#> tau              0.87 0.20  0.55  0.73  0.84  0.98  1.33     1215     1784    1
+#> tau_consistency  0.84 0.18  0.55  0.71  0.82  0.94  1.26     1380     2233    1
+#> 
+#> Residual deviance: 53.7 (on 50 data points)
+#>                pD: 44.1
+#>               DIC: 97.8
+#> 
+#> Bayesian p-value: 0.47
+#> 
+#> ------------------------------------ Node-split Self-help vs. Group counselling ---- 
+#> 
+#>                  mean   sd  2.5%   25%   50%   75% 97.5% Bulk_ESS Tail_ESS Rhat
+#> d_net           -0.59 0.48 -1.54 -0.89 -0.59 -0.29  0.35     3241     2871    1
+#> d_dir           -0.60 0.64 -1.85 -1.03 -0.62 -0.18  0.66     4527     3486    1
+#> d_ind           -0.62 0.69 -2.02 -1.06 -0.59 -0.17  0.70     2189     2282    1
+#> omega            0.01 0.90 -1.68 -0.59 -0.01  0.57  1.82     2583     2545    1
+#> tau              0.88 0.20  0.58  0.74  0.85  0.99  1.33     1305     2176    1
+#> tau_consistency  0.84 0.18  0.55  0.71  0.82  0.94  1.26     1380     2233    1
+#> 
+#> Residual deviance: 54 (on 50 data points)
+#>                pD: 44.4
+#>               DIC: 98.4
+#> 
+#> Bayesian p-value: 1
+#> 
+#> ------------------------------- Node-split Self-help vs. Individual counselling ---- 
+#> 
+#>                  mean   sd  2.5%   25%   50%   75% 97.5% Bulk_ESS Tail_ESS Rhat
+#> d_net           -0.34 0.41 -1.17 -0.60 -0.33 -0.08  0.47     2744     2564    1
+#> d_dir            0.06 0.66 -1.28 -0.36  0.06  0.49  1.35     3032     2784    1
+#> d_ind           -0.61 0.53 -1.66 -0.96 -0.60 -0.27  0.42     1647     2275    1
+#> omega            0.67 0.82 -0.93  0.14  0.66  1.19  2.31     2064     2246    1
+#> tau              0.85 0.19  0.55  0.72  0.82  0.96  1.30     1018     1891    1
+#> tau_consistency  0.84 0.18  0.55  0.71  0.82  0.94  1.26     1380     2233    1
+#> 
+#> Residual deviance: 54 (on 50 data points)
+#>                pD: 44.1
+#>               DIC: 98.1
+#> 
+#> Bayesian p-value: 0.38

The DIC of each inconsistency model is unchanged from the consistency model, no node-splits result in reduced heterogeneity standard deviation \(\tau\) compared to the consistency @@ -724,59 +724,59 @@

Node-splitting

information for the Individual counselling vs. No intervention comparison, so the network (consistency) estimate is very similar to the direct estimate for this comparison.

-
plot(smk_nodesplit) +
-  ggplot2::theme(legend.position = "bottom", legend.direct = "horizontal")
-

+
plot(smk_nodesplit) +
+  ggplot2::theme(legend.position = "bottom", legend.direction = "horizontal")
+

Further results

Pairwise relative effects, for all pairwise contrasts with all_contrasts = TRUE.

-
(smk_releff <- relative_effects(smkfit, all_contrasts = TRUE))
-#>                                                  mean   sd  2.5%   25%   50%   75% 97.5% Bulk_ESS
-#> d[Group counselling vs. No intervention]         1.11 0.44  0.27  0.81  1.09  1.40  2.01     2031
-#> d[Individual counselling vs. No intervention]    0.84 0.25  0.39  0.67  0.83  0.99  1.36     1221
-#> d[Self-help vs. No intervention]                 0.50 0.40 -0.27  0.24  0.48  0.75  1.31     1930
-#> d[Individual counselling vs. Group counselling] -0.27 0.42 -1.13 -0.54 -0.27  0.01  0.54     2773
-#> d[Self-help vs. Group counselling]              -0.61 0.50 -1.58 -0.92 -0.61 -0.29  0.37     2816
-#> d[Self-help vs. Individual counselling]         -0.34 0.41 -1.15 -0.61 -0.35 -0.08  0.47     2229
-#>                                                 Tail_ESS Rhat
-#> d[Group counselling vs. No intervention]            2307    1
-#> d[Individual counselling vs. No intervention]       1875    1
-#> d[Self-help vs. No intervention]                    2362    1
-#> d[Individual counselling vs. Group counselling]     2754    1
-#> d[Self-help vs. Group counselling]                  2787    1
-#> d[Self-help vs. Individual counselling]             2425    1
-plot(smk_releff, ref_line = 0)
-

+
(smk_releff <- relative_effects(smkfit, all_contrasts = TRUE))
+#>                                                  mean   sd  2.5%   25%   50%   75% 97.5% Bulk_ESS
+#> d[Group counselling vs. No intervention]         1.11 0.43  0.32  0.81  1.09  1.38  2.00     1750
+#> d[Individual counselling vs. No intervention]    0.84 0.23  0.41  0.68  0.83  0.99  1.33     1257
+#> d[Self-help vs. No intervention]                 0.49 0.40 -0.29  0.24  0.49  0.74  1.28     1639
+#> d[Individual counselling vs. Group counselling] -0.26 0.41 -1.07 -0.53 -0.26  0.01  0.52     2784
+#> d[Self-help vs. Group counselling]              -0.61 0.48 -1.58 -0.92 -0.60 -0.29  0.32     2982
+#> d[Self-help vs. Individual counselling]         -0.35 0.41 -1.19 -0.61 -0.34 -0.08  0.45     2219
+#>                                                 Tail_ESS Rhat
+#> d[Group counselling vs. No intervention]            2675    1
+#> d[Individual counselling vs. No intervention]       1773    1
+#> d[Self-help vs. No intervention]                    2303    1
+#> d[Individual counselling vs. Group counselling]     2698    1
+#> d[Self-help vs. Group counselling]                  2767    1
+#> d[Self-help vs. Individual counselling]             2406    1
+plot(smk_releff, ref_line = 0)
+

Treatment rankings, rank probabilities, and cumulative rank probabilities. We set lower_better = FALSE since a higher log odds of cessation is better (the outcome is positive).

-
(smk_ranks <- posterior_ranks(smkfit, lower_better = FALSE))
-#>                              mean   sd 2.5% 25% 50% 75% 97.5% Bulk_ESS Tail_ESS Rhat
-#> rank[No intervention]        3.90 0.31    3   4   4   4     4     2465       NA    1
-#> rank[Group counselling]      1.36 0.63    1   1   1   2     3     2834     2832    1
-#> rank[Individual counselling] 1.94 0.63    1   2   2   2     3     2409     2647    1
-#> rank[Self-help]              2.80 0.69    1   3   3   3     4     2581       NA    1
-plot(smk_ranks)
+
(smk_ranks <- posterior_ranks(smkfit, lower_better = FALSE))
+#>                              mean   sd 2.5% 25% 50% 75% 97.5% Bulk_ESS Tail_ESS Rhat
+#> rank[No intervention]        3.90 0.30    3   4   4   4     4     2443       NA    1
+#> rank[Group counselling]      1.35 0.60    1   1   1   2     3     3345     3426    1
+#> rank[Individual counselling] 1.93 0.63    1   2   2   2     3     2611     3205    1
+#> rank[Self-help]              2.82 0.66    1   3   3   3     4     2414       NA    1
+plot(smk_ranks)

-
(smk_rankprobs <- posterior_rank_probs(smkfit, lower_better = FALSE))
-#>                           p_rank[1] p_rank[2] p_rank[3] p_rank[4]
-#> d[No intervention]             0.00      0.00      0.10       0.9
-#> d[Group counselling]           0.71      0.22      0.07       0.0
-#> d[Individual counselling]      0.23      0.60      0.17       0.0
-#> d[Self-help]                   0.06      0.18      0.67       0.1
-plot(smk_rankprobs)
-

-
(smk_cumrankprobs <- posterior_rank_probs(smkfit, lower_better = FALSE, cumulative = TRUE))
-#>                           p_rank[1] p_rank[2] p_rank[3] p_rank[4]
-#> d[No intervention]             0.00      0.00       0.1         1
-#> d[Group counselling]           0.71      0.93       1.0         1
-#> d[Individual counselling]      0.23      0.83       1.0         1
-#> d[Self-help]                   0.06      0.24       0.9         1
-plot(smk_cumrankprobs)
-

+
(smk_rankprobs <- posterior_rank_probs(smkfit, lower_better = FALSE))
+#>                           p_rank[1] p_rank[2] p_rank[3] p_rank[4]
+#> d[No intervention]             0.00      0.00      0.10       0.9
+#> d[Group counselling]           0.72      0.22      0.06       0.0
+#> d[Individual counselling]      0.23      0.60      0.16       0.0
+#> d[Self-help]                   0.05      0.17      0.68       0.1
+plot(smk_rankprobs)
+

+
(smk_cumrankprobs <- posterior_rank_probs(smkfit, lower_better = FALSE, cumulative = TRUE))
+#>                           p_rank[1] p_rank[2] p_rank[3] p_rank[4]
+#> d[No intervention]             0.00      0.00       0.1         1
+#> d[Group counselling]           0.72      0.94       1.0         1
+#> d[Individual counselling]      0.23      0.84       1.0         1
+#> d[Self-help]                   0.05      0.22       0.9         1
+plot(smk_cumrankprobs)
+

References