 Research article
 Open Access
 Open Peer Review
 Published:
Network metaanalysis of multiple outcome measures accounting for borrowing of information across outcomes
BMC Medical Research Methodology volume 14, Article number: 92 (2014)
Abstract
Background
Network metaanalysis (NMA) enables simultaneous comparison of multiple treatments while preserving randomisation. When summarising evidence to inform an economic evaluation, it is important that the analysis accurately reflects the dependency structure within the data, as correlations between outcomes may have implication for estimating the net benefit associated with treatment. A multivariate NMA offers a framework for evaluating multiple treatments across multiple outcome measures while accounting for the correlation structure between outcomes.
Methods
The standard NMA model is extended to multiple outcome settings in two stages. In the first stage, information is borrowed across outcomes as well across studies through modelling the withinstudy and betweenstudy correlation structure. In the second stage, we make use of the additional assumption that intervention effects are exchangeable between outcomes to predict effect estimates for all outcomes, including effect estimates on outcomes where evidence is either sparse or the treatment had not been considered by any one of the studies included in the analysis. We apply the methods to binary outcome data from a systematic review evaluating the effectiveness of nine home safety interventions on uptake of three poisoning prevention practices (safe storage of medicines, safe storage of other household products, and possession of poison centre control telephone number) in households with children. Analyses are conducted in WinBUGS using Markov Chain Monte Carlo (MCMC) simulations.
Results
Univariate and the first stage multivariate models produced broadly similar point estimates of intervention effects but the uncertainty around the multivariate estimates varied depending on the prior distribution specified for the betweenstudy covariance structure. The second stage multivariate analyses produced more precise effect estimates while enabling intervention effects to be predicted for all outcomes, including intervention effects on outcomes not directly considered by the studies included in the analysis.
Conclusions
Accounting for the dependency between outcomes in a multivariate metaanalysis may or may not improve the precision of effect estimates from a network metaanalysis compared to analysing each outcome separately.
Background
Metaanalysis or the quantitative synthesis of evidence, usually from systematic reviews, has become a popular tool in healthcare evaluations [1, 2]. Largely driven by a desire for more realistic synthesis of complex healthcare evidence, increasingly sophisticated methodology has been developed. One area of metaanalysis that has seen significant methodological development is the application of multivariate statistical methods for the comparison of treatments on two or more endpoints (usually known as multivariate metaanalysis) [3–8]. These methods are appealing because many studies and systematic reviews focus on broad health effects and therefore typically report several outcome measures [4, 6, 9]. In such instances, the multivariate approach offers some advantages over separate univariate analyses including the ability to account for the interrelationship between outcomes and borrow strength across studies as well as across outcomes [10] through modelling the correlation structure [7, 11]. This can potentially reduce outcome reporting bias [12] and the uncertainty with which intervention effects are estimated. Additionally, in a decision making context where the synthesis is meant to inform a health economic evaluation, accounting for the correlations between effect estimates on different outcomes is important as the dependence between outcomes may have implication for estimating quality of life or economic consequences associated with treatment [13]. An example is the situation where a particularly effective treatment for a disease condition is associated with a large side effect profile. Ignoring information about the interrelationships between beneficial and ‘side effect’ endpoints in such instances may have implications for quantifying the benefits associated with treatment.
When summarising effectiveness evidence, correlations between the effectiveness estimates typically arise at either withinstudy and or betweenstudy levels. At the withinstudy level, correlations arise mainly due to differences in patientlevel characteristics. They are rarely reported in the published literature and usually have to be estimated from external sources such as individual patient level data if available or elicited from expert opinion [8, 11, 14, 15]. At the betweenstudy level, correlations arise from i) differences in the distribution of patientlevel characteristics across studies, in which case they will be related to the withinstudy correlations and/or ii) differences in the distribution of other studylevel characteristics such as study design, population and baseline disease severity [16]. The withinstudy correlations thus give an indication of the association between multiple endpoints within a study while the betweenstudy correlations indicate how the underlying true studyspecific effects on different outcomes vary jointly across studies.
A second area of rapid methodological development is network metaanalysis (NMA) [17], also known as mixed treatment comparison metaanalysis [18–20] or multiple treatment metaanalysis [21–23]. NMA methods extend standard pairwise metaanalysis to enable simultaneous comparison of multiple treatments while maintaining randomisation of individual studies [18]. The method enables ‘direct’ evidence (i.e. evidence from studies directly comparing two interventions of interest) and ‘indirect’ evidence (i.e. evidence from studies that do not compare the two interventions directly) to be pooled under the assumption of evidence consistency [24]. Estimates of intervention effects can then be obtained, including effects between treatments not directly compared within any one individual study [19]. NMA methods thus provide a coherent framework for appraising all available evidence relevant to a specific decision problem. The results from such analyses are increasingly being used to inform economic evaluations in healthcare decision making where coherent decisions (about judicious use of scarce resource) need to be made based on sound appraisal of all available evidence.
Approaches to extend NMA methodology to multiple outcome settings have been proposed in the literature [13, 25–27], initially focusing on mutually exclusive competing risk outcomes [13] or a single outcome measured at multiple time points [26, 28]. More recently, Efthimiou et al.[14] proposed a method for modelling multiple correlated outcomes in networks of evidence with binary outcome measures. The proposed method accounts for both the withinstudy and betweenstudy correlation structure and includes a strategy for eliciting expert opinion to inform the withinstudy correlations. This paper contributes to the growing literature on the simultaneous evaluation of correlated outcomes. We do this in two stages. In the first stage (labelled as model 2 in the remainder of the paper), information is borrowed across studies as well as across outcomes through modelling the correlations between effectiveness estimates on different outcomes. In the second stage (labelled as model 3 in the remainder of the paper), additional information is borrowed across outcomes based on ideas for combining evidence across human and animal studies originally proposed by DuMouchel and Harris [29] and also revisited by Jones et al.[30]. The proposed second stage analysis methods allows: i) disconnected treatments to be incorporated as nodes in a network of evidence and ii) prediction of intervention effects for outcomes where evidence from primary studies is either sparse or not directly available from any one study included in the analysis. The motivating application area is injury prevention in children where a broad array of outcomes and intervention packages have been evaluated with the aim of increasing safety practices around the home (to ultimately reduce household injuries).
The remainder of this paper is structured as follows: the example dataset is first described followed by a Methods section describing the statistical models developed and implementation of the models. These are followed by sections presenting the results of applying the methods to the motivating dataset and a discussion.
Dataset
The example data comes from a recently updated Cochrane systematic review of home safety education and provision of safety equipment for injury prevention in children [31]. The models developed in this paper are applied to a subset of the review evidence relating to the prevention of poisoning injuries. Table 1 presents the data from 22 studies for the following outcomes:

a)
Safe storage of medicines

b)
Safe storage of other household products (e.g. cleaning products) and

c)
Possession of a poison control centre (PCC) telephone number.
Thirteen of the 22 studies considered at least two of the three outcomes. Of these, 8 considered storage of medicines and storage of other household products, 2 considered storage of other household products and possession of a PCC telephone number, and 3 considered all three outcome measures. Individual patient data (IPD) were available for 8 of the 13 studies, of which 7 were in a format suitable for the analysis reported here as explained by the footnotes in Table 1.
We classified the interventions trialled in the 22 studies into 9 relatively homogenous treatment packages:

1)
Usual care (UC)

2)
Education (E)

3)
Education + provision of free/low cost equipment (E + FE)

4)
Education + provision of free/low cost equipment + home safety inspection (E + FE + HSI)

5)
Education + provision of free/low cost equipment + fitting of equipment (E + FE + F)

6)
Education + home safety inspection (E + HSI)

7)
Education + provision of free/low cost equipment + home safety inspection + fitting of equipment (E + FE + HSI + F)

8)
Education + home visit (E + HV)

9)
Provision of free/low cost equipment (FE).
Figure 1 shows the comparisons between the interventions that were made by individual studies and the number of comparisons in each network. All studies compared 2 intervention strategies, except Babul et al. (2007) [33] which compared 3 strategies. Data on each outcome was not available for all interventions; i.e. for the storage of medicines and other household products outcomes, interventions E + HSI and E + HV were not investigated in any of the included studies, and for possession of a PCC number interventions, E + FE + F and FE were not available.
Methods
In this section, we first present the NMA statistical model for one binary outcome measure and then extend it to compare multiple interventions across multiple outcomes. Throughout the paper, we refer to these single and multiple outcome models as univariate and multivariate NMAs respectively. Where studies report multiple outcomes, these will not be independent as each household provides information on the different outcome measures within intervention arms. The multivariate model takes this correlation structure into account by allowing the intervention effects measured by one outcome to be correlated with the intervention effects measured by other outcomes.
Model 1: Univariate NMA
Given armlevel binary data of the form presented in Table 1, a random effects NMA may be specified using the method of Lu and Ades [20]. It is assumed that the occurrence of r _{ ik } events out of a total of n _{ ik } households in the ktharm (k = A,B,C,….,) of the ithstudy follow a binomial distribution with underlying event probability p _{ ik }:
where μ _{ ib } is a studyspecific baseline effect (i.e. the logodds for the control group in study i with baseline treatment b), _{ i(bk) } is a studyspecific logodds ratio, d _{ (bk) } is the pooled effect of treatment k relative to treatment b (a quantity usually of interest in a metaanalysis) and ${\sigma}_{\left(\mathit{bk}\right)}^{2}$ is the betweenstudy variance or heterogeneity parameter. Random effects NMA assumes that intervention effects are exchangeable across the network regardless of whether or not treatments b and k are included in study i[18]. This assumption implies that the pooled effects d _{ (bk) }, can be expressed as functions of basic parameters with reference to a common comparator or baseline treatment (i.e. d _{(bk)} = d _{(Ak)} − d _{(Ab)}) [24]. Throughout this paper, we take usual care (UC) intervention to be the reference or ‘baseline’ treatment (i.e. UC is taken as treatment A of equation (2) above). Multiarm studies (i.e. studies with more than 2 treatment groups) present a special problem in network metaanalysis because they produce evidence on multiple treatment effects that are correlated through sharing a common reference or ‘baseline’ treatment. Under a homogenous variance assumption (${\sigma}_{\left(\mathit{bk}\right)}^{2}={\sigma}^{2}$), the covariance between any two effects that share a common reference treatment is $\frac{{\sigma}^{2}}{2}$[20]. The homogeneous variance assumption allows for the distribution of effects (in a study with an arbitrary number of arms) to be expressed as a univariate marginal distribution and a series of univariate conditional distributions. Specifically, for the ithstudy with p + 1 arms and p treatment effect estimates relative to the reference treatment, if
then the marginal and conditional univariate distributions for arm j, given the previous 1, ⋯, (j − 1) arms are:
The analysis is conducted within a Bayesian framework requiring prior distributions to be specified for all model parameters. Accordingly, we specified minimally informative prior distributions corresponding to a Normal (0,10^{3}) prior distribution for the pooled mean effects relative to usual care, d _{(Ak)} and the studyspecific baseline effects, μ _{ ib } and a Uniform (0,2) prior distribution for the betweenstudy standard deviation on the log odds ratio scale σ[34].
Model 2: Multivariate NMA
We extend the univariate NMA model defined above to the multiple outcomes settings in order to account for correlations between intervention effects on different outcomes. The method presented here follows from Ades et al. (2010) NMA with competing risks model [13] where only the withinstudy correlations are taken into account. We extend their method to account for the betweenstudy correlation as well.
Note that in Ades et al. (2010), a multinomial likelihood was appropriate as the three binary outcomes (relapse during treatment for Schizophrenia, discontinuation because of intolerable side effects, and discontinuation for other reasons) are mutually exclusive and event probabilities sum to 1 across outcomes. A multinomial likelihood will not be appropriate for our example dataset because each household can have one, two or all three outcomes simultaneously so that the event probabilities do not sum to 1 across outcomes. Instead, we use a Normal distribution on the logodds scale to take account of the withinstudy correlations between outcomes. We assume that in each study i and for each kth arm, the estimates ${\widehat{\theta}}_{\mathit{ikm}}$ of the observed logodds of an event on the mth outcome (m = 1, 2, ⋯, M) jointly follow a multivariate normal distribution:
where (μ _{ ib1}, μ _{ ib2}, ⋯, μ _{ ibM }) and (δ _{ i(bk)1}, δ _{ i(bk)2}, ⋯, δ _{ i(bk)M }) represent vectors of ‘true’ baseline and studyspecific effects in study i with baseline treatment b respectively. The quantities $\left({\widehat{\theta}}_{\mathit{ik}1},\phantom{\rule{0.25em}{0ex}}{\widehat{\theta}}_{\mathit{ik}2},\cdots ,\phantom{\rule{0.5em}{0ex}}{\theta}_{\mathit{ikM}}\right)$ and (θ _{ ik1}, θ _{ ik2}, ⋯, θ _{ ikM }) represent vectors of observed and ‘true’ logodds of response in arm k of study i and S _{ ik } is the associated withinstudy covariance matrix usually assumed known but estimated in practice from the data [35]. We calculated $\left({\widehat{\theta}}_{\mathit{ik}1},\phantom{\rule{0.5em}{0ex}}{\widehat{\theta}}_{\mathit{ik}2}\phantom{\rule{0.25em}{0ex}},\cdots ,\phantom{\rule{0.25em}{0ex}}{\widehat{\theta}}_{\mathit{ikM}}\right)$ and the diagonal elements of S _{ ik } using standard formulae for logodds and variance of the logodds [2]. We applied a continuity correction by adding 0.5 to the numerators and 1 to the denominators of studies with 0% or 100% event rate in one of the treatment arms [36, 37]. The offdiagonals of S _{ ik } were calculated from estimates of withinstudy correlations ${r}_{\mathit{ik}}^{\mathit{mn}}$ between outcomes m and n (m ≠ n) in arm k of study i obtained from studies with IPD (see Additional file 1: Table S1). The method used to estimate the correlations from the IPD is described in the implementation section below.
When summarising evidence across multiple endpoints, it is common to encounter instances where some studies do not report information for all outcomes of interest leading to incomplete vectors with missing studyspecific effects for the outcomes not reported [5, 10]. Such studies can be included in our model under the assumption that the effects for outcomes not reported are missing at random. When implemented using the WinBUGS software, the missing study effects and standard errors are coded as NA in the data, a strategy previously outlined in Bujkiewicz et al.[10] and Dakin et al.[28]. This enables WinBUGS to automatically ‘impute’ values for the missing information under missing at random assumption with predicted distributions.
We refer to equation (5) as the withinstudy model and the model describing the distribution of the ‘true’ effects across studies (presented below) as the betweenstudy model following standard terminology in multivariate metaanalysis [5, 6, 10, 11, 38, 39]. For the network of twoarm trials, the betweenstudy model for the ith study is thus given by:
where the ‘true’ effects δ _{ i(bk)m } (m = 1, 2, ⋯, M) jointly follow a Normal distribution with mean effects d _{(bk)m }. The parameters in equation (6) have the same interpretation as in equation (2) except that they are now specific to each outcome. The covariance matrix Σ _{(bk)} contains terms for the betweenstudy variances, ${\sigma}_{\left(\mathit{bk}\right)m}^{2}$ for each outcome m and the betweenstudy correlations ${\rho}_{\mathit{bk}}^{\mathit{mn}}$ between effects measured by outcome m and n (m ≠ n) specific to each k versus b comparison. Fitting the full model would thus require a large number of possibly multiarm studies in order to make Σ _{(bk)} identifiable [5, 13]. The number of parameters in Σ _{(bk)}, can however be reduced if reasonable assumptions can be made about the covariance structure. In particular, most practical applications of NMA methods involve the assumption of a common betweenstudy variance across treatment arms, often referred to as a homogenous variance assumption [18, 40, 41]. Therefore, to simplify Σ _{(bk)} we make the additional assumption in this context of a common betweenstudy correlation (${\rho}_{\mathit{bk}}^{\mathit{mn}}\phantom{\rule{0.25em}{0ex}}=\phantom{\rule{0.25em}{0ex}}{\rho}^{\mathit{mn}}$) leading to the following simplified betweenstudy covariance structure for twoarm studies:
where, as in the univariate case, σ _{ m } represent the common betweenstudy standard deviation or heterogeneity parameter specific to outcome m.
To include multiarm studies in our model, we extend equations (3) and (4) to the multiple outcome setting. We show in Appendix A, that under evidence consistency and the homogenous betweenstudy covariance structure (${\sigma}_{\left(\mathit{bk}\right)m}^{2}={\sigma}_{m}^{2}\phantom{\rule{0.25em}{0ex}}\mathrm{and}\phantom{\rule{0.5em}{0ex}}{\rho}_{\mathit{bk}}^{\mathit{mn}}\phantom{\rule{0.25em}{0ex}}=\phantom{\rule{0.25em}{0ex}}{\rho}^{\mathit{mn}}$), equation (3) can be extended to the multiple outcome settings by formulating the distribution of effects in a multiarm study i with p + 1 arms reporting on m = 1, 2, ⋯, M outcomes as follows:
where p is the number of treatment effect estimates. The corresponding marginal and conditional distributions for arm j, given the previous 1, ⋯, (j − 1) arms are:
To complete model 2, μ _{ ibm } and d _{(1k)m } are given minimally informative prior distributions:
Prior distributions also need to be specified for Σ _{(M × M)} which, in general, is nontrivial because of the positive definite constraint. Initially we specified an InverseWishart distribution [42]:
where K is M × M scale matrix and M is the total number of outcomes. Specifying minimally informative InverseWishart prior distribution is, however, problematic, especially when the amount of data is small relative to the dimensions of Σ _{(M × M)} as is the case for our example data. Therefore, to allow for flexibility in formulating a prior distribution for Σ _{(M × M)}, we also followed a strategy outlined by Lu and Ades (2009) [43] and more recently by Wei and Higgins (2013) [39] to express Σ _{(M × M)} in terms of a diagonal matrix of standard deviations V^{1/2} and squared positive semidefinite matrix of correlations R based on a separation strategy (Barnard et al.[44]):
where the offdiagonal elements of R contain correlation terms and diagonal elements equal 1. Lu and Ades [43] and also Wei and Higgins [39] showed that R can be written as R = L ^{T} L using Cholesky decomposition where L is an upper triangular matrix. The spherical parameterization technique [39, 43] can be used to express R in terms of sine and cosine functions of the elements in L. Using this later technique, we specified Uniform (0,π) prior distributions for the spherical coordinate φ _{ mn } in our model to ensure that elements of the correlation matrix R lie in the interval (−1,1). Finally, the elements of V ^{1/2} correspond to the betweenstudy standard deviation terms in Σ _{(M × M)} and are given independent Uniform (0,2) prior distributions as in the univariate case (model 1).
Model 3: Borrowing strength across interventions and outcomes
From Table 1, it can be seen that none of the studies had considered the interventions E + HSI and E + HV for storage of medicines and storage of other household products. Similarly, interventions E + FE + F and FE were not trialled by any of the included studies on possession of a PCC number. To estimate the full set of 24 basic intervention effects relative to usual care from 9 interventions on 3 outcomes, we applied methods originally proposed by DuMouchel and Harris [29] and revisited by DuMouchel and Groer [45] and Jones et al.[30]. We assume that the pooled effects of treatment k relative to usual care intervention d _{(Ak)m }, can be expressed as a sum of treatmentspecific effect α _{ k } and outcomespecific effect γ _{ m }. This assumption replaces the minimally informative prior distribution N(0, 10^{3}) specified for d _{(Ak)m } in model 2 with:
where K is the total number of treatments being evaluated across M outcomes, and for k = 1 (i.e. reference treatment A), d _{(Ak)m } equal to zero. Note that on the logarithmic scale, this would imply that the ratio of any intervention effects is constant across outcomes as the γ _{ m } cancel, i.e.
Equation (10) thus embodies an assumption of equal or constant relative potency of treatments across outcomes which imply exchangeability of the relative effects between the nonreference/baseline treatments indicated by equation (11). For our example dataset, this implies that missing intervention effects for comparisons with the usual care intervention can be predicted directly from equation (10) as a linear combination of γ _{ m } and α _{ k } assuming that each treatment effect relative to usual care is reported on at least one outcome. The missing intervention effects between non reference/baseline treatments if required can similarly be predicted directly from the model as linear combinations of the intervention effects relative to usual care.
The parameter τ controls the accuracy of the constant relative potency assumption. Values of τ close to zero would thus indicate a high degree of confidence (and support from the data) in the parallelism of effect profiles across outcomes and the constant relative potency assumption. Conversely, larger values of τ would indicate otherwise.
Multiarm studies are included in model 3 based on equations (8) and (9) in the same way as in model 2. To complete model 3, the parameters α _{ k }, γ _{ m } and τ are given minimally informative prior distributions. For the mean effects, this is a normal distribution with zero mean and large variance:
We give τ a Uniform (0, 2) prior distribution, considered to be minimally informative on the logodds ratio scale. Sensitivity analyses were conducted to assess the impact of specifying alternative prior distributions for τ that are also considered minimally informative [46]:

i.
Normal prior distribution centred on 0 with large variance and constrained to be positive, τ ~ N(0, 10^{2}), τ ≥ 0

ii.
Gamma prior distribution placed on the precision: τ ^{− 2} ~ Gamma(0.001, 0.001).
The results of the sensitivity analyses are presented in Additional file 1: Figure S1.
There is a limitation to the number of data (i.e. intervention effects relative to the usual care) on outcomes allowed to be missing for the model hyperparameters to be identifiable. For K interventions and M outcomes, there will be (K − 1) × M equation (10) that are used to estimate a total of (K − 1) + M hyperparameters (i.e. (K − 1) of α _{ k } and M of γ _{ m } hyperparameters). Therefore no more than ((K − 1) × M) − ((K − 1) + M) missing values in total are allowed. For example, for K = 3 treatments and M = 2 outcomes, data has to be available on both outcomes for both treatment comparisons with the baseline when the prior distributions are noninformative. When large number of data on outcomes is missing, placing informative prior distributions on the hyperparameters can improve convergence.
Implementation of the models
We fitted a total of four models, models 1 and 3 as described above and two versions of model 2. In model 2a, we specified an inverseWishart prior distribution for the betweenstudy covariance matrix Σ _{(M × M)} whilst in model 2b, we specified a prior distribution for Σ _{(M × M)} based on the separation strategy. All four models allowed for multiarm trials to be included in the analysis. To fit the multivariate NMA models, the quantities $\left({\hat{\theta}}_{\mathit{ik}1},\phantom{\rule{0.5em}{0ex}}{\hat{\theta}}_{\mathit{ik}2},\phantom{\rule{0.25em}{0ex}}{\hat{\theta}}_{\mathit{ik}3}\right)$ and the diagonals of S _{ ik } were estimated using standard 2×2–table formulae [2]. Next, we obtained estimates of the withinstudy correlations from the IPD studies using the following three methods: i) Pearson correlation coefficient between the observed outcome events ii) Bootstrapping as described in Daniel and Hughes (1998), and iii) Generalised Estimating Equations (See details in Additional file 1). All three methods produced identical estimates of the correlations between pairs of outcome specific logodds of event from each IPD study (Additional file 1: Table S1). Therefore, we formulated informative prior distributions for the correlation terms in S _{ ik } of equation (5) using the estimates of the correlations between the observed outcome events (Pearson correlation) as follows:
where ${r}_{\mathit{ik}}^{\mathit{mn}}$ is the withinstudy correlation between the outcomes m and n effects measured on the logodds scale in arm k of study i, and $\overline{{r}^{\mathit{mn}}}$ and var (r ^{mn})are the mean and variance of the withinstudy correlation between outcomes m and n effects measured by IPD respectively.
We also assessed consistency of the evidence within each network using a method based on node splitting [24]. We found no evidence of conflict between the direct and indirect sources (Additional file 1: Table S2) in all three outcome networks.
We fitted all models described above in WinBUGS [47] using Markov Chain Monte Carlo (MCMC) simulations. The univariate models were fitted separately for each outcome using WinBUGS code available from Dias et al.[48]. The WinBUGS code for the multivariate models is provided in the appendices 1 and 2 of the Additional file 1. Convergence was assessed by examination of the trace and autocorrelation plots and the RubinGelman statistic after running 400 000 simulations and discarding the first 200 000 samples as ‘burn in samples’.
Results
Univariate and multivariate analyses
Parameters of interest were the posterior median estimate (and 95% credible intervals) of the pooled intervention effects relative to the usual care intervention and the posterior median estimate (and 95% credible intervals) of the betweenstudy standard deviation and correlation terms. Summary forest plots displaying effectiveness estimates relative to usual care on the odds ratio (OR) scale are presented in Figure 2. It can be seen that, all 4 models produced broadly similar estimates when the treatment effect is not extreme compared to the other effect estimates for the same outcome. Compared to the univariate analysis, the multivariate models produced noticeably less extreme estimates of intervention effects. This can be seen in the effect of FE + HSI (3) on possession of PCC number being shifted towards the line of no effect from an OR of 39.35 (95% CrI 2.37 to 732.30) in model 1 to 23.55 (95% CrI 1.39, to 456.80) in model 2a, 20.37 (95% CrI 0.72, to 706.00) in model 2b and 4.20 (95% CrI 1.59 to 13.16) in model 3. Similarly, the OR for FE (9) on safe storage of other household products shifted from 0.37 (95% CrI 0.00 to 15.10) in model 1 to 1.81 (95% CrI 0.63, to 5.37) in model 3.
Posterior median and 95% credible intervals of the betweenstudy standard deviations and correlations are presented in Table 2. The posterior medians of the betweenstudy correlations from the multivariate models were small and estimated with considerable uncertainty (i.e. all had large variances). Estimates of the betweenstudy standard deviations were broadly similar for the univariate NMA (model 1) and the multivariate NMA using the separation strategy (model 2b), and relatively high for multivariate NMA using the inverseWishart prior distribution (model 2a).
Borrowing strength across outcomes
It can be seen from Figure 2 that the effect of E + HSI and E + HV relative to usual care intervention on safe storage of medicines and safe storage of other household products, and E + FE + F and FE on possession of a PCC telephone number were only estimated in model 3 as none of the studies had trialled these interventions on the respective outcomes. In this model, estimates of relative effects between non reference/baseline treatments were assumed to be exchangeable across outcomes, which enabled estimates to be obtained for all outcomes by predicting effects where the interventions have not been considered for the particular outcome of interest. For the intervention/outcome pair where data from trials were available, this extrapolation step had the additional effect of producing more precise estimates of the treatment effect in comparison to the models that do not assume exchangeability effects across outcomes.
The posterior median and 95% credible intervals of intervention effects relative to usual care were unaffected by placing alternative minimally informative prior distributions on τ (Additional file 1: Appendix 2). The posterior median and credible intervals for τ (Table 2) were similarly not sensitive to the choice of prior distribution placed on τ in the primary and sensitivity analyses. The posterior median estimates of τ were all close to zero, which suggest that assumptions about the parallelism of effect profiles across outcomes is supported by the data. The estimates of τ would thus suggest with 95% probability, that based on the information in our example dataset, the extrapolation model could be accurate to within a factor of about e ^{(2 × 0.10)} = 1.24 (95% CrI 1.01 to 2.87).
Discussion
We have presented methods for simultaneous comparison of multiple treatments across multiple outcome measures while preserving the internal randomisation of individual studies. Our method may be viewed as an extension of Ades et al.’s (2010) NMA with competing risks paper [13] wherein only the withinstudy correlation is taken into account. We have extended their method to account for the dependency between outcome effects across studies as well as withinstudies.
In this particular application of the multivariate approach, accounting for the correlation between outcomes alone (models 2a and 2b) did not reduce the uncertainty around estimates of intervention effects compared to analysing each outcome separately (model 1). Assuming that intervention effects are exchangeable across outcome did however lead to a modest reduction in uncertainty around effectiveness estimates (model 3).
The betweenstudy correlations were estimated with considerable uncertainty (Table 2) and appear to have little impact on overall effect estimates. This may be because the betweenstudy correlation arises due to, among other things, differences in studylevel characteristics that also give rise to betweenstudy heterogeneity in a metaanalysis. Based on a criterion outlined in Spiegelhalter et al.[49] the posterior median estimates of the betweenstudy standard deviations, σ _{1} and σ _{2} on the log odds ratio scale (Table 2) could be interpreted as indicating evidence of low to moderate heterogeneity for storage of medicines and storage of other household products outcomes. Only the estimates for possession of poison control centre number exhibited a considerable degree of heterogeneity. Consequently, the posterior medians of the betweenstudy correlations were small. There was therefore very little gain (in terms of increasing the precision of estimates) from formulating the betweenstudy covariance structure described for the analysis presented here. Accounting for the betweenstudy correlation is likely to be beneficial in situations where the betweenstudy variance (heterogeneity) is large relative to withinstudy variances.
We opted to incorporate the withinstudy correlation through the armspecific effects (logodds) rather than the studyspecific treatment difference (logodds ratio) as is often done in multivariate metaanalysis [3, 4, 15, 38]. This approach greatly simplifies the likelihood for multiarm studies because treatment arms can be considered independent as a consequence of randomisation. Hence, there is no requirement to account for the additional correlations between effect estimates which share a common comparator treatment in the model likelihood [50]. The armbased approach is also likely to be useful when (as is typical with many practical application of multivariate metaanalysis) the withinstudy correlations are not available [10, 12, 15, 51] and have to be obtained from an external source such as expert opinion [14]. In such situations, formulating questions about correlations between outcomespecific event probabilities (which can be used directly in an armbased approach) is more likely to be intuitive and easily understood by nonstatistician healthcare experts than questions about correlations between intervention effects. It is acknowledged however, that the correlations between the intervention effects if required can easily be obtained from the correlations between the outcomes [14, 51].
At the betweenstudy level, we assumed a common correlation structure across treatments in addition to the common variance assumption underlying most practical application of NMA methods. The common correlation assumption implies that if several separate multivariate metaanalyses were conducted with the same outcomes, each with a different set of k versus b comparison, the assumption is that the betweenstudy correlations would be the same across the different sets of bk comparisons. We suggested this structure to simplify the covariance structure and reduce the number of parameters in the model. Appropriateness of such modelling assumptions would need to be considered carefully and assessed when it is feasible to do so.
Initially we specified an inverseWishart distribution for the betweenstudy covariance matrix Σ _{(M × M)}. However, we believe this prior distribution to be influential due to the small number of studies in our example dataset relative to the number of outcomes. Under these conditions, the inverseWishart prior distribution produced upwardlybiased estimates of σ _{1} and σ _{2} and downward bias in the estimate for σ _{3} when compared to the corresponding estimates obtained from the univariate model (Table 2). These findings are consistent with observations in the univariate case where the use of a Gamma prior distribution (which is the univariate analogue of the InverseWishart prior distribution) can lead to an overestimation of the heterogeneity parameter when the true value is close to 0 [46, 52]. As an alternative to an inverseWishart prior distribution therefore, we followed the spherical decomposition technique suggested by Lu and Ades [43]. This parameterization offered greater flexibility in formulating independent prior distributions for the standard deviation and correlation terms in Σ _{(M × M)}.
An obvious limitation to implementation of the multivariate models presented in this paper is the limited availability of data including i) the problem of missing withinstudy correlations and ii) the requirement for a relatively large number of studies to estimate all model parameters. The problem of missing withinstudy correlations has traditionally hampered the widespread application of multivariate metaanalysis [7, 10, 15]. In our example, IPD was available from a proportion of the included studies and we have used correlations estimated from the IPD to formulate informative prior distributions for the withinstudy model. Alternative approaches to dealing with missing withinstudy correlations when IPD is not available include: i) using the observed correlation from the summary studyspecific effects [12], ii) eliciting information about the correlations from external sources such as clinical experts [14] and iii) specifying ‘vague’ prior distributions for analysis conducted within a Bayesian framework [6].
The second data issue concerns the number of studies needed to estimate the full unstructured betweenstudy covariance matrix presented in equation (6). We anticipate a large number of multiarm studies reporting across the three outcomes will be needed to identify Σ _{(bk)} and estimate all model parameters. This can be problematic considering the fact that most applications of network metaanalysis typically include mostly twoarm studies with very small numbers of multiarm studies. Even with the simplification of the betweenstudy covariance matrix given in equation (5), a relatively large number of studies in comparison to the total number of outcomes being considered may still be needed. We are unable to answer the question of how many studies should be considered large enough for a NMA with multiple outcomes. As a guide, Wei and Higgins [39] recently estimated 15, 27 and 42 studies as a minimum for multivariate pairwise metaanalysis with two, three and fouroutcomes respectively. Hence, we believe an even larger number of studies will be required for the NMA with multiple outcomes.
Another limitation of the multivariate models presented here is that they rely on the normal approximation to binomial distribution to incorporate the withinstudy correlations in the model. The normal approximation frequently fails and may not provide adequate fit to the data in the presence of studies with zero or a small number of events, necessitating use of continuity corrections. We were unable to use the exact binomial distribution as our primary interest was to develop models for summary binary data where outcomes are not mutually exclusive, and where it is not reasonable to assume that withinstudy correlations are zero so that the likelihood factorises easily as in Arends et al.[3]. Further methodological investigations into modelling multivariate summary data that is not normally distributed will therefore be useful. An example is provided in Chu et al.[53] where parameterization of the withinstudy model enabled the special case of diagnostic sensitivity and specificity to be jointly modelled with disease prevalence using a trivariate binomial likelihood. In the interim, an alternative formulation which bypasses the need for approximating normal distributions is to directly model the IPD where this is available. This will require extending Saramago et al.’s[54] NMA model with aggregate and individual participant level data from single outcome to multiple outcome settings.
We assessed the consistency of each outcome network separately using the method of node splitting [24]. We found no evidence of conflict between the direct and indirect sources on pairwise contrasts that have both sources of evidence in model 1. We did not assess the consistency of the multivariate estimates partly because we are unaware of current methods for carrying out this type of assessment. We are investigating extensions of the nodesplit method to multiple outcome networks and investigate the effect of jointly synthesising evidence across multiple endpoints on evidence consistency in a simulation study.
Our initial motivation for a multiple outcome NMA was to estimate intervention effects for all the outcomes, including effects of interventions on outcomes not considered by any of the studies included in the analysis. This requires the correlation structure between effects on multiple outcomes to be appropriately modelled and also ensuring the mechanism of “borrow strength” across outcomes through the assumption of exchangeability of the random effect across outcomes. This implies a priori assumption that outcomes are related but different and that there is no way of knowing the order of magnitude of effects on outcomes. If this assumption does not hold, it may potentially lead to worse or more biased effectiveness estimates. In our example, the outcomes are similar and measured on the same scale. It would be clearly inappropriate to assume that intervention effects are exchangeable across outcomes that are different in some important respects such as being measured on different scales (e.g. where one outcome reports a weighted mean difference and another outcome reports a logodds ratio) as such estimates will differ in terms of the precision with which they are estimated.
Conclusion
Our aim in this paper was to present methods for simultaneous comparison of multiple treatments across multiple outcome measures while preserving the internal randomisation of individual studies. Application of the method to the poison prevention data yielded similar point estimates of treatment effect to those obtained from a univariate NMA but the uncertainty around the multivariate estimates increased or decreased depending on the prior distribution specified for the betweenstudy covariance structure. The proposed method followed the usual hierarchical approach to multivariate metaanalysis where correlations between outcomes are modelled at the withinstudy and or betweenstudy levels.
Betweenstudy covariance for multiarm studies reporting multiple outcomes
For a multiarm study i with K treatments labelled A, B, C, …, K reporting a total of M outcomes labelled 1, 2, …, M. A random effects betweenstudy model can be represented as:
Where δ _{ i(Ak)m } and d _{(Ak)m } are studyspecific and mean effect of treatment k relative A (reference treatment) on outcome m in study i respectively and Σ _{ FULL } is the full (K1) × (K1) blocks of M × M withintreatment betweenoutcome covariance matrix. The parameters in Σ _{ FULL } have the following interpretation:(${\sigma}_{\left(\mathit{Ak}\right)m}^{2}$) indicate the variance of the effect of treatment k (k = B, C, ⋯, K) relative to A on outcome m across studies.
${\rho}_{\left(\mathit{Ak},\mathit{Ak}\right)}^{\mathit{mn}}$indicate the correlation between δ _{ i(Ak)m } and δ _{ i(Ak)n } (i.e. the correlation between the effect of treatment k relative to A on outcome m and the effect of treatment k relative to A on outcome n (m ≠ n)) specific to the Ak comparison.
${\rho}_{\left(\mathit{Ah},\mathit{Ak}\right)}^{\mathit{mm}}$indicate the correlation between δ _{ i(Ah)m } and δ _{ i(Ak)m } (i.e. the correlation between the effect of treatment h relative to A on outcome m and the effect of treatment k relative to A (h ≠ k) on outcome m because they share a common comparator A).The diagonal block matrices in Σ _{ FULL } thus carry terms for the betweenstudy variance (${\sigma}_{\left(\mathit{Ak}\right)m}^{2}$) while the offdiagonal blocks carry terms for the betweenstudy covariance. We make two assumptions to simplify and reduce the number of parameters in Σ _{ FULL }. First, we assume homogenous variances for intervention effects within outcomes [20]. This implies ${\sigma}_{\left(\mathit{Ak}\right)m}^{2}={\sigma}_{m}^{2}$ and ${\rho}_{\left(\mathit{Ah},\mathit{Ak}\right)}^{\mathit{mm}}=\phantom{\rule{0.25em}{0ex}}\frac{1}{2}$ as in the single outcome network metaanalysis case [20, 34]. Second, we make the assumption of homogenous betweenstudy correlations for the intervention effects from different outcomes. Under this assumption we can express ${\rho}_{\left(\mathit{Ah},\mathit{Ah}\right)}^{\mathit{mn}}$ and ${\rho}_{\left(\mathit{Ak},\mathit{Ak}\right)}^{\mathit{mn}}$ in terms of a common correlation parameter ρ ^{mn} by noting that for any 3treatment ( A, h, k ) configuration, the covariance between outcome m and n effects across studies can be expressed as a covariance between two sums under evidence consistency:
The homogenous betweenstudy correlation assumption implies ${\rho}_{\left(\mathit{Ah},\mathit{Ah}\right)}^{\mathit{mn}}=\phantom{\rule{0.25em}{0ex}}{\rho}_{\left(\mathit{Ak},\mathit{Ak}\right)}^{\mathit{mn}}={\rho}^{\mathit{mn}}$ and ${\rho}_{\left(\mathit{Ak},\mathit{Ah}\right)}^{\mathit{mn}}=\phantom{\rule{0.25em}{0ex}}\frac{1}{2}{\rho}^{\mathit{mn}}$ for the inequality $1\le \left({\rho}_{\left(\mathit{Ak},\mathit{Ak}\right)}^{\mathit{mn}}+{\rho}_{\left(\mathit{Ah},\mathit{Ah}\right)}^{\mathit{mn}}2{\rho}_{\left(\mathit{Ak},\mathit{Ah}\right)}^{\mathit{mn}}\right)\phantom{\rule{0.25em}{0ex}}\le 1$ to hold. Substituting these expressions into equation (A1), we see that the betweenstudy correlation terms equal ρ ^{mn} in the diagonal block of matrices and $\frac{1}{2}{\rho}^{\mathit{mn}}$ in the offdiagonal block of matrices of in Σ _{ FULL } leading to the following simplificationfollowing simplification of the betweenstudy covariance matrix:
Finally by relabeling the reference treatment A as b, (δ _{ i(AB)1}, ⋯, δ _{ i(AK)m }) as $\left({\delta}_{i\left(b{k}_{1}\right)1},\phantom{\rule{0.5em}{0ex}}\cdots ,{\delta}_{i\left(b{k}_{j}\right)M}\phantom{\rule{0.25em}{0ex}}\right)$ and (d _{(AB)1}, ⋯, d _{(AK)M }) as $\left({d}_{\left(b{k}_{1}\right)M},\phantom{\rule{0.5em}{0ex}}\cdots ,{d}_{\left(b{k}_{j}\right)M}\phantom{\rule{0.25em}{0ex}}\right)$, equation (A3) can be rewritten as equation (A4).
Abbreviations
 CrI:

Credible interval
 IPD:

Individual patient data
 MCMC:

Markov Chain Monte Carlo
 NMA:

Network metaanalysis
 PCC:

Poison centre control number.
References
 1.
Cooper NJ, Sutton AJ, Ades AE, Paisley S, Jones DR, W. G. O. U. o. E. i. E. D: Models, Use of evidence in economic decision models: Practical issues and methodological challenges. Health Economics. 2007, 16 (12): 12771286.
 2.
Sutton A, Abrams K, Jones DR, Sheldon TA, Song F: Methods for Metaanalysis in Medical Research. 2000, Chichester: Wiley
 3.
Arends LR, Voko Z, Stijnen T: Combining multiple outcome measures in a metaanalysis: an application. Statistics in Medicine. 2003, 22 (8): 13351353.
 4.
Berkey CS, Hoaglin DC, AntczakBouckoms A, Mosteller F, Colditz GA: Metaanalysis of multiple outcomes by regression with random effects. Statistics in Medicine. 1998, 17 (22): 25372550.
 5.
Jackson D, Riley R, White IR: Multivariate metaanalysis: Potential and promise. Statistics in Medicine. 2011, 30 (20): 24812498.
 6.
Nam IS, Mengersen K, Garthwaite P: Multivariate metaanalysis. Statistics in Medicine. 2003, 22 (14): 23092333.
 7.
Riley RD, Abrams KR, Sutton AJ, Lambert PC, Thompson JR: Bivariate randomeffects metaanalysis and the estimation of betweenstudy correlation. BMC Medical Research Methodology. 2007, 7:
 8.
Riley RD, Thompson JR, Abrams KR: An alternative model for bivariate randomeffects metaanalysis when the withinstudy correlations are unknown. Biostatistics. 2008, 9 (1): 172186.
 9.
Kendrick D, Coupland C, Mulvaney C, Simpson J, Smith SJ, Sutton A, Watson M, Woods A: Home safety education and provision of safety equipment for injury prevention. Cochrane Database Systematic reviews. 2007, 1: CD005014
 10.
Bujkiewicz S, Thompson JR, Sutton AJ, Cooper NJ, Harrison MJ, Symmons DPM, Abrams KR: Multivariate metaanalysis of mixed outcomes: a Bayesian approach. Statistics in Medicine. 2013, 32 (22): 39263943.
 11.
Riley RD, Abrams KR, Lambert PC, Sutton AJ, Thompson JR: An evaluation of bivariate randomeffects metaanalysis for the joint synthesis of two correlated outcomes. Statistics in Medicine. 2007, 26 (1): 7897.
 12.
Kirkham JJ, Riley RD, Williamson PR: A multivariate metaanalysis approach for reducing the impact of outcome reporting bias in systematic reviews. Statistics in Medicine. 2012, 31 (20): 21792195.
 13.
Ades AE, Mavranezouli I, Dias S, Welton NJ, Whittington C, Kendall T: Network metaanalysis with competing risk outcomes. Value Health. 2010, 13 (8): 976983.
 14.
Efthimiou O, Mavridis D, Cipriani A, Leucht S, Bagos P, Salanti G: An approach for modelling multiple correlated outcomes in a network of interventions using odds ratios. Statistics in Medicine. 2014, 33 (13): 22752287.
 15.
Riley RD: Multivariate metaanalysis: the effect of ignoring withinstudy correlation. Journal of the Royal Statistical Society Stat Soc. 2009, 172: 789811.
 16.
Rodgers M, Sowden A, Petticrew M, Arai L, Roberts H, Britten N, Popay J: Testing Methodological Guidance on the Conduct of Narrative Synthesis in Systematic Reviews. Evaluation. 2009, 15 (1): 4973.
 17.
Lumley T: Network metaanalysis for indirect treatment comparisons. Statistics in Medicine. 2002, 21 (16): 23132324.
 18.
Caldwell DM, Ades AE, Higgins JP: Simultaneous comparison of multiple treatments: combining direct and indirect evidence. British Medical Journal. 2005, 331 (7521): 897900.
 19.
Caldwell DM, Welton NJ, Ades AE: Mixed treatment comparison analysis provides internally coherent treatment effect estimates based on overviews of reviews and can reveal inconsistency. Journal of Clinical Epidemiology. 2010, 63 (8): 875882.
 20.
Lu G, Ades AE: Combination of direct and indirect evidence in mixed treatment comparisons. Statistics in Medicine. 2004, 23 (20): 31053124.
 21.
Salanti G, Ades AE, Ioannidis JP: Graphical methods and numerical summaries for presenting results from multipletreatment metaanalysis: an overview and tutorial. Journal of Clinical Epidemiology. 2010, 64 (2): 163171.
 22.
Salanti G, Higgins JP, Ades AE, Ioannidis JP: Evaluation of networks of randomized trials. Statistical Methods in Medical Research. 2008, 17 (3): 279301.
 23.
Salanti G, Marinho V, Higgins JP: A case study of multipletreatments metaanalysis demonstrates that covariates should be considered. Journal of Clinical Epidemiology. 2009, 62 (8): 857864.
 24.
Dias S, Welton NJ, Caldwell DM, Ades AE: Checking consistency in mixed treatment comparison metaanalysis. Statistics in Medicine. 2010, 29 (7–8): 932944.
 25.
Hong H, Carlin BP, Chu H, Shamliyan TA, Wang S, Kane RL: Methods Research Report: A Bayesian Missing Data Framework for Multiple Continuous Outcome Mixed Treatment Comparisons. 2013, Avaialable online from http://www.ncbi.nlm.nih.gov/books/NBK116689/pdf/TOC.pdf
 26.
Lu G, Ades AE, Sutton AJ, Cooper NJ, Briggs AH, Caldwell DM: Metaanalysis of mixed treatment comparisons at multiple followup times. Statistics in Medicine. 2007, 26 (20): 36813699.
 27.
Welton NJ, Cooper NJ, Ades AE, Lu G, Sutton AJ: Mixed treatment comparison with multiple outcomes reported inconsistently across trials: evaluation of antivirals for treatment of influenza A and B. Statistics in Medicine. 2008, 27 (27): 56205639.
 28.
Dakin HA, Welton NJ, Ades AE, Collins S, Orme M, Kelly S: Mixed treatment comparison of repeated measurements of a continuous endpoint: an example using topical treatments for primary openangle glaucoma and ocular hypertension. Statistics in Medicine. 2011, 30 (20): 25112535.
 29.
DuMouchel WH, Harris JE: Bayes Methods for Combining the Results of Cancer Studies in Humans and Other Species. Journal of the American Statistical Association. 1983, 78 (382): 293308.
 30.
Jones DR, Peters JL, Rushton L, Sutton AJ, Abrams KR: Interspecies extrapolation in environmental exposure standard setting: A Bayesian synthesis approach. Regulatory Toxicology Pharmacology. 2009, 53 (3): 217225.
 31.
Kendrick D, Young B, MasonJones AJ, Ilyas N, Achana FA, Cooper NJ, Hubbard SJ, Sutton AJ, Smith S, Wynn P, Mulvaney C, Watson MC, Coupland C: Home safety education and provision of safety equipment for injury prevention. Cochrane Database Systematic reviews. 2012, 9
 32.
Johnston BD, Huebner CE, Anderson ML, Tyll LT, Thompson RS: Healthy steps in an integrated delivery system  Child and parent outcomes at 30 months. Archives of Pediatrics and Adolescent Medicine. 2006, 160 (8): 793800.
 33.
Babul S, Olsen L, Janssen P, McIntee P, Raina P: A randomized trial to assess the effectiveness of an infant home safety programme. International Journal of Control and Safety Promotion. 2007, 14 (2): 109117.
 34.
Dias S, Sutton AJ, Ades AE, Welton NJ: A Generalized Linear Modeling Framework for Pairwise and Network Metaanalysis of Randomized Controlled Trials. Medical Decision Making. 2012, 33 (5): 60717.
 35.
van Houwelingen HC, Arends LR, Stijnen T: Advanced methods in metaanalysis: multivariate approach and metaregression. Statistics in Medicine. 2002, 21 (4): 589624.
 36.
Dershewitz RA, Williamson JW: Prevention of Childhood Household Injuries  Controlled ClinicalTrial. American Journal of Public Health. 1977, 67 (12): 11481153.
 37.
Kelly B, Sein C, McCarthy PL: Safety education in a pediatric primary care setting. Pediatrics. 1987, 79 (5): 818824.
 38.
Mavridis D, Salanti G: A practical introduction to multivariate metaanalysis. Statistical Methods in Medical Research. 2013, 22 (2): 133158.
 39.
Wei Y, Higgins JPT: Bayesian multivariate metaanalysis with multiple outcomes. Statistics in Medicine. 2013, 32 (17): 29112934.
 40.
Carter P, Achana F, Troughton J, Gray LJ, Khunti K, Davies MJ: A Mediterranean diet improves HbA1c but not fasting blood glucose compared to alternative dietary strategies: a network metaanalysis. Journal of Human Nutrition and Dietetics. 2014, 27 (3): 280297.
 41.
Cooper NJ, Sutton AJ, Lu G, Khunti K: Mixed comparison of stroke prevention treatments in individuals with nonrheumatic atrial fibrillation. Archives Internal Medicine. 2006, 166 (12): 12691275.
 42.
Spiegelhalter D, Thomas A, Best N, Lunn D: WinBUGS User Manual. Available from http://www.politicalbubbles.org/bayes_beach/manual14.pdf, 2003. Accessed 10^{th} December 2013
 43.
Lu GB, Ades A: Modeling betweentrial variance structure in mixed treatment comparisons. Biostatistics. 2009, 10 (4): 792805.
 44.
Barnard J, McCulloch R, Meng XL: Modeling covariance matrices in terms of standard deviations and correlations, with application to shrinkage. Statistica Sinica. 2000, 10 (4): 12811311.
 45.
DuMouchel W, Groer PG: A Bayesian Methodology for Scaling Radiation Studies from Animals to Man. Health Physics. 1989, 57: 411418.
 46.
Lambert PC, Sutton AJ, Burton PR, Abrams KR, Jones DR: How vague is vague? A simulation study of the impact of the use of vague prior distributions in MCMC using WinBUGS. Statistics in Medicine. 2005, 24 (15): 24012428.
 47.
Lunn DJ, Thomas A, Best N, Spiegelhalter D: WinBUGS  A Bayesian modelling framework: Concepts, structure, and extensibility. Statistical Computing. 2000, 10 (4): 325337.
 48.
Dias S, Welton NJ, Sutton AJ, Ades AE: NICE DSU Technical Support Document 2: A Generalised Linear Modelling Framework for Pairwise and Network Metaanalysis of Randomised Controlled Trials. 2011, Available from http://www.nicedsu.org.uk
 49.
Spiegelhalter DJ, Abrams KR, Myles JP: Bayesian Approaches to Clinical Trials and HealthCare Evaluation. Statistics in Practice. 2004, Chichester: John Wiley & Sons
 50.
Franchini AJ, Dias S, Ades AE, Jansen JP, Welton NJ: Accounting for correlation in network metaanalysis with multiarm trials. Research Synthesis Methodology. 2012, 3 (2): 142160.
 51.
Wei Y, Higgins JPT: Estimating withinstudy covariances in multivariate metaanalysis with multiple outcomes. Statistics in Medicine. 2013, 32 (7): 11911205.
 52.
Gelman A: Prior distributions for variance parameters in hierarchical models(Comment on an Article by Browne and Draper). Bayesian Analysis. 2006, 1 (3): 515533.
 53.
Chu HT, Nie L, Cole SR, Poole C: Metaanalysis of diagnostic accuracy studies accounting for disease prevalence: Alternative parameterizations and model selection. Statistics in Medicine. 2009, 28 (18): 23842399.
 54.
Saramago P, Sutton AJ, Cooper NJ, Manca A: Mixed treatment comparisons using aggregate and individual participant level data. Statistics in Medicine. 2012, 31 (28): 35163536.
Prepublication history
The prepublication history for this paper can be accessed here:http://www.biomedcentral.com/14712288/14/92/prepub
Acknowledgements
FAA is funded by the National Institute for Health Research (NIHR) under its Programme Grants for Applied Research programme (Reference Number RPPG 0407–10231) to pursue a PhD in Public Health Evaluation of Interventions.
SB is supported by the MRC Methodology Research Programme (grant number MR/L009854/1).
Author information
Additional information
Competing interests
The authors declare that they have no competing interests.
Author’s contribution
FAA drafted the initial manuscript, conducted the analysis and coordinated contributions from all authors in drafting the final manuscript. NJC designed the study and revised the manuscript to improve general readability. SB contributed to developing the statistical analysis and interpretation of the results. SJH contributed to preparing the data for analysis, drafting the network diagrams and revision of the initial and final drafts. DK provided the data for analysis, defined the outcomes and interventions and revision of the final drafts.DRJ participated in designing the study and revision of the final manuscript. AJS conceived the original idea and participated developing the statistical model. All authors read and approved the final manuscript.
Authors’ original submitted files for images
Below are the links to the authors’ original submitted files for images.
Rights and permissions
This article is published under license to BioMed Central Ltd. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly credited. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.
About this article
Received
Accepted
Published
DOI
Keywords
 Network metaanalysis
 Mixed treatment comparisons
 Multiple outcomes
 Multivariate
 WinBUGS