add_integration {multinma}

Add numerical integration points to aggregate data

Description

The add_integration() generic creates Quasi-Monte Carlo numerical integration points using a Gaussian copula and Sobol' sequences, as described in Phillippo et al. (2020). Methods are available for networks stored in nma_data objects, and for data frames. The function unnest_integration() unnests integration points stored in a data frame, to aid plotting or other exploration.

Usage

add_integration(x, ...)

## Default S3 method:
add_integration(x, ...)

## S3 method for class 'data.frame'
add_integration(
  x,
  ...,
  cor = NULL,
  cor_adjust = NULL,
  n_int = 64L,
  int_args = list()
)

## S3 method for class 'nma_data'
add_integration(
  x,
  ...,
  cor = NULL,
  cor_adjust = NULL,
  n_int = 64L,
  int_args = list()
)

unnest_integration(data)

Arguments

x

An nma_data object, as created by the ⁠set_*()⁠ functions or combine_network(), or data frame

...

Distributions for covariates, see "Details"

cor

Correlation matrix to use for generating the integration points. By default, this takes a weighted correlation matrix from all IPD studies. Rows and columns should match the order of covariates specified in ....

cor_adjust

Adjustment to apply to the correlation matrix given by cor (or computed from the IPD if cor = NULL) to obtain the Gaussian copula correlations, either "spearman", "pearson", or "none", see "Details". The default when cor = NULL is "spearman", otherwise the default is "pearson".

n_int

Number of integration points to generate, default 64. Powers of 2 are recommended, as these are expected to be particularly efficient for QMC integration.

int_args

A named list of arguments to pass to sobol()

data

Data frame with nested integration points, stored in list columns as ⁠.int_<variable name>⁠

Details

The arguments passed to ... specify distributions for the covariates. Argument names specify the name of the covariate, which should match a covariate name in the IPD (if IPD are present). The required marginal distribution is then specified using the function distr().

The argument cor_adjust specifies how the correlation matrix given by cor (or computed from the IPD if cor = NULL) is adjusted to obtain the correlation matrix for the Gaussian copula, using the formulae in Xiao and Zhou (2018).

• cor_adjust = "spearman" should be used when the correlations cor have been computed using Spearman's rank correlation. Correlations between continuous covariates will be reproduced exactly by the integration points. Correlations between discrete covariates will be reproduced approximately. This is the default when cor = NULL and correlations are calculated from the IPD studies.

• cor_adjust = "pearson" should be used when the correlations cor have been computed using Pearson's product-moment correlation. Correlations between Normal covariates will be reproduced exactly by the integration points, all others will be reproduced approximately. Correlations between discrete covariates will be reproduced approximately (and identically to cor_adjust = "spearman"). This is the default when cor is provided by the user, since cor() defaults to method = "pearson" and Pearson correlations are most likely reported in published data. However, we recommend providing Spearman correlations (e.g. from cor(., method = "spearman")) and using cor_adjust = "spearman" where possible.

• cor_adjust = "none" allows the user to specify the correlation matrix for the Gaussian copula directly; no adjustment is applied.

• cor_adjust = "legacy" is also available, which reproduces exactly the behaviour from version 0.3.0 and earlier. This is similar to cor_adjust = "none", but unadjusted Spearman correlations are used if cor = NULL.

When adding integration points to a network object the correlation matrix used is stored in ⁠$int_cor⁠, and the copula correlation matrix and adjustment used are stored as attributes of ⁠$int_cor⁠. If this correlation matrix is passed again to add_integration() (e.g. to reuse the correlations for an external target population) this will be detected, and the correct setting for cor_adjust will automatically be applied.

Value

For the nma_data method, an object of class nma_data. For the data.frame method, the input data frame is returned (as a tibble) with an added column for each covariate (prefixed with ".int_"), containing the numerical integration points nested as length-n_int vectors within each row. For unnest_integration(), a data frame with integration points unnested.

References

Phillippo DM, Dias S, Ades AE, Belger M, Brnabic A, Schacht A, Saure D, Kadziola Z, Welton NJ (2020). “Multilevel Network Meta-Regression for population-adjusted treatment comparisons.” Journal of the Royal Statistical Society: Series A (Statistics in Society), 183(3), 1189–1210. doi:10.1111/rssa.12579.

Xiao Q, Zhou S (2018). “Matching a correlation coefficient by a Gaussian copula.” Communications in Statistics - Theory and Methods, 48(7), 1728–1747. doi:10.1080/03610926.2018.1439962.

Examples

## Plaque psoriasis ML-NMR - network setup and adding integration points
# Set up plaque psoriasis network combining IPD and AgD
library(dplyr)
pso_ipd <- filter(plaque_psoriasis_ipd,
                  studyc %in% c("UNCOVER-1", "UNCOVER-2", "UNCOVER-3"))

pso_agd <- filter(plaque_psoriasis_agd,
                  studyc == "FIXTURE")

head(pso_ipd)
head(pso_agd)

pso_ipd <- pso_ipd %>%
  mutate(# Variable transformations
    bsa = bsa / 100,
    prevsys = as.numeric(prevsys),
    psa = as.numeric(psa),
    weight = weight / 10,
    durnpso = durnpso / 10,
    # Treatment classes
    trtclass = case_when(trtn == 1 ~ "Placebo",
                         trtn %in% c(2, 3, 5, 6) ~ "IL blocker",
                         trtn == 4 ~ "TNFa blocker"),
    # Check complete cases for covariates of interest
    complete = complete.cases(durnpso, prevsys, bsa, weight, psa)
  )

pso_agd <- pso_agd %>%
  mutate(
    # Variable transformations
    bsa_mean = bsa_mean / 100,
    bsa_sd = bsa_sd / 100,
    prevsys = prevsys / 100,
    psa = psa / 100,
    weight_mean = weight_mean / 10,
    weight_sd = weight_sd / 10,
    durnpso_mean = durnpso_mean / 10,
    durnpso_sd = durnpso_sd / 10,
    # Treatment classes
    trtclass = case_when(trtn == 1 ~ "Placebo",
                         trtn %in% c(2, 3, 5, 6) ~ "IL blocker",
                         trtn == 4 ~ "TNFa blocker")
  )

# Exclude small number of individuals with missing covariates
pso_ipd <- filter(pso_ipd, complete)

pso_net <- combine_network(
  set_ipd(pso_ipd,
          study = studyc,
          trt = trtc,
          r = pasi75,
          trt_class = trtclass),
  set_agd_arm(pso_agd,
              study = studyc,
              trt = trtc,
              r = pasi75_r,
              n = pasi75_n,
              trt_class = trtclass)
)

# Print network details
pso_net

# Add integration points to the network
pso_net <- add_integration(pso_net,
  durnpso = distr(qgamma, mean = durnpso_mean, sd = durnpso_sd),
  prevsys = distr(qbern, prob = prevsys),
  bsa = distr(qlogitnorm, mean = bsa_mean, sd = bsa_sd),
  weight = distr(qgamma, mean = weight_mean, sd = weight_sd),
  psa = distr(qbern, prob = psa),
  n_int = 64)


## Adding integration points to a data frame, e.g. for prediction
# Define a data frame of covariate summaries
new_agd_int <- data.frame(
  bsa_mean = 0.6,
  bsa_sd = 0.3,
  prevsys = 0.1,
  psa = 0.2,
  weight_mean = 10,
  weight_sd = 1,
  durnpso_mean = 3,
  durnpso_sd = 1)

# Adding integration points, using the weighted average correlation matrix
# computed for the plaque psoriasis network
new_agd_int <- add_integration(new_agd_int,
  durnpso = distr(qgamma, mean = durnpso_mean, sd = durnpso_sd),
  prevsys = distr(qbern, prob = prevsys),
  bsa = distr(qlogitnorm, mean = bsa_mean, sd = bsa_sd),
  weight = distr(qgamma, mean = weight_mean, sd = weight_sd),
  psa = distr(qbern, prob = psa),
  cor = pso_net$int_cor,
  n_int = 64)

# Here, since we reused the correlation matrix pso_net$int_cor from the
# network, the correct setting of cor_adjust = "spearman" is automatically
# applied

new_agd_int