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 |
Adds follow-up time (fups
, fupcount
) and arm (arms
, narms
) indices to a dataset.
add_index(data.ab, reference = 1)
add_index(data.ab, reference = 1)
data.ab |
A data frame of arm-level data in "long" format containing the columns:
|
reference |
A number or character (depending on the format of |
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).
# 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)
# 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)
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.
alog_pcfb
alog_pcfb
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
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
.
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
binplot( network, overlay.nma = c(0, stats::quantile(network$data.ab$time)), method = "common", link = "identity", lim = "cred", plot.bins = TRUE, legend = TRUE, ... )
binplot( network, overlay.nma = c(0, stats::quantile(network$data.ab$time)), method = "common", link = "identity", lim = "cred", plot.bins = TRUE, legend = TRUE, ... )
network |
An object of class |
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 |
method |
Can take |
link |
Can take either |
lim |
Specifies calculation of either 95% credible intervals ( |
plot.bins |
Plot time bin boundaries as vertical dashed lines. Setting |
legend |
|
... |
Arguments to be sent to R2jags. |
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.
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
.
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.
# 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")
# 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")
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.
copd
copd
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
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
.
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
cumrank(x, sucra = TRUE, ...)
cumrank(x, sucra = TRUE, ...)
x |
An object of class |
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 |
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
.
# 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)
# 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)
This function creates JAGS code snippets for default MBNMA model priors.
default.priors(fun = tloglin())
default.priors(fun = tloglin())
fun |
An object of class |
A list, each element of which is a named JAGS snippet corresponding to a prior in the MBNMA JAGS code.
default.priors(fun=temax()) default.priors(fun=titp(p.expon=TRUE))
default.priors(fun=temax()) default.priors(fun=titp(p.expon=TRUE))
Plot deviance contributions from an MBNMA model
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, ... )
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, ... )
mbnma |
An S3 object of class |
dev.type |
Deviances to plot - can be either residual deviances ( |
plot.type |
Deviances can be plotted either as scatter points ( |
xaxis |
A character object that indicates whether deviance contributions should be plotted
by time ( |
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 |
n.thin |
The thinning rate. Must be a positive integer. Default is the value used in |
... |
Arguments to be sent to |
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.
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.
# 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)
# 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)
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.
diabetes
diabetes
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
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
.
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
fitplot( mbnma, treat.labs = NULL, disp.obs = TRUE, n.iter = round(mbnma$BUGSoutput$n.iter/4), n.thin = mbnma$BUGSoutput$n.thin, ... )
fitplot( mbnma, treat.labs = NULL, disp.obs = TRUE, n.iter = round(mbnma$BUGSoutput$n.iter/4), n.thin = mbnma$BUGSoutput$n.thin, ... )
mbnma |
An S3 object of class |
treat.labs |
A character vector of treatment labels with which to name graph panels.
Can use |
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
|
... |
Arguments to be sent to |
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.
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.
# 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)
# 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
gen.parameters.to.save(fun, model)
gen.parameters.to.save(fun, model)
fun |
An object of class |
model |
A JAGS model written as a character object |
A character vector of parameter names that should be monitored in the model
Get large vector of distinct colours using Rcolorbrewer
genmaxcols()
genmaxcols()
Generates spline basis matrices for fitting to time-course function
genspline( x, spline = "bs", knots = 1, degree = 1, max.time = max(x), boundaries = NULL )
genspline( x, spline = "bs", knots = 1, degree = 1, max.time = max(x), boundaries = NULL )
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 ( |
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, |
degree |
a positive integer giving the degree of the polynomial from which the spline function is composed
(e.g. |
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 |
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.
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)
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)
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.
get.closest.time(network, t = stats::median(network$data.ab$time))
get.closest.time(network, t = stats::median(network$data.ab$time))
network |
An object of class |
t |
The time point at which |
A data frame in long format of responses at the closest time point to t in each arm of each study.
# 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)
# 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)
Takes the earliest time point from each arm in each study within an
mb.network
object. Useful for network plots.
get.earliest.time(network)
get.earliest.time(network)
network |
An object of class |
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
)
# 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)
# 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)
Takes the latest time point from each arm in each study within an
mb.network
object. Useful for network plots.
get.latest.time(network)
get.latest.time(network)
network |
An object of class |
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
)
# 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)
# 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)
Extracts specific information required for prediction from a time-course MBNMA model
get.model.vals( mbnma, E0 = 0, level = "treatments", lim = "cred", link = "identity" )
get.model.vals( mbnma, E0 = 0, level = "treatments", lim = "cred", link = "identity" )
mbnma |
An S3 object of class |
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.
|
level |
Can take either |
lim |
Specifies calculation of either 95% credible intervals ( |
link |
Can take either |
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
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.
get.prior(model)
get.prior(model)
model |
A character object of JAGS MBNMA model code |
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.
A character vector, each element of which is a line of JAGS code corresponding to a prior in the JAGS code.
# 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)
# 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)
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.
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" )
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" )
mbnma |
An S3 object of class |
mbnma.add |
An S3 object of |
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 |
classes |
A character vector of class names for which to calculate relative effects/mean
differences. Must be a subset of |
lim |
Specifies calculation of either 95% credible intervals ( |
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.
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.
# 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)
# 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)
Converts MBNMA data frame to a list for use in JAGS model
getjagsdata( data.ab, fun = NULL, class = FALSE, rho = NULL, covstruct = "CS", link = "identity", sdscale = FALSE, cfb = NULL )
getjagsdata( data.ab, fun = NULL, class = FALSE, rho = NULL, covstruct = "CS", link = "identity", sdscale = FALSE, cfb = NULL )
data.ab |
A data frame of arm-level data in "long" format containing the columns:
|
fun |
An object of class |
class |
A boolean object indicating whether or not |
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. |
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 |
link |
Can take either |
sdscale |
Logical object to indicate whether to write a model that specifies a reference SD
for standardising when modelling using Standardised Mean Differences. Specifying |
cfb |
A logical vector whose length is equal to the unique number of studies in |
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.
# 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")
# 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")
Converts data frame to a list for use in JAGS NMA model
getnmadata(data.ab, link = "identity", sdscale = FALSE)
getnmadata(data.ab, link = "identity", sdscale = FALSE)
data.ab |
A data frame of arm-level data in "long" format containing the columns:
|
link |
Can take either |
sdscale |
Logical object to indicate whether to write a model that specifies a reference SD
for standardising when modelling using Standardised Mean Differences. Specifying |
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
# 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)
# 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)
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.
goutSUA_CFB
goutSUA_CFB
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.
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
.
Pfizer Ltd.
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.
goutSUA_CFBcomb
goutSUA_CFBcomb
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.
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
.
Pfizer Ltd.
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.
hyalarthritis
hyalarthritis
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
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
.
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 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).
inconsistency.loops(data)
inconsistency.loops(data)
data |
A data frame containing variables |
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
.
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.
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/.
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)
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 contrasts within a network that make up all the head-to-head comparisons. Repetitions of the same treatment comparison are grouped together.
mb.comparisons(data)
mb.comparisons(data)
data |
A data frame containing variables |
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.
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)
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)
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.
mb.make.contrast(network, datatype = NULL, format = "wide")
mb.make.contrast(network, datatype = NULL, format = "wide")
network |
An object of class |
datatype |
A string indicating the data type. Can be |
format |
A string indicating the data format. Can be |
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
# 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")
# 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 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).
mb.nodesplit.comparisons(network)
mb.nodesplit.comparisons(network)
network |
An object of class |
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
.
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.
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/.
# Create mb.network object network <- mb.network(osteopain) # Identify comparisons informed by direct and indirect evidence mb.nodesplit.comparisons(network)
# Create mb.network object network <- mb.network(osteopain) # Identify comparisons informed by direct and indirect evidence mb.nodesplit.comparisons(network)
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).
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, ... )
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, ... )
network |
An object of class |
fun |
An object of class |
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 |
link |
Can take either |
sdscale |
Logical object to indicate whether to write a model that specifies a reference SD
for standardising when modelling using Standardised Mean Differences. Specifying |
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. |
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:
|
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). |
corparam |
A boolean object that indicates whether correlation should be modeled
between relative effect time-course parameters. Default is |
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
( |
UME |
Can take either |
pd |
Can take either:
|
parallel |
A boolean value that indicates whether JAGS should be run 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.thin |
thinning rate. Must be a positive integer. Set |
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 |
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 |
... |
Arguments to be sent to R2jags. |
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.
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
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()
When modelling correlation between observations using rho
, values for rho
must imply a
positive semidefinite covariance matrix.
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.
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.
# 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")
# 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
mb.update( mbnma, param = "theta", n.iter = mbnma$BUGSoutput$n.iter, n.thin = mbnma$BUGSoutput$n.thin )
mb.update( mbnma, param = "theta", n.iter = mbnma$BUGSoutput$n.iter, n.thin = mbnma$BUGSoutput$n.thin )
mbnma |
An S3 object of class |
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 |
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 |
n.thin |
The thinning rate. Must be a positive integer. Default is the value used in |
A data frame containing posterior means for the specified param
at each observation, arm and study.
# 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)
# 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
mb.validate.data(data.ab, single.arm = FALSE, CFB = TRUE)
mb.validate.data(data.ab, single.arm = FALSE, CFB = TRUE)
data.ab |
A data frame of arm-level data in "long" format containing the columns:
|
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 ( |
CFB |
A boolean object to indicate if the dataset is composed of studies measuring change from
baseline ( |
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
An error or warnings if checks are not passed. Runs silently if checks are passed
Writes JAGS code for a Bayesian time-course model for model-based network meta-analysis (MBNMA).
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 )
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 )
fun |
An object of class |
link |
Can take either |
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 |
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. |
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:
|
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). |
corparam |
A boolean object that indicates whether correlation should be modeled
between relative effect time-course parameters. Default is |
sdscale |
Logical object to indicate whether to write a model that specifies a reference SD
for standardising when modelling using Standardised Mean Differences. Specifying |
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
( |
UME |
Can take either |
A single long character string containing the JAGS model generated based on the arguments passed to the function.
# 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"))
# 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
nma.run( data.ab, treatments = NULL, method = "common", link = "identity", sdscale = FALSE, ... )
nma.run( data.ab, treatments = NULL, method = "common", link = "identity", sdscale = FALSE, ... )
data.ab |
A data frame of arm-level data in "long" format containing the columns:
|
treatments |
A vector of treatment names. If left as |
method |
Can take |
link |
Can take either |
sdscale |
Logical object to indicate whether to write a model that specifies a reference SD
for standardising when modelling using Standardised Mean Differences. Specifying |
... |
Options for plotting in |
Returns an object of class("nma", "rjags")
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")
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")
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.
obesityBW_CFB
obesityBW_CFB
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
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
.
Pfizer Ltd.
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.
osteopain
osteopain
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
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
.
Pfizer Ltd.
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).
pDcalc( obs1, obs2, fups = NULL, narm, NS, theta.result, resdev.result, likelihood = "normal", type = "time" )
pDcalc( obs1, obs2, fups = NULL, narm, NS, theta.result, resdev.result, likelihood = "normal", type = "time" )
obs1 |
A matrix (study x arm) or array (study x arm x time point) containing
observed data for |
obs2 |
A matrix (study x arm) or array (study x arm x time point) containing
observed data for |
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 |
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:
|
type |
The type of MBNMA model fitted. Can be either |
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).
A numeric value for the effective number of parameters, pD, calculated via the plugin method
TO ADD pV REF
# 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 )
# 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 )
mb.network
objectCreates an object of class("mb.network")
. Various MBNMA functions can subsequently be applied
to this object.
## 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")
## 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")
x |
An object of class |
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 |
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 |
level |
A string indicating whether nodes/facets should represent |
remove.loops |
A boolean value indicating whether to include loops that indicate comparisons within a node. |
v.color |
Can take either |
v.scale |
A number with which to scale the size of the nodes. If the variable |
layout |
An igraph layout specification. This is a function specifying an igraph
layout that determines the arrangement of the vertices (nodes). The default
|
legend |
A boolean value indicating whether or not to plot a legend with class names if |
legend.x |
Can be either a string or a numerical x-coordinate indicating where the legend should be
plotted (see |
legend.y |
A numerical y-coordinate indicating where the legend should be plotted - only required if |
... |
Options for plotting in |
data.ab |
A data frame of arm-level data in "long" format containing the columns:
|
reference |
A number or character (depending on the format of |
cfb |
A logical vector whose length is equal to the unique number of studies in |
description |
Optional. Short description of the network. |
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.
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.
plot(mb.network)
: Generate a network plot
# 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)
# 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
## 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, ... )
## 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, ... )
x |
An object of class |
disp.obs |
A boolean object to indicate whether to show shaded sections
of the plot for where there is observed data ( |
overlay.ref |
A boolean object indicating whether to overlay a line
showing the median network reference treatment response over time on the
plot (
|
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 |
method |
Can take |
col |
A character indicating the colour to use for shading if |
max.col.scale |
Rarely requires adjustment. The maximum count of
observations (therefore the darkest shaded color) only used if |
treat.labs |
A vector of treatment labels in the same order as treatment codes.
Easiest to use treatment labels stored by |
plot.bins |
Plot time bin boundaries as vertical dashed lines. Setting |
... |
Arguments for |
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
.
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.
# 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
# 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
## S3 method for class 'mb.rank' plot(x, treat.labs = NULL, ...)
## S3 method for class 'mb.rank' plot(x, treat.labs = NULL, ...)
x |
An object of class |
treat.labs |
A vector of treatment labels in the same order as treatment codes.
Easiest to use treatment labels stored by |
... |
Arguments to be sent to |
A histogram that shows rankings for each treatment/agent/prediction.
The object returned is an object of class c("gg", "ggplot")
.
# 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)
# 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)
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()
## S3 method for class 'mbnma' plot(x, params = NULL, treat.labs = NULL, class.labs = NULL, ...)
## S3 method for class 'mbnma' plot(x, params = NULL, treat.labs = NULL, class.labs = NULL, ...)
x |
An S3 object of class |
params |
A character vector of time-course parameters to plot.
Parameters must be given the same name as monitored nodes in |
treat.labs |
A character vector of treatment labels. If left as |
class.labs |
A character vector of class labels if |
... |
Arguments to be sent to |
A forest plot of class c("gg", "ggplot")
that has separate panels for different time-course parameters
# 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")
# 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")
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).
## 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", ... )
## 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", ... )
x |
An object of |
plot.type |
A character string that can take the value of |
params |
A character vector corresponding to a time-course parameter(s) for which to plot results.
If left as |
... |
Arguments to be sent to |
network |
An object of class |
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: |
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 |
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 |
lim |
Specifies calculation of either 95% credible intervals ( |
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.
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
plot(nodesplit)
: Plot outputs from nodesplit models
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/.
# 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)
# 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)
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).
## 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", ... )
## 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", ... )
object |
An S3 object of |
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 |
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.
|
treats |
A character vector of treatment/class names or a numeric vector of treatment/class codes (as coded
in |
level |
Can take either |
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
|
synth |
A character object that can take the value |
lim |
Specifies calculation of either 95% credible intervals ( |
... |
Arguments to be sent to R2jags for synthesis of the network
reference treatment effect (using |
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.
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
# 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)
# 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
## S3 method for class 'mb.network' print(x, ...)
## S3 method for class 'mb.network' print(x, ...)
x |
An object of class |
... |
further arguments passed to or from other methods |
Prints the contents of x
to the console.
Print summary information from an mb.predict object
## S3 method for class 'mb.predict' print(x, ...)
## S3 method for class 'mb.predict' print(x, ...)
x |
An object of |
... |
further arguments passed to or from other methods |
Prints a summary of rankings for each parameter
## S3 method for class 'mb.rank' print(x, ...)
## S3 method for class 'mb.rank' print(x, ...)
x |
An object of class |
... |
further arguments passed to or from other methods |
Prints summary details of treatment rankings to the console
Prints basic results from a node-split to the console
## S3 method for class 'nodesplit' print(x, groupby = "time.param", ...)
## S3 method for class 'nodesplit' print(x, groupby = "time.param", ...)
x |
An object of class |
groupby |
A character object that can take the value |
... |
arguments to be sent to |
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
## S3 method for class 'relative.array' print(x, digits = 2, ...)
## S3 method for class 'relative.array' print(x, digits = 2, ...)
x |
An object of class |
digits |
An integer indicating the number of significant digits to be used. |
... |
further arguments passed to |
Prints a league table of treatment effects to the console
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.
radian.rescale(x, start = 0, direction = 1)
radian.rescale(x, start = 0, direction = 1)
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 |
A numeric vector of rescaled values
https://gist.github.com/kjhealy/834774/a4e677401fd6e4c319135dabeaf9894393f9392c
Set rank as a method
rank(x, ...)
rank(x, ...)
x |
An object on which to apply the rank method |
... |
Arguments to be passed to methods |
Uses the rank method
Rank predictions at a specific time point
## S3 method for class 'mb.predict' rank( x, time = max(x$summary[[1]]$time), lower_better = FALSE, treats = names(x$summary), ... )
## S3 method for class 'mb.predict' rank( x, time = max(x$summary[[1]]$time), lower_better = FALSE, treats = names(x$summary), ... )
x |
an object of |
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 |
lower_better |
Indicates whether negative responses are better ( |
treats |
A character vector of treatment/class names for which responses have been predicted
in |
... |
Arguments to be passed to methods |
Returns an object of class("mb.rank")
containing ranked predictions
# 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)
# 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)
Ranks desired parameters saved from a time-course MBNMA model from "best" to "worst".
## S3 method for class 'mbnma' rank( x, param = "auc", lower_better = FALSE, treats = NULL, int.range = NULL, n.iter = x$BUGSoutput$n.sims, ... )
## S3 method for class 'mbnma' rank( x, param = "auc", lower_better = FALSE, treats = NULL, int.range = NULL, n.iter = x$BUGSoutput$n.sims, ... )
x |
An S3 object of |
param |
A character object containing any model parameter monitored
in |
lower_better |
Indicates whether negative responses are better ( |
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 |
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 |
n.iter |
The number of iterations for which to calculate AUC (if |
... |
Arguments to be sent to |
"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.
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.
# 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))
# 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
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, ... )
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, ... )
mbnma |
An S3 object of class |
lower_better |
Indicates whether negative responses are better ( |
treats |
A character vector of treatment/class names (depending on the value of |
level |
Can take either |
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 |
n.iter |
The number of iterations for which to calculate AUC (if |
subdivisions |
The number of subdivisions over which to integrate (see |
... |
Arguments to be sent to R2jags for synthesis of the network
reference treatment effect (using |
"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.
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 contrasts relative to each study reference within a network. Repetitions of the same treatment comparison are grouped together.
ref.comparisons(data)
ref.comparisons(data)
data |
A data frame containing variables |
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.
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.
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, ... )
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, ... )
data.ab |
A data frame of arm-level data in "long" format containing the columns:
|
mbnma |
An S3 object of class |
synth |
A character object that can take the value |
link |
Can take either |
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.thin |
thinning rate. Must be a positive integer. Set
|
n.chains |
number of Markov chains (default: 3) |
... |
Arguments to be sent to R2jags for synthesis of the network
reference treatment effect (using |
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.
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.
# 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")
# 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")
Ensures ref.resp
takes the correct form to allow for synthesis of network
reference treatment effect if data is provided for meta-analysis
ref.validate(data.ab)
ref.validate(data.ab)
data.ab |
A data frame of arm-level data in "long" format containing the columns:
|
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
remove.loops(model)
remove.loops(model)
model |
A character object of JAGS MBNMA model code |
A character vector of JAGS MBNMA model code that has had empty loops removed from it
Identical to replace.prior()
in MBNMAdose.
replace.prior(priors, model = NULL, mbnma = NULL)
replace.prior(priors, model = NULL, mbnma = NULL)
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 |
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).
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
## S3 method for class 'mb.network' summary(object, ...)
## S3 method for class 'mb.network' summary(object, ...)
object |
An object of class |
... |
further arguments passed to or from other methods |
Prints summary details of x
to the console.
Prints a summary table of the mean of MCMC iterations at each time point for each treatment
## S3 method for class 'mb.predict' summary(object, ...)
## S3 method for class 'mb.predict' summary(object, ...)
object |
An object of class |
... |
further arguments passed to or from other methods |
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.
# 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)
# 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
## S3 method for class 'mbnma' summary(object, ...)
## S3 method for class 'mbnma' summary(object, ...)
object |
An S3 object of |
... |
further arguments passed to |
Prints summary details of the model results to the console
Takes node-split results and produces summary data frame
## S3 method for class 'nodesplit' summary(object, ...)
## S3 method for class 'nodesplit' summary(object, ...)
object |
An object of class |
... |
further arguments passed to or from other methods |
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)
** 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 **
temax( pool.emax = "rel", method.emax = "common", pool.et50 = "rel", method.et50 = "common", pool.hill = NULL, method.hill = NULL, p.expon = FALSE )
temax( pool.emax = "rel", method.emax = "common", pool.et50 = "rel", method.et50 = "common", pool.hill = NULL, method.hill = NULL, p.expon = FALSE )
pool.emax |
Pooling for Emax parameter. Can take |
method.emax |
Method for synthesis of Emax parameter. Can take |
pool.et50 |
Pooling for ET50 parameter. Can take |
method.et50 |
Method for synthesis of ET50 parameter. Can take |
pool.hill |
Pooling for Hill parameter. Can take |
method.hill |
Method for synthesis of Hill parameter. Can take |
p.expon |
Should parameters that can only take positive values be modeled on the exponential scale ( |
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:
With Hill parameter:
An object of class("timefun")
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()
.
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/.
# 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")
# 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")
As first described for use in Network Meta-Analysis by Jansen et al. (2015).
tfpoly( degree = 1, pool.1 = "rel", method.1 = "common", pool.2 = "rel", method.2 = "common", method.power1 = 0, method.power2 = 0 )
tfpoly( degree = 1, pool.1 = "rel", method.1 = "common", pool.2 = "rel", method.2 = "common", method.power1 = 0, method.power2 = 0 )
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 |
method.1 |
Method for synthesis of the 1st fractional polynomial coefficient. Can take |
pool.2 |
Pooling for the 2nd fractional polynomial coefficient. Can take |
method.2 |
Method for synthesis of the 2nd fractional polynomial coefficient. Can take |
method.power1 |
Value for the 1st fractional polynomial power. Must take any numeric value in the set |
method.power2 |
Value for the 2nd fractional polynomial power. Must take any numeric value in the set |
represents the 1st coefficient.
represents the 2nd coefficient.
represents the 1st power
represents the 2nd power
For a polynomial of degree=1
:
For a polynomial of degree=2
:
is a regular power except where
, where
.
If a fractional polynomial power
repeats within the function it is multiplied by another
.
An object of class("timefun")
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()
.
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.
# 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)
# 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
timeplot(network, level = "treatment", plotby = "arm", link = "identity", ...)
timeplot(network, level = "treatment", plotby = "arm", link = "identity", ...)
network |
An object of class |
level |
A string indicating whether nodes/facets should represent |
plotby |
A character object that can take either |
link |
Can take either |
... |
Arguments to be sent to |
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.
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()
.
# 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")
# 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")
Similar parameterisation to the Emax model but with non-asymptotic maximal effect (Emax). Proposed by proposed by Fu and Manner (2010)
titp( pool.emax = "rel", method.emax = "common", pool.rate = "rel", method.rate = "common", p.expon = FALSE )
titp( pool.emax = "rel", method.emax = "common", pool.rate = "rel", method.rate = "common", p.expon = FALSE )
pool.emax |
Pooling for exponential Emax parameter. Can take |
method.emax |
Method for synthesis of exponential Emax parameter. Can take |
pool.rate |
Pooling for parameter controlling rate of onset. Default is |
method.rate |
Method for synthesis of parameter controlling rate of onset. Can take |
p.expon |
Should parameters that can only take positive values be modeled on the exponential scale ( |
An object of class("timefun")
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. |
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/.
titp(pool.emax="rel", method.emax="random") titp(pool.emax="abs")
titp(pool.emax="rel", method.emax="random") titp(pool.emax="abs")
tloglin(pool.rate = "rel", method.rate = "common")
tloglin(pool.rate = "rel", method.rate = "common")
pool.rate |
Pooling for rate parameter. Can take |
method.rate |
Method for synthesis of rate parameter. Can take |
An object of class("timefun")
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. |
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/.
tloglin(pool.rate="rel", method.rate="random") tloglin(pool.rate="abs")
tloglin(pool.rate="rel", method.rate="random") tloglin(pool.rate="abs")
Polynomial time-course function
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" )
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" )
degree |
The degree of the polynomial - e.g. |
pool.1 |
Pooling for the 1st polynomial coefficient. Can take |
method.1 |
Method for synthesis of the 1st polynomial coefficient.Can take |
pool.2 |
Pooling for the 2nd polynomial coefficient. Can take |
method.2 |
Method for synthesis of the 2nd polynomial coefficient. Can take |
pool.3 |
Pooling for the 3rd polynomial coefficient. Can take |
method.3 |
Method for synthesis of the 3rd polynomial coefficient. Can take |
pool.4 |
Pooling for the 4th polynomial coefficient. Can take |
method.4 |
Method for synthesis of the 4th polynomial coefficient. Can take |
represents the 1st coefficient.
represents the 2nd coefficient.
represents the 3rd coefficient.
represents the 4th coefficient.
Linear model:
Quadratic model:
Cubic model:
Quartic model:
An object of class("timefun")
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()
.
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/.
# 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")
# 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")
Used to fit B-splines, natural cubic splines, and piecewise linear splines(Perperoglu et al. 2019).
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" )
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" )
type |
The type of spline. Can take |
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. |
pool.1 |
Pooling for the 1st coefficient. Can take |
method.1 |
Method for synthesis of the 1st coefficient. Can take |
pool.2 |
Pooling for the 2nd coefficient. Can take |
method.2 |
Method for synthesis of the 2nd coefficient. Can take |
pool.3 |
Pooling for the 3rd coefficient. Can take |
method.3 |
Method for synthesis of the 3rd coefficient. Can take |
pool.4 |
Pooling for the 4th coefficient. Can take |
method.4 |
Method for synthesis of the 4th coefficient. Can take |
An object of class("timefun")
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()
.
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.
# 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")
# 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
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" )
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" )
fun |
A formula specifying any relationship including |
pool.1 |
Pooling for |
method.1 |
Method for synthesis of |
pool.2 |
Pooling for |
method.2 |
Method for synthesis of |
pool.3 |
Pooling for |
method.3 |
Method for synthesis of |
pool.4 |
Pooling for |
method.4 |
Method for synthesis of |
An object of class("timefun")
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()
.
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/.
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")
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
write.beta(model, timecourse, fun, UME, class.effect)
write.beta(model, timecourse, fun, UME, class.effect)
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 |
UME |
Can take either |
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
( |
A character vector of JAGS MBNMA model code that includes beta parameter components of the model
Checks validity of arguments for mb.write
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() )
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() )
fun |
An object of class |
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 |
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. |
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:
|
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). |
link |
Can take either |
sdscale |
Logical object to indicate whether to write a model that specifies a reference SD
for standardising when modelling using Standardised Mean Differences. Specifying |
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
( |
UME |
Can take either |
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.
A boolean object that indicates whether the arguments imply modelling correlation between time points.
This uses a Wishart prior as default for modelling the correlation
write.cor(model, fun, omega = NULL, class.effect = list())
write.cor(model, fun, omega = NULL, class.effect = list())
model |
A character object of JAGS MBNMA model code |
fun |
An object of class |
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). |
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
( |
Adds sections of JAGS code for an MBNMA model that correspond to the likelihood
write.likelihood( model, timecourse, rho = 0, covar = "varadj", link = "identity", sdscale = FALSE, fun )
write.likelihood( model, timecourse, rho = 0, covar = "varadj", link = "identity", sdscale = FALSE, fun )
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. |
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:
|
link |
Can take either |
sdscale |
Logical object to indicate whether to write a model that specifies a reference SD
for standardising when modelling using Standardised Mean Differences. Specifying |
fun |
An object of class |
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
write.model()
write.model()
A character vector of JAGS model code
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.
write.ref.synth( fun = tpoly(degree = 1), link = "identity", positive.scale = TRUE, intercept = TRUE, rho = 0, covar = "varadj", mu.synth = "random", priors = NULL )
write.ref.synth( fun = tpoly(degree = 1), link = "identity", positive.scale = TRUE, intercept = TRUE, rho = 0, covar = "varadj", mu.synth = "random", priors = NULL )
fun |
An object of class |
link |
Can take either |
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 |
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. |
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:
|
mu.synth |
A string that takes the value |
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)). |
A character object of JAGS MBNMA model code that includes beta parameter components of the model
# 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
# 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
write.timecourse(model, fun, intercept, positive.scale)
write.timecourse(model, fun, intercept, positive.scale)
model |
A character string representing the MBNMA model in JAGS code |
fun |
An object of class |
intercept |
A boolean object that indicates whether an intercept (written
as |
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 |
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