estimate_lucid: Fit LUCID models with one or multiple omics layers

In LUCIDus: LUCID with Multiple Omics Data

Fit LUCID models with one or multiple omics layers

Description

EM algorithm to estimate LUCID with one or multiple omics layers

Usage

estimate_lucid(
  lucid_model = c("early", "parallel", "serial"),
  G,
  Z,
  Y,
  CoG = NULL,
  CoY = NULL,
  K,
  init_omic.data.model = "EEV",
  useY = TRUE,
  tol = 0.001,
  max_itr = 1000,
  max_tot.itr = 10000,
  Rho_G = 0,
  Rho_Z_Mu = 0,
  Rho_Z_Cov = 0,
  family = c("normal", "binary"),
  seed = 123,
  init_impute = c("lod", "mix"),
  init_par = c("mclust", "random"),
  verbose = FALSE,
  n_starts = 1L
)

Arguments

lucid_model	Specifying LUCID model, "early" for early integration, "parallel" for lucid in parallel, "serial" for lucid in serial
G	an N by P matrix representing exposures
Z	Omics data, if "early", an N by M matrix; If "parallel", a list, each element i is a matrix with N rows and P_i features; If "serial", a list, each element i is a matrix with N rows and p_i features or a list with two or more matrices with N rows and a certain number of features
Y	a length N vector
CoG	an N by V matrix representing covariates to be adjusted for G -> X
CoY	an N by K matrix representing covariates to be adjusted for X -> Y
K	Number of latent clusters. If "early", an integer greater or equal to 2; If "parallel", an integer vector, same length as Z, with each element being an integer greater or equal to 2; If "serial", a list, each element is either an integer like that for "early" or an list of integers like that for "parallel", same length as Z
init_omic.data.model	a vector of strings specifies the geometric model of omics data. If NULL, See more in ?mclust::mclustModelNames
useY	logical, if TRUE, EM algorithm fits a supervised LUCID; otherwise unsupervised LUCID.
tol	stopping criterion for the EM algorithm
max_itr	Maximum iterations of the EM algorithm. If the EM algorithm iterates more than max_itr without converging, the EM algorithm is forced to stop.
max_tot.itr	Max number of total iterations for estimate_lucid function. estimate_lucid may conduct EM algorithm for multiple times if the algorithm fails to converge.
Rho_G	A scalar. This parameter is the LASSO penalty to regularize exposure coefficients in the G-to-X model. CoG adjustment covariates are included unpenalized. If user wants to tune the penalty, use the wrapper function lucid. Penalty tuning is supported for "early" and "parallel". For "serial", only scalar penalty inputs are supported.
Rho_Z_Mu	A scalar. This parameter is the LASSO penalty to regularize cluster-specific means for omics data (Z). If user wants to tune the penalty, use the wrapper function lucid. Penalty tuning is supported for "early" and "parallel". For "serial", only scalar penalty inputs are supported.
Rho_Z_Cov	A scalar. This parameter is the graphical LASSO penalty to estimate sparse cluster-specific variance-covariance matrices for omics data (Z). If user wants to tune the penalty, use the wrapper function lucid. Penalty tuning is supported for "early" and "parallel". For "serial", only scalar penalty inputs are supported.
family	The distribution of the outcome
seed	Random seed to initialize the EM algorithm
init_impute	Method to initialize the imputation of missing values in LUCID. lod (the default) initializes the imputation via replacing missing values by LOD / sqrt(2), where LOD is determined by the minimum of each variable in omics data; mix uses mclust::imputeData to implement EM Algorithm for Unrestricted General Location Model via the mix package to impute the missing values in omics data. mix is archived on CRAN and must be installed manually (e.g. from the CRAN Archive) to use this option; a request for init_impute = "mix" without mix installed raises an informative error.
init_par	For "early", an interface to initialize EM algorithm, if mclust, initiate the parameters using the mclust package, if random, initiate the parameters by drawing from a uniform distribution; For "parallel", mclust is the default for quick convergence; For "serial", each sub-model follows the above depending on it is a "early" or "parallel"
verbose	Logging level for fitting progress. If FALSE, concise start/finish status lines are printed. If TRUE, detailed iteration-level traces (including log-likelihood updates) are printed.
n_starts	Number of independent random starts for the EM algorithm (default 1). The EM algorithm converges only to a local optimum, so with n_starts > 1 the model is fitted from that many starting points and the fit with the highest observed-data log-likelihood is returned. Per-start log-likelihoods are recorded in em_control$start_loglik, which is worth inspecting: a wide spread indicates the likelihood surface is multi-modal and that a single start would have been unreliable.

Value

An object of class early_lucid, lucid_parallel or lucid_serial according to lucid_model. All three are lists; the components common to every fit are:

res_Beta

Estimates of the exposure-to-cluster (G -> X) association. For "early", a K by (1 + P + V) matrix of multinomial logistic coefficients, cluster 1 as reference. For "parallel" and "serial", a list holding the fitted object and the coefficient matrix per layer or stage.

res_Mu

Cluster-specific omics means (the mu of X -> Z). A K by M matrix for "early"; a list by layer or stage otherwise.

res_Sigma

Cluster-specific omics variance-covariance matrices (the sigma of X -> Z). A list of K matrices for "early"; a list by layer or stage otherwise.

res_Gamma

Estimates of the cluster-to-outcome (X -> Y) association, holding beta (absolute cluster levels), the reference-coded cluster_effect contrasts printed by summary(), any covariate coefficients, the residual sigma for a normal outcome, and the parameterization used.

inclusion.p

Posterior probability of cluster membership for each observation, r_{ij} of Eq 3. An N by K matrix for "early"; a list by layer or stage otherwise.

K

Number of latent clusters: an integer for "early", a list of integers for "parallel" and "serial".

var.names

Names of the G, Z and Y variables, as list(Gnames, Znames, Ynames).

init_omic.data.model

The mclust geometric model used for the omics covariances.

family

Outcome distribution, "normal" or "binary".

useY

Whether the outcome was used in fitting, i.e. whether the model is supervised.

Z

The omics data. For "early" and "parallel" this is the data the model was fitted to, so sporadically missing cells hold their imputed values and listwise-missing rows remain NA. For "serial" it is the omics data as supplied, still containing every missing value: imputation happens inside each stage, so the imputed omics for stage i are in submodel[[i]]$Z.

init_impute

The imputation method used to initialize missing omics values.

init_par

The parameter-initialization method used.

Rho

The penalties actually applied, as list(Rho_G, Rho_Z_Mu, Rho_Z_Cov).

missing_summary

How much omics data was missing and in what pattern, using the taxonomy of the incomplete-omics extension: complete_rows, listwise_rows (a whole omics layer missing for that subject) and sporadic_rows (some features missing), the corresponding proportions, and cell-level counts total_missing_cells, sporadic_missing_cells and prop_total_missing_cells. A list by layer for "parallel"; for "serial", n_stages and a per-stage breakdown.

em_control

The stopping controls used (tol, max_itr, max_tot.itr), which bootstrap refits reuse, together with convergence diagnostics. For "early" and "parallel" these are converged (FALSE, with a warning, if the fit exhausted max_itr without meeting tol), n_iter, n_restart, loglik_trace, n_starts and n_starts_ok. A "serial" fit runs no EM loop of its own, so it reports converged (TRUE only if every sub-model converged), n_iter (the total across sub-models), the per-sub-model submodel_converged and submodel_n_iter, and submodel_loglik_trace (each sub-model's own loglik_trace, by stage); the full diagnostics for a stage are in submodel[[i]]$em_control. For all three model types, the log-likelihood trace is checked for monotonicity as it is recorded: a decrease beyond a small, majorization-step-aware slack triggers a warning() naming the iteration (or stage) and the two values, since that indicates a numerical problem rather than expected EM behaviour.

The remaining components appear only for some model types:

likelihood

The observed-data log-likelihood at the returned estimates, present for all three model types. For "early" and "parallel" this is the single joint log-likelihood from the EM fit. A "serial" fit runs no joint EM loop – it is a sequence of conditionally fitted stages – so its likelihood is instead the sum of each stage's own log-likelihood (matching cal_loglik_serial()); the individual per-stage values are in submodel[[i]]$likelihood.

select

Feature-selection indicators, present for all three model types. select$selectG and select$selectZ are logical vectors of retained exposures and omics features. For "parallel", select$selectG is the exposure-wise union across layers (selected in at least one layer), select$selectG_layer holds the per-layer exposure selection, and select$selectZ holds the per-layer omics selection (a list with at least one layer's selection reported). For "serial", select is stage 1's own selection only (submodel[[1]]$select) – stage 1 is the only stage whose G is the user's actual exposures; every later stage's G is the previous stage's posterior cluster-membership probabilities, so Rho_G is always 0 there and that stage's selectG is not a meaningful exposure-selection result (its selectZ still is). The complete per-stage selection record remains available at submodel[[i]]$select for every i.

N

Number of observations. Present for "parallel" and "serial"; for "early", use nrow(fit$inclusion.p).

z

The E-step responsibilities over the joint cluster configuration across layers, before they are marginalised into inclusion.p. Present for "parallel" only – nothing in the package reads this field, and there is no plan to extend it to "serial".

res_Delta

Estimates of the between-stage cluster transition associations, one element per transition (length n_stages - 1). Element i is the coefficient object of stage i + 1 fitted with stage i's cluster assignment in place of the exposures, so it has the same structure as res_Beta for that stage. Present for "serial" only.

submodel

The fitted sub-models, one per stage, each itself an early_lucid or lucid_parallel object. Present for "serial" only.

Examples

i <- 1008
set.seed(i)
G <- matrix(rnorm(500), nrow = 100)
Z1 <- matrix(rnorm(1000), nrow = 100)
Z2 <- matrix(rnorm(1000), nrow = 100)
Z3 <- matrix(rnorm(1000), nrow = 100)
Z <- list(Z1 = Z1, Z2 = Z2, Z3 = Z3)
Y <- rnorm(100)
CoY <- matrix(rnorm(200), nrow = 100)
CoG <- matrix(rnorm(200), nrow = 100)
fit1 <- estimate_lucid(G = G, Z = Z, Y = Y, K = list(2, 2, 2),
lucid_model = "serial",
family = "normal",
seed = i,
CoG = CoG, CoY = CoY,
useY = TRUE,
max_itr = 20, max_tot.itr = 50)

LUCIDus documentation built on Sept. 3, 2026, 1:06 a.m.