Package 'MBNMAtime'

Title: Run Time-Course Model-Based Network Meta-Analysis (MBNMA) Models
Description: Fits Bayesian time-course models for model-based network meta-analysis (MBNMA) that allows inclusion of multiple time-points from studies. Repeated measures over time are accounted for within studies by applying different time-course functions, following the method of Pedder et al. (2019) <doi:10.1002/jrsm.1351>. The method allows synthesis of studies with multiple follow-up measurements that can account for time-course for a single or multiple treatment comparisons. Several general time-course functions are provided; others may be added by the user. Various characteristics can be flexibly added to the models, such as correlation between time points and shared class effects. The consistency of direct and indirect evidence in the network can be assessed using unrelated mean effects models and/or by node-splitting.
Authors: Hugo Pedder [aut, cre]
Maintainer: Hugo Pedder <[email protected]>
License: GPL-3
Version: 0.2.4
Built: 2024-11-20 09:27:05 UTC
Source: https://github.com/cranhaven/cranhaven.r-universe.dev

Help Index


Add follow-up time and arm indices to a dataset

Description

Adds follow-up time (fups, fupcount) and arm (arms, narms) indices to a dataset.

Usage

add_index(data.ab, reference = 1)

Arguments

data.ab

A data frame of arm-level data in "long" format containing the columns:

  • studyID Study identifiers

  • time Numeric data indicating follow-up times

  • treatment Treatment identifiers (can be numeric, factor or character)

  • class An optional column indicating a particular class code. Treatments with the same identifier must also have the same class code.

reference

A number or character (depending on the format of treatment within data.ab) indicating the reference treatment in the network (i.e. those for which estimated relative treatment effects estimated by the model will be compared to).

Value

A data frame similar to data.ab but with additional columns:

  • arm Arm identifiers coded for each study

  • fupcount Follow-up identifiers coded for each study

  • fups The total number of follow-up measurements in each study

  • narm The total number of arms in each study

If treatment or class are non-numeric or non-sequential (i.e. with missing numeric codes), treatments/classes in the returned data frame will be numbered and recoded to enforce sequential numbering (a warning will be shown stating this).

Examples

# Add indices to osteoarthritis pain dataset
data.ab <- add_index(osteopain)

# Add indices to dataset using different network reference treatment
data.ab <- add_index(osteopain, reference=3)

Studies of alogliptin for lowering blood glucose concentration in patients with type II diabetes

Description

A dataset from a systematic review of Randomised-Controlled Trials (RCTs) comparing different doses of alogliptin with placebo (Langford et al. 2016). The systematic review was simply performed and was intended to provide data to illustrate a statistical methodology rather than for clinical inference. Alogliptin is a treatment aimed at reducing blood glucose concentration in type II diabetes. The outcome is continuous, and aggregate data responses correspond to the mean change in HbA1c from baseline to follow-up. The dataset includes 14 Randomised-Controlled Trials (RCTs), comparing 5 different doses of alogliptin with placebo, leading to 6 different treatments (combination of dose and agent) within the network.

Usage

alog_pcfb

Format

A data frame in long format (one row per arm and study), with 46 rows and 9 variables:

  • studyID Study identifiers

  • clinicaltrialGov_ID The clinicaltrial.gov ID code

  • agent Character data indicating the agent to which participants were randomised

  • dose Numeric data indicating the standardised dose received

  • treatment Character data indicating the treatment (combination of agent and dose) to which participants were randomised

  • time Numeric data indicating the time at which the observation was measured (given in weeks)

  • y Numeric data indicating the mean change from baseline in blood glucose concentration (mg/dL) in a study arm

  • se Numeric data indicating the standard error for the mean change from baseline in blood glucose concentration (mg/dL) in a study arm

  • n Numeric data indicating the number in each arm at each follow-up time

Details

alog_pcfb is a data frame in long format (one row per observation, arm and study), with the variables studyID, clinicaltrialGov_ID, agent, dose, treatment, time, y, se, and n.

References

Langford O, Aronson JK, van Valkenhoef G, Stevens RJ (2016). “Methods for meta-analysis of pharmacodynamic dose-response data with application to multi-arm studies of alogliptin.” Stat Methods Med Res. ISSN 1477-0334 (Electronic) 0962-2802 (Linking), doi:10.1177/0962280216637093, https://pubmed.ncbi.nlm.nih.gov/26994216/.


Plot relative effects from NMAs performed at multiple time-bins

Description

Plot relative effects from NMAs performed at multiple time-bins

Usage

binplot(
  network,
  overlay.nma = c(0, stats::quantile(network$data.ab$time)),
  method = "common",
  link = "identity",
  lim = "cred",
  plot.bins = TRUE,
  legend = TRUE,
  ...
)

Arguments

network

An object of class "mb.network".

overlay.nma

Numeric vector used to overlay the results from a standard NMA model that "lumps" time-points together within the time bin ranges specified in overlay.nma. The numbers in overlay.nma define the boundaries of the time bins within which to perform a standard NMA. Length must be >=2, or can be left as NULL (the default) to indicate that no NMA should be perfomed. overlay.nma can only be specified if overlay.ref==TRUE. See Details for further information.

method

Can take "common" or "random" to indicate the type of NMA model used to synthesise data points given in overlay.nma. The default is "random" since this assumes different time-points in overlay.nma have been lumped together to estimate the NMA.

link

Can take either "identity" (the default), "log" (for modelling Ratios of Means (Friedrich et al. 2011)) or "smd" (for modelling Standardised Mean Differences - although this also corresponds to an identity link function).

lim

Specifies calculation of either 95% credible intervals (lim="cred") or 95% prediction intervals (lim="pred").

plot.bins

Plot time bin boundaries as vertical dashed lines. Setting plot.bins=TRUE if overlay.nma is specified also sets x-axis ticks to time bin boundaries automatically.

legend

TRUE/FALSE to indicate whether a legend should be plotted.

...

Arguments to be sent to R2jags.

Details

Performs several standard NMAs at different time "bins", time periods within which treatment effects are assumed to be constant over time. Separate NMAs are then performed within each time bin on data points from studies that fall within the time bin (only a single follow-up time is taken from each study to avoid double counting).

Note that the wider the time bin boundaries specified by the user, the larger the potential range of included follow-up times and this can introduce heterogeneity or inconsistency.

Results are plotted versus the network reference and are plotted on the specified link scale. Each time bin window is marked on the plot by vertical dashed lines. The NMA estimates within each time bin are plotted as a horizontal solid black line (the posterior median) with a shaded region indicating the 95% credible interval (prediction intervals can instead be plotted). The width of these shaded regions is equal to the range of study time-points included in the NMA performed within that timebin, which may therefore be more narrow than the time bin specified in the binplot() command due to the follow-up times at which data is available in included studies.

Value

Plots treatment effects from NMAs performed within discrete time bins. The object returned is a list containing the plot and a sublist of NMA results and predictions from each time bin specified in overlay.nma.

Overlaying NMA results

overlay.nma indicates regions of the data (defined as "time bins") over which it may be reasonable to "lump" different follow-up times from different studies together and assume a standard NMA model. For example:

  • overlay.nma=c(5,10) indicates a single NMA of studies with follow-up times ⁠>5⁠ and ⁠<=10⁠

  • overlay.nma=c(5,10,15) indicates two NMAs should be performed of studies with follow-up times ⁠>5⁠ and ⁠<=10⁠ of studies with follow-up times ⁠>10⁠ and ⁠<=15⁠

When used with MBNMA (via predict.mbnma()) this allows comparison to MBNMA results over a specific range of time within each time bin. It can be useful to assess which time-course function might be suitable when using binplot(), or to to assess if the MBNMA predictions are in agreement with predictions from an NMA model when using plot.mb.predict() for a specific range of time-points. This can be a general indicator of the fit of the time-course model.

However, it is important to note that the wider the range specified in overlay.nma, the more likely it is that different time-points are included, and therefore that there is greater heterogeneity/inconsistency in the NMA model. If overlay.nma includes several follow-up times for any study then only a single time-point will be taken (the one closest to mean(overlay.nma)). The NMA predictions are plotted over the range specified in overlay.nma as a horizontal line, with the 95%CrI shown by a grey rectangle. The NMA predictions represent those for any time-points within this range since they lump together data at all these time-points. Predictions for treatments that are disconnected from the network reference treatment at data points specified within overlay.nma cannot be estimated so are not included.

It is important to note that the NMA model is not necessarily the "correct" model, since it "lumps" different time-points together and ignores potential differences in treatment effects that may arise from this. The wider the range specified in overlay.nma, the greater the effect of "lumping" and the stronger the assumption of similarity between studies.

For an NMA model to be estimated and a corresponding prediction to be made from it, each time bin must include the network reference treatment (treatment=1) evaluated in at least 1 connected study in the time bin. If a given time bin does not meet this criteria then an NMA will not be calculated for it.

Examples

# Create an mb.network object from a dataset
alognet <- mb.network(alog_pcfb)

# Plot relative effects from NMAs calculated for a single time-bins
# Do not plot time-bin boundaries
binplot(alognet, overlay.nma=c(0,5), plot.bins=FALSE)

# Plot relative effects from NMAs at multiple time-bins
# With random treatment effects
binplot(alognet, overlay.nma=c(5,10,15,20),
  method="random")

Studies comparing Tiotropium, Aclidinium and Placebo for maintenance treatment of moderate to severe chronic obstructive pulmonary disease

Description

A dataset from a systematic review of Randomised-Controlled Trials (RCTs) for maintenance treatment of moderate to severe chronic obstructive pulmonary disease (COPD) (Karabis et al. 2013). Data are extracted from (Tallarita et al. 2019). SEs were imputed for three studies, and number of patients randomised were imputed for one study (LAS 39) in which they were missing, using the median standard deviation calculated from other studies in the dataset. The outcome is trough Forced Expiratory Volume in 1 second (FEV1), measured in litres and reported in each study arm as mean change from baseline to follow-up. The dataset includes 13 Randomised-Controlled Trials (RCTs), comparing 2 treatments (Tiotropium and Aclidinium) and placebo.

Usage

copd

Format

A data frame in long format (one row per arm and study), with 80 rows and 6 variables:

  • studyID Study identifiers

  • time Numeric data indicating the time at which the observation was measured (given in weeks)

  • y Numeric data indicating the mean change from baseline in FEV1 (litres) in a study arm

  • se Numeric data indicating the standard error for the mean change from baseline in FEV1 in a study arm

  • treatment Factor data indicating the treatment to which participants were randomised

  • n Numeric data indicating the number of participants randomised to each arm

Details

copd is a data frame in long format (one row per observation, arm and study), with the variables studyID, time, y, se, treatment, and n.

References

Karabis A, Lindner L, Mocarski M, Huisman E, Greening A (2013). “Comparative efficacy of aclidinium versus glycopyrronium and tiotropium, as maintenance treatment of moderate to severe COPD patients: a systematic review and network meta-analysis.” Int J Chron Obstruct Pulmon Dis, 8, 405-423. doi:10.2147/COPD.S48967, https://pubmed.ncbi.nlm.nih.gov/24043936/.

Tallarita M, De lorio M, Baio G (2019). “A comparative review of network meta-analysis models in longitudinal randomized controlled trial.” Statistics in Medicine, 38(16), 3053-3072. doi:10.1002/sim.8169.


Plot cumulative ranking curves from MBNMA models

Description

Plot cumulative ranking curves from MBNMA models

Usage

cumrank(x, sucra = TRUE, ...)

Arguments

x

An object of class "mb.rank" generated by rank.mbnma()

sucra

A logical object to indicate whether Surface Under Cumulative Ranking Curve (SUCRA) values should be calculated and returned as a data frame. Areas calculated using trapezoid approach.

...

Arguments to be sent to ggplot::geom_line()

Value

Line plots showing the cumulative ranking probabilities for each agent/class for the ranked dose response paramtere in x. The object returned is a list which contains the plot (an object of ⁠class(c("gg", "ggplot")⁠) and a data frame of SUCRA values if sucra = TRUE.

Examples

# Using the alogliptin data
network <- mb.network(alog_pcfb)

# Estimate rankings  from an Emax dose-response MBNMA
emax <- mb.run(network, fun=temax())
ranks <- rank(emax, param=c("emax"))

# Plot cumulative rankings for both dose-response parameters simultaneously
# Note that SUCRA values are also returned
cumrank(ranks)

Sets default priors for JAGS model code

Description

This function creates JAGS code snippets for default MBNMA model priors.

Usage

default.priors(fun = tloglin())

Arguments

fun

An object of class "timefun" generated (see Details) using any of tloglin(), tpoly(), titp(), temax(), tfpoly(), tspline() or tuser()

Value

A list, each element of which is a named JAGS snippet corresponding to a prior in the MBNMA JAGS code.

Examples

default.priors(fun=temax())

default.priors(fun=titp(p.expon=TRUE))

Plot deviance contributions from an MBNMA model

Description

Plot deviance contributions from an MBNMA model

Usage

devplot(
  mbnma,
  dev.type = "dev",
  plot.type = "box",
  xaxis = "time",
  facet = TRUE,
  n.iter = round(mbnma$BUGSoutput$n.iter/4),
  n.thin = mbnma$BUGSoutput$n.thin,
  ...
)

Arguments

mbnma

An S3 object of class "mbnma" generated by running a time-course MBNMA model

dev.type

Deviances to plot - can be either residual deviances ("resdev") or deviances ("dev", the default)

plot.type

Deviances can be plotted either as scatter points ("scatter" - using geom_point()) or as boxplots ("box", the default)

xaxis

A character object that indicates whether deviance contributions should be plotted by time ("time") or by follow-up count ("fup")

facet

A boolean object that indicates whether or not to facet by treatment

n.iter

The number of iterations to update the model whilst monitoring additional parameters (if necessary). Must be a positive integer. Default is the value used in mbnma.

n.thin

The thinning rate. Must be a positive integer. Default is the value used in mbnma.

...

Arguments to be sent to ggplot2::ggplot()

Details

Deviances should only be plotted for models that have converged successfully. If deviance contributions have not been monitored in mbnma$parameters.to.save then additional iterations will have to be run to get results for these.

Deviance contributions cannot be calculated for models with a multivariate likelihood (i.e. those that account for correlation between observations) because the covariance matrix in these models is treated as unknown (if rho="estimate") and deviance contributions will be correlated.

Value

Generates a plot of deviance contributions and returns a list containing the plot (as an object of class c("gg", "ggplot")), and a data.frame of posterior mean deviance/residual deviance contributions for each observation.

Examples

# Make network
alognet <- mb.network(alog_pcfb)

# Run MBNMA
mbnma <- mb.run(alognet, fun=tpoly(degree=2), intercept=FALSE)

# Plot residual deviance contributions in a scatterplot
devplot(mbnma)

# Plot deviance contributions in boxplots at each follow-up measurement
# Monitor for 500 additional iterations
devplot(mbnma, dev.type="dev", plot.type="box", xaxis="fup", n.iter=500)

Studies comparing treatments for type 2 diabetes

Description

A dataset from of trials for reduction of haemoglobin A1c (HbA1c) in patients with type 2 diabetes(Ding and Fu 2013). Data are reported in each study arm as mean change from baseline to follow-up. The dataset includes 4 Randomised-Controlled Trials (RCTs), comparing 4 treatments.

Usage

diabetes

Format

A data frame in long format (one row per arm and study), with 28 rows and 7 variables:

  • studyID Study identifiers

  • treatment Numeric data indicating the treatment to which participants were randomised

  • time Numeric data indicating the time at which the observation was measured (given in weeks)

  • y Numeric data indicating the mean change from baseline in HbA1c in a study arm

  • se Numeric data indicating the standard error for the mean change from baseline in HbA1c in a study arm

  • sd Numeric data indicating the standard deviation for the mean change from baseline in HbA1c in a study arm

  • n Numeric data indicating the number of participants in each arm at each time-point

Details

diabetes is a data frame in long format (one row per observation, arm and study), with the variables studyID, treatment, time, y, se, sd, and n.

References

Ding Y, Fu H (2013). “Bayesian indirect and mixed treatment comparisons across longitudinal time points.” Statistics in Medicine, 32, 2613-2628. doi:10.1002/sim.5688.


Plot fitted values from MBNMA model

Description

Plot fitted values from MBNMA model

Usage

fitplot(
  mbnma,
  treat.labs = NULL,
  disp.obs = TRUE,
  n.iter = round(mbnma$BUGSoutput$n.iter/4),
  n.thin = mbnma$BUGSoutput$n.thin,
  ...
)

Arguments

mbnma

An S3 object of class "mbnma" generated by running a time-course MBNMA model

treat.labs

A character vector of treatment labels with which to name graph panels. Can use mb.network()[["treatments"]] with original dataset if in doubt.

disp.obs

A boolean object to indicate whether raw data responses should be plotted as points on the graph

n.iter

number of total iterations per chain (including burn in; default: 2000)

n.thin

thinning rate. Must be a positive integer. Set n.thin > 1 to save memory and computation time if n.iter is large. Default is max(1, floor(n.chains * (n.iter-n.burnin) / 1000)) which will only thin if there are at least 2000 simulations.

...

Arguments to be sent to ggplot2::ggplot()

Details

Fitted values should only be plotted for models that have converged successfully. If fitted values (theta) have not been monitored in mbnma$parameters.to.save then additional iterations will have to be run to get results for these.

Value

Generates a plot of fitted values from the MBNMA model and returns a list containing the plot (as an object of class c("gg", "ggplot")), and a data.frame of posterior mean fitted values for each observation.

Examples

# Make network
painnet <- mb.network(osteopain)

# Run MBNMA
mbnma <- mb.run(painnet,
  fun=temax(pool.emax="rel", method.emax="common",
    pool.et50="abs", method.et50="random"))

# Plot fitted values from the model
# Monitor fitted values for 500 additional iterations
fitplot(mbnma, n.iter=500)

Automatically generate parameters to save for a time-course MBNMA model

Description

Automatically generate parameters to save for a time-course MBNMA model

Usage

gen.parameters.to.save(fun, model)

Arguments

fun

An object of class "timefun" generated (see Details) using any of tloglin(), tpoly(), titp(), temax(), tfpoly(), tspline() or tuser()

model

A JAGS model written as a character object

Value

A character vector of parameter names that should be monitored in the model


Get large vector of distinct colours using Rcolorbrewer

Description

Get large vector of distinct colours using Rcolorbrewer

Usage

genmaxcols()

Generates spline basis matrices for fitting to time-course function

Description

Generates spline basis matrices for fitting to time-course function

Usage

genspline(
  x,
  spline = "bs",
  knots = 1,
  degree = 1,
  max.time = max(x),
  boundaries = NULL
)

Arguments

x

A numeric vector indicating all time points available in the dataset

spline

Indicates the type of spline function. Can be either a piecewise linear spline ("ls"), natural cubic spline ("ns") or B-spline ("bs").

knots

The number/location of knots. If a single integer is given it indicates the number of knots (they will be equally spaced across the range of time-points). If a numeric vector is given it indicates the quantiles of the knots as a proportion of the maximum study follow-up in the dataset. For example, if the maximum follow-up time in the dataset is 10 months, knots=c(0.1,0.5) would indicate knots should be fitted at 1 and 5 months follow-up.

degree

a positive integer giving the degree of the polynomial from which the spline function is composed (e.g. degree=3 represents a cubic spline).

max.time

A number indicating the maximum time between which to calculate the spline function.

boundaries

A positive numeric vector of length 2 that represents the time-points at which to anchor the B-spline or natural cubic spline basis matrix. This allows data to extend beyond the boundary knots, or for the basis parameters to not depend on x. The default (boundaries=NULL)is the range of x.

Value

A spline basis matrix with number of rows equal to length(x) and the number of columns equal to the number of coefficients in the spline.

Examples

x <- 0:100

genspline(x)

# Generate a quadratic B-spline with 1 equally spaced internal knot
genspline(x, spline="bs", knots=2, degree=2)

# Generate a natural spline with 2 knots at selected quantiles
genspline(x, spline="ns", knots=c(0.1, 0.5))

# Generate a piecewise linear spline with 3 equally spaced knots
genspline(x, spline="ls", knots=3)

Create a dataset with a single time point from each study closest to specified time

Description

Takes the closest time point from each arm in each study to a specified time (t) within an mb.network object. Useful for network plots or exploring standard NMA.

Usage

get.closest.time(network, t = stats::median(network$data.ab$time))

Arguments

network

An object of class "mb.network".

t

The time point at which

Value

A data frame in long format of responses at the closest time point to t in each arm of each study.

Examples

# Using the alogliptin dataset
network <- mb.network(alog_pcfb)

# Take a single follow-up time from each study...
# ...closest to 7
get.closest.time(network, t=7)

# ...closest to 20
get.closest.time(network, t=7)

# ...closest to the median follow-up across all studies
get.closest.time(network, t=26)

Create a dataset with the earliest time point only

Description

Takes the earliest time point from each arm in each study within an mb.network object. Useful for network plots.

Usage

get.earliest.time(network)

Arguments

network

An object of class "mb.network".

Value

A list containing:

  • a data frame in long format of responses at the earliest time point in each arm of each study

  • a vector of studyIDs

  • a vector of treatment names

  • a vector of class names (if included in network)

  • a data frame of treatment -> class codings (if included in network)

Examples

# Using the alogliptin dataset
network <- mb.network(alog_pcfb)

# Generate a data set with only the earliest time point included in each study
get.earliest.time(network)

Create a dataset with the latest time point only

Description

Takes the latest time point from each arm in each study within an mb.network object. Useful for network plots.

Usage

get.latest.time(network)

Arguments

network

An object of class "mb.network".

Value

A list containing:

  • a data frame in long format of responses at the latest time point in each arm of each study

  • a vector of studyIDs

  • a vector of treatment names

  • a vector of class names (if included in network)

  • a data frame of treatment -> class codings (if included in network)

Examples

# Using the alogliptin dataset
network <- mb.network(alog_pcfb)

# Generate a data frame with only the latest time point included in each study
get.latest.time(network)

Get MBNMA model values

Description

Extracts specific information required for prediction from a time-course MBNMA model

Usage

get.model.vals(
  mbnma,
  E0 = 0,
  level = "treatments",
  lim = "cred",
  link = "identity"
)

Arguments

mbnma

An S3 object of class "mbnma" generated by running a time-course MBNMA model

E0

An object to indicate the value(s) to use for the response at time = 0 in the prediction. This can take a number of different formats depending on how it will be used/calculated. The default is 0 but this may lead to non-sensical predictions if Ratio of Means are modeled.

  • numeric() A single numeric value representing the deterministic response at time = 0

  • formula() A formula representing a stochastic distribution for the response at time = 0. This is specified as a random number generator (RNG) given as a string, and can take any RNG distribution for which a function exists in R. For example: ~rnorm(n, 7, 0.5).

level

Can take either "treatment" to make predictions for treatments, or "class" to make predictions for classes (in which case object must be a class effect model).

lim

Specifies calculation of either 95% credible intervals (lim="cred") or 95% prediction intervals (lim="pred").

link

Can take either "identity" (the default), "log" (for modelling Ratios of Means (Friedrich et al. 2011)) or "smd" (for modelling Standardised Mean Differences - although this also corresponds to an identity link function).

Value

A list containing named elements that correspond to different time-course parameters in mbnma. These elements contain MCMC results either taken directly from mbnma or (in the case of random time-course parameters specified as method="random") randomly generated using parameter values estimated in mbnma.

Additional elements contain the following values:

  • timecourse A character object that specifies the time-course used in mbnma in terms of alpha, beta, mu, d and time. Consistency relative time-course parameters are specified in terms of mu and d.

  • time.params A character vector that indicates the different time-course parameters that are required for the prediction

@noRd


Get current priors from JAGS model code

Description

Identical to get.prior() in MBNMAdose. This function takes JAGS model presented as a string and identifies what prior values have been used for calculation.

Usage

get.prior(model)

Arguments

model

A character object of JAGS MBNMA model code

Details

Even if an MBNMA model that has not initialised successfully and results have not been calculated, the JAGS model for it is saved in MBNMA$model.arg$jagscode and therefore priors can still be obtained. This allows for priors to be changed even in failing models, which may help solve issues with initialisation.

Value

A character vector, each element of which is a line of JAGS code corresponding to a prior in the JAGS code.

Examples

# Create mb.network object using an MBNMAtime dataset
network <- mb.network(osteopain)

# Create mb.network object using an MBNMAdose dataset

# Run linear MBNMA
result <- mb.run(network, fun=tpoly(degree=1,
    pool.1="rel", method.1="random"))

# Obtain model prior values
get.prior(result$model.arg$jagscode)

# ...also equivalent to
print(result$model.arg$priors)

Calculates relative effects/mean differences at a particular time-point

Description

Uses mbnma time-course parameter estimates to calculate treatment differences between treatments or classes at a particular time-point. Can be used to compare treatments evaluated in studies at different follow-up times, or even to compare treatments in different MBNMA models via a common comparator.

Usage

get.relative(
  mbnma,
  mbnma.add = NULL,
  time = max(mbnma$model.arg$jagsdata$time, na.rm = TRUE),
  treats = unique(c(mbnma$network$treatments, mbnma.add$network$treatments)),
  classes = NULL,
  lim = "cred"
)

Arguments

mbnma

An S3 object of class "mbnma" generated by running a time-course MBNMA model

mbnma.add

An S3 object of class("mbnma") generated by running a time-course MBNMA model. This should only be specified if results from two different MBNMA models are to be combined to perform a 2-stage MBNMA (see Details).

time

A numeric value for the time at which to estimate relative effects/mean differences.

treats

A character vector of treatment names for which to calculate relative effects/mean differences. Must be a subset of mbnma$network$treatments.

classes

A character vector of class names for which to calculate relative effects/mean differences. Must be a subset of mbnma$network$classes. Only works for class effect models.

lim

Specifies calculation of either 95% credible intervals (lim="cred") or 95% prediction intervals (lim="pred").

Details

get.relative() can also be used to perform a 2-stage MBNMA that allows synthesis of results from two different MBNMA models via a single common comparator. In an MBNMA model, all treatments must share the same time-course function. However, a 2-stage approach can enable fitting of different time-course functions to different sets ("subnetworks") of treatments. For example, some treatments may have rich time-course information, allowing for a more complex time-course function to be used, whereas others may be sparse, requiring a simpler time-course function.

Relative comparisons between treatments in the two datasets at specific follow-up times can then be estimated from MBNMA predicted effects versus a common comparator using the Bucher method and assuming consistency. See the MBNMAtime vignette for further details.

Value

An object of class "relative.array" list containing:

  • The time-point for which results are estimated

  • Matrices of posterior means, medians, SDs and upper and lower 95% credible intervals for the differences between each treatment

  • An array containing MCMC results for the differences between all treatments specified in treats or all classes specified in classes.

Results are reported in tables as the row-defined treatment minus the column-defined treatment.

Examples

# Create an mb.network object from a dataset
alognet <- mb.network(alog_pcfb)

# Run a quadratic time-course MBNMA using the alogliptin dataset
mbnma <- mb.run(alognet,
  fun=tpoly(degree=2,
  pool.1="rel", method.1="random",
  pool.2="rel", method.2="common"
  )
)

# Calculate differences between all treatments at 20 weeks follow-up
allres <- get.relative(mbnma, time=20)

# Calculate difference between a subset of treatments at 10 weeks follow-up
subres <- get.relative(mbnma, time=10,
  treats=c("alog_50", "alog_25", "placebo"))



###########################
##### 2-stage MBNMA #####
###########################

# Using the osteoarthritis dataset
# With placebo (Pl_0) as common comparator between subnetworks

#### Sparse model ####

# Treatments on which time-course data is limited
sparse.trt <- c("Ce_100", "Ce_400", "Du_90", "Lu_200", "Lu_400",
  "Lu_NA", "Et_5", "Ox_44")

# Create a subnetwork of studies comparing these treatments
sparse.df <- osteopain %>% dplyr::group_by(studyID) %>%
  dplyr::filter(any(treatment %in% sparse.trt)) %>%
  dplyr::ungroup() %>%
  subset(treatment %in% c("Pl_0", sparse.trt))

sparse.net <- mb.network(sparse.df)

# Run a ITP MBNMA with a known rate
sparse.mbnma <- mb.run(sparse.net, fun=titp(method.rate=0.8, pool.rate="abs"))


#### Complex model ####

# Treatments on which time-course data is rich
rich.trt <- levels(osteopain$treatment)[!levels(osteopain$treatment) %in%
  c("Pl_0", "Ce_100", "Ce_400", "Du_90", "Lu_200",
    "Lu_400", "Lu_NA", "Et_5", "Ox_44")]

# Create a subnetwork of studies comparing these treatments
rich.df <- osteopain %>% dplyr::group_by(studyID) %>%
  dplyr::filter(any(treatment %in% rich.trt)) %>%
  dplyr::ungroup() %>%
  subset(treatment %in% c("Pl_0", rich.trt))

rich.net <- mb.network(rich.df)

# Run a Emax MBNMA
rich.mbnma <- mb.run(rich.net, temax(p.expon = FALSE))


#### Calculate relative effects between models ####

# At 10 weeks follow-up
rels.sparse <- get.relative(sparse.mbnma, time=10)
rels.rich <- get.relative(rich.mbnma, time=10)

rels.all <- get.relative(mbnma=rich.mbnma,
  mbnma.add=sparse.mbnma, time=10)

print(rels.all$median)

Prepares data for JAGS

Description

Converts MBNMA data frame to a list for use in JAGS model

Usage

getjagsdata(
  data.ab,
  fun = NULL,
  class = FALSE,
  rho = NULL,
  covstruct = "CS",
  link = "identity",
  sdscale = FALSE,
  cfb = NULL
)

Arguments

data.ab

A data frame of arm-level data in "long" format containing the columns:

  • studyID Study identifiers

  • time Numeric data indicating follow-up times

  • y Numeric data indicating the aggregate response for a given observation (e.g. mean)

  • se Numeric data indicating the standard error for a given observation

  • treatment Treatment identifiers (can be numeric, factor or character)

  • class An optional column indicating a particular class identifier. Observations with the same treatment identifier must also have the same class identifier.

  • n An optional column indicating the number of participants used to calculate the response at a given observation (required if modelling using Standardised Mean Differences)

  • standsd An optional column of numeric data indicating reference SDs used to standardise treatment effects when modelling using Standardised Mean Differences (SMD).

fun

An object of class "timefun" generated (see Details) using any of tloglin(), tpoly(), titp(), temax(), tfpoly(), tspline() or tuser()

class

A boolean object indicating whether or not data.ab contains information on different classes of treatments

rho

The correlation coefficient when modelling within-study correlation between time points. The default is a string representing a prior distribution in JAGS, indicating that it be estimated from the data (e.g. rho="dunif(0,1)"). rho also be assigned a numeric value (e.g. rho=0.7), which fixes rho in the model to this value (e.g. for use in a deterministic sensitivity analysis). If set to rho=0 (the default) then this implies modelling no correlation between time points.

covstruct

A character to indicate the covariance structure required for modelling correlation between time points (if any), since this determines some of the data. Can be either "CS" (compound symmetry), "AR1" (autoregressive AR1) or "varadj" (variance-adjustment).

link

Can take either "identity" (the default), "log" (for modelling Ratios of Means (Friedrich et al. 2011)) or "smd" (for modelling Standardised Mean Differences - although this also corresponds to an identity link function).

sdscale

Logical object to indicate whether to write a model that specifies a reference SD for standardising when modelling using Standardised Mean Differences. Specifying sdscale=TRUE will therefore only modify the model if link function is set to SMD (link="smd").

cfb

A logical vector whose length is equal to the unique number of studies in data.ab, where each element is TRUE if the study data reported is change-from-baseline and FALSE otherwise. If left as NULL (the default) then this will be identified from the data by assuming any study for which there is no data at time=0 reports change-from-baseline.

Value

A named list of numbers, vector, matrices and arrays to be sent to JAGS. List elements are:

  • y An array of mean responses for each observation in each arm within each study

  • se An array of standard errors for each observation in each arm within each study

  • time A matrix of follow-up times within each study

  • fups A numeric vector with the number of follow-up measurements per study

  • narm A numeric vector with the number of arms per study

  • NS The total number of studies in the dataset

  • NT The total number of treatments in the dataset

  • treat A matrix of treatment codes within each study

  • Nclass Optional. The total number of classes in the dataset

  • class Optional. A matrix of class codes within each study

  • classkey Optional. A vector of class codes that correspond to treatment codes. Same length as the number of treatment codes.

  • mat.triangle Optional. A matrix with number indicating how to fill covariance matrices within the JAGS code.

  • mat.order Optional. A matrix with number indicating what order to fill covariance matrices within the JAGS code.

  • timedif.0 Optional. A vector of the difference in times between the first and second follow-up time in each study.

Examples

# Using the alogliptin dataset
network <- mb.network(alog_pcfb)
jagsdat <- getjagsdata(network$data.ab)


# Get JAGS data with class
netclass <- mb.network(goutSUA_CFBcomb)
jagsdat <- getjagsdata(netclass$data.ab, class=TRUE)


# Get JAGS data that allows for modelling correlation between time points
painnet <- mb.network(osteopain)
jagsdat <- getjagsdata(painnet$data.ab, rho="dunif(0,1)", covstruct="AR1")

Prepares NMA data for JAGS

Description

Converts data frame to a list for use in JAGS NMA model

Usage

getnmadata(data.ab, link = "identity", sdscale = FALSE)

Arguments

data.ab

A data frame of arm-level data in "long" format containing the columns:

  • studyID Study identifiers

  • time Numeric data indicating follow-up times

  • y Numeric data indicating the aggregate response for a given observation (e.g. mean)

  • se Numeric data indicating the standard error for a given observation

  • treatment Treatment identifiers (can be numeric, factor or character)

  • class An optional column indicating a particular class identifier. Observations with the same treatment identifier must also have the same class identifier.

  • n An optional column indicating the number of participants used to calculate the response at a given observation (required if modelling using Standardised Mean Differences)

  • standsd An optional column of numeric data indicating reference SDs used to standardise treatment effects when modelling using Standardised Mean Differences (SMD).

link

Can take either "identity" (the default), "log" (for modelling Ratios of Means (Friedrich et al. 2011)) or "smd" (for modelling Standardised Mean Differences - although this also corresponds to an identity link function).

sdscale

Logical object to indicate whether to write a model that specifies a reference SD for standardising when modelling using Standardised Mean Differences. Specifying sdscale=TRUE will therefore only modify the model if link function is set to SMD (link="smd").

Value

A named list of numbers, vector, matrices and arrays to be sent to JAGS. List elements are:

  • y An array of mean responses for each observation in each arm within each study

  • se An array of standard errors for each observation in each arm within each study

  • narm A numeric vector with the number of arms per study

  • NS The total number of studies in the dataset

  • NT The total number of treatments in the dataset

  • treat A matrix of treatment codes within each study

Examples

# Using the alogliptin dataset
network <- mb.network(alog_pcfb)

# Construct a dataset with the latest time point in each study
data.ab <- get.latest.time(network)$data.ab
getnmadata(data.ab)

Studies of treatments for reducing serum uric acid in patients with gout

Description

A dataset from a systematic review of interventions for lowering Serum Uric Acid (SUA) concentration in patients with gout (not published previously). The outcome is continuous, and aggregate data responses correspond to the mean change from baseline in SUA in mg/dL. Overall there are 41 treatments of 8 agents in the network. Standard deviations have been imputed for 181 observations.

Usage

goutSUA_CFB

Format

A data frame with 224 rows and 7 variables:

  • studyID Study identifiers

  • time Numeric data indicating follow-up times

  • y Numeric data indicating the mean response for a given observation

  • se Numeric data indicating the standard error for a given observation

  • treatment Treatment identifiers as factors. Labels are shortened treatment names.

  • treatname Character data giving the full names of each treatment in the format agent_dose

  • class Shortened agent names stored as factors.

Details

goutSUA_CFB is a data frame in long format (one row per observation, arm and study), with the variables studyID, time, y, se, treatment, treatname and class.

Source

Pfizer Ltd.


Studies of combined treatments for reducing serum uric acid in patients with gout

Description

A dataset from a systematic review of interventions for lowering Serum Uric Acid (SUA) concentration in patients with gout (not published previously). The outcome is continuous, and aggregate data responses correspond to the mean change from baseline in SUA in mg/dL. Treatments with similar doses have been pooled together to improve network connectivity and facilitate evidence synthesis, resulting in 19 treatments of 7 agents included in the network. Standard deviations have been imputed for 181 observations.

Usage

goutSUA_CFBcomb

Format

A data frame with 224 rows and 7 variables:

  • studyID Study identifiers

  • time Numeric data indicating follow-up times

  • y Numeric data indicating the mean response for a given observation

  • se Numeric data indicating the standard error for a given observation

  • treatment Treatment identifiers as factors. Labels are shortened treatment names.

  • treatname Character data giving the full names of each treatment in the format agent_dose

  • class Shortened agent names stored as factors.

Details

goutSUA_CFBcomb is a data frame in long format (one row per observation, arm and study), with the variables studyID, time, y, se, treatment, treatname and class.

Source

Pfizer Ltd.


Studies comparing hyaluronan (HA)–based viscosupplements for osteoarthritis

Description

A dataset from of trials for pain reduction for patients with osteoarthritis treated with HA-based viscosupplements(Jansen et al. 2015). Data are reported in each study arm as mean change from baseline to follow-up on a visual analogue scale (0-100). The dataset includes 16 Randomised-Controlled Trials (RCTs), comparing 6 treatments and placebo.

Usage

hyalarthritis

Format

A data frame in long format (one row per arm and study), with 150 rows and 6 variables:

  • studyID Study identifiers

  • time Numeric data indicating the time at which the observation was measured (given in weeks)

  • treatment Factor data indicating the treatment to which participants were randomised

  • n Numeric data indicating the number of participants randomised to each arm

  • y Numeric data indicating the mean change from baseline in HbA1c in a study arm

  • se Numeric data indicating the standard error for the mean change from baseline in HbA1c in a study arm

Details

hyalarthritis is a data frame in long format (one row per observation, arm and study), with the variables studyID, time, treatment, n, y, and se.

References

Jansen JP, Vieira MC, Cope S (2015). “Network meta-analysis of longitudinal data using fractional polynomials.” Stat Med, 34(15), 2294-311. ISSN 1097-0258 (Electronic) 0277-6715 (Linking), doi:10.1002/sim.6492, https://pubmed.ncbi.nlm.nih.gov/25877808/.


Identify comparisons in loops that fulfil criteria for node-splitting

Description

Identify comparisons informed by both direct and indirect evidence from independent sources, which therefore fulfil the criteria for testing for inconsistency via node-splitting. Follows the method of van Valkenhoef van Valkenhoef et al. (2016).

Usage

inconsistency.loops(data)

Arguments

data

A data frame containing variables studyID and treatment (as numeric codes) that indicate which treatments are used in which studies.

Details

Similar to gemtc::mtc.nodesplit() but uses a fixed reference treatment and therefore suggests fewer loops in which to test for inconsistency. Heterogeneity can also be parameterised as inconsistency and so testing for inconsistency in additional loops whilst changing the reference treatment would also be identifying heterogeneity. Depends on igraph.

Value

A data frame of comparisons that are informed by direct and indirect evidence from independent sources. Each row of the data frame is a different treatment comparison. Numerical codes in t1 and t2 correspond to treatment codes.

References

van Valkenhoef G, Dias S, Ades AE, Welton NJ (2016). “Automated generation of node-splitting models for assessment of inconsistency in network meta-analysis.” Res Synth Methods, 7(1), 80-93. ISSN 1759-2887 (Electronic) 1759-2879 (Linking), doi:10.1002/jrsm.1167, https://pubmed.ncbi.nlm.nih.gov/26461181/.

Examples

data <- data.frame(studyID=c(1,1,2,2,3,3,4,4,5,5,5),
  treatment=c(1,2,1,3,2,3,3,4,1,2,4)
  )

# Identify comparisons informed by direct and indirect evidence
inconsistency.loops(data)

Identify unique comparisons within a network (identical to MBNMAdose)

Description

Identify unique contrasts within a network that make up all the head-to-head comparisons. Repetitions of the same treatment comparison are grouped together.

Usage

mb.comparisons(data)

Arguments

data

A data frame containing variables studyID and treatment (as numeric codes) that indicate which treatments are used in which studies.

Value

A data frame of unique comparisons in which each row represents a different comparison. t1 and t2 indicate the treatment codes that make up the comparison. nr indicates the number of times the given comparison is made within the network.

If there is only a single observation for each study within the dataset (i.e. as for standard network meta-analysis) nr will represent the number of studies that compare treatments t1 and t2.

If there are multiple observations for each study within the dataset (as in time-course MBNMA) nr will represent the number of time points in the dataset in which treatments t1 and t2 are compared.

Examples

data <- data.frame(studyID=c(1,1,2,2,3,3,4,4,5,5,5),
  treatment=c(1,2,1,3,2,3,3,4,1,2,4)
  )

# Identify comparisons informed by direct and indirect evidence
mb.comparisons(data)

Convert arm-based MBNMA data to contrast data

Description

Converts an object of class mb.network from arm-based long MBNMA data to a data frame with contrast data (a separate contrast for each treatment comparison at each time point within each study). Data can be either long or wide.

Usage

mb.make.contrast(network, datatype = NULL, format = "wide")

Arguments

network

An object of class mb.network

datatype

A string indicating the data type. Can be binomial or normal

format

A string indicating the data format. Can be wide (two additional columns for each variable - contrast arms) or long.

Value

A data frame with the following columns. In wide format, some columns are given the indices 1 and 2 to indicate each arm in a given treatment comparison.:

  • t The treatment in each arm

  • TE The treatment effect (mean difference, log-odds) for the treatment in arm 1 versus the treatment in arm 2

  • seTE The standard error for the treatment effect (mean difference, log-odds) for the treatment in arm 1 versus the treatment in arm 2

  • y The mean response in each arm

  • se The standard error of the mean in each arm

  • r The number of responders in each arm

  • n The total number of participants in each arm

  • fupcount Follow-up identifier

  • time The time the data are reported

  • studyID Study identifier

Examples

# Create mb.network
network <- mb.network(osteopain)

# Convert to wide contrast data
mb.make.contrast(network, format="wide")

# Convert to long contrast data
mb.make.contrast(network, format="long")

Identify comparisons in time-course MBNMA datasets that fulfil criteria for node-splitting

Description

Identify comparisons informed by both direct and indirect evidence from independent sources in MBNMA datasets with repeated measurements in each study. These comparisons are therefore those which fulfil the criteria for testing for inconsistency via node-splitting, following the method of van Valkenhoef van Valkenhoef et al. (2016).

Usage

mb.nodesplit.comparisons(network)

Arguments

network

An object of class "mb.network".

Details

Similar to gemtc::mtc.nodesplit() but uses a fixed reference treatment and therefore suggests fewer loops in which to test for inconsistency. Heterogeneity can also be parameterised as inconsistency and so testing for inconsistency in additional loops whilst changing the reference treatment would also be identifying heterogeneity. Depends on igraph.

Value

A data frame of comparisons that are informed by direct and indirect evidence from independent sources. Each row of the data frame is a different treatment comparison. Numerical codes in t1 and t2 correspond to treatment codes.

References

van Valkenhoef G, Dias S, Ades AE, Welton NJ (2016). “Automated generation of node-splitting models for assessment of inconsistency in network meta-analysis.” Res Synth Methods, 7(1), 80-93. ISSN 1759-2887 (Electronic) 1759-2879 (Linking), doi:10.1002/jrsm.1167, https://pubmed.ncbi.nlm.nih.gov/26461181/.

Examples

# Create mb.network object
network <- mb.network(osteopain)

# Identify comparisons informed by direct and indirect evidence
mb.nodesplit.comparisons(network)

Run MBNMA time-course models

Description

Fits a Bayesian time-course model for model-based network meta-analysis (MBNMA) that can account for repeated measures over time within studies by applying a desired time-course function. Follows the methods of Pedder et al. (2019).

Usage

mb.run(
  network,
  fun = tpoly(degree = 1),
  positive.scale = FALSE,
  intercept = NULL,
  link = "identity",
  sdscale = FALSE,
  parameters.to.save = NULL,
  rho = 0,
  covar = "varadj",
  omega = NULL,
  corparam = FALSE,
  class.effect = list(),
  UME = FALSE,
  pd = "pv",
  parallel = FALSE,
  priors = NULL,
  n.iter = 20000,
  n.chains = 3,
  n.burnin = floor(n.iter/2),
  n.thin = max(1, floor((n.iter - n.burnin)/1000)),
  model.file = NULL,
  jagsdata = NULL,
  ...
)

Arguments

network

An object of class "mb.network".

fun

An object of class "timefun" generated (see Details) using any of tloglin(), tpoly(), titp(), temax(), tfpoly(), tspline() or tuser()

positive.scale

A boolean object that indicates whether all continuous mean responses (y) are positive and therefore whether the baseline response should be given a prior that constrains it to be positive (e.g. for scales that cannot be <0).

intercept

A boolean object that indicates whether an intercept (written as alpha in the model) is to be included. If left as NULL (the default), an intercept will be included only for studies reporting absolute means, and will be excluded for studies reporting change from baseline (as indicated in network$cfb).

link

Can take either "identity" (the default), "log" (for modelling Ratios of Means (Friedrich et al. 2011)) or "smd" (for modelling Standardised Mean Differences - although this also corresponds to an identity link function).

sdscale

Logical object to indicate whether to write a model that specifies a reference SD for standardising when modelling using Standardised Mean Differences. Specifying sdscale=TRUE will therefore only modify the model if link function is set to SMD (link="smd").

parameters.to.save

A character vector containing names of parameters to monitor in JAGS

rho

The correlation coefficient when modelling within-study correlation between time points. The default is a string representing a prior distribution in JAGS, indicating that it be estimated from the data (e.g. rho="dunif(0,1)"). rho also be assigned a numeric value (e.g. rho=0.7), which fixes rho in the model to this value (e.g. for use in a deterministic sensitivity analysis). If set to rho=0 (the default) then this implies modelling no correlation between time points.

covar

A character specifying the covariance structure to use for modelling within-study correlation between time-points. This can be done by specifying one of the following:

  • "varadj" - a univariate likelihood with a variance adjustment to assume a constant correlation between subsequent time points (Jansen et al. 2015). This is the default.

  • "CS" - a multivariate normal likelihood with a compound symmetry structure

  • "AR1" - a multivariate normal likelihood with an autoregressive AR1 structure

omega

DEPRECATED IN VERSION 0.2.3 ONWARDS (~uniform(-1,1) now used for correlation between parameters rather than a Wishart prior). A scale matrix for the inverse-Wishart prior for the covariance matrix used to model the correlation between time-course parameters (see Details for time-course functions). omega must be a symmetric positive definite matrix with dimensions equal to the number of time-course parameters modelled using relative effects (pool="rel"). If left as NULL (the default) a diagonal matrix with elements equal to 1 is used.

corparam

A boolean object that indicates whether correlation should be modeled between relative effect time-course parameters. Default is FALSE and this is automatically set to FALSE if class effects are modeled. Setting it to TRUE models correlation between time-course parameters. This can help identify parameters that are estimated poorly for some treatments by allowing sharing of information between parameters for different treatments in the network, but may also cause some shrinkage.

class.effect

A list of named strings that determines which time-course parameters to model with a class effect and what that effect should be ("common" or "random"). For example: list(emax="common", et50="random").

UME

Can take either TRUE or FALSE (for an unrelated mean effects model on all or no time-course parameters respectively) or can be a vector of parameter name strings to model as UME. For example: c("beta.1", "beta.2").

pd

Can take either:

  • pv only pV will be reported (as automatically outputted by R2jags).

  • plugin calculates pD by the plug-in method (Spiegelhalter et al. 2002). It is faster, but may output negative non-sensical values, due to skewed deviances that can arise with non-linear models.

  • pd.kl (the default) calculates pD by the Kullback–Leibler divergence (Plummer 2008). This will require running the model for additional iterations but will always produce a sensical result.

  • popt calculates pD using an optimism adjustment which allows for calculation of the penalized expected deviance (Plummer 2008)

parallel

A boolean value that indicates whether JAGS should be run in parallel (TRUE) or not (FALSE). If TRUE then the number of cores to use is automatically calculated. Functions that involve updating the model (e.g. devplot(), fitplot()) cannot be used with models implemented in parallel.

priors

A named list of parameter values (without indices) and replacement prior distribution values given as strings using distributions as specified in JAGS syntax (see Plummer (2017)).

n.iter

number of total iterations per chain (including burn in; default: 20000)

n.chains

number of Markov chains (default: 3)

n.burnin

length of burn in, i.e. number of iterations to discard at the beginning. Default is ⁠n.iter/2``, that is, discarding the first half of the simulations. If ⁠n.burnin' is 0, jags() will run 100 iterations for adaption.

n.thin

thinning rate. Must be a positive integer. Set ⁠n.thin > 1`` to save memory and computation time if ⁠n.iter⁠is large. Default is⁠max(1, floor(n.chains * (n.iter-n.burnin) / 1000))“ which will only thin if there are at least 2000 simulations.

model.file

The file path to a JAGS model (.jags file) that can be used to overwrite the JAGS model that is automatically written based on the specified options in MBNMAtime. Useful for adding further model flexibility.

jagsdata

A named list of the data objects to be used in the JAGS model. Only required if users are defining their own JAGS model using model.file. Format should match that of standard models fitted in MBNMAtime (see mbnma$model.arg$jagsdata)

...

Arguments to be sent to R2jags.

Value

An object of S3 class ⁠c("mbnma", "rjags")`` containing parameter results from the model. Can be summarized by ⁠print()⁠and can check traceplots using⁠R2jags::traceplot()⁠or various functions from the package⁠mcmcplots‘.#’

If there are errors in the JAGS model code then the object will be a list consisting of two elements - an error message from JAGS that can help with debugging and model.arg, a list of arguments provided to mb.run() which includes jagscode, the JAGS code for the model that can help users identify the source of the error.

Time-course parameters

Nodes that are automatically monitored (if present in the model) have the same name as in the time-course function for named time-course parameters (e.g. emax). However, for named only as beta.1, beta.2, beta.3 or beta.4 parameters may have an alternative interpretation.

Details of the interpretation and model specification of different parameters can be shown by using the summary() method on an "mbnma" object generated by mb.run().

Parameters modelled using relative effects

  • If pooling is relative (e.g. pool.1="rel") for a given parameter then the named parameter (e.g. emax) or a numbered d parameter (e.g. d.1) corresponds to the pooled relative effect for a given treatment compared to the network reference treatment for this time-course parameter.

  • sd. followed by a named (e.g. emax, beta.1) is the between-study SD (heterogeneity) for relative effects, reported if pooling for a time-course parameter is relative (e.g. pool.1="rel") and the method for synthesis is random (e.g. ⁠method.1="random⁠).

  • If class effects are modelled, parameters for classes are represented by the upper case name of the time-course parameter they correspond to. For example if class.effect=list(emax="random"), relative class effects will be represented by EMAX. The SD of the class effect (e.g. sd.EMAX, sd.BETA.1) is the SD of treatments within a class for the time-course parameter they correspond to.

Parameters modelled using absolute effects

  • If pooling is absolute (e.g. pool.1="abs") for a given parameter then the named parameter (e.g. emax) or a numbered beta parameter (e.g. beta.1) corresponds to the estimated absolute effect for this time-course parameter.

  • For an absolute time-course parameter if the corresponding method is common (e.g. method.1="common") the parameter corresponds to a single common parameter estimated across all studies and treatments. If the corresponding method is random (e.g. method.1="random") then parameter is a mean effect around which the study-level absolute effects vary with SD corresponding to sd. followed by the named parameter (e.g. sd.emax, sd.beta.1) .

Other model parameters

  • rho The correlation coefficient for correlation between time-points. Its interpretation will differ depending on the covariance structure specified in covar

  • totresdev The residual deviance of the model

  • deviance The deviance of the model

Time-course function

Several general time-course functions with up to 4 time-course parameters are provided, but a user-defined time-course relationship can instead be used. Details can be found in the respective help files for each function.

Available time-course functions are:

  • Log-linear: tloglin()

  • Polynomial: tpoly()

  • Integrated Two-Component Prediction (ITP): titp()

  • Emax: temax()

  • Fractional polynomial: tfpoly()

  • Splines (various spline types can be used): tspline()

  • User-defined: tuser()

Correlation between observations

When modelling correlation between observations using rho, values for rho must imply a positive semidefinite covariance matrix.

Advanced options

model.file and jagsdata can be used to run an edited JAGS model and dataset. This allows users considerably more modelling flexibility than is possible using the basic MBNMAtime syntax, though requires strong understanding of JAGS and the MBNMA modelling framework. Treatment-specific priors, meta-regression and bias-adjustment are all possible in this way, and it allows users to make use of the subsequent functions in MBNMAtime (plotting, prediction, ranking) whilst fitting these more complex models.

References

Friedrich JO, Adhikari NKJ, Beyene J (2011). “Ratio of means for analyzing continuous outcomes in meta-analysis performed as well as mean difference methods.” Journal of Clinical Epidemiology, 64(5), 556-564. doi:10.1016/j.jclinepi.2010.09.016.

Jansen JP, Vieira MC, Cope S (2015). “Network meta-analysis of longitudinal data using fractional polynomials.” Stat Med, 34(15), 2294-311. ISSN 1097-0258 (Electronic) 0277-6715 (Linking), doi:10.1002/sim.6492, https://pubmed.ncbi.nlm.nih.gov/25877808/.

Pedder H, Dias S, Bennetts M, Boucher M, Welton NJ (2019). “Modelling time-course relationships with multiple treatments: Model-Based Network Meta-Analysis for continuous summary outcomes.” Res Synth Methods, 10(2), 267-286.

Plummer M (2008). “Penalized loss functions for Bayesian model comparison.” Biostatistics, 9(3), 523-39. ISSN 1468-4357 (Electronic) 1465-4644 (Linking), https://pubmed.ncbi.nlm.nih.gov/18209015/.

Plummer M (2017). JAGS user manual. https://people.stat.sc.edu/hansont/stat740/jags_user_manual.pdf.

Spiegelhalter DJ, Best NG, Carlin BP, van der Linde A (2002). “Bayesian measures of model complexity and fit.” J R Statistic Soc B, 64(4), 583-639.

Examples

# Create mb.network object
network <- mb.network(osteopain)

# Fit a linear time-course MBNMA with:
# random relative treatment effects on the slope
mb.run(network, fun=tpoly(degree=1, pool.1="rel", method.1="random"))

# Fit an emax time-course MBNMA with:
# fixed relative treatment effects on emax
# a common parameter estimated independently of treatment
# a common Hill parameter estimated independently of treatment
# a prior for the Hill parameter (normal with mean 0 and precision 0.1)
# data reported as change from baseline
result <- mb.run(network, fun=temax(pool.emax="rel", method.emax="common",
                                    pool.et50="abs", method.et50="common",
                                    pool.hill="abs", method.hill="common"),
                 priors=list(hill="dunif(0.5, 2)"),
                 intercept=TRUE)


#### commented out to prevent errors from JAGS version in github actions build ####
# Fit a log-linear MBNMA with:
# random relative treatment effects on the rate
# an autoregressive AR1 covariance structure
# modelled as standardised mean differences
# copdnet <- mb.network(copd)
# result <- mb.run(copdnet, fun=tloglin(pool.rate="rel", method.rate="random"),
#                covar="AR1", rho="dunif(0,1)", link="smd")



####### Examine MCMC diagnostics (using mcmcplots package) #######

# Traceplots
# mcmcplots::traplot(result)

# Plots for assessing convergence
# mcmcplots::mcmcplot(result, c("rate", "sd.rate", "rho"))

########## Output ###########

# Print R2jags output and summary
print(result)
summary(result)

# Plot forest plot of results
plot(result)


###### Additional model arguments ######

# Use gout dataset
goutnet <- mb.network(goutSUA_CFBcomb)

# Define a user-defined time-course relationship for use in mb.run
timecourse <- ~ exp(beta.1 * time) + (time^beta.2)

# Run model with:
# user-defined time-course function
# random relative effects on beta.1
# default common effects on beta.2
# default relative pooling on beta.1 and beta.2
# common class effect on beta.2
mb.run(goutnet, fun=tuser(fun=timecourse, method.1="random"),
       class.effect=list(beta.1="common"))

# Fit a log-linear MBNMA
# with variance adjustment for correlation between time-points
result <- mb.run(network, fun=tloglin(),
                 rho="dunif(0,1)", covar="varadj")

Update MBNMA to obtain deviance contributions or fitted values

Description

Update MBNMA to obtain deviance contributions or fitted values

Usage

mb.update(
  mbnma,
  param = "theta",
  n.iter = mbnma$BUGSoutput$n.iter,
  n.thin = mbnma$BUGSoutput$n.thin
)

Arguments

mbnma

An S3 object of class "mbnma" generated by running a time-course MBNMA model

param

A character object that represents the parameter within the model to monitor when updating. Can currently only be used for monitoring fitted values and deviance contributions and so can take either "dev" (for deviance contributions), "resdev" (for residual deviance contributions) or "theta" (for fitted values).

n.iter

The number of iterations to update the model whilst monitoring additional parameters (if necessary). Must be a positive integer. Default is the value used in mbnma.

n.thin

The thinning rate. Must be a positive integer. Default is the value used in mbnma.

Value

A data frame containing posterior means for the specified param at each observation, arm and study.

Examples

# Using the alogliptin dataset
network <- mb.network(alog_pcfb)

# Run Emax model
emax <- mb.run(network, fun=temax())

# Update model for 500 iterations to monitor fitted values
mb.update(emax, param="theta", n.iter=500)

# Update model for 500 iterations to monitor residual deviance contributions
mb.update(emax, param="resdev", n.iter=500)

# Update model for 500 iterations to monitor deviance contributions
mb.update(emax, param="dev", n.iter=500)

Validates that a dataset fulfils requirements for MBNMA

Description

Validates that a dataset fulfils requirements for MBNMA

Usage

mb.validate.data(data.ab, single.arm = FALSE, CFB = TRUE)

Arguments

data.ab

A data frame of arm-level data in "long" format containing the columns:

  • studyID Study identifiers

  • time Numeric data indicating follow-up times

  • y Numeric data indicating the aggregate response for a given observation (e.g. mean)

  • se Numeric data indicating the standard error for a given observation

  • treatment Treatment identifiers (can be numeric, factor or character)

  • class An optional column indicating a particular class identifier. Observations with the same treatment identifier must also have the same class identifier.

  • n An optional column indicating the number of participants used to calculate the response at a given observation (required if modelling using Standardised Mean Differences)

  • standsd An optional column of numeric data indicating reference SDs used to standardise treatment effects when modelling using Standardised Mean Differences (SMD).

single.arm

A boolean object to indicate whether or not function should allow singe arm studies to be allowed in the network without returning an error. Default is not to allow their inclusion (single.arm=FALSE)

CFB

A boolean object to indicate if the dataset is composed of studies measuring change from baseline (TRUE) or not (FALSE). It is not essential to specify this correctly but failing to do so may lead to warnings.

Details

Checks done within the validation:

  • Checks data.ab has required column names

  • Checks there are no NAs

  • Checks that all SEs are positive

  • Checks that studies have baseline measurement (unless change from baseline data is being used)

  • Checks that arms are balanced at each time point

  • Checks that class codes are consistent within each treatment

  • Checks that treatment codes are consistent across different time points within a study

  • Checks that studies have at least two arms (if single.arm = FALSE)

  • Checks that standsd values are consistent within a study

Value

An error or warnings if checks are not passed. Runs silently if checks are passed


Write MBNMA time-course models JAGS code

Description

Writes JAGS code for a Bayesian time-course model for model-based network meta-analysis (MBNMA).

Usage

mb.write(
  fun = tpoly(degree = 1),
  link = "identity",
  positive.scale = TRUE,
  intercept = NULL,
  rho = 0,
  covar = "varadj",
  omega = NULL,
  corparam = TRUE,
  sdscale = FALSE,
  class.effect = list(),
  UME = FALSE
)

Arguments

fun

An object of class "timefun" generated (see Details) using any of tloglin(), tpoly(), titp(), temax(), tfpoly(), tspline() or tuser()

link

Can take either "identity" (the default), "log" (for modelling Ratios of Means (Friedrich et al. 2011)) or "smd" (for modelling Standardised Mean Differences - although this also corresponds to an identity link function).

positive.scale

A boolean object that indicates whether all continuous mean responses (y) are positive and therefore whether the baseline response should be given a prior that constrains it to be positive (e.g. for scales that cannot be <0).

intercept

A boolean object that indicates whether an intercept (written as alpha in the model) is to be included. If left as NULL (the default), an intercept will be included only for studies reporting absolute means, and will be excluded for studies reporting change from baseline (as indicated in network$cfb).

rho

The correlation coefficient when modelling within-study correlation between time points. The default is a string representing a prior distribution in JAGS, indicating that it be estimated from the data (e.g. rho="dunif(0,1)"). rho also be assigned a numeric value (e.g. rho=0.7), which fixes rho in the model to this value (e.g. for use in a deterministic sensitivity analysis). If set to rho=0 (the default) then this implies modelling no correlation between time points.

covar

A character specifying the covariance structure to use for modelling within-study correlation between time-points. This can be done by specifying one of the following:

  • "varadj" - a univariate likelihood with a variance adjustment to assume a constant correlation between subsequent time points (Jansen et al. 2015). This is the default.

  • "CS" - a multivariate normal likelihood with a compound symmetry structure

  • "AR1" - a multivariate normal likelihood with an autoregressive AR1 structure

omega

DEPRECATED IN VERSION 0.2.3 ONWARDS (~uniform(-1,1) now used for correlation between parameters rather than a Wishart prior). A scale matrix for the inverse-Wishart prior for the covariance matrix used to model the correlation between time-course parameters (see Details for time-course functions). omega must be a symmetric positive definite matrix with dimensions equal to the number of time-course parameters modelled using relative effects (pool="rel"). If left as NULL (the default) a diagonal matrix with elements equal to 1 is used.

corparam

A boolean object that indicates whether correlation should be modeled between relative effect time-course parameters. Default is FALSE and this is automatically set to FALSE if class effects are modeled. Setting it to TRUE models correlation between time-course parameters. This can help identify parameters that are estimated poorly for some treatments by allowing sharing of information between parameters for different treatments in the network, but may also cause some shrinkage.

sdscale

Logical object to indicate whether to write a model that specifies a reference SD for standardising when modelling using Standardised Mean Differences. Specifying sdscale=TRUE will therefore only modify the model if link function is set to SMD (link="smd").

class.effect

A list of named strings that determines which time-course parameters to model with a class effect and what that effect should be ("common" or "random"). For example: list(emax="common", et50="random").

UME

Can take either TRUE or FALSE (for an unrelated mean effects model on all or no time-course parameters respectively) or can be a vector of parameter name strings to model as UME. For example: c("beta.1", "beta.2").

Value

A single long character string containing the JAGS model generated based on the arguments passed to the function.

Examples

# Write a linear time-course MBNMA:
# random treatment effects on beta.1
# equal baselines in study arms
model <- mb.write(fun=tpoly(degree=1, pool.1="rel", method.1="random"))

# Write an emax time-course MBNMA with:
# a Hill parameter
# no intercept
model <- mb.write(fun=temax(pool.emax="rel", method.emax="common",
    pool.et50="abs", method.et50="common", pool.hill="abs", method.hill="common"),
  intercept=TRUE)

# Write a log-linear time-course MBNMA with:
# AR1 correlation between time points
model <- mb.write(fun=tloglin(),
  rho="dunif(0,1)", covar="AR1")

# Define a user-defined time-course relationship for the MBNMA JAGS model
userfun <- ~ (exp(beta.1 * time) / (beta.2 * time))
model <- mb.write(fun=tuser(fun=userfun,
    pool.1="rel", method.1="random",
    pool.2="rel", method.2="common"))

Run an NMA model

Description

Run an NMA model

Usage

nma.run(
  data.ab,
  treatments = NULL,
  method = "common",
  link = "identity",
  sdscale = FALSE,
  ...
)

Arguments

data.ab

A data frame of arm-level data in "long" format containing the columns:

  • studyID Study identifiers

  • time Numeric data indicating follow-up times

  • y Numeric data indicating the aggregate response for a given observation (e.g. mean)

  • se Numeric data indicating the standard error for a given observation

  • treatment Treatment identifiers (can be numeric, factor or character)

  • class An optional column indicating a particular class identifier. Observations with the same treatment identifier must also have the same class identifier.

  • n An optional column indicating the number of participants used to calculate the response at a given observation (required if modelling using Standardised Mean Differences)

  • standsd An optional column of numeric data indicating reference SDs used to standardise treatment effects when modelling using Standardised Mean Differences (SMD).

treatments

A vector of treatment names. If left as NULL it will use the treatment coding given in data.ab

method

Can take "common" or "random" to indicate the type of NMA model used to synthesise data points given in overlay.nma. The default is "random" since this assumes different time-points in overlay.nma have been lumped together to estimate the NMA.

link

Can take either "identity" (the default), "log" (for modelling Ratios of Means (Friedrich et al. 2011)) or "smd" (for modelling Standardised Mean Differences - although this also corresponds to an identity link function).

sdscale

Logical object to indicate whether to write a model that specifies a reference SD for standardising when modelling using Standardised Mean Differences. Specifying sdscale=TRUE will therefore only modify the model if link function is set to SMD (link="smd").

...

Options for plotting in igraph.

Value

Returns an object of class("nma", "rjags")

Examples

network <- mb.network(osteopain)

# Get the latest time point
late.time <- get.latest.time(network)

# Get the closest time point to a given value (t)
early.time <- get.closest.time(network, t=7)

# Run NMA on the data
nma.run(late.time$data.ab, treatments=late.time$treatments,
  method="random")

Studies of treatments for reducing body weight in patients with obesity

Description

A dataset from a systematic review of pharmacological treatments for reducing body weight in patients with obesity. The outcome is continuous, and aggregate data responses are given as mean change from baseline in body weight (KG). Overall there are 35 RCTs investigating 26 treatments of 16 agents (/combinations of agents) in the network. Standard deviations have been imputed for 421 observations.

Usage

obesityBW_CFB

Format

A data frame with 710 rows and 7 variables:

  • studyID Study identifiers

  • time Numeric data indicating follow-up times

  • y Numeric data indicating the mean response for a given observation

  • se Numeric data indicating the standard error for a given observation

  • n Numeric data indicating the number of participants used to calculate means for each observation

  • treatment Treatment identifiers as factors. Labels are shortened treatment names.

  • treatname Character data giving the full names of each treatment in the format agent_dose

  • agent Agent (drug) names stored as characters

  • class The drug class of the agent (a broader category than agent) stored as characters

Details

obesityBW_CFB is a data frame in long format (one row per observation, arm and study), with the variables studyID, time, y, se, n, treatment, treatname, agent and class.

Source

Pfizer Ltd.


Studies of pain relief medications for osteoarthritis

Description

A dataset containing results on the WOMAC pain scale (0-10) over time for studies investigating 29 treatments for pain relief in patients with osteoarthritis. Standard deviations have been imputed for 269 observations.

Usage

osteopain

Format

A data frame with 417 rows and 7 variables:

  • studyID Study identifiers

  • time Numeric data indicating follow-up times

  • y Numeric data indicating the mean response for a given observation

  • se Numeric data indicating the standard error for a given observation

  • treatment Treatment identifiers as factors. Labels are shortened treatment names.

  • arm Arm identifiers coded for each study

  • treatname Character data giving the full names of each treatment

Details

osteopain is a data frame in long format (one row per observation, arm and study), with the variables studyID, time, y, se, treatment, arm and treatname.

Source

Pfizer Ltd.


Calculate plugin pD from a JAGS model with univariate likelihood for studies with repeated measurements

Description

Uses results from MBNMA JAGS models to calculate pD via the plugin method (Spiegelhalter et al. 2002). Can only be used for models with known standard errors or covariance matrices (typically univariate).

Usage

pDcalc(
  obs1,
  obs2,
  fups = NULL,
  narm,
  NS,
  theta.result,
  resdev.result,
  likelihood = "normal",
  type = "time"
)

Arguments

obs1

A matrix (study x arm) or array (study x arm x time point) containing observed data for y (normal likelihood) or r (binomial or Poisson likelihood) in each arm of each study. This will be the same array used as data for the JAGS model.

obs2

A matrix (study x arm) or array (study x arm x time point) containing observed data for se (normal likelihood), N (binomial likelihood) or E (Poisson likelihood) in each arm of each study. This will be the same array used as data for the JAGS model.

fups

A numeric vector of length equal to the number of studies, containing the number of follow-up mean responses reported in each study. Required for time-course MBNMA models (if type="time")

narm

A numeric vector of length equal to the number of studies, containing the number of arms in each study.

NS

A single number equal to the number of studies in the dataset.

theta.result

A matrix (study x arm) or array (study x arm x time point) containing the posterior mean predicted means/probabilities/rate in each arm of each study. This will be estimated by the JAGS model.

resdev.result

A matrix (study x arm) or array (study x arm x time point) containing the posterior mean residual deviance contributions in each arm of each study. This will be estimated by the JAGS model.

likelihood

A character object of any of the following likelihoods:

  • univariate

  • binomial (does not work with time-course MBNMA models)

  • multivar.normal (does not work with time-course MBNMA models)

type

The type of MBNMA model fitted. Can be either "time" or "dose"

Details

Method for calculating pD via the plugin method proposed by (Spiegelhalter et al. 2002). Standard errors / covariance matrices must be assumed to be known. To obtain values for theta.result and resdev.result these parameters must be monitored when running the JAGS model.

For non-linear time-course MBNMA models residual deviance contributions may be skewed, which can lead to non-sensical results when calculating pD via the plugin method. Alternative approaches are to use pV (pv) as an approximation (Plummer 2008) or pD calculated by Kullback–Leibler divergence (pd.kl) or using an optimism adjustment (popt) (Plummer 2008).

Value

A numeric value for the effective number of parameters, pD, calculated via the plugin method

References

TO ADD pV REF

Examples

# Using the alogliptin dataset
network <- mb.network(alog_pcfb)

# Run Emax model saving predicted means and residual deviance contributions
emax <- mb.run(network, fun=temax(),
  parameters.to.save=c("theta", "resdev"), intercept=FALSE)

# Get matrices of observed data
jagsdat <- getjagsdata(network$data.ab)

# Plugin estimation of pD is problematic with non-linear models as it often leads to
#negative values, hence use of pV, pd.kl and popt as other measures for the effective
#number of parameters
pDcalc(obs1=jagsdat$y, obs2=jagsdat$se,
  fups=jagsdat$fups, narm=jagsdat$narm, NS=jagsdat$NS,
  theta.result = emax$BUGSoutput$mean$theta,
  resdev.result = emax$BUGSoutput$mean$resdev
  )

Create an mb.network object

Description

Creates an object of class("mb.network"). Various MBNMA functions can subsequently be applied to this object.

Usage

## S3 method for class 'mb.network'
plot(
  x,
  edge.scale = 1,
  label.distance = 0,
  level = "treatment",
  remove.loops = FALSE,
  v.color = "connect",
  v.scale = NULL,
  layout = igraph::in_circle(),
  legend = TRUE,
  legend.x = "bottomleft",
  legend.y = NULL,
  ...
)

mb.network(data.ab, reference = 1, cfb = NULL, description = "Network")

Arguments

x

An object of class mb.network.

edge.scale

A number to scale the thickness of connecting lines (edges). Line thickness is proportional to the number of studies for a given comparison. Set to 0 to make thickness equal for all comparisons.

label.distance

A number scaling the distance of labels from the nodes to improve readability. The labels will be directly on top of the nodes if the default of 0 is used. Option only applicable if layout_in_circle is set to TRUE.

level

A string indicating whether nodes/facets should represent treatment or class in the plot. Can be used to examine the expected impact of modelling class/agent effects.

remove.loops

A boolean value indicating whether to include loops that indicate comparisons within a node.

v.color

Can take either "connect" (the default) to indicate that nodes should only be coloured if they are connected to the network reference treatment (indicates network connectivity) or "class" to colour nodes by class (this requires that the variable class be included in the dataset).

v.scale

A number with which to scale the size of the nodes. If the variable N (to indicate the numbers of participants at each observation) is included in the dataset then the size of the nodes will be proportional to the number of participants within a treatment/class in the network at the earliest time point reported in each study.

layout

An igraph layout specification. This is a function specifying an igraph layout that determines the arrangement of the vertices (nodes). The default igraph::as_circle() arranged vertices in a circle. Two other useful layouts for network plots are: igraph::as_star(), igraph::with_fr(). Others can be found in layout_

legend

A boolean value indicating whether or not to plot a legend with class names if v.color="class"

legend.x

Can be either a string or a numerical x-coordinate indicating where the legend should be plotted (see legend).

legend.y

A numerical y-coordinate indicating where the legend should be plotted - only required if legend.x is also a numeric co-ordinate.

...

Options for plotting in igraph.

data.ab

A data frame of arm-level data in "long" format containing the columns:

  • studyID Study identifiers

  • time Numeric data indicating follow-up times

  • y Numeric data indicating the aggregate response for a given observation (e.g. mean)

  • se Numeric data indicating the standard error for a given observation

  • treatment Treatment identifiers (can be numeric, factor or character)

  • class An optional column indicating a particular class identifier. Observations with the same treatment identifier must also have the same class identifier.

  • n An optional column indicating the number of participants used to calculate the response at a given observation (required if modelling using Standardised Mean Differences)

  • standsd An optional column of numeric data indicating reference SDs used to standardise treatment effects when modelling using Standardised Mean Differences (SMD).

reference

A number or character (depending on the format of treatment within data.ab) indicating the reference treatment in the network (i.e. those for which estimated relative treatment effects estimated by the model will be compared to).

cfb

A logical vector whose length is equal to the unique number of studies in data.ab, where each element is TRUE if the study data reported is change-from-baseline and FALSE otherwise. If left as NULL (the default) then this will be identified from the data by assuming any study for which there is no data at time=0 reports change-from-baseline.

description

Optional. Short description of the network.

Details

The S3 method plot() on an mb.network object generates a network plot that shows how different treatments are connected within the network via study comparisons. This can be used to identify how direct and indirect evidence are informing different treatment comparisons. Depends on igraph.

Missing values (NA) cannot be included in the dataset. Studies must have a baseline measurement and more than a single follow-up time (unless change from baseline data are being used). Data must be present for all arms within a study at each follow-up time.

Value

Returns an object of class "igraph", which can be modified by other functions within the igraph package.

An object of class("mb.network") which is a list containing:

  • description A short description of the network

  • data.ab A data frame containing the arm-level network data (treatment identifiers will have been recoded to a sequential numeric code)

  • studyID A character vector with the IDs of included studies.

  • cfb A logical vector indicating which studies report change from baseline data

  • treatments A character vector indicating the treatment identifiers that correspond to the new treatment codes.

  • classes A character vector indicating the class identifiers (if included in the original data) that correspond to the new class codes.

Methods (by generic)

  • plot(mb.network): Generate a network plot

Examples

# Create an mb.network object from the data
network <- mb.network(osteopain)

# Arrange network plot in a star with the reference treatment in the centre
plot(network, layout=igraph::as_star())

# Generate a network plot at the class level that removes loops indicating comparisons
#within a node
goutnet <- mb.network(goutSUA_CFB)
plot(goutnet, level="class", remove.loops=TRUE)

# Generate a network plot at the treatment level that colours nodes by class
plot(goutnet, v.color="class", remove.loops=TRUE)

# Plot network in which node size is proportional to number of participants
alognet <- mb.network(alog_pcfb)
plot(alognet, v.scale=2)



# Using the osteoarthritis dataset
print(osteopain)

# Define network
network <- mb.network(osteopain, description="Osteoarthritis Dataset")

# Define network with different network reference treatment
network <- mb.network(osteopain, reference="Ce_200")


# Using the alogliptin dataset
network <- mb.network(alog_pcfb, description="Alogliptin Dataset")

# Examine networks
print(network)
plot(network)

Plots predicted responses from a time-course MBNMA model

Description

Plots predicted responses from a time-course MBNMA model

Usage

## S3 method for class 'mb.predict'
plot(
  x,
  disp.obs = FALSE,
  overlay.ref = TRUE,
  overlay.nma = NULL,
  method = "random",
  col = "blue",
  max.col.scale = NULL,
  treat.labs = NULL,
  plot.bins = TRUE,
  ...
)

Arguments

x

An object of class "mb.predict" generated by predict("mbnma")

disp.obs

A boolean object to indicate whether to show shaded sections of the plot for where there is observed data (TRUE) or not (FALSE)

overlay.ref

A boolean object indicating whether to overlay a line showing the median network reference treatment response over time on the plot (TRUE) or not (FALSE). The network reference treatment (treatment

  1. must be included in predict for this to display the network reference treatment properly.

overlay.nma

Numeric vector used to overlay the results from a standard NMA model that "lumps" time-points together within the time bin ranges specified in overlay.nma. The numbers in overlay.nma define the boundaries of the time bins within which to perform a standard NMA. Length must be >=2, or can be left as NULL (the default) to indicate that no NMA should be perfomed. overlay.nma can only be specified if overlay.ref==TRUE. See Details for further information.

method

Can take "common" or "random" to indicate the type of NMA model used to synthesise data points given in overlay.nma. The default is "random" since this assumes different time-points in overlay.nma have been lumped together to estimate the NMA.

col

A character indicating the colour to use for shading if disp.obs is set to TRUE. Can be either "blue", "green", or "red"

max.col.scale

Rarely requires adjustment. The maximum count of observations (therefore the darkest shaded color) only used if disp.obs is used. This allows consistency of shading between multiple plotted graphs. It should always be at least as high as the maximum count of observations plotted

treat.labs

A vector of treatment labels in the same order as treatment codes. Easiest to use treatment labels stored by mb.network()

plot.bins

Plot time bin boundaries as vertical dashed lines. Setting plot.bins=TRUE if overlay.nma is specified also sets x-axis ticks to time bin boundaries automatically.

...

Arguments for ggplot() or R2jags()

Details

For the S3 method plot(), if disp.obs is set to TRUE it is advisable to ensure predictions in predict are estimated using an even sequence of time points to avoid misrepresentation of shaded densities. Shaded counts of observations will be relative to the treatment plotted in each panel rather than to the network reference treatment if disp.obs is set to TRUE.

Overlaying NMA results

overlay.nma indicates regions of the data (defined as "time bins") over which it may be reasonable to "lump" different follow-up times from different studies together and assume a standard NMA model. For example:

  • overlay.nma=c(5,10) indicates a single NMA of studies with follow-up times ⁠>5⁠ and ⁠<=10⁠

  • overlay.nma=c(5,10,15) indicates two NMAs should be performed of studies with follow-up times ⁠>5⁠ and ⁠<=10⁠ of studies with follow-up times ⁠>10⁠ and ⁠<=15⁠

When used with MBNMA (via predict.mbnma()) this allows comparison to MBNMA results over a specific range of time within each time bin. It can be useful to assess which time-course function might be suitable when using binplot(), or to to assess if the MBNMA predictions are in agreement with predictions from an NMA model when using plot.mb.predict() for a specific range of time-points. This can be a general indicator of the fit of the time-course model.

However, it is important to note that the wider the range specified in overlay.nma, the more likely it is that different time-points are included, and therefore that there is greater heterogeneity/inconsistency in the NMA model. If overlay.nma includes several follow-up times for any study then only a single time-point will be taken (the one closest to mean(overlay.nma)). The NMA predictions are plotted over the range specified in overlay.nma as a horizontal line, with the 95%CrI shown by a grey rectangle. The NMA predictions represent those for any time-points within this range since they lump together data at all these time-points. Predictions for treatments that are disconnected from the network reference treatment at data points specified within overlay.nma cannot be estimated so are not included.

It is important to note that the NMA model is not necessarily the "correct" model, since it "lumps" different time-points together and ignores potential differences in treatment effects that may arise from this. The wider the range specified in overlay.nma, the greater the effect of "lumping" and the stronger the assumption of similarity between studies.

For an NMA model to be estimated and a corresponding prediction to be made from it, each time bin must include the network reference treatment (treatment=1) evaluated in at least 1 connected study in the time bin. If a given time bin does not meet this criteria then an NMA will not be calculated for it.

Examples

# Create an mb.network object from a dataset
copdnet <- mb.network(copd)

# Run an MBNMA model with a log-linear time-course
loglin <- mb.run(copdnet,
  fun=tloglin(pool.rate="rel", method.rate="common"),
  rho="dunif(0,1)", covar="varadj")

# Predict responses using the original dataset to estimate the network reference
#treatment response
df.ref <- copd[copd$treatment=="Placebo",]
predict <- predict(loglin, times=c(0:20), E0=0, ref.resp=df.ref)

# Plot the predicted responses with observations displayed on plot as green shading
plot(predict, disp.obs=TRUE, overlay.ref=FALSE, col="green")

# Plot the predicted responses with the median network reference treatment response overlayed
#on the plot
plot(predict, disp.obs=FALSE, overlay.ref=TRUE)

# Plot predictions from NMAs calculated between different time-points
plot(predict, overlay.nma=c(5,10), n.iter=20000)
plot(predict, overlay.nma=c(5,10,15,20), n.iter=20000)
# Time-course fit may be less well at 15-20 weeks follow-up

Plot histograms of rankings from MBNMA models

Description

Plot histograms of rankings from MBNMA models

Usage

## S3 method for class 'mb.rank'
plot(x, treat.labs = NULL, ...)

Arguments

x

An object of class "mb.rank" generated by rank.mbnma()

treat.labs

A vector of treatment labels in the same order as treatment codes. Easiest to use treatment labels stored by mb.network()

...

Arguments to be sent to ggplot2::ggplot()

Value

A histogram that shows rankings for each treatment/agent/prediction. The object returned is an object of class c("gg", "ggplot").

Examples

# Create an mb.network object from a dataset
painnet <- mb.network(osteopain)

# Run an MBNMA model with an Emax time-course
emax <- mb.run(painnet,
  fun=temax(pool.emax="rel", method.emax="common",
    pool.et50="abs", method.et50="random"),
  positive.scale=TRUE)

# Calculate treatment rankings for AUC and emax
ranks <- rank(emax,
  param=c("auc"),
  int.range=c(0,15), n.iter=500)

# Plot histograms for ranking by AUC
plot(ranks)

Forest plot for results from time-course MBNMA models

Description

Generates a forest plot for time-course parameters of interest from results from time-course MBNMA models. Posterior densities are plotted above each result using ggdist:stat_:halfeye()

Usage

## S3 method for class 'mbnma'
plot(x, params = NULL, treat.labs = NULL, class.labs = NULL, ...)

Arguments

x

An S3 object of class "mbnma" generated by running a time-course MBNMA model

params

A character vector of time-course parameters to plot. Parameters must be given the same name as monitored nodes in mbnma and must vary by treatment or class. Can be set to NULL to include all available time-course parameters estimated by mbnma.

treat.labs

A character vector of treatment labels. If left as NULL (the default) then labels will be used as defined in the data.

class.labs

A character vector of class labels if mbnma was modelled using class effects If left as NULL (the default) then labels will be used as defined in the data.

...

Arguments to be sent to ggdist::stat_halfeye()

Value

A forest plot of class c("gg", "ggplot") that has separate panels for different time-course parameters

Examples

# Create an mb.network object from a dataset
alognet <- mb.network(alog_pcfb)

# Run an MBNMA model with an Emax time-course
emax <- mb.run(alognet,
  fun=temax(pool.emax="rel", method.emax="common",
    pool.et50="rel", method.et50="common"),
  intercept=FALSE)

# Generate forest plot
plot(emax)

# Plot results for only one time-course parameter
plot(emax, params="emax")

Perform node-splitting on a MBNMA time-course network

Description

Within a MBNMA time-course network, split contributions into direct and indirect evidence and test for consistency between them. Closed loops of treatments in which it is possible to test for consistency are those in which direct and indirect evidence are available from independent sources van Valkenhoef van Valkenhoef et al. (2016).

Usage

## S3 method for class 'nodesplit'
plot(x, plot.type = NULL, params = NULL, ...)

mb.nodesplit(
  network,
  comparisons = mb.nodesplit.comparisons(network),
  nodesplit.parameters = "all",
  fun = tpoly(degree = 1),
  times = NULL,
  lim = "cred",
  ...
)

Arguments

x

An object of class("nodesplit")

plot.type

A character string that can take the value of "forest" to plot only forest plots, "density" to plot only density plots, or left as NULL (the default) to plot both types of plot.

params

A character vector corresponding to a time-course parameter(s) for which to plot results. If left as NULL (the default), nodes-split results for all time-course parameters will be plotted.

...

Arguments to be sent to mb.run()

network

An object of class "mb.network".

comparisons

A data frame specifying the comparisons to be split (one row per comparison). The frame has two columns indicating each treatment for each comparison: t1 and t2.

nodesplit.parameters

A character vector of named time-course parameters on which to node-split (e.g. c("beta.1", "beta.2")). Can use "all" to split on all time-course parameters.

fun

An object of class "timefun" generated (see Details) using any of tloglin(), tpoly(), titp(), temax(), tfpoly(), tspline() or tuser()

times

A sequence of positive numbers indicating which time points to predict mean responses for (or at which to conduct a node-split if used with mb.nodesplit())

lim

Specifies calculation of either 95% credible intervals (lim="cred") or 95% prediction intervals (lim="pred").

Details

The S3 method plot() on an mb.nodesplit object generates either forest plots of posterior medians and 95\% credible intervals, or density plots of posterior densities for direct and indirect evidence.

Note that by specifying the times argument a user can perform a node-split of treatment effects at a specific time-point. This will give the treatment effect for both direct, indirect, and MBNMA estimates at this time point.

Value

Plots the desired graph(s) and returns an object (or list of objects if plot.type=NULL) of class(c("gg", "ggplot")), which can be edited using ggplot commands.

A an object of class("mb.nodesplit") that is a list containing elements d.X.Y (treatment 1 = X, treatment 2 = Y). Each element (corresponding to each comparison) contains additional numbered elements corresponding to each parameter in the time-course function on which node splitting was performed. These elements then contain:

  • ⁠overlap matrix⁠ MCMC results for the difference between direct and indirect evidence

  • p.values Bayesian p-value for the test of consistency between direct and indirect evidence

  • quantiles

  • forest.plot

  • density.plot

  • direct MCMC results for the direct evidence

  • indirect MCMC results for the indirect evidence

Functions

  • plot(nodesplit): Plot outputs from nodesplit models

References

van Valkenhoef G, Dias S, Ades AE, Welton NJ (2016). “Automated generation of node-splitting models for assessment of inconsistency in network meta-analysis.” Res Synth Methods, 7(1), 80-93. ISSN 1759-2887 (Electronic) 1759-2879 (Linking), doi:10.1002/jrsm.1167, https://pubmed.ncbi.nlm.nih.gov/26461181/.

Examples

# Create mb.network object
painnet <- mb.network(osteopain)

# Identify comparisons informed by direct and indirect evidence
splits <- mb.nodesplit.comparisons(painnet)

# Fit a log-linear time-course MBNMA (takes a while to run)
result <- mb.nodesplit(painnet, comparisons=splits, nodesplit.parameters="all",
  fun=tloglin(pool.rate="rel", method.rate="common"),
  rho="dunif(0,1)", covar="varadj"
  )

# Fit an emax time-course MBNMA with a node-split on emax parameters only
result <- mb.nodesplit(painnet, comparisons=splits, nodesplit.parameters="emax",
  fun=temax(pool.emax="rel", method.emax="common",
    pool.et50="rel", method.et50="common"))

# Inspect results
print(result)
summary(result)

# Plot results
plot(result)

Predict effects over time in a given population based on MBNMA time-course models

Description

Used to predict effects over time for different treatments or to predict the results of a new study. For MBNMA models that include consistency relative effects on time-course parameters, this is calculated by combining relative treatment effects with a given reference treatment response (specific to the population of interest).

Usage

## S3 method for class 'mbnma'
predict(
  object,
  times = seq(0, max(object$model.arg$jagsdata$time, na.rm = TRUE), length.out = 30),
  E0 = 0,
  treats = NULL,
  level = "treatment",
  ref.resp = NULL,
  synth = "common",
  lim = "cred",
  ...
)

Arguments

object

An S3 object of class("mbnma") generated by running a time-course MBNMA model

times

A sequence of positive numbers indicating which time points to predict mean responses for (or at which to conduct a node-split if used with mb.nodesplit())

E0

An object to indicate the value(s) to use for the response at time = 0 in the prediction. This can take a number of different formats depending on how it will be used/calculated. The default is 0 but this may lead to non-sensical predictions if Ratio of Means are modeled.

  • numeric() A single numeric value representing the deterministic response at time = 0

  • formula() A formula representing a stochastic distribution for the response at time = 0. This is specified as a random number generator (RNG) given as a string, and can take any RNG distribution for which a function exists in R. For example: ~rnorm(n, 7, 0.5).

treats

A character vector of treatment/class names or a numeric vector of treatment/class codes (as coded in mbnma) that indicates which treatments/classes to calculate predictions for. If left as NULL then predictions will be calculated for all treatments/classes. Whether the vector should correspond to treatments or classes depends on the value of level.

level

Can take either "treatment" to make predictions for treatments, or "class" to make predictions for classes (in which case object must be a class effect model).

ref.resp

An object to indicate the value(s) to use for the reference treatment response in MBNMA models in which the reference treatment response is not estimated within the model (i.e. those that model any time- course parameters using pool="rel"). This can take a number of different formats depending on how it will be used/calculated. There are two approaches for this:

  1. The reference response can be estimated from a dataset of studies investigating the reference treatment using meta-analysis. This dataset could be a set of observational studies that are specific to the population on which to make predictions, or it could be a subset of the study arms within the MBNMA dataset that investigate the reference treatment. The data should be provided to ref.resp as a data.frame() containing the data in long format (one row per observation). See ref.synth()

  2. Values for the reference treatment response can be assigned to different time-course parameters within the model that have been modelled using consistency relative effects (pool="rel"). These are given as a list, in which each named element corresponds to a time-course parameter modelled in mbnma, specified on the corresponding scale (i.e. specified on the log scale if modelled on the log scale using Ratios of Means). Their values can be either of the following:

  • numeric() A numeric value representing the deterministic value of the time-course parameter in question in individuals given the reference treatment. 0 is used as the default, which assumes no effect of time on the reference treatment (i.e. mean differences / relative effects versus the reference treatment are modeled).

  • formula() A formula representing a stochastic distribution for the value of the time-course parameter in question. This is specified as a random number generator (RNG) given as a formula, and can take any RNG distribution for which a function exists in R. For example: ~rnorm(n, -3, 0.2).

synth

A character object that can take the value "common" or "random" that specifies the the type of pooling to use for synthesis of ref.resp. Using "random" rather than "common" for synth will result in wider 95\% CrI for predictions.

lim

Specifies calculation of either 95% credible intervals (lim="cred") or 95% prediction intervals (lim="pred").

...

Arguments to be sent to R2jags for synthesis of the network reference treatment effect (using ref.synth())

Details

By default the network reference treatment baseline (E0) and time-course parameter values are set to zero so that predict() estimates mean differences (/relative treatment effects) over time versus the network reference treatment.

ref.resp only needs to be specified if mbnma has been estimated using consistency relative effects (pool="rel") for any time-course parameters, as these inform the absolute values of the network reference treatment parameters which can then be added to the relative effects to calculate specific predictions.

Value

An S3 object of class mb.predict that contains the following elements:

  • summary A named list of data frames. Each data frame contains a summary of predicted responses at follow-up times specified in times for each treatment specified in treats

  • pred.mat A named list of matrices. Each matrix contains the MCMC results of predicted responses at follow-up times specified in times for each treatment specified in treats

Examples

# Create an mb.network object from a dataset
network <- mb.network(osteopain)

# Run an MBNMA model with an Emax time-course
emax <- mb.run(network,
  fun=temax(pool.emax="rel", method.emax="common",
    pool.et50="abs", method.et50="common"))

# Predict responses using a stochastic baseline (E0) and a distribution for the
#network reference treatment
preds <- predict(emax, times=c(0:10),
  E0=~rnorm(n, 7, 0.5),
  ref.resp=list(emax=~rnorm(n, -0.5, 0.05)))
summary(preds)

# Predict responses using the original dataset to estimate the network reference
#treatment response
paindata.ref <- osteopain[osteopain$treatname=="Placebo_0",]
preds <- predict(emax, times=c(5:15),
  E0=10,
  ref.resp=paindata.ref)
summary(preds)

# Repeat the above prediction but using a random effects meta-analysis of the
#network reference treatment response
preds <- predict(emax, times=c(5:15),
  E0=10,
  ref.resp=paindata.ref,
  synth="random")
summary(preds)

Print mb.network information to the console

Description

Print mb.network information to the console

Usage

## S3 method for class 'mb.network'
print(x, ...)

Arguments

x

An object of class mb.network.

...

further arguments passed to or from other methods

Value

Prints the contents of x to the console.


Print summary information from an mb.predict object

Description

Print summary information from an mb.predict object

Usage

## S3 method for class 'mb.predict'
print(x, ...)

Arguments

x

An object of class("mb.predict") generated by predict.mbnma()

...

further arguments passed to or from other methods


Prints a summary of rankings for each parameter

Description

Prints a summary of rankings for each parameter

Usage

## S3 method for class 'mb.rank'
print(x, ...)

Arguments

x

An object of class "mb.rank" generated by rank.mbnma()

...

further arguments passed to or from other methods

Value

Prints summary details of treatment rankings to the console


Prints basic results from a node-split to the console

Description

Prints basic results from a node-split to the console

Usage

## S3 method for class 'nodesplit'
print(x, groupby = "time.param", ...)

Arguments

x

An object of class "nodesplit" generated by mb.nodeplit()

groupby

A character object that can take the value "time.param" to present results grouped by time-course parameter (the default) or "comparison" to present results grouped by treatment comparison.

...

arguments to be sent to knitr::kable()

Value

Prints summary details of nodesplit results to the console


Print posterior medians (95% credible intervals) for table of relative effects/mean differences between treatments/classes

Description

Print posterior medians (95% credible intervals) for table of relative effects/mean differences between treatments/classes

Usage

## S3 method for class 'relative.array'
print(x, digits = 2, ...)

Arguments

x

An object of class "relative.array" generated by get.relative()

digits

An integer indicating the number of significant digits to be used.

...

further arguments passed to knitr::kable

Value

Prints a league table of treatment effects to the console


Calculate position of label with respect to vertex location within a circle

Description

Useful for graphs drawn using igraph to reposition labels relative to vertices when vertices are laid out in a circle (as is common in network plots). igraph interprets position within vertex.label.degree as radians, so it is necessary to convert locations into radian values. This is the main role of this function.

Usage

radian.rescale(x, start = 0, direction = 1)

Arguments

x

A numeric vector of positions around a circle, typically sequentially numbered.

start

A number giving the offset from 12 o'clock in radians for the label locations.

direction

Either 1 for clockwise numbering (based on the order of x) or -1 for anti-clockwise.

Value

A numeric vector of rescaled values

References

https://gist.github.com/kjhealy/834774/a4e677401fd6e4c319135dabeaf9894393f9392c


Set rank as a method

Description

Set rank as a method

Usage

rank(x, ...)

Arguments

x

An object on which to apply the rank method

...

Arguments to be passed to methods

Value

Uses the rank method


Rank predictions at a specific time point

Description

Rank predictions at a specific time point

Usage

## S3 method for class 'mb.predict'
rank(
  x,
  time = max(x$summary[[1]]$time),
  lower_better = FALSE,
  treats = names(x$summary),
  ...
)

Arguments

x

an object of class("mb.predict") that contains predictions from an MBNMA model

time

a number indicating the time point at which predictions should be ranked. It must be one of the time points for which predictions in x are available.

lower_better

Indicates whether negative responses are better (lower_better=TRUE) or positive responses are better (lower_better=FALSE)

treats

A character vector of treatment/class names for which responses have been predicted in x As default, rankings will be calculated for all treatments/classes in x.

...

Arguments to be passed to methods

Value

Returns an object of class("mb.rank") containing ranked predictions

Examples

# Create an mb.network object from a dataset
network <- mb.network(osteopain)

# Run an MBNMA model with an Emax time-course
emax <- mb.run(network,
  fun=temax(pool.emax="rel", method.emax="common",
    pool.et50="abs", method.et50="common"))

# Predict responses using a stochastic baseline (E0) and a distribution for the
#network reference treatment
preds <- predict(emax, E0=7,
  ref.resp=list(emax=~rnorm(n, -0.5, 0.05)))

# Rank predictions at latest predicted time-point
rank(preds, lower_better=TRUE)


#### Rank predictions at 5 weeks follow-up ####

# First ensure responses are predicted at 5 weeks
preds <- predict(emax, E0=7,
  ref.resp=list(emax=~rnorm(n, -0.5, 0.05)),
  times=c(0,5,10))

# Rank predictions at 5 weeks follow-up
ranks <- rank(preds, lower_better=TRUE, time=5)

# Plot ranks
plot(ranks)

Rank parameters from a time-course MBNMA

Description

Ranks desired parameters saved from a time-course MBNMA model from "best" to "worst".

Usage

## S3 method for class 'mbnma'
rank(
  x,
  param = "auc",
  lower_better = FALSE,
  treats = NULL,
  int.range = NULL,
  n.iter = x$BUGSoutput$n.sims,
  ...
)

Arguments

x

An S3 object of class("mbnma") generated by running a time-course MBNMA model

param

A character object containing any model parameter monitored in mbnma for which ranking is desired (e.g. "beta.1", "emax"). Parameters must vary by treatment for ranking to be possible. Can also be specified as "auc" (see details).

lower_better

Indicates whether negative responses are better (lower_better=TRUE) or positive responses are better (lower_better=FALSE)

treats

A character vector of treatment/class names (depending on the parameter to be ranked) or a numeric vector of treatment/class codes (as coded in mbnma) that indicate which treatments/classes to calculate rankings for. If left 'NULL“ then rankings will be calculated for all treatments/classes.

int.range

A numeric vector with two elements that indicates the range over which to calculate AUC. Takes the form c(lower bound, upper bound). If left as NULL (the default) then the range will be between zero and the maximum follow-up time in the dataset.

n.iter

The number of iterations for which to calculate AUC (if "auc" is included in params). Must be a positive integer. Default is the value used in mbnma.

...

Arguments to be sent to integrate()

Details

"auc" can be specified in param to rank treatments based on Area Under the Curve (AUC). This accounts for the effect of multiple time-course parameters simultaneously on the treatment response, but will be impacted by the range of time over which AUC is calculated (int.range). This requires integration over int.range and can take some time to run (particularly) for spline functions as this uses the trapezoid method rather than adaptive quadrature). Note that "auc" can only be calculated at the treatment-level in class effect models.

As with other post-estimation functions, rank() should only be performed on models which have successfully converged. Note that rankings can be very sensitive to even small changes in treatment effects and therefore failure to converge in only one parameter may have substantial impact on rankings.

Value

A named list whose elements include:

  • summary.rank A data frame containing mean, sd, and quantiles for the ranks of each treatment given in treats

  • prob.matrix A matrix of the proportions of MCMC results for which each treatment/class in treats ranked in which position for the given parameter

  • rank.matrix A matrix of the ranks of MCMC results for each treatment/class in treats for the given parameter.

Examples

# Create an mb.network object from a dataset
network <- mb.network(alog_pcfb)

# Run an MBNMA model with an Emax time-course
emax <- mb.run(network,
  fun=temax(pool.emax="rel", method.emax="common",
            pool.et50="rel", method.et50="random"),
  intercept=FALSE)

# Rank treatments by time-course parameter from the model with lower scores being better
rank(emax, param=c("emax"), lower_better=TRUE)

# Rank treatments 1-3 by AUC
rank(emax, param="auc", treats=c(1:3), lower_better=TRUE,
  int.range=c(0,20))

Calculates ranking probabilities for AUC from a time-course MBNMA

Description

Calculates ranking probabilities for AUC from a time-course MBNMA

Usage

rankauc(
  mbnma,
  lower_better = FALSE,
  treats = NULL,
  level = "treatments",
  int.range = c(0, max(mbnma$network$data.ab$time)),
  n.iter = mbnma$BUGSoutput$n.sims,
  subdivisions = 100,
  ...
)

Arguments

mbnma

An S3 object of class "mbnma" generated by running a time-course MBNMA model

lower_better

Indicates whether negative responses are better (lower_better=TRUE) or positive responses are better (lower_better=FALSE)

treats

A character vector of treatment/class names (depending on the value of level). If left ⁠NULL`` then rankings will be calculated for all treatments/classes. Note that unlike ⁠rank.mbnma()' this argument cannot take a numeric vector.

level

Can take either "treatment" to make predictions for treatments, or "class" to make predictions for classes (in which case object must be a class effect model).

int.range

A numeric vector with two elements that indicates the range over which to calculate AUC. Takes the form c(lower bound, upper bound). If left as NULL (the default) then the range will be between zero and the maximum follow-up time in the dataset.

n.iter

The number of iterations for which to calculate AUC (if "auc" is included in params). Must be a positive integer. Default is the value used in mbnma.

subdivisions

The number of subdivisions over which to integrate (see integrate)

...

Arguments to be sent to R2jags for synthesis of the network reference treatment effect (using ref.synth())

Details

"auc" can be specified in param to rank treatments based on Area Under the Curve (AUC). This accounts for the effect of multiple time-course parameters simultaneously on the treatment response, but will be impacted by the range of time over which AUC is calculated (int.range). This requires integration over int.range and can take some time to run (particularly) for spline functions as this uses the trapezoid method rather than adaptive quadrature). Note that "auc" can only be calculated at the treatment-level in class effect models.

As with other post-estimation functions, rank() should only be performed on models which have successfully converged. Note that rankings can be very sensitive to even small changes in treatment effects and therefore failure to converge in only one parameter may have substantial impact on rankings.

Value

A named list whose elements include:

  • summary.rank A data frame containing mean, sd, and quantiles for the ranks of each treatment given in treats

  • prob.matrix A matrix of the proportions of MCMC results for which each treatment/class in treats ranked in which position for the given parameter

  • rank.matrix A matrix of the ranks of MCMC results for each treatment/class in treats for the given parameter.


Identify unique comparisons relative to study reference treatment within a network

Description

Identify unique contrasts relative to each study reference within a network. Repetitions of the same treatment comparison are grouped together.

Usage

ref.comparisons(data)

Arguments

data

A data frame containing variables studyID and treatment (as numeric codes) that indicate which treatments are used in which studies.

Value

A data frame of unique comparisons in which each row represents a different comparison. t1 and t2 indicate the treatment codes that make up the comparison. nr indicates the number of times the given comparison is made within the network.

If there is only a single observation for each study within the dataset (i.e. as for standard network meta-analysis) nr will represent the number of studies that compare treatments t1 and t2.

If there are multiple observations for each study within the dataset (as in MBNMAtime) nr will represent the number of time points in the dataset in which treatments t1 and t2 are compared.


Synthesise single arm studies with repeated observations of the same treatment over time

Description

Synthesises single arm studies with repeated measures by applying a particular time-course function. Used in predicting mean responses from a time-course MBNMA. The same parameterisation of the time course must be used as in the MBNMA.

Usage

ref.synth(
  data.ab,
  mbnma,
  synth = "random",
  link = mbnma$model.arg$link,
  n.iter = mbnma$BUGSoutput$n.iter,
  n.burnin = mbnma$BUGSoutput$n.burnin,
  n.thin = mbnma$BUGSoutput$n.thin,
  n.chains = mbnma$BUGSoutput$n.chains,
  ...
)

Arguments

data.ab

A data frame of arm-level data in "long" format containing the columns:

  • studyID Study identifiers

  • time Numeric data indicating follow-up times

  • y Numeric data indicating the mean response for a given observation

  • se Numeric data indicating the standard error for a given observation

mbnma

An S3 object of class "mbnma" generated by running a time-course MBNMA model

synth

A character object that can take the value "common" or "random" that specifies the the type of pooling to use for synthesis of ref.resp. Using "random" rather than "common" for synth will result in wider 95\% CrI for predictions.

link

Can take either "identity" (the default), "log" (for modelling Ratios of Means (Friedrich et al. 2011)) or "smd" (for modelling Standardised Mean Differences - although this also corresponds to an identity link function).

n.iter

number of total iterations per chain (including burn in; default: 2000)

n.burnin

length of burn in, i.e. number of iterations to discard at the beginning. Default is n.iter/2, that is, discarding the first half of the simulations. If n.burnin is 0, jags() will run 100 iterations for adaption.

n.thin

thinning rate. Must be a positive integer. Set n.thin > 1 to save memory and computation time if n.iter is large. Default is max(1, floor(n.chains * (n.iter-n.burnin) / 1000)) which will only thin if there are at least 2000 simulations.

n.chains

number of Markov chains (default: 3)

...

Arguments to be sent to R2jags for synthesis of the network reference treatment effect (using ref.synth())

Details

data.ab can be a collection of studies that closely resemble the population of interest intended for the prediction, which could be different to those used to estimate the MBNMA model, and could be include single arms of RCTs or observational studies. If other data is not available, the data used to estimate the MBNMA model can be used by selecting only the studies and arms that specify the network reference treatment responses.

Value

A list of named elements corresponding to each time-course parameter within an MBNMA model that contain the median posterior value for the network reference treatment response.

Examples

# Create an mb.network object from a dataset
network <- mb.network(osteopain)

# Run an MBNMA model with an Emax time-course
emax <- mb.run(network,
  fun=temax(pool.emax="rel", method.emax="common",
    pool.et50="abs", method.et50="random"))

# Generate a set of studies with which to estimate the network reference treatment response
paindata.ref <- osteopain[osteopain$treatname=="Placebo_0",]

# Estimate the network reference treatment effect using common effects meta-analysis
ref.synth(data.ab=paindata.ref, mbnma=emax, synth="common")

# Estimate the network reference treatment effect using random effects meta-analysis
ref.synth(data.ab=paindata.ref, mbnma=emax, synth="random")

Checks the validity of ref.resp if given as data frame

Description

Ensures ref.resp takes the correct form to allow for synthesis of network reference treatment effect if data is provided for meta-analysis

Usage

ref.validate(data.ab)

Arguments

data.ab

A data frame of arm-level data in "long" format containing the columns:

  • studyID Study identifiers

  • time Numeric data indicating follow-up times

  • y Numeric data indicating the mean response for a given observation

  • se Numeric data indicating the standard error for a given observation

Value

Returns data.ab, the data frame used to estimate the network reference treatment time-course function, with additional required indices added.


Removes any loops from MBNMA model JAGS code that do not contain any expressions

Description

Removes any loops from MBNMA model JAGS code that do not contain any expressions

Usage

remove.loops(model)

Arguments

model

A character object of JAGS MBNMA model code

Value

A character vector of JAGS MBNMA model code that has had empty loops removed from it


Replace original priors in an MBNMA model with new priors

Description

Identical to replace.prior() in MBNMAdose.

Usage

replace.prior(priors, model = NULL, mbnma = NULL)

Arguments

priors

A named list of parameter values (without indices) and replacement prior distribution values given as strings using distributions as specified in JAGS syntax (see Plummer (2017)).

model

A character object of JAGS MBNMA model code

mbnma

An S3 object of class c("mbnma", "rjags") generated by running a time-course MBNMA model.

Details

This function takes new priors, as specified by the user, and adds them to the JAGS code from an MBNMA model. New priors replace old priors in the JAGS model.

Values in priors can include any JAGS functions/distributions (e.g. censoring/truncation).

Value

A character object of JAGS MBNMA model code that includes the new priors in place of original priors


Print summary mb.network information to the console

Description

Print summary mb.network information to the console

Usage

## S3 method for class 'mb.network'
summary(object, ...)

Arguments

object

An object of class mb.network.

...

further arguments passed to or from other methods

Value

Prints summary details of x to the console.


Prints summary of mb.predict object

Description

Prints a summary table of the mean of MCMC iterations at each time point for each treatment

Usage

## S3 method for class 'mb.predict'
summary(object, ...)

Arguments

object

An object of class "mb.predict"

...

further arguments passed to or from other methods

Value

A matrix containing times at which responses have been predicted (time) and an additional column for each treatment for which responses have been predicted. Each row represents mean MCMC predicted responses for each treatment at a particular time.

Examples

# Define network
network <- mb.network(obesityBW_CFB, reference="plac")

# Run an MBNMA with a quadratic time-course function
quad <- mb.run(network,
  fun=tpoly(degree=2, pool.1="rel", method.1="common",
    pool.2="rel", method.2="common"),
  intercept=TRUE)

# Predict responses
pred <- predict(quad, times=c(0:50), treats=c(1:5),
  ref.resp = network$data.ab[network$data.ab$treatment==1,],
  E0=10)

# Generate summary of predictions
summary(pred)

Print summary MBNMA results to the console

Description

Print summary MBNMA results to the console

Usage

## S3 method for class 'mbnma'
summary(object, ...)

Arguments

object

An S3 object of class("mbnma") generated by running a time-course MBNMA model

...

further arguments passed to knitr::kable

Value

Prints summary details of the model results to the console


Takes node-split results and produces summary data frame

Description

Takes node-split results and produces summary data frame

Usage

## S3 method for class 'nodesplit'
summary(object, ...)

Arguments

object

An object of class "nodesplit" generated by mb.nodeplit()

...

further arguments passed to or from other methods

Value

A data frame of summary node-split results with the following variables:

  • Comparison The treatment comparison on which a node-split has been performed

  • Time.Param The time-course parameter on which a node-split has been performed

  • Evidence The evidence contribution for the given comparison (either "Direct" or "Indirect")

  • Median The posterior median

  • ⁠2.5%⁠ The lower 95% credible interval limit

  • ⁠97.5%⁠ The upper 95% credible interval limit

  • p.value The Bayesian p-value for the overlap between direct and indirect evidence for the given comparison (it will therefore have an identical value for direct and indirect evidence within a particular comparison and time-course parameter)


Emax time-course function

Description

** For version 0.2.3: to ensure positive posterior values, et50 and hill parameters are now modeled on the natural scale using a half-normal prior rather than a symmetrical prior on the exponential scale to improve model stability **

Usage

temax(
  pool.emax = "rel",
  method.emax = "common",
  pool.et50 = "rel",
  method.et50 = "common",
  pool.hill = NULL,
  method.hill = NULL,
  p.expon = FALSE
)

Arguments

pool.emax

Pooling for Emax parameter. Can take "rel" or "abs" (see details).

method.emax

Method for synthesis of Emax parameter. Can take ⁠"common⁠, "random", or be assigned a numeric value (see details).

pool.et50

Pooling for ET50 parameter. Can take "rel" or "abs" (see details).

method.et50

Method for synthesis of ET50 parameter. Can take ⁠"common⁠, "random", or be assigned a numeric value (see details).

pool.hill

Pooling for Hill parameter. Can take "rel" or "abs" (see details).

method.hill

Method for synthesis of Hill parameter. Can take ⁠"common⁠, "random", or be assigned a numeric value (see details).

p.expon

Should parameters that can only take positive values be modeled on the exponential scale (TRUE) or should they be assigned a prior that restricts the posterior to positive values (FALSE)

Details

  • Emax represents the maximum response.

  • ET50 represents the time at which 50% of the maximum response is achieved. This can only take positive values and so is modeled on the exponential scale and assigned a symmetrical normal prior Alternatively it can be assigned a normal prior truncated at zero (half-normal) (this will be the default in MBNMAtime version >=0.2.3).

  • Hill is the Hill parameter, which allows for a sigmoidal function. This can only take positive values and so is modeled on the exponential scale and assigned a symmetrical normal prior Alternatively it can be assigned a normal prior truncated at zero (half-normal) (this will be the default in MBNMAtime version >=0.2.3).

Without Hill parameter:

Emax×xET50+x\frac{E_{max}\times{x}}{ET_{50}+x}

With Hill parameter:

Emax×xhillET50×hill+xhill\frac{E_{max}\times{x^{hill}}}{ET_{50}\times{hill}+x^{hill}}

Value

An object of class("timefun")

Time-course parameters

Time-course parameters in the model must be specified using a pool and a method prefix.

pool is used to define the approach used for pooling of a given time-course parameter and can take any of:

Argument Model specification
"rel" Indicates that relative effects should be pooled for this time-course parameter. Relative effects preserve randomisation within included studies, are likely to vary less between studies (only due to effect modification), and allow for testing of consistency between direct and indirect evidence. Pooling follows the general approach for Network Meta-Analysis proposed by Lu and Ades (2004).
"abs" Indicates that study arms should be pooled across the whole network for this time-course parameter independently of assigned treatment to estimate an absolute effect. This implies estimating a single value across the network for this time-course parameter, and may therefore be making strong assumptions of similarity.

method is used to define the model used for meta-analysis for a given time-course parameter and can take any of the following values:

Argument Model specification
"common" Implies that all studies estimate the same true effect (often called a "fixed effect" meta-analysis)
"random" Implies that all studies estimate a separate true effect, but that each of these true effects vary randomly around a true mean effect. This approach allows for modelling of between-study heterogeneity.
numeric() Assigned a numeric value, indicating that this time-course parameter should not be estimated from the data but should be assigned the numeric value determined by the user. This can be useful for fixing specific time-course parameters (e.g. Hill parameters in Emax functions, power parameters in fractional polynomials) to a single value.

When relative effects are modelled on more than one time-course parameter, correlation between them is automatically estimated using a vague inverse-Wishart prior. This prior can be made slightly more informative by specifying the scale matrix omega and by changing the degrees of freedom of the inverse-Wishart prior using the priors argument in mb.run().

References

Lu G, Ades AE (2004). “Combination of direct and indirect evidence in mixed treatment comparisons.” Stat Med, 23(20), 3105-24. ISSN 0277-6715 (Print) 0277-6715 (Linking), doi:10.1002/sim.1875, https://pubmed.ncbi.nlm.nih.gov/15449338/.

Examples

# Model without a Hill parameter
temax(pool.emax="rel", method.emax="random", pool.et50="abs", method.et50="common")

# Model including a Hill parameter and defaults for Emax and ET50 parameters
temax(pool.hill="abs", method.hill="common")

Fractional polynomial time-course function

Description

As first described for use in Network Meta-Analysis by Jansen et al. (2015).

Usage

tfpoly(
  degree = 1,
  pool.1 = "rel",
  method.1 = "common",
  pool.2 = "rel",
  method.2 = "common",
  method.power1 = 0,
  method.power2 = 0
)

Arguments

degree

The degree of the fractional polynomial as defined in Royston and Altman (1994)

pool.1

Pooling for the 1st fractional polynomial coefficient. Can take "rel" or "abs" (see details).

method.1

Method for synthesis of the 1st fractional polynomial coefficient. Can take ⁠"common⁠, "random", or be assigned a numeric value (see details).

pool.2

Pooling for the 2nd fractional polynomial coefficient. Can take "rel" or "abs" (see details).

method.2

Method for synthesis of the 2nd fractional polynomial coefficient. Can take ⁠"common⁠, "random", or be assigned a numeric value (see details).

method.power1

Value for the 1st fractional polynomial power. Must take any numeric value in the set ⁠-2, -1, -0.5, 0, 0.5, 1, 2, 3⁠. pool for this parameter is set to "abs".

method.power2

Value for the 2nd fractional polynomial power. Must take any numeric value in the set ⁠-2, -1, -0.5, 0, 0.5, 1, 2, 3⁠. pool for this parameter is set to "abs".

Details

  • β1\beta_1 represents the 1st coefficient.

  • β2\beta_2 represents the 2nd coefficient.

  • p1p_1 represents the 1st power

  • p2p_2 represents the 2nd power

For a polynomial of degree=1:

β1xp1{\beta_1}x^{p_1}

For a polynomial of degree=2:

β1xp1+β2xp2{\beta_1}x^{p_1}+{\beta_2}x^{p_2}

x(p)x^{(p)} is a regular power except where p=0p=0, where x(0)=ln(x)x^{(0)}=ln(x). If a fractional polynomial power pmp_m repeats within the function it is multiplied by another ln(x)ln(x).

Value

An object of class("timefun")

Time-course parameters

Time-course parameters in the model must be specified using a pool and a method prefix.

pool is used to define the approach used for pooling of a given time-course parameter and can take any of:

Argument Model specification
"rel" Indicates that relative effects should be pooled for this time-course parameter. Relative effects preserve randomisation within included studies, are likely to vary less between studies (only due to effect modification), and allow for testing of consistency between direct and indirect evidence. Pooling follows the general approach for Network Meta-Analysis proposed by Lu and Ades (2004).
"abs" Indicates that study arms should be pooled across the whole network for this time-course parameter independently of assigned treatment to estimate an absolute effect. This implies estimating a single value across the network for this time-course parameter, and may therefore be making strong assumptions of similarity.

method is used to define the model used for meta-analysis for a given time-course parameter and can take any of the following values:

Argument Model specification
"common" Implies that all studies estimate the same true effect (often called a "fixed effect" meta-analysis)
"random" Implies that all studies estimate a separate true effect, but that each of these true effects vary randomly around a true mean effect. This approach allows for modelling of between-study heterogeneity.
numeric() Assigned a numeric value, indicating that this time-course parameter should not be estimated from the data but should be assigned the numeric value determined by the user. This can be useful for fixing specific time-course parameters (e.g. Hill parameters in Emax functions, power parameters in fractional polynomials) to a single value.

When relative effects are modelled on more than one time-course parameter, correlation between them is automatically estimated using a vague inverse-Wishart prior. This prior can be made slightly more informative by specifying the scale matrix omega and by changing the degrees of freedom of the inverse-Wishart prior using the priors argument in mb.run().

References

Jansen JP, Vieira MC, Cope S (2015). “Network meta-analysis of longitudinal data using fractional polynomials.” Stat Med, 34(15), 2294-311. ISSN 1097-0258 (Electronic) 0277-6715 (Linking), doi:10.1002/sim.6492, https://pubmed.ncbi.nlm.nih.gov/25877808/.

Lu G, Ades AE (2004). “Combination of direct and indirect evidence in mixed treatment comparisons.” Stat Med, 23(20), 3105-24. ISSN 0277-6715 (Print) 0277-6715 (Linking), doi:10.1002/sim.1875, https://pubmed.ncbi.nlm.nih.gov/15449338/.

Royston P, Altman D (1994). “Regression Using Fractional Polynomials of Continuous Covariates: Parsimonious Parametric Modelling.” Journal of the Royal Statistical Society: Series C, 43(3), 429-467.

Examples

# 1st order fractional polynomial with random effects
tfpoly(pool.1="rel", method.1="random")

# 2nd order fractional polynomial
# with a single absolute parameter estimated for the 2nd coefficient
# 1st power equal to zero
tfpoly(degree=2, pool.1="rel", method.1="common",
  pool.2="abs", method.2="random",
  method.power1=0)

Plot raw responses over time by treatment or class

Description

Plot raw responses over time by treatment or class

Usage

timeplot(network, level = "treatment", plotby = "arm", link = "identity", ...)

Arguments

network

An object of class "mb.network".

level

A string indicating whether nodes/facets should represent treatment or class in the plot. Can be used to examine the expected impact of modelling class/agent effects.

plotby

A character object that can take either "arm" to indicate that raw responses should be plotted separately for each study arm, or "rel" to indicate that the within-study relative effects/treatment differences should be plotted. In this way the time-course of both the absolute effects and the relative effects can be examined.

link

Can take either "identity" (the default), "log" (for modelling Ratios of Means (Friedrich et al. 2011)) or "smd" (for modelling Standardised Mean Differences - although this also corresponds to an identity link function).

...

Arguments to be sent to ggplot()

Details

Plots can be faceted by either treatment (level="treatment") or class (level="class") to investigate similarity of treatment responses within classes/agents. Points represent observed responses and lines connect between observations within the same study and arm.

Value

The function returns an object of ⁠class(c("gg", "ggplot")⁠. Characteristics of the object can therefore be amended as with other plots generated by ggplot().

Examples

# Make network
goutnet <- mb.network(goutSUA_CFB)

# Use timeplot to plot responses grouped by treatment
timeplot(goutnet)

# Use timeplot ot plot resposes grouped by class
timeplot(goutnet, level="class")

# Plot matrix of relative effects
timeplot(goutnet, level="class", plotby="rel")

# Plot using Standardised Mean Differences
copdnet <- mb.network(copd)
timeplot(copdnet, plotby="rel", link="smd")

Integrated Two-Component Prediction (ITP) function

Description

Similar parameterisation to the Emax model but with non-asymptotic maximal effect (Emax). Proposed by proposed by Fu and Manner (2010)

Usage

titp(
  pool.emax = "rel",
  method.emax = "common",
  pool.rate = "rel",
  method.rate = "common",
  p.expon = FALSE
)

Arguments

pool.emax

Pooling for exponential Emax parameter. Can take "rel" or "abs" (see details).

method.emax

Method for synthesis of exponential Emax parameter. Can take ⁠"common⁠, "random", or be assigned a numeric value (see details).

pool.rate

Pooling for parameter controlling rate of onset. Default is NULL which avoids including this parameter (i.e. fixes it to 1 for all treatments). Can take "rel" or "abs" (see details).

method.rate

Method for synthesis of parameter controlling rate of onset. Can take ⁠"common⁠, "random", or be assigned a numeric value (see details).

p.expon

Should parameters that can only take positive values be modeled on the exponential scale (TRUE) or should they be assigned a prior that restricts the posterior to positive values (FALSE)

Details

Emax×(1exp(rate×x))(1exp(rate×max(x))){E_{max}}\times\frac{(1-exp(-{rate}\times{x}))}{(1-exp(-{rate}\times{max(x)}))}

Value

An object of class("timefun")

Time-course parameters

Time-course parameters in the model must be specified using a pool and a method prefix.

pool is used to define the approach used for pooling of a given time-course parameter and can take any of:

Argument Model specification
"rel" Indicates that relative effects should be pooled for this time-course parameter. Relative effects preserve randomisation within included studies, are likely to vary less between studies (only due to effect modification), and allow for testing of consistency between direct and indirect evidence. Pooling follows the general approach for Network Meta-Analysis proposed by Lu and Ades (2004).
"abs" Indicates that study arms should be pooled across the whole network for this time-course parameter independently of assigned treatment to estimate an absolute effect. This implies estimating a single value across the network for this time-course parameter, and may therefore be making strong assumptions of similarity.

method is used to define the model used for meta-analysis for a given time-course parameter and can take any of the following values:

Argument Model specification
"common" Implies that all studies estimate the same true effect (often called a "fixed effect" meta-analysis)
"random" Implies that all studies estimate a separate true effect, but that each of these true effects vary randomly around a true mean effect. This approach allows for modelling of between-study heterogeneity.
numeric() Assigned a numeric value, indicating that this time-course parameter should not be estimated from the data but should be assigned the numeric value determined by the user. This can be useful for fixing specific time-course parameters (e.g. Hill parameters in Emax functions, power parameters in fractional polynomials) to a single value.

References

Fu H, Manner D (2010). “Bayesian adaptive dose-finding studies with delayed responses.” J Biopharm Stat, 20(5), 1055-1070. doi:10.1080/10543400903315740.

Lu G, Ades AE (2004). “Combination of direct and indirect evidence in mixed treatment comparisons.” Stat Med, 23(20), 3105-24. ISSN 0277-6715 (Print) 0277-6715 (Linking), doi:10.1002/sim.1875, https://pubmed.ncbi.nlm.nih.gov/15449338/.

Examples

titp(pool.emax="rel", method.emax="random")
titp(pool.emax="abs")

Log-linear (exponential) time-course function

Description

rate×log(x+1)rate\times{log(x + 1)}

Usage

tloglin(pool.rate = "rel", method.rate = "common")

Arguments

pool.rate

Pooling for rate parameter. Can take "rel" or "abs" (see details).

method.rate

Method for synthesis of rate parameter. Can take ⁠"common⁠, "random", or be assigned a numeric value (see details).

Value

An object of class("timefun")

Time-course parameters

Time-course parameters in the model must be specified using a pool and a method prefix.

pool is used to define the approach used for pooling of a given time-course parameter and can take any of:

Argument Model specification
"rel" Indicates that relative effects should be pooled for this time-course parameter. Relative effects preserve randomisation within included studies, are likely to vary less between studies (only due to effect modification), and allow for testing of consistency between direct and indirect evidence. Pooling follows the general approach for Network Meta-Analysis proposed by Lu and Ades (2004).
"abs" Indicates that study arms should be pooled across the whole network for this time-course parameter independently of assigned treatment to estimate an absolute effect. This implies estimating a single value across the network for this time-course parameter, and may therefore be making strong assumptions of similarity.

method is used to define the model used for meta-analysis for a given time-course parameter and can take any of the following values:

Argument Model specification
"common" Implies that all studies estimate the same true effect (often called a "fixed effect" meta-analysis)
"random" Implies that all studies estimate a separate true effect, but that each of these true effects vary randomly around a true mean effect. This approach allows for modelling of between-study heterogeneity.
numeric() Assigned a numeric value, indicating that this time-course parameter should not be estimated from the data but should be assigned the numeric value determined by the user. This can be useful for fixing specific time-course parameters (e.g. Hill parameters in Emax functions, power parameters in fractional polynomials) to a single value.

References

Lu G, Ades AE (2004). “Combination of direct and indirect evidence in mixed treatment comparisons.” Stat Med, 23(20), 3105-24. ISSN 0277-6715 (Print) 0277-6715 (Linking), doi:10.1002/sim.1875, https://pubmed.ncbi.nlm.nih.gov/15449338/.

Examples

tloglin(pool.rate="rel", method.rate="random")
tloglin(pool.rate="abs")

Polynomial time-course function

Description

Polynomial time-course function

Usage

tpoly(
  degree = 1,
  pool.1 = "rel",
  method.1 = "common",
  pool.2 = "rel",
  method.2 = "common",
  pool.3 = "rel",
  method.3 = "common",
  pool.4 = "rel",
  method.4 = "common"
)

Arguments

degree

The degree of the polynomial - e.g. degree=1 for linear, degree=2 for quadratic, degree=3 for cubic.

pool.1

Pooling for the 1st polynomial coefficient. Can take "rel" or "abs" (see details).

method.1

Method for synthesis of the 1st polynomial coefficient.Can take ⁠"common⁠, "random", or be assigned a numeric value (see details).

pool.2

Pooling for the 2nd polynomial coefficient. Can take "rel" or "abs" (see details).

method.2

Method for synthesis of the 2nd polynomial coefficient. Can take ⁠"common⁠, "random", or be assigned a numeric value (see details).

pool.3

Pooling for the 3rd polynomial coefficient. Can take "rel" or "abs" (see details).

method.3

Method for synthesis of the 3rd polynomial coefficient. Can take ⁠"common⁠, "random", or be assigned a numeric value (see details).

pool.4

Pooling for the 4th polynomial coefficient. Can take "rel" or "abs" (see details).

method.4

Method for synthesis of the 4th polynomial coefficient. Can take ⁠"common⁠, "random", or be assigned a numeric value (see details).

Details

  • β1\beta_1 represents the 1st coefficient.

  • β2\beta_2 represents the 2nd coefficient.

  • β3\beta_3 represents the 3rd coefficient.

  • β4\beta_4 represents the 4th coefficient.

Linear model:

β1x\beta_1{x}

Quadratic model:

β1x+β2x2\beta_1{x} + \beta_2{x^2}

Cubic model:

β1x+β2x2+β3x3\beta_1{x} + \beta_2{x^2} + \beta_3{x^3}

Quartic model:

β1x+β2x2+β3x3+β4x4\beta_1{x} + \beta_2{x^2} + \beta_3{x^3} + \beta_4{x^4}

Value

An object of class("timefun")

Time-course parameters

Time-course parameters in the model must be specified using a pool and a method prefix.

pool is used to define the approach used for pooling of a given time-course parameter and can take any of:

Argument Model specification
"rel" Indicates that relative effects should be pooled for this time-course parameter. Relative effects preserve randomisation within included studies, are likely to vary less between studies (only due to effect modification), and allow for testing of consistency between direct and indirect evidence. Pooling follows the general approach for Network Meta-Analysis proposed by Lu and Ades (2004).
"abs" Indicates that study arms should be pooled across the whole network for this time-course parameter independently of assigned treatment to estimate an absolute effect. This implies estimating a single value across the network for this time-course parameter, and may therefore be making strong assumptions of similarity.

method is used to define the model used for meta-analysis for a given time-course parameter and can take any of the following values:

Argument Model specification
"common" Implies that all studies estimate the same true effect (often called a "fixed effect" meta-analysis)
"random" Implies that all studies estimate a separate true effect, but that each of these true effects vary randomly around a true mean effect. This approach allows for modelling of between-study heterogeneity.
numeric() Assigned a numeric value, indicating that this time-course parameter should not be estimated from the data but should be assigned the numeric value determined by the user. This can be useful for fixing specific time-course parameters (e.g. Hill parameters in Emax functions, power parameters in fractional polynomials) to a single value.

When relative effects are modelled on more than one time-course parameter, correlation between them is automatically estimated using a vague inverse-Wishart prior. This prior can be made slightly more informative by specifying the scale matrix omega and by changing the degrees of freedom of the inverse-Wishart prior using the priors argument in mb.run().

References

Lu G, Ades AE (2004). “Combination of direct and indirect evidence in mixed treatment comparisons.” Stat Med, 23(20), 3105-24. ISSN 0277-6715 (Print) 0277-6715 (Linking), doi:10.1002/sim.1875, https://pubmed.ncbi.nlm.nih.gov/15449338/.

Examples

# Linear model with random effects
tpoly(pool.1="rel", method.1="random")

# Quadratic model with a single absolute parameter estimated for the 2nd coefficient
tpoly(pool.1="rel", method.1="common", pool.2="abs", method.2="random")

Spline time-course functions

Description

Used to fit B-splines, natural cubic splines, and piecewise linear splines(Perperoglu et al. 2019).

Usage

tspline(
  type = "bs",
  knots = 1,
  degree = 1,
  pool.1 = "rel",
  method.1 = "common",
  pool.2 = "rel",
  method.2 = "common",
  pool.3 = "rel",
  method.3 = "common",
  pool.4 = "rel",
  method.4 = "common"
)

Arguments

type

The type of spline. Can take "bs" (B-spline), "ns" (natural cubic spline) or "ls" (piecewise linear spline)

knots

The number/location of spline internal knots. If a single number is given it indicates the number of knots (they will be equally spaced across the range of time points). If a numeric vector is given it indicates the location of the knots.

degree

The degree of the piecewise B-spline polynomial - e.g. degree=1 for linear, degree=2 for quadratic, degree=3 for cubic.

pool.1

Pooling for the 1st coefficient. Can take "rel" or "abs" (see details).

method.1

Method for synthesis of the 1st coefficient. Can take ⁠"common⁠, "random", or be assigned a numeric value (see details).

pool.2

Pooling for the 2nd coefficient. Can take "rel" or "abs" (see details).

method.2

Method for synthesis of the 2nd coefficient. Can take ⁠"common⁠, "random", or be assigned a numeric value (see details).

pool.3

Pooling for the 3rd coefficient. Can take "rel" or "abs" (see details).

method.3

Method for synthesis of the 3rd coefficient. Can take ⁠"common⁠, "random", or be assigned a numeric value (see details).

pool.4

Pooling for the 4th coefficient. Can take "rel" or "abs" (see details).

method.4

Method for synthesis of the 4th coefficient. Can take ⁠"common⁠, "random", or be assigned a numeric value (see details).

Value

An object of class("timefun")

Time-course parameters

Time-course parameters in the model must be specified using a pool and a method prefix.

pool is used to define the approach used for pooling of a given time-course parameter and can take any of:

Argument Model specification
"rel" Indicates that relative effects should be pooled for this time-course parameter. Relative effects preserve randomisation within included studies, are likely to vary less between studies (only due to effect modification), and allow for testing of consistency between direct and indirect evidence. Pooling follows the general approach for Network Meta-Analysis proposed by Lu and Ades (2004).
"abs" Indicates that study arms should be pooled across the whole network for this time-course parameter independently of assigned treatment to estimate an absolute effect. This implies estimating a single value across the network for this time-course parameter, and may therefore be making strong assumptions of similarity.

method is used to define the model used for meta-analysis for a given time-course parameter and can take any of the following values:

Argument Model specification
"common" Implies that all studies estimate the same true effect (often called a "fixed effect" meta-analysis)
"random" Implies that all studies estimate a separate true effect, but that each of these true effects vary randomly around a true mean effect. This approach allows for modelling of between-study heterogeneity.
numeric() Assigned a numeric value, indicating that this time-course parameter should not be estimated from the data but should be assigned the numeric value determined by the user. This can be useful for fixing specific time-course parameters (e.g. Hill parameters in Emax functions, power parameters in fractional polynomials) to a single value.

When relative effects are modelled on more than one time-course parameter, correlation between them is automatically estimated using a vague inverse-Wishart prior. This prior can be made slightly more informative by specifying the scale matrix omega and by changing the degrees of freedom of the inverse-Wishart prior using the priors argument in mb.run().

References

Lu G, Ades AE (2004). “Combination of direct and indirect evidence in mixed treatment comparisons.” Stat Med, 23(20), 3105-24. ISSN 0277-6715 (Print) 0277-6715 (Linking), doi:10.1002/sim.1875, https://pubmed.ncbi.nlm.nih.gov/15449338/.

Perperoglu A, Sauerbrei W, Abrahamowicz M, Schmid M (2019). “A review of spline function procedures in R.” BMC Medical Research Methodology, 19(46), 1-16. doi:10.1186/s12874-019-0666-3.

Examples

# Second order B spline with 2 knots and random effects on the 2nd coefficient
tspline(type="bs", knots=2, degree=2,
  pool.1="rel", method.1="common",
  pool.2="rel", method.2="random")

# Piecewise linear spline with knots at 0.1 and 0.5 quantiles
# Single parameter independent of treatment estimated for 1st coefficient
#with random effects
tspline(type="ls", knots=c(0.1,0.5),
  pool.1="abs", method.1="random",
  pool.2="rel", method.2="common")

User-defined time-course function

Description

User-defined time-course function

Usage

tuser(
  fun,
  pool.1 = "rel",
  method.1 = "common",
  pool.2 = "rel",
  method.2 = "common",
  pool.3 = "rel",
  method.3 = "common",
  pool.4 = "rel",
  method.4 = "common"
)

Arguments

fun

A formula specifying any relationship including time and one/several of: beta.1, beta.2, beta.3, beta.4.

pool.1

Pooling for beta.1. Can take ⁠"common⁠, "random", or be assigned a numeric value (see details).

method.1

Method for synthesis of beta.1. Can take ⁠"common⁠ or "random" (see details).

pool.2

Pooling for beta.2. Can take ⁠"common⁠, "random", or be assigned a numeric value (see details).

method.2

Method for synthesis of ⁠beta.2. Can take ⁠"commonor"random"' (see details).

pool.3

Pooling for beta.3. Can take ⁠"common⁠, "random", or be assigned a numeric value (see details).

method.3

Method for synthesis of ⁠beta.3. Can take ⁠"commonor"random"' (see details).

pool.4

Pooling for beta.4. Can take ⁠"common⁠, "random", or be assigned a numeric value (see details).

method.4

Method for synthesis of beta.4. Can take ⁠"common⁠, "random", or be assigned a numeric value (see details).

Value

An object of class("timefun")

Time-course parameters

Time-course parameters in the model must be specified using a pool and a method prefix.

pool is used to define the approach used for pooling of a given time-course parameter and can take any of:

Argument Model specification
"rel" Indicates that relative effects should be pooled for this time-course parameter. Relative effects preserve randomisation within included studies, are likely to vary less between studies (only due to effect modification), and allow for testing of consistency between direct and indirect evidence. Pooling follows the general approach for Network Meta-Analysis proposed by Lu and Ades (2004).
"abs" Indicates that study arms should be pooled across the whole network for this time-course parameter independently of assigned treatment to estimate an absolute effect. This implies estimating a single value across the network for this time-course parameter, and may therefore be making strong assumptions of similarity.

method is used to define the model used for meta-analysis for a given time-course parameter and can take any of the following values:

Argument Model specification
"common" Implies that all studies estimate the same true effect (often called a "fixed effect" meta-analysis)
"random" Implies that all studies estimate a separate true effect, but that each of these true effects vary randomly around a true mean effect. This approach allows for modelling of between-study heterogeneity.
numeric() Assigned a numeric value, indicating that this time-course parameter should not be estimated from the data but should be assigned the numeric value determined by the user. This can be useful for fixing specific time-course parameters (e.g. Hill parameters in Emax functions, power parameters in fractional polynomials) to a single value.

When relative effects are modelled on more than one time-course parameter, correlation between them is automatically estimated using a vague inverse-Wishart prior. This prior can be made slightly more informative by specifying the scale matrix omega and by changing the degrees of freedom of the inverse-Wishart prior using the priors argument in mb.run().

References

Lu G, Ades AE (2004). “Combination of direct and indirect evidence in mixed treatment comparisons.” Stat Med, 23(20), 3105-24. ISSN 0277-6715 (Print) 0277-6715 (Linking), doi:10.1002/sim.1875, https://pubmed.ncbi.nlm.nih.gov/15449338/.

Examples

timecourse <- ~ beta.1 * (1/(time+1)) + beta.2 * time^2
tuser(fun=timecourse,
  pool.1="abs", method.1="common",
  pool.2="rel", method.2="common")

Adds sections of JAGS code for an MBNMA model that correspond to beta parameters

Description

Adds sections of JAGS code for an MBNMA model that correspond to beta parameters

Usage

write.beta(model, timecourse, fun, UME, class.effect)

Arguments

model

A character object of JAGS MBNMA model code

timecourse

A character object representing the time-course used in the MBNMA model

fun

An object of class "timefun" generated (see Details) using any of tloglin(), tpoly(), titp(), temax(), tfpoly(), tspline() or tuser()

UME

Can take either TRUE or FALSE (for an unrelated mean effects model on all or no time-course parameters respectively) or can be a vector of parameter name strings to model as UME. For example: c("beta.1", "beta.2").

class.effect

A list of named strings that determines which time-course parameters to model with a class effect and what that effect should be ("common" or "random"). For example: list(emax="common", et50="random").

Value

A character vector of JAGS MBNMA model code that includes beta parameter components of the model


Checks validity of arguments for mb.write

Description

Checks validity of arguments for mb.write

Usage

write.check(
  fun = tpoly(degree = 1),
  positive.scale = TRUE,
  intercept = NULL,
  rho = 0,
  covar = NULL,
  omega = NULL,
  link = "identity",
  sdscale = FALSE,
  class.effect = list(),
  UME = c()
)

Arguments

fun

An object of class "timefun" generated (see Details) using any of tloglin(), tpoly(), titp(), temax(), tfpoly(), tspline() or tuser()

positive.scale

A boolean object that indicates whether all continuous mean responses (y) are positive and therefore whether the baseline response should be given a prior that constrains it to be positive (e.g. for scales that cannot be <0).

intercept

A boolean object that indicates whether an intercept (written as alpha in the model) is to be included. If left as NULL (the default), an intercept will be included only for studies reporting absolute means, and will be excluded for studies reporting change from baseline (as indicated in network$cfb).

rho

The correlation coefficient when modelling within-study correlation between time points. The default is a string representing a prior distribution in JAGS, indicating that it be estimated from the data (e.g. rho="dunif(0,1)"). rho also be assigned a numeric value (e.g. rho=0.7), which fixes rho in the model to this value (e.g. for use in a deterministic sensitivity analysis). If set to rho=0 (the default) then this implies modelling no correlation between time points.

covar

A character specifying the covariance structure to use for modelling within-study correlation between time-points. This can be done by specifying one of the following:

  • "varadj" - a univariate likelihood with a variance adjustment to assume a constant correlation between subsequent time points (Jansen et al. 2015). This is the default.

  • "CS" - a multivariate normal likelihood with a compound symmetry structure

  • "AR1" - a multivariate normal likelihood with an autoregressive AR1 structure

omega

DEPRECATED IN VERSION 0.2.3 ONWARDS (~uniform(-1,1) now used for correlation between parameters rather than a Wishart prior). A scale matrix for the inverse-Wishart prior for the covariance matrix used to model the correlation between time-course parameters (see Details for time-course functions). omega must be a symmetric positive definite matrix with dimensions equal to the number of time-course parameters modelled using relative effects (pool="rel"). If left as NULL (the default) a diagonal matrix with elements equal to 1 is used.

link

Can take either "identity" (the default), "log" (for modelling Ratios of Means (Friedrich et al. 2011)) or "smd" (for modelling Standardised Mean Differences - although this also corresponds to an identity link function).

sdscale

Logical object to indicate whether to write a model that specifies a reference SD for standardising when modelling using Standardised Mean Differences. Specifying sdscale=TRUE will therefore only modify the model if link function is set to SMD (link="smd").

class.effect

A list of named strings that determines which time-course parameters to model with a class effect and what that effect should be ("common" or "random"). For example: list(emax="common", et50="random").

UME

Can take either TRUE or FALSE (for an unrelated mean effects model on all or no time-course parameters respectively) or can be a vector of parameter name strings to model as UME. For example: c("beta.1", "beta.2").

Details

Used to check if the arguments given to mb.write are valid. The function will return informative errors if arguments are mispecified and will return an object that indicates whether the arguments imply modelling a correlation between time points if it passes.

Value

A boolean object that indicates whether the arguments imply modelling correlation between time points.


Adds correlation between time-course relative effects

Description

This uses a Wishart prior as default for modelling the correlation

Usage

write.cor(model, fun, omega = NULL, class.effect = list())

Arguments

model

A character object of JAGS MBNMA model code

fun

An object of class "timefun" generated (see Details) using any of tloglin(), tpoly(), titp(), temax(), tfpoly(), tspline() or tuser()

omega

DEPRECATED IN VERSION 0.2.3 ONWARDS (~uniform(-1,1) now used for correlation between parameters rather than a Wishart prior). A scale matrix for the inverse-Wishart prior for the covariance matrix used to model the correlation between time-course parameters (see Details for time-course functions). omega must be a symmetric positive definite matrix with dimensions equal to the number of time-course parameters modelled using relative effects (pool="rel"). If left as NULL (the default) a diagonal matrix with elements equal to 1 is used.

class.effect

A list of named strings that determines which time-course parameters to model with a class effect and what that effect should be ("common" or "random"). For example: list(emax="common", et50="random").


Adds sections of JAGS code for an MBNMA model that correspond to the likelihood

Description

Adds sections of JAGS code for an MBNMA model that correspond to the likelihood

Usage

write.likelihood(
  model,
  timecourse,
  rho = 0,
  covar = "varadj",
  link = "identity",
  sdscale = FALSE,
  fun
)

Arguments

model

A character object of JAGS MBNMA model code

timecourse

A character object representing the time-course used in the MBNMA model

rho

The correlation coefficient when modelling within-study correlation between time points. The default is a string representing a prior distribution in JAGS, indicating that it be estimated from the data (e.g. rho="dunif(0,1)"). rho also be assigned a numeric value (e.g. rho=0.7), which fixes rho in the model to this value (e.g. for use in a deterministic sensitivity analysis). If set to rho=0 (the default) then this implies modelling no correlation between time points.

covar

A character specifying the covariance structure to use for modelling within-study correlation between time-points. This can be done by specifying one of the following:

  • "varadj" - a univariate likelihood with a variance adjustment to assume a constant correlation between subsequent time points (Jansen et al. 2015). This is the default.

  • "CS" - a multivariate normal likelihood with a compound symmetry structure

  • "AR1" - a multivariate normal likelihood with an autoregressive AR1 structure

link

Can take either "identity" (the default), "log" (for modelling Ratios of Means (Friedrich et al. 2011)) or "smd" (for modelling Standardised Mean Differences - although this also corresponds to an identity link function).

sdscale

Logical object to indicate whether to write a model that specifies a reference SD for standardising when modelling using Standardised Mean Differences. Specifying sdscale=TRUE will therefore only modify the model if link function is set to SMD (link="smd").

fun

An object of class "timefun" generated (see Details) using any of tloglin(), tpoly(), titp(), temax(), tfpoly(), tspline() or tuser()

Value

A character vector of JAGS MBNMA model code that includes likelihood components of the model


Write the basic JAGS model code for MBNMA to which other lines of model code can be added

Description

Write the basic JAGS model code for MBNMA to which other lines of model code can be added

Usage

write.model()

Value

A character vector of JAGS model code


Write MBNMA time-course models JAGS code for synthesis of studies investigating reference treatment

Description

Writes JAGS code for a Bayesian time-course model for model-based network meta-analysis (MBNMA) that pools reference treatment effects from different studies. This model only pools single study arms and therefore does not pool relative effects.

Usage

write.ref.synth(
  fun = tpoly(degree = 1),
  link = "identity",
  positive.scale = TRUE,
  intercept = TRUE,
  rho = 0,
  covar = "varadj",
  mu.synth = "random",
  priors = NULL
)

Arguments

fun

An object of class "timefun" generated (see Details) using any of tloglin(), tpoly(), titp(), temax(), tfpoly(), tspline() or tuser()

link

Can take either "identity" (the default), "log" (for modelling Ratios of Means (Friedrich et al. 2011)) or "smd" (for modelling Standardised Mean Differences - although this also corresponds to an identity link function).

positive.scale

A boolean object that indicates whether all continuous mean responses (y) are positive and therefore whether the baseline response should be given a prior that constrains it to be positive (e.g. for scales that cannot be <0).

intercept

A boolean object that indicates whether an intercept (written as alpha in the model) is to be included. If left as NULL (the default), an intercept will be included only for studies reporting absolute means, and will be excluded for studies reporting change from baseline (as indicated in network$cfb).

rho

The correlation coefficient when modelling within-study correlation between time points. The default is a string representing a prior distribution in JAGS, indicating that it be estimated from the data (e.g. rho="dunif(0,1)"). rho also be assigned a numeric value (e.g. rho=0.7), which fixes rho in the model to this value (e.g. for use in a deterministic sensitivity analysis). If set to rho=0 (the default) then this implies modelling no correlation between time points.

covar

A character specifying the covariance structure to use for modelling within-study correlation between time-points. This can be done by specifying one of the following:

  • "varadj" - a univariate likelihood with a variance adjustment to assume a constant correlation between subsequent time points (Jansen et al. 2015). This is the default.

  • "CS" - a multivariate normal likelihood with a compound symmetry structure

  • "AR1" - a multivariate normal likelihood with an autoregressive AR1 structure

mu.synth

A string that takes the value fixed or random, indicating the type of synthesis model to use

priors

A named list of parameter values (without indices) and replacement prior distribution values given as strings using distributions as specified in JAGS syntax (see Plummer (2017)).

Value

A character object of JAGS MBNMA model code that includes beta parameter components of the model

Examples

# Write a log-linear time-course MBNMA synthesis model with:
# Common effects for synthesis of mu
# Modelled as ratio of means
model <- write.ref.synth(fun=tloglin(pool.rate="rel", method.rate="common"),
  mu.synth="common", link="log")

cat(model) # Concatenates model representations making code more easily readable

Adds sections of JAGS code for an MBNMA model that correspond to alpha parameters

Description

Adds sections of JAGS code for an MBNMA model that correspond to alpha parameters

Usage

write.timecourse(model, fun, intercept, positive.scale)

Arguments

model

A character string representing the MBNMA model in JAGS code

fun

An object of class "timefun" generated (see Details) using any of tloglin(), tpoly(), titp(), temax(), tfpoly(), tspline() or tuser()

intercept

A boolean object that indicates whether an intercept (written as alpha in the model) is to be included. If left as NULL (the default), an intercept will be included only for studies reporting absolute means, and will be excluded for studies reporting change from baseline (as indicated in network$cfb).

positive.scale

A boolean object that indicates whether all continuous mean responses (y) are positive and therefore whether the baseline response should be given a prior that constrains it to be positive (e.g. for scales that cannot be <0).

timecourse

A character object that contains JAGS code for the time-course component of the model

Value

A list of named elements: model is a character vector of JAGS MBNMA model code that includes alpha parameter components of the model timecourse is a character object that contains JAGS code for the time-course component of the model, for which alpha will be indexed correctly