Title: | Unified Cross-Omics Deconvolution |
---|---|
Description: | UNIfied Cross-Omics deconvolution (Unico) deconvolves standard 2-dimensional bulk matrices of samples by features into a 3-dimensional tensors representing samples by features by cell types. Unico stands out as the first principled model-based deconvolution method that is theoretically justified for any heterogeneous genomic data. For more details see Chen and Rahmani et al. (2024) <doi:10.1101/2024.01.27.577588>. |
Authors: | Zeyuan Chen [aut, cre], Elior Rahmani [aut] |
Maintainer: | Zeyuan Chen <[email protected]> |
License: | GPL-3 |
Version: | 0.1.0 |
Built: | 2025-02-09 05:00:52 UTC |
Source: | https://github.com/cozygene/unico |
Performs asymptotic statistical testing on (1) the marginal effect of each covariate in C1
at source-specific level (2) non-source-specific effect for each covariate in C2
. In the context of bulk genomic data containing a mixture of cell types, these correspond to the marginal effect of each covariate in C1
(potentially including the phenotype of interest) at each cell type and tissue-level effect for each covariate in C2
.
association_asymptotic( X, Unico.mdl, slot_name = "asymptotic", diag_only = FALSE, intercept = TRUE, X_max_stds = 2, Q_max_stds = Inf, V_min_qlt = 0.05, parallel = TRUE, num_cores = NULL, log_file = "Unico.log", verbose = FALSE, debug = FALSE )
association_asymptotic( X, Unico.mdl, slot_name = "asymptotic", diag_only = FALSE, intercept = TRUE, X_max_stds = 2, Q_max_stds = Inf, V_min_qlt = 0.05, parallel = TRUE, num_cores = NULL, log_file = "Unico.log", verbose = FALSE, debug = FALSE )
X |
An |
Unico.mdl |
The entire set of model parameters estimated by Unico on the 2D mixture matrix (i.e. the list returned by applying function |
slot_name |
A string indicating the key for storing the results under |
diag_only |
A logical value indicating whether to only use the estimated source-level variances (and thus ignoring the estimate covariance) for controlling the heterogeneity in the observed mixture. if set to FALSE, Unico instead estimates the observation- and feature-specific variance in the mixture by leveraging the entire |
intercept |
A logical value indicating whether to fit the intercept term when performing the statistical testing. |
X_max_stds |
A non-negative numeric value indicating, for each feature, the portions of data that are considered as outliers due to the observed mixture value. Only samples whose observed mixture value fall within |
Q_max_stds |
A non-negative numeric value indicating, for each feature, the portions of data that are considered as outliers due to the estimated mixture variance. Only samples whose estimated mixture variance fall within |
V_min_qlt |
A non-negative numeric value indicating, for each feature, the portions of data that are considered as outliers due to the estimated moment condition variance. This value should be between 0 and 1. Only samples whose estimated moment condition variance fall outside the bottom |
parallel |
A logical value indicating whether to use parallel computing (possible when using a multi-core machine). |
num_cores |
A numeric value indicating the number of cores to use (activated only if |
log_file |
A path to an output log file. Note that if the file |
verbose |
A logical value indicating whether to print logs. |
debug |
A logical value indicating whether to set the logger to a more detailed debug level; set |
Under no distribution assumption, we can solve for the following weighted least square problem, which is similar to the heteroskedastic regression view described in association_parametric.
denotes the design matrix formed by stacking samples in the rows and dependent variables
on the columns.
denotes the corresponding effect sizes on the dependent variables.
denotes the feature-specific weighting scheme. Similar to the parametric counterpart,
, where for each sample
, its corresponding weight will be the inverse of the estimated variance in the mixture:
.
Marginal testing can thus be carried out on each dependent variable via the asymptotic distribution of the estimator
.
An updated Unico.mdl
object with the the following list of effect size and p-value estimates stored in an additional key specified by slot_name
gammas_hat |
An |
betas_hat |
An |
gammas_hat_pvals |
An |
betas_hat_pvals |
An |
Q |
An |
masks |
An |
fphi_hat |
An |
phi_hat |
An |
phi_se |
An |
phi_hat_pvals |
An |
data = simulate_data(n=100, m=2, k=3, p1=1, p2=1, taus_std=0, log_file=NULL) res = list() res$params.hat = Unico(data$X, data$W, data$C1, data$C2, parallel=FALSE, log_file=NULL) res$params.hat = association_asymptotic(data$X, res$params.hat, parallel=FALSE, log_file=NULL)
data = simulate_data(n=100, m=2, k=3, p1=1, p2=1, taus_std=0, log_file=NULL) res = list() res$params.hat = Unico(data$X, data$W, data$C1, data$C2, parallel=FALSE, log_file=NULL) res$params.hat = association_asymptotic(data$X, res$params.hat, parallel=FALSE, log_file=NULL)
Performs parametric statistical testing (T-test) on (1) the marginal effect of each covariate in C1
at source-specific level (2) the joint effect across all sources for each covariate in C1
(3) non-source-specific effect for each covariate in C2
. In the context of bulk genomic data containing a mixture of cell types, these correspond to the marginal effect of each covariate in C1
(potentially including the phenotype of interest) at each cell type, joint tissue-level effect for each covariate in C1
, and tissue-level effect for each covariate in C2
.
association_parametric( X, Unico.mdl, slot_name = "parametric", diag_only = FALSE, intercept = TRUE, X_max_stds = 2, Q_max_stds = Inf, XQ_max_stds = Inf, parallel = TRUE, num_cores = NULL, log_file = "Unico.log", verbose = FALSE, debug = FALSE )
association_parametric( X, Unico.mdl, slot_name = "parametric", diag_only = FALSE, intercept = TRUE, X_max_stds = 2, Q_max_stds = Inf, XQ_max_stds = Inf, parallel = TRUE, num_cores = NULL, log_file = "Unico.log", verbose = FALSE, debug = FALSE )
X |
An |
Unico.mdl |
The entire set of model parameters estimated by Unico on the 2D mixture matrix (i.e. the list returned by applying function |
slot_name |
A string indicating the key for storing the results under |
diag_only |
A logical value indicating whether to only use the estimated source-level variances (and thus ignoring the estimate covariance) for controlling the heterogeneity in the observed mixture. if set to FALSE, Unico instead estimates the observation- and feature-specific variance in the mixture by leveraging the entire |
intercept |
A logical value indicating whether to fit the intercept term when performing the statistical testing. |
X_max_stds |
A non-negative numeric value indicating, for each feature, the portions of data that are considered as outliers due to the observed mixture value. Only samples whose observed mixture value fall within |
Q_max_stds |
A non-negative numeric value indicating, for each feature, the portions of data that are considered as outliers due to the estimated mixture variance. Only samples whose estimated mixture variance fall within |
XQ_max_stds |
A non-negative numeric value indicating, for each feature, the portions of data that are considered as outliers due to the weighted mixture value. Only samples whose weighted mixture value fall within |
parallel |
A logical value indicating whether to use parallel computing (possible when using a multi-core machine). |
num_cores |
A numeric value indicating the number of cores to use (activated only if |
log_file |
A path to an output log file. Note that if the file |
verbose |
A logical value indicating whether to print logs. |
debug |
A logical value indicating whether to set the logger to a more detailed debug level; set |
If we assume that source-specific values are normally distributed, under the Unico model, we have the following:
For a given feature under test, the above equation corresponds to a heteroskedastic regression problem with
as the dependent variable and
as the set of independent variables.
This view allows us to perform parametric statistical testing (T-test for marginal effects and partial F-test for joint effects) by solving a generalized least squares problem with sample
scaled by the inverse of its estimated standard deviation.
An updated Unico.mdl
object with the the following list of effect size and p-value estimates stored in an additional key specified by slot_name
gammas_hat |
An |
betas_hat |
An |
gammas_hat_pvals |
An |
betas_hat_pvals |
An |
gammas_hat_pvals.joint |
An |
Q |
An |
masks |
An |
phi_hat |
An |
phi_se |
An |
phi_hat_pvals |
An |
data = simulate_data(n=100, m=2, k=3, p1=1, p2=1, taus_std=0, log_file=NULL) res = list() res$params.hat = Unico(data$X, data$W, data$C1, data$C2, parallel=FALSE, log_file=NULL) res$params.hat = association_parametric(data$X, res$params.hat, parallel=FALSE, log_file=NULL)
data = simulate_data(n=100, m=2, k=3, p1=1, p2=1, taus_std=0, log_file=NULL) res = list() res$params.hat = Unico(data$X, data$W, data$C1, data$C2, parallel=FALSE, log_file=NULL) res$params.hat = association_parametric(data$X, res$params.hat, parallel=FALSE, log_file=NULL)
Simulate all model parameters and sample source specific data from multivariate gaussian with full covariance structure.
simulate_data( n, m, k, p1, p2, mus_mean = 10, mus_std = 2, gammas_mean = 1, gammas_std = 0.1, betas_mean = 1, betas_std = 0.1, sigmas_lb = 0, sigmas_ub = 1, taus_std = 0.1, log_file = "Unico.log", verbose = FALSE )
simulate_data( n, m, k, p1, p2, mus_mean = 10, mus_std = 2, gammas_mean = 1, gammas_std = 0.1, betas_mean = 1, betas_std = 0.1, sigmas_lb = 0, sigmas_ub = 1, taus_std = 0.1, log_file = "Unico.log", verbose = FALSE )
n |
A positive integer indicating the number of observations to simulate. |
m |
A positive integer indicating the number of features to simulate. |
k |
A positive integer indicating the number of sources to simulate. |
p1 |
A non-negative integer indicating the number of source-specific covariates to simulate. |
p2 |
A non-negative indicating the number of non-source-specific covariates to simulate. |
mus_mean |
A numerical value indicating the average of the source specific means. |
mus_std |
A positive value indicating the variation of the source specific means across difference sources. |
gammas_mean |
A numerical value indicating the average effect sizes of the source-specific covariates. |
gammas_std |
A non-negative numerical value indicating the variation of the effect sizes of the source-specific covariates. |
betas_mean |
A numerical value indicating the average effect sizes of the non-source-specific covariates. |
betas_std |
A non-negative numerical value indicating the variation of the effect sizes of the non-source-specific covariates. |
sigmas_lb |
A numerical value indicating the lower bound of a uniform distribution from which we sample entries of matrix |
sigmas_ub |
A numerical value indicating the upper bound of a uniform distribution from which we sample entries of matrix |
taus_std |
non-negative numerical value indicating the variation of the measurement noise across difference features. |
log_file |
A path to an output log file. Note that if the file |
verbose |
A logical value indicating whether to print logs. |
Simulate data based on the generative model described in function Unico.
A list of simulated model parameters, covariates, observed mixture, and source-specific data.
X |
An |
W |
An |
C1 |
An |
C2 |
An |
Z |
A |
mus |
An |
gammas |
An |
betas |
An |
sigmas |
An |
taus |
An |
data = simulate_data(n=100, m=2, k=3, p1=1, p2=1, taus_std=0, log_file=NULL)
data = simulate_data(n=100, m=2, k=3, p1=1, p2=1, taus_std=0, log_file=NULL)
Infers the underlying (sources by features by observations) 3D tensor from the observed (features by observations) 2D mixture, under the assumption of the Unico model that each observation is a mixture of unique source-specific values (in each feature in the data). In the context of bulk genomics containing a mixture of cell types (i.e. the input could be CpG sites by individuals for DNA methylation and genes by individuals for RNA expression), tensor
allows to estimate the cell-type-specific levels for each individual in each CpG site/gene (i.e. a tensor of CpG sites/genes by individuals by cell types).
tensor( X, W, C1, C2, Unico.mdl, parallel = TRUE, num_cores = NULL, log_file = "Unico.log", verbose = FALSE, debug = FALSE )
tensor( X, W, C1, C2, Unico.mdl, parallel = TRUE, num_cores = NULL, log_file = "Unico.log", verbose = FALSE, debug = FALSE )
X |
An |
W |
An |
C1 |
An |
C2 |
An |
Unico.mdl |
The entire set of model parameters estimated by Unico on the 2D mixture matrix (i.e. the list returned by applying function |
parallel |
A logical value indicating whether to use parallel computing (possible when using a multi-core machine). |
num_cores |
A numeric value indicating the number of cores to use (activated only if |
log_file |
A path to an output log file. Note that if the file |
verbose |
A logical value indicating whether to print logs. |
debug |
A logical value indicating whether to set the logger to a more detailed debug level; set |
After obtaining all the estimated parameters in the Unico model (by calling Unico), tensor
uses the conditional distribution for estimating the
source-specific levels of each sample
at each feature
.
A k
by m
by n
array with the estimated source-specific values. The first axis/dimension in the array corresponds to the different sources.
data = simulate_data(n=100, m=2, k=3, p1=1, p2=1, taus_std=0, log_file=NULL) res = list() res$params.hat = Unico(data$X, data$W, data$C1, data$C2, parallel=FALSE, log_file=NULL) res$Z = tensor(data$X, data$W, data$C1, data$C2, res$params.hat, parallel=FALSE, log_file=NULL)
data = simulate_data(n=100, m=2, k=3, p1=1, p2=1, taus_std=0, log_file=NULL) res = list() res$params.hat = Unico(data$X, data$W, data$C1, data$C2, parallel=FALSE, log_file=NULL) res$Z = tensor(data$X, data$W, data$C1, data$C2, res$params.hat, parallel=FALSE, log_file=NULL)
Fits the Unico model for an input matrix of features by observations that are coming from a mixture of k
sources, under the assumption that each observation is a mixture of unique (unobserved) source-specific values (in each feature in the data). Specifically, for each feature, it standardizes the data and learns the source-specific mean and full k
by k
variance-covariance matrix.
Unico( X, W, C1, C2, fit_tau = FALSE, mean_penalty = 0, var_penalty = 0.01, covar_penalty = 0.01, mean_max_iterations = 2, var_max_iterations = 3, nloptr_opts_algorithm = "NLOPT_LN_COBYLA", max_stds = 2, init_weight = "default", max_u = 1, max_v = 1, parallel = TRUE, num_cores = NULL, log_file = "Unico.log", verbose = FALSE, debug = FALSE )
Unico( X, W, C1, C2, fit_tau = FALSE, mean_penalty = 0, var_penalty = 0.01, covar_penalty = 0.01, mean_max_iterations = 2, var_max_iterations = 3, nloptr_opts_algorithm = "NLOPT_LN_COBYLA", max_stds = 2, init_weight = "default", max_u = 1, max_v = 1, parallel = TRUE, num_cores = NULL, log_file = "Unico.log", verbose = FALSE, debug = FALSE )
X |
An |
W |
An |
C1 |
An |
C2 |
An |
fit_tau |
A logical value indicating whether to fit the standard deviation of the measurement noise (i.e. the i.i.d. component of variation in the model denoted as |
mean_penalty |
A non-negative numeric value indicating the regularization strength on the source-specific mean estimates. |
var_penalty |
A non-negative numeric value indicating the regularization strength on the diagonal entries of the full |
covar_penalty |
A non-negative numeric value indicating the regularization strength on the off diagonal entries of the full |
mean_max_iterations |
A non-negative numeric value indicating the number of iterative updates performed on the mean estimates. |
var_max_iterations |
A non-negative numeric value indicating the number of iterative updates performed on the variance-covariance matrix. |
nloptr_opts_algorithm |
A string indicating the optimization algorithm to use. |
max_stds |
A non-negative numeric value indicating, for each feature, the portions of data that are considered as outliers. Only samples within |
init_weight |
A string indicating the initial weights on the samples to start the iterative optimization. |
max_u |
A non-negative numeric value indicating the maximum weights/influence a sample can have on mean estimates. |
max_v |
A non-negative numeric value indicating the maximum weights/influence a sample can have on variance-covariance estimates. |
parallel |
A logical value indicating whether to use parallel computing (possible when using a multi-core machine). |
num_cores |
A numeric value indicating the number of cores to use (activated only if |
log_file |
A path to an output log file. Note that if the file |
verbose |
A logical value indicating whether to print logs. |
debug |
A logical value indicating whether to set the logger to a more detailed debug level; set |
Unico assumes the following model:
The mixture value at sample feature
:
is modeled as a weighted linear combination, specified by weights
, of a total of
source-specific levels, specified by
.
In addition, we also consider global-level covariates
that systematically affect the observed mixture values and their effect sizes
.
denotes the i.i.d measurement noise with variance
across all samples.
Weights have be to non-negative and sum up to
across all sources for each sample. In practice, we assume that the weights are fixed and estimated by external methods.
Source specific profiles are further modeled as:
denotes the population level mean of feature
at source
.
We also consider covariates
that systematically affect the source-specific values and their effect sizes
on each source.
Finally, we actively model the
by
covariance structure of a given feature
across all
sources
.
A list with the estimated parameters of the model. This list can be then used as the input to other functions such as tensor.
W |
An |
C1 |
An |
C2 |
An |
mus_hat |
An |
gammas_hat |
An |
betas_hat |
An |
sigmas_hat |
An |
taus_hat |
An |
scale.factor |
An |
config |
A list with hyper-parameters used for fitting the model and configurations for in the optimization algorithm. |
Us_hat_list |
A list tracking, for each feature, the sample weights used for each iteration of the mean optimization (activated only if |
Vs_hat_list |
A list tracking, for each feature, the sample weights used for each iteration of the variance-covariance optimization (activated only if |
Ls_hat_list |
A list tracking, for each feature, the computed estimates of the upper triangular cholesky decomposition of variance-covariance matrix at each iteration of the variance-covariance optimization (activated only if |
sigmas_hat_list |
A list tracking, for each feature, the computed estimates of the variance-covariance matrix at each iteration of the variance-covariance optimization (activated only if |
data = simulate_data(n=100, m=2, k=3, p1=1, p2=1, taus_std=0, log_file=NULL) res = list() res$params.hat = Unico(data$X, data$W, data$C1, data$C2, parallel=FALSE, log_file=NULL)
data = simulate_data(n=100, m=2, k=3, p1=1, p2=1, taus_std=0, log_file=NULL) res = list() res$params.hat = Unico(data$X, data$W, data$C1, data$C2, parallel=FALSE, log_file=NULL)