Integrative Bayesian Multi-Platform
Regression with IntegMultiReg
Sinian Zhang, Thierry Chekouo
Overview
IntegMultiReg implements the integrative
multi-regression (IMR) model of Chekouo, Stingo, Doecke and Do
(2017) and extends it from time-to-event outcomes to continuous
(Gaussian) and binary (probit) outcomes. The central idea is to
integrate all available subjects across several
molecular platforms, even when many subjects are measured on only a
subset of the platforms.
Given \(K\) platforms, subjects are
partitioned into the \(2^K-1\)
non-empty availability subgroups of a \(K\)-set Venn diagram. One regression model
is built per subgroup, and information is shared across availability
subgroups through
- non-local (product moment) priors on the regression
coefficients, which sharply penalise small effects and yield
parsimonious models, and
- a Markov random field (MRF) prior on the
variable-selection indicators, which encourages the same biomarkers to
be selected across availability subgroups.
library(IntegMultiReg)
data("simIMR", package = "IntegMultiReg")
sapply(simIMR$platforms, dim)
#> genomic proteomic metabolomic
#> [1,] 240 180 180
#> [2,] 21 11 9
The bundled simIMR data set has three platforms
(genomic, proteomic, metabolomic)
measured on overlapping but partially missing sets of subjects, clinical
covariates for everyone, and three outcome types generated from the same
latent signal. The truly associated features are recorded in
simIMR$truth.
Real TCGA-KIRC
example
The package also includes kircIMR, a reduced real-data
example derived from public UCSC Xena TCGA kidney renal clear cell
carcinoma (KIRC) sampleMap files, not from controlled-access TCGA/GDC
files. It mirrors the case-study structure of Chekouo et al. (2017):
mRNA expression, miRNA expression and DNA methylation measured on
overlapping patient sets, clinical covariates, and a right-censored
survival outcome. Patient IDs are package-internal labels such as
KIRC001; the package does not distribute TCGA barcodes, and
users should not attempt participant re-identification or linkage to
external resources.
data("kircIMR", package = "IntegMultiReg")
sapply(kircIMR$platforms, dim)
#> mrna mirna methylation
#> [1,] 521 235 312
#> [2,] 51 31 51
kircIMR$model_subgroup_sizes
#> bitstrings
#> 011 101 001 111
#> 63 139 147 172
The original Biometrics analysis used a much larger screened panel
(776 mRNA, 91 miRNA and 729 methylation features) and long MCMC chains.
The package object keeps a smaller Cox-screened panel (50, 30 and 50
features) so that examples remain lightweight while retaining the same
scientific data structure. A short survival fit follows the same call
pattern:
kirc_fit <- imr(
platform_data_list = kircIMR$platforms,
outcome = kircIMR$outcome.survival,
cov = kircIMR$covariates,
type_outcome = "right.censored",
nu = c(-4, -3, -4),
sample_mcmc = c(4000, 1000),
ssize = 30,
seed = 1
)
summary(kirc_fit)
plot_top_features(kirc_fit, top = 12)
Fitting the model
The workhorse is imr(). Here we fit a binary (probit)
outcome with a short MCMC run for illustration; in practice use longer
chains.
fit <- imr(
platform_data_list = simIMR$platforms,
outcome = simIMR$outcome.binary,
cov = simIMR$covariates,
type_outcome = "binary",
nu = c(-4, -3, -4),
sample_mcmc = c(1500, 500),
ssize = 30,
seed = 1
)
fit
#> Integrative Bayesian Multi-Platform Regression (IMR)
#> ----------------------------------------------------
#> Call:
#> imr(platform_data_list = simIMR$platforms, outcome = simIMR$outcome.binary,
#> cov = simIMR$covariates, type_outcome = "binary", ssize = 30,
#> nu = c(-4, -3, -4), sample_mcmc = c(1500, 500), seed = 1)
#>
#> Outcome type : binary
#> Method : IMR
#> Platforms : 3 (genomic, proteomic, metabolomic)
#> MCMC : 1500 retained draws after 500 burn-in
#>
#> Platform key:
#> P1 = genomic
#> P2 = proteomic
#> P3 = metabolomic
#>
#> Availability subgroups modelled (bitstring : platforms : size):
#> 011 : P1 + P2 : 120
#> 100 : P3 : 60
#> 101 : P1 + P3 : 60
#> 111 : P1 + P2 + P3 : 60
#>
#> Features with mPIP > 0.50 (in any subgroup):
#> genomic : 2 of 20
#> proteomic : 2 of 10
#> metabolomic : 2 of 8
The four modelled availability subgroups correspond to the
availability patterns 011 (genomic + proteomic),
111 (all three), 101 (genomic + metabolomic)
and 100 (metabolomic only); the bitstring digits run from
the first platform (right) to the last (left).
plot_subgroup_sizes() shows how the sample splits across
them.
Inspecting selected
biomarkers
summary() ranks, per platform, the features whose
marginal posterior inclusion probability (mPIP) exceeds a threshold.
summary(fit, threshold = 0.5)
#> Integrative Bayesian Multi-Platform Regression (IMR) -- summary
#> --------------------------------------------------------------
#> Outcome type : binary Method: IMR
#> Selection threshold (mPIP) : 0.50
#>
#> Platform 'genomic': 2 selected feature(s)
#> feature max_mpip subgroup
#> G02 1.000 011
#> G01 0.892 111
#>
#> Platform 'proteomic': 2 selected feature(s)
#> feature max_mpip subgroup
#> P01 0.956 011
#> P03 0.698 111
#>
#> Platform 'metabolomic': 2 selected feature(s)
#> feature max_mpip subgroup
#> M03 0.855 101
#> M01 0.711 100
The posterior inclusion probabilities themselves are available
through coef() (a list of subgroup \(\times\) feature matrices) and can be
visualised either as per-platform heatmaps or as a single ranked bar
chart of the most strongly supported features across all platforms.
plot(fit, type = "selection", platform = 1)
plot_top_features(fit, top = 8)
The estimated MRF interaction parameters, which measure how strongly
each pair of subgroups shares biomarkers, are shown with
type = "theta", and the log-posterior trace with
type = "trace".
plot(fit, type = "trace")
Predicting new
subjects
predict() routes new subjects to the appropriate
availability subgroup, standardises their features with the training
scaling, and performs Bayesian model averaging. When the new platforms
are supplied in the same order as at training,
platform_names can be omitted. Binary outcomes are returned
on the probability scale.
new_x <- simIMR$platforms$genomic[1:20, ]
new_p <- simIMR$platforms$proteomic[1:20, ]
pred <- predict(fit, newdata = list(new_x, new_p),
covariates = simIMR$covariates)
head(pred[["model:011"]])
#> id predict
#> 1 1 0.380
#> 2 2 0.388
#> 3 3 0.703
#> 4 4 0.327
#> 5 5 0.737
#> 6 6 0.250
Other outcome
types
The same interface handles continuous and right-censored outcomes;
only type_outcome and the outcome data frame
change.
# continuous (Gaussian) outcome
fit_c <- imr(
simIMR$platforms, simIMR$outcome.continuous, cov = simIMR$covariates,
type_outcome = "continuous", nu = c(-4, -3, -4),
sample_mcmc = c(1500, 500), ssize = 30, seed = 1)
# right-censored survival outcome
fit_s <- imr(
simIMR$platforms, simIMR$outcome.survival, cov = simIMR$covariates,
type_outcome = "right.censored", nu = c(-4, -3, -4),
sample_mcmc = c(1500, 500), ssize = 30, seed = 1)
References
Chekouo T, Stingo FC, Doecke JD, Do K-A (2017). “A Bayesian
Integrative Approach for Multi-Platform Genomic Data: A Kidney Cancer
Case Study.” Biometrics, 73(2), 615–624. doi:10.1111/biom.12587