Skip to content
Projects

kstawiski/OmicSelector has 36 GitHub stars and 4 forks. GitHot measured 0 stars in 24 hours and 0 in 7 days. Last checked 2026-11-09 16:37 UTC. Last push 2026-10-30.

How these counts are measured

OmicSelector - Environment, docker-based application and R package for biomarker signiture selection (feature selection) & deep learning diagnostic tool development from high-throughput high-throughput omics experiments and other multidimensional datasets. Initially developed for miRNA-seq, RNA-seq and qPCR.

What you can do

  • FilterGOF_KS - Kolmogorov-Smirnov filter for sparse data
  • FilterHurdle - Two-part hurdle model filter
  • FilterZeroProp - Zero-proportion difference filter
  • make_autotuner_glmnet() - Bayesian-tuned elastic net

Measured star growth

in 24h
0
in 7 days
0
in 15 days
Update delayed
in 1 month
Update delayed

Net star changes over overlapping periods, not daily totals. Missing measurements are left blank. Stars measure attention, not software quality.

README.md

Original

OmicSelector 2.6.10

OmicSelector provides 73 registered scoring methods and audits whether a frozen model can score one specimen independently. All existing methods, validation routes and outputs remain available.

Quick start: core R installation

install.packages("remotes")
remotes::install_github("kstawiski/[email protected]", dependencies = NA)
library(OmicSelector)

set.seed(2610)
X <- matrix(rexp(41 * 30), 41, 30,
            dimnames = list(NULL, paste0("marker", 1:30)))
y <- rep(0:1, each = 20)
X[21:40, 1:5] <- X[21:40, 1:5] * 3
X_train <- X[1:40, , drop = FALSE]
x_new <- X[41, , drop = FALSE]
frozen <- deploy_singlesample(X_train, y, method = "ws-rclr-trimmed")
saveRDS(frozen, "reference.rds")
frozen <- readRDS("reference.rds")
score_specimen(frozen, x_new, strict = TRUE)
is_singlesample_deployable(frozen, X_train[1:8, , drop = FALSE])

This simulated example needs no Python or suggested R packages. The reference returns a signed panel score, not a disease probability. The audit checks score independence from companion specimens and frozen model state. Clinical performance and decision thresholds require validation for the intended use. See Deploying a frozen single-sample classifier.

Choosing a method

Start with the reference ws-rclr-trimmed. In the prespecified benchmark, none of 48 estimable methods in the 30-unit primary synthesis improved on the reference after BY multiplicity correction. TabPFN had a pooled AUC difference of +0.0245 (95% CI +0.0060 to +0.0431), TabICL +0.0228 (-0.0071 to +0.0526), and selbal +0.0116 (-0.0227 to +0.0460). Lack of a significant difference does not establish equivalence. An exploratory input-quality exclusion favored selbal (+0.0302; BY P = 2.68e-7), with 98.63% of the estimate's weight from Toray. These are distinct analyses. Source: the R12a primary synthesis and sensitivity tables in the completed paper reproduction bundle.

requirements <- method_requirements()
requirements[, c("method_id", "family", "r_packages", "python_modules",
                 "model_weights", "gpu", "fit_seconds",
                 "score_seconds_per_specimen", "peak_memory_bytes",
                 "frozen_model_bytes", "measuring_node", "cost_status")]
available_methods()

PB1 costs summarize completed CPU measurements. The measurements attribute retains dataset-specific results, hardware and source information. Missing measurements are shown as NA. Transfer methods may require kit, technology, biofluid or other training metadata. Dependency availability does not replace these input requirements or the deployment audit.

Method Fit (s) Score one specimen (s) Peak fit memory (MiB) Frozen model (MiB) CPU hardware
ws-rclr-trimmed 0.3649 0.02435 170.4 0.06907 Intel(R) Xeon(R) Platinum 8380 CPU @ 2.30GHz
reo-ktsp 6.612 0.02383 189.8 0.1144 Intel(R) Xeon(R) Platinum 8380 CPU @ 2.30GHz
ai-tabpfn 751.2 7.855 1.426e+04 2836 Intel(R) Xeon(R) Platinum 8380 CPU @ 2.30GHz
tab-tabdpt 15.32 0.5271 1400 248.6 Intel(R) Xeon(R) Platinum 8380 CPU @ 2.30GHz
tab-tabicl 331.9 4.974 1.955e+04 5896 Intel(R) Xeon(R) Platinum 8380 CPU @ 2.30GHz

These are medians over completed CPU measurement cells; dataset sizes and hardware are retained in the shipped tables.

Full installation

# Preview the required packages and checkpoints without changing the environment.
install_method_dependencies(c("ai-tabpfn", "tab-tabdpt", "tab-tabicl"), dry_run = TRUE)
# In a fresh R session, install every registered method's required dependencies.
install_method_dependencies()
# Or select methods or roster families.
install_method_dependencies(families = "A", python = FALSE)
available_methods()

The installer uses the omicselector reticulate virtual environment by default and reuses packages and weights already present. In later R sessions, set RETICULATE_PYTHON to that environment's interpreter before loading OmicSelector. Optional alternative implementations are reported separately. Other package features, such as the web interface and plotting, retain their documented suggested packages; install those with remotes::install_github(..., dependencies = TRUE) where their repositories are available.

The release workflow builds kstawiski/omicselector:2.6.10. Its Dockerfile provides the RStudio and web-interface environment. Add selected Python-backed methods with the same installer inside a persistent container and retain the virtual environment and weight caches. See Dockerfile.


Testing and continuous integration

Correctness is enforced by a testthat suite that runs on every push through GitHub Actions (R-CMD-check, badge above). The suite covers the pipeline, the nested-CV leakage guarantees, feature-stability metrics, calibration, batch correction, multi-omics integration, and the single-sample scoring roster. Run it locally with:

devtools::test()

Overview

OmicSelector is an R package for biomarker discovery that enforces methodologically sound machine learning practices. Built on the mlr3 ecosystem, it guarantees:

  • Zero Data Leakage: Feature selection occurs strictly inside cross-validation folds
  • Proper Nested CV: Separation of outer (evaluation) and inner (selection) loops
  • Feature Stability: Nogueira Stability Index for reproducible biomarker sets
  • Multi-Objective Selection: Balance performance, stability, and parsimony

Single-sample scoring methods

The within-sample, reference-free scoring layer, introduced in 2.6.0, is exported as documented R functions:

  1. Five within-sample compositional scoring methods: ws_rclr_trimmed, ws_balance_ilr, ws_alr_pivot, ws_mad_logratio, and ws_dominance_score.
  2. Frozen-reference denoising add-ons: fit_frozen_ruv / apply_frozen_ruv, fit_robust_pca_residual / apply_robust_pca_residual, and hemolysis_index_blondal with fit_hemolysis_rr / apply_hemolysis_rr.
  3. qPCR non-detect handling: qpcr_nondetect_impute wraps Bioconductor nondetects, with qpcr_nondetect_lod_fallback as the documented Ct-limit fallback.
  4. Matched-null benchmark: singlesample_matched_null_benchmark and singlesample_matched_null_benchmark_cv implement stratified random-panel benchmarking with default K = 10000.
  5. Provenance pre-flight: os_provenance_preflight audits a manifest TSV and returns NO_OVERLAP, KNOWN_OVERLAP, or UNKNOWN_ACCESSION.
  6. Dependent-p combination: singlesample_bh_fdr_correct_blocked reports both conservative and textbook Brown 1975 / Kost-McDermott 2002 reference distributions.
  7. Block-aware FDR: singlesample_bh_fdr_correct_blocked and singlesample_bh_fdr_correct_matched_null provide the block-aware and matched-null BH correction helpers.
  8. AUC confidence intervals: singlesample_hanley_mcneil_auc_ci exposes the Hanley-McNeil 1982 95% AUC CI calculation used by the matched-null benchmark outputs.
  9. ILR coverage gate: ws_balance_ilr exposes min_balance_coverage so panels that preserve the fixed balance dictionary can activate ILR scoring at the documented coverage threshold.
  10. Full single-sample method-bank registry: singlesample_method_bank maps every method family (A-J) to exported package functions, including the kit-aware, provenance-aware, anchor/reference, technology-aware, biofluid-aware, learned/adversarial, Group-DRO, and Sinkhorn-OT scorer families.
  11. External-validation split helpers: singlesample_make_loco_splits and singlesample_make_locto_splits construct leave-one-cohort-out and leave-one-cancer-type-out splits with same provenance-block exclusion.

Reproducible scores across machines

Set CPU dispatch before starting R or Rscript. Some numerical libraries read these settings when loaded. Use the same package, Python environment, model weights, seeds, thread counts and OPENBLAS_CORETYPE on every machine. The validated Linux x86-64 environment used this recipe:

export ATEN_CPU_CAPABILITY=avx2
export MKL_ENABLE_INSTRUCTIONS=AVX2
export ONEDNN_MAX_CPU_ISA=AVX2
export DNNL_MAX_CPU_ISA=AVX2
export OPENBLAS_CORETYPE=Zen
export OMP_NUM_THREADS=4 OPENBLAS_NUM_THREADS=4 MKL_NUM_THREADS=4
export NUMEXPR_NUM_THREADS=4 NUMBA_NUM_THREADS=4
Rscript analysis.R

OPENBLAS_CORETYPE=Zen must be supported by the OpenBLAS build in use; the measured environment used the same build on all hosts. Record these variables, the host CPU, library versions, device and thread counts with the results.

Under these matched startup settings, measured TabPFN and TabDPT scores on x86-64 CPU workstations, a CPU server and AVX-512 cluster nodes agreed within 8.05e-7, with zero changed classifications and zero AUC differences. Without matched dispatch, AVX-512 scores differed by up to 0.030 for TabDPT and 0.0003 for TabPFN. Repeated CPU TabICL runs differed by up to 5.96e-7. A fixed seed alone does not guarantee bitwise equality across machines. These are computational checks of the specified environment and datasets; they do not establish clinical validity.

TabPFN also bounds query batches using a 4 GiB float32 attention-workspace estimate for three dense query-by-context matrices. CPU and CUDA use the same bound. It excludes fitted context, model parameters and other process memory. The default training context is capped at 4,096 rows. Allow additional memory for fitting and scoring; the workspace bound is not a total memory limit.

Within-Sample Biomarker Panel Classification

Since 2.6.0, OmicSelector has provided within-sample classification methods for reference-cohort-free biomarker panels and provenance-aware validation helpers:

  • Four image encodings: encode_simple_grid(), encode_corr_grid(), encode_deepinsight(), and encode_ratio_image()
  • Two dense CoDA learners: train_ratio_cnn() for pairwise log-ratio images and train_clr_mlp() for CLR-transformed vectors
  • A sparse balance interface: train_codacore(), which uses the official codacore backend when available and otherwise falls back to a deterministic sparse-balance model
  • mlr3 learners registered on load: classif.ratio_cnn, classif.clr_mlp, and classif.codacore
  • Group-aware resampling, provenance-floor, identifiability, matched-null, and operating-point gates for public-omics biomarker-panel audits
library(OmicSelector)

X <- matrix(rnorm(14 * 40), nrow = 40, ncol = 14)
y <- factor(rep(c("control", "case"), each = 20), levels = c("control", "case"))

clr_fit <- train_clr_mlp(X, y, backend = "glm", verbose = FALSE)
ratio_fit <- train_ratio_cnn(
  X_train = X,
  y_train = y,
  X_test = X[1:5, , drop = FALSE],
  epochs = 5,
  verbose = FALSE
)

mlr3::lrn("classif.clr_mlp", backend = "glm", verbose = FALSE)

Installation

# Install the release described in the accompanying article
remotes::install_github("kstawiski/[email protected]")

# Or install the current development version
remotes::install_github("kstawiski/OmicSelector")

# Load the package
library(OmicSelector)

Version 2.6.9 changes every TabDPT fit with more than 128 features. It removes undetermined directions and fixes arbitrary projection signs using exact float64 CPU SVD and a positive largest-absolute-loading convention. Only numerically determined directions are retained, up to 128. The fitted classifier records projection_rank, projection_q, projection_rank_tolerance, projection_min_retained_ratio and projection_sign_min_margin, the smallest gap between the two largest absolute loadings in any retained direction. A zero margin marks a tie. Fits with at most 128 features, TabPFN and TabICL keep their existing paths.


Data Format

OmicSelector expects a data.frame with:

  • Feature columns: Numeric values (gene expression, miRNA counts, etc.)
  • Target column: Factor or character for classification, numeric for regression
  • Optional metadata: patient_id, batch, etc.

Example Data Structure

# Your data should look like this:
head(my_data)
#>   gene_A  gene_B  gene_C  gene_D  gene_E  outcome
#> 1   2.34    1.56    3.21    0.89    1.45     Case
#> 2   1.12    2.89    0.45    2.11    0.78  Control
#> 3   3.45    0.23    2.67    1.34    2.90     Case
#> ...

# Required structure:
# - All feature columns must be numeric
# - Target column can be factor, character (classification) or numeric (regression)
# - No row names required (use a patient_id column if needed)

Creating Test Data

# Simulate expression data for 100 samples, 50 genes
set.seed(42)
n_samples <- 100
n_genes <- 50

# Create feature matrix
features <- matrix(rnorm(n_samples * n_genes), nrow = n_samples)
colnames(features) <- paste0("gene_", 1:n_genes)

# Add some signal to first 5 genes
outcome <- factor(rep(c("Case", "Control"), each = n_samples/2))
features[outcome == "Case", 1:5] <- features[outcome == "Case", 1:5] + 1.5

# Combine into data.frame
my_data <- data.frame(features, outcome = outcome)

Multi-Omics Data

For multi-omics, pass a named list of data.frames (same samples, same row order):

# RNA-seq data (100 samples x 1000 genes)
rna_data <- data.frame(matrix(rnorm(100 * 1000), nrow = 100))
names(rna_data) <- paste0("gene_", 1:1000)

# miRNA data (100 samples x 200 miRNAs)
mirna_data <- data.frame(matrix(rnorm(100 * 200), nrow = 100))
names(mirna_data) <- paste0("miR_", 1:200)

# Target vector
outcome <- factor(rep(c("Case", "Control"), each = 50))

# Create multi-omics pipeline
pipeline <- OmicPipeline$new(
  data = list(rna = rna_data, mirna = mirna_data),
  target = outcome,
  positive = "Case"
)
# Features are automatically namespaced: rna::gene_1, mirna::miR_1, etc.

Quick Start

Step 1: Create Pipeline

library(OmicSelector)

# Create pipeline from data
pipeline <- OmicPipeline$new(
  data = my_data,          # Your data.frame
  target = "outcome",      # Name of the target column
  positive = "Case"        # Positive class for AUC calculation
)
#> OmicPipeline created: 100 samples, 50 features, target='outcome'

Step 2: Create GraphLearner

# Create a learner with embedded feature selection
learner <- pipeline$create_graph_learner(
  filter = "anova",        # Feature selection method
  model = "ranger",        # Classification model
  n_features = 10          # Select top 10 features
)

Step 3: Run Nested Cross-Validation

# Create benchmark service
benchmark <- BenchmarkService$new(
  task = pipeline,         # Your pipeline
  outer_folds = 5,         # Outer CV folds (evaluation)
  inner_folds = 3,         # Inner CV folds (selection)
  seed = 42                # For reproducibility
)

# Add learner and run
benchmark$add_learner(learner)
result <- benchmark$run()
#> Running nested cross-validation...
#> Completed in 12.3 seconds

# View results
print(result)

Web UI (Shiny)

Run the interactive WebUI locally:

shiny::runApp("shiny/omicselector", host = "0.0.0.0", port = 3838)

Run via Docker (Shiny WebUI + RStudio Server):

docker compose up

Then open:

  • WebUI: http://localhost:3838
  • RStudio Server: http://localhost:8787 (user: rstudio, password: omicselector)

Uploads and results are saved under /OmicSelector/analyses/<analysis_id> inside the container and persisted to ./analyses on the host. Set OMICSELECTOR_ANALYSES_DIR to override the storage location.


Recipes

Binary Classification (Single-Omics)

pipeline <- OmicPipeline$new(
  data = my_data,
  target = "outcome",
  positive = "Case"
)

learner <- pipeline$create_graph_learner(
  filter = "anova",
  model = "rpart",
  n_features = 20
)

fit <- pipeline$fit(learner, seed = 1)
fit$selected_features

Multi-Omics (Named List Input)

pipeline <- OmicPipeline$new(
  data = list(rna = rna_data, mirna = mirna_data),
  target = "outcome",
  positive = "Case"
)

learner <- pipeline$create_graph_learner(
  filter = "mrmr",
  model = "rpart",
  n_features = 30
)

Benchmark (Nested CV + Screening + Caching)

learner <- pipeline$create_graph_learner(
  filter = "anova",
  model = "rpart",
  n_features = 20,
  screening = TRUE,
  screening_frac = 0.2
)

result <- pipeline$benchmark(
  learners = learner,
  outer_folds = 5,
  inner_folds = 3,
  seed = 42,
  cache_dir = "cache",
  parallel = TRUE,
  threads = 1
)

Deep Learning (mlr3torch + Autoencoder)

# Requires: torch + mlr3torch
learner <- pipeline$create_graph_learner(
  filter = "anova",
  model = "mlp",
  n_features = 50,
  autoencoder = list(
    latent_dim = 32,
    hidden_layers = c(128, 64),
    epochs = 50,
    batch_size = 64,
    early_stopping = TRUE,
    patience = 10
  )
)

result <- pipeline$benchmark(
  learners = learner,
  outer_folds = 3,
  inner_folds = 2,
  seed = 1
)

Transfer learning via autoencoder pretraining:

ae <- autoencoder_fit(x_train, latent_dim = 32, epochs = 50)
autoencoder_save(ae, "pretrained_ae.pt")

learner <- pipeline$create_graph_learner(
  filter = "anova",
  model = "mlp",
  n_features = 50,
  autoencoder = list(
    pretrained = "pretrained_ae.pt",
    freeze_encoder = TRUE
  )
)

Checkpoint export/import for fine‑tuning:

fit <- pipeline$fit(learner, seed = 1)
export_omicfit_checkpoint(fit, "mlp_checkpoint.pt")

finetuned <- finetune_mlr3torch_checkpoint(
  learner = learner,
  task = pipeline$get_task(),
  path = "mlp_checkpoint.pt",
  epochs = 30
)

# Or stay at the OmicFit level
finetuned_fit <- finetune_omicfit_checkpoint(
  fit = fit,
  task = pipeline$get_task(),
  path = "mlp_checkpoint.pt",
  epochs = 30
)

Configuration Options

Feature Selection Methods

Use the filter parameter in create_graph_learner():

Univariate Statistical Tests

Method Code Description Best For
ANOVA F-test "anova" Tests mean differences between classes Default choice, continuous features
Kruskal-Wallis "kruskal" Non-parametric rank-based test Non-normal distributions
Chi-Squared "chi_squared" Chi-squared test Discrete/categorical features

Variance-Based

Method Code Description Best For
Variance "variance" Removes low-variance features Pre-filtering, removing constants

Correlation-Based

Method Code Description Best For
Correlation "correlation" Correlation with target Quick univariate filter

Information-Theoretic

Method Code Description Best For
Information Gain "information_gain" Entropy-based importance Mixed feature types
Gain Ratio "gain_ratio" Normalized information gain Avoiding bias toward high-cardinality
mRMR "mrmr" Min Redundancy Max Relevance Reducing feature redundancy
CMIM "cmim" Conditional Mutual Info Max Complex feature dependencies
JMIM "jmim" Joint Mutual Info Max Capturing feature interactions
JMI "jmi" Joint Mutual Information Similar to JMIM

Model-Based Importance

Method Code Description Best For
AUC "auc" AUC of univariate models Classification performance
Relief "relief" Instance-based algorithm Detecting interactions
RF Importance "importance" Random Forest importance Non-linear relationships
Permutation "permutation" Permutation importance Model-agnostic importance
# Example: Try different filters
learner_anova <- pipeline$create_graph_learner(filter = "anova", model = "ranger", n_features = 20)
learner_mrmr <- pipeline$create_graph_learner(filter = "mrmr", model = "ranger", n_features = 20)
learner_relief <- pipeline$create_graph_learner(filter = "relief", model = "xgboost", n_features = 15)

Classification Models

Use the model parameter in create_graph_learner():

Tree-Based Ensembles

Model Code Description Strengths
Random Forest "ranger" Fast RF implementation Handles interactions, robust, default choice
XGBoost "xgboost" Gradient boosting High performance, handles missing values
LightGBM "lightgbm" Fast gradient boosting Very fast, memory efficient
Decision Tree "rpart" Single CART tree Interpretable, fast

Linear Models

Model Code Description Strengths
Elastic Net "glmnet" L1/L2 regularized regression Interpretable coefficients, feature selection
Logistic Regression "log_reg" Simple logistic regression Baseline, highly interpretable

Distance-Based Methods

Model Code Description Strengths
SVM "svm" Support Vector Machine (RBF) High-dimensional data
k-NN "kknn" K-Nearest Neighbors Non-parametric, simple

Probabilistic Classifiers

Model Code Description Strengths
Naive Bayes "naive_bayes" Probabilistic classifier Fast, works with small data
LDA "lda" Linear Discriminant Analysis Dimensionality reduction
QDA "qda" Quadratic Discriminant Analysis Non-linear decision boundaries

Neural Networks

Model Code Description Strengths
Neural Net "nnet" Single-layer neural network Non-linear relationships
# Example: Compare models
learner_rf <- pipeline$create_graph_learner(filter = "anova", model = "ranger", n_features = 15)
learner_xgb <- pipeline$create_graph_learner(filter = "anova", model = "xgboost", n_features = 15)
learner_glmnet <- pipeline$create_graph_learner(filter = "anova", model = "glmnet", n_features = 15)

Using Any mlr3 Learner

OmicSelector is built on the mlr3 ecosystem. You can use any mlr3 learner directly with BenchmarkService:

library(mlr3learners)

# Get the task from the pipeline
task <- pipeline$get_task()

# Create any mlr3 learner with probability predictions
xgboost_learner <- lrn("classif.xgboost", predict_type = "prob")
nnet_learner <- lrn("classif.nnet", predict_type = "prob")

# Wrap in a graph with feature selection (optional)
library(mlr3pipelines)
library(mlr3filters)

xgb_graph <- po("filter", filter = flt("anova"), filter.nfeat = 20) %>>%
  po("learner", xgboost_learner)

# Convert to GraphLearner
xgb_graph_learner <- as_learner(xgb_graph)
xgb_graph_learner$id <- "xgboost_anova_20"

# Use in benchmark
benchmark <- BenchmarkService$new(task = pipeline, outer_folds = 5, seed = 42)
benchmark$add_learner(xgb_graph_learner)
result <- benchmark$run()

See the mlr3learners documentation for all available learners.

Full Parameter Reference

learner <- pipeline$create_graph_learner(
  # Feature Selection
  filter = "anova",           # Filter method (see table above)
  n_features = 20,            # Number of features to select
                              # Use decimal < 1 for proportion (e.g., 0.1 = top 10%)

  # Preprocessing
  impute_method = "median",   # Handle missing values: "median", "mean", "sample", "pmm"
  scale = TRUE,               # Scale features to zero mean, unit variance

  # Class Imbalance (applied INSIDE each CV fold)
  oversample = NULL,          # NULL (none), "smote", or "rose"

  # Batch Effect Correction
  batch_correct = FALSE       # TRUE to apply FrozenComBat (requires batch column)
)

Advanced Imputation: PMM (Predictive Mean Matching)

For complex missing data patterns, use PMM imputation via the mice package:

# PMM uses predictive models to impute missing values
# Best for: MAR (Missing At Random) data, preserving distributions
learner <- pipeline$create_graph_learner(
  filter = "anova",
  model = "ranger",
  n_features = 20,
  impute_method = "pmm"       # Predictive Mean Matching via mice
)

# PMM advantages:
# - Preserves data distribution (unlike mean imputation)
# - Handles complex correlation structures
# - Produces plausible values from observed data
# - Applied properly inside each CV fold to prevent leakage

Comparing Multiple Configurations

# Create benchmark service
benchmark <- BenchmarkService$new(pipeline, outer_folds = 5, inner_folds = 3, seed = 42)

# Add multiple learners with different configurations
benchmark$add_learner(
  pipeline$create_graph_learner(filter = "anova", model = "ranger", n_features = 10),
  id = "anova_rf_10"
)
benchmark$add_learner(
  pipeline$create_graph_learner(filter = "anova", model = "ranger", n_features = 20),
  id = "anova_rf_20"
)
benchmark$add_learner(
  pipeline$create_graph_learner(filter = "mrmr", model = "glmnet", n_features = 15),
  id = "mrmr_glmnet_15"
)
benchmark$add_learner(
  pipeline$create_graph_learner(filter = "variance", model = "svm", n_features = 25),
  id = "var_svm_25"
)

# Run all configurations
result <- benchmark$run()

# Compare performance
print(result$performance)
#>           learner_id classif.auc classif.acc classif.bbrier
#> 1       anova_rf_10       0.892       0.840          0.112
#> 2       anova_rf_20       0.901       0.850          0.098
#> 3   mrmr_glmnet_15        0.885       0.830          0.121
#> 4       var_svm_25        0.878       0.820          0.134

Advanced Features

Phase 5: Cutting-Edge Biomarker Discovery

OmicSelector includes state-of-the-art methods validated on TCGA data:

GOF Filters for Sparse/Zero-Inflated Data

Traditional filters (ANOVA, variance) miss biologically important features that are "off" in one condition but expressed in another. GOF filters detect these patterns:

library(OmicSelector)

# Use KS filter - captures distributional differences including zero-inflation
learner <- pipeline$create_graph_learner(
  filter = "gof_ks",      # Kolmogorov-Smirnov GOF filter
  model = "ranger",
  n_features = 30
)

# Other GOF options:
# filter = "hurdle"      # Two-part model (zero-frequency + magnitude)
# filter = "zero_prop"   # Simple zero-proportion difference

# GOF filters discovered miR-200c/miR-141 family in TCGA kidney data
# that ANOVA missed - these are established EMT regulators in cancer

Bayesian Hyperparameter Optimization

Replace grid search with intelligent Bayesian optimization:

# Bayesian-tuned glmnet with omics-optimized search space
autotuner <- make_autotuner_glmnet(
  task,
  n_evals = 20,           # Budget of evaluations
  inner_folds = 3         # Inner CV for tuning
)
autotuner$train(task)

# Get optimal parameters
params <- get_optimal_params(autotuner)
# Returns: alpha, s (lambda) optimized for your data

AutoXAI: Interpretability with Correlation Warnings

Get SHAP-based interpretability with automatic warnings for correlated features:

# Train a model
learner <- lrn("classif.ranger", predict_type = "prob", importance = "impurity")
learner$train(task)

# Run XAI pipeline
xai <- xai_pipeline(
  learner = learner,
  task = task,
  top_k = 20,
  cor_threshold = 0.7     # Warn if features correlated > 0.7
)

# View results
print(xai$top_features)          # Top features by permutation importance
print(xai$correlations$clusters) # Correlated feature clusters
plot_xai_importance(xai)         # Visualization

# Warnings prevent SHAP misinterpretation when features are collinear

Bootstrap Stability Ensemble

Find reproducible biomarkers across resamples:

# Create stability ensemble with multiple filters
ensemble <- create_stability_ensemble(
  preset = "default",     # Uses mRMR, AUC, Information Gain
  n_bootstrap = 100,      # 100 bootstrap iterations
  n_features = 30
)

# Fit and get stable features
ensemble$fit(task, seed = 42)
stable <- ensemble$get_feature_importance(30)

# Features with >90% selection frequency are highly reproducible
# TCGA validation found 6 biomarkers with 100% stability

SMOTE for Imbalanced Data

Proper SMOTE application inside CV folds:

# Balance classes with SMOTE (applied inside each fold)
task_balanced <- smote_augment(task, ratio = 1.0, k = 5)

# Validate synthetic data quality
validation <- validate_synthetic(real_data, synthetic_data)
print(validation$quality_score)  # 0-1, higher is better

Feature Stability Analysis

# Compute Nogueira Stability Index with Bootstrap Confidence Intervals
stability <- compute_stability_from_resample(result$benchmark_result)
print(stability)
#> Nogueira Stability Index: 0.823
#> Interpretation: Good stability

# Bootstrap CI for uncertainty quantification
stability_ci <- compute_nogueira_stability(
  feature_matrix,
  n_bootstrap = 1000,      # Number of bootstrap iterations
  conf_level = 0.95        # 95% confidence interval
)
print(stability_ci)
#> Nogueira Stability: 0.823
#> 95% Bootstrap CI: [0.756, 0.889]
#> Interpretation: Good stability (CI does not include 0.7 threshold)

# Get features selected in ALL folds (consensus)
consensus_features <- stability$consensus_features
print(consensus_features)
#> [1] "gene_1" "gene_3" "gene_5"

# Get selection frequency per feature
head(stability$selection_frequency)
#>   gene_1   gene_3   gene_5   gene_2   gene_4
#>     1.00     1.00     1.00     0.80     0.60

Optimal Signature Selection

# Select best signature balancing multiple objectives
best <- select_best_signature(
  result,
  mode = "weighted",         # Selection mode
  stability_weight = 0.3,    # Weight for stability (vs performance)
  parsimony_weight = 0.1     # Weight for fewer features
)

# Available modes:
# - "best_auc": Highest AUC only
# - "best_stability": Highest Nogueira index only
# - "weighted": Weighted combination (default)
# - "pareto": Pareto-optimal with distance-to-utopia

Model Calibration

# After getting predictions from your model
probabilities <- predict(model, newdata = test_data, type = "prob")
true_labels <- test_data$outcome

# Compute calibration metrics
ece <- compute_ece(probabilities, true_labels)
brier <- decompose_brier(probabilities, true_labels)

# Fit calibration (Platt scaling or isotonic regression)
calibrator <- fit_platt_scaling(probabilities, true_labels)
# or: calibrator <- fit_isotonic_calibration(probabilities, true_labels)

# Apply to new predictions (calibrator is a function)
calibrated_probs <- calibrator(new_probabilities)

Batch Correction (FrozenComBat)

# Include batch information in pipeline
pipeline <- OmicPipeline$new(
  data = my_data,
  target = "outcome",
  batch = "batch_id"         # Column containing batch labels
)

# Enable batch correction (parameters fit on training data only)
learner <- pipeline$create_graph_learner(
  filter = "anova",
  model = "ranger",
  n_features = 20,
  batch_correct = TRUE       # Applies FrozenComBat inside CV
)

SHAP Interpretability

# Create explainer for a trained model
explainer <- create_explainer(trained_learner, task)

# Compute SHAP values
shap <- shap_values(explainer, new_data)

# Check for problematic correlations
# (SHAP can be misleading when features are highly correlated)
correlations <- check_feature_correlations(
  data = my_data,
  features = selected_features,
  threshold = 0.7            # Warn if correlation > 0.7
)

Multi-Omics Late Integration

# Late-integrate modalities with a stacked ensemble: one base learner per
# modality plus a meta-learner. Samples must be aligned across modalities.
ensemble <- stack_omics(
  data_list    = list(rna = rna_data, mirna = mirna_data),
  learner_list = list(rna = "classif.ranger", mirna = "classif.glmnet"),
  y            = outcome,
  meta_learner = "classif.log_reg"
)

# Predict on new, aligned modalities
predictions <- ensemble$predict(list(rna = new_rna, mirna = new_mirna))

For early integration, pass a named list of modalities directly to OmicPipeline$new(); features are namespaced automatically as modality::feature.

TRIPOD+AI Reporting

Generate a TRIPOD+AI-oriented reporting object from a fitted model and result for transparent model documentation:

report <- create_report_data(result)
generate_tripod_report(report, output_file = "tripod_report.html")

YAML Configuration

For reproducible batch processing:

# config.yaml
data:
  file: "data/expression.csv"
  target: "outcome"
  positive: "Case"

pipeline:
  filter: "anova"
  model: "ranger"
  n_features: 20
  scale: true
  impute: "median"

benchmark:
  outer_folds: 5
  inner_folds: 3
  seed: 42
  stratify: true

output:
  dir: "./results"
  report_format: "html"

Run via CLI:

Rscript inst/bin/omicselector run --config=config.yaml

Docker

# Build image
docker build -f Dockerfile -t omicselector .

# Run interactive R session
docker run -it --rm -v $(pwd):/workspace omicselector R

# Run analysis with config
docker run -it --rm -v $(pwd):/workspace omicselector \
  Rscript inst/bin/omicselector run --config=/workspace/config.yaml

Module Reference

Core Pipeline

Module Description
OmicPipeline Build mlr3 graphs with zero-leakage preprocessing
BenchmarkService Nested CV with proper inner/outer loop separation
select_best_signature() Multi-objective signature selection
compute_nogueira_stability() Feature selection stability metrics

Batch & Calibration

Module Description
FrozenComBat Batch correction with frozen parameters
fit_platt_scaling() Platt calibration
fit_isotonic_calibration() Isotonic regression calibration

Phase 5: Advanced Features

Module Description
FilterGOF_KS Kolmogorov-Smirnov filter for sparse data
FilterHurdle Two-part hurdle model filter
FilterZeroProp Zero-proportion difference filter
make_autotuner_glmnet() Bayesian-tuned elastic net
make_autotuner_xgboost() Bayesian-tuned XGBoost
make_autotuner_ranger() Bayesian-tuned Random Forest
xai_pipeline() DALEX-based interpretability with correlation warnings
plot_xai_importance() Feature importance visualization
create_stability_ensemble() Bootstrap stability feature selection
smote_augment() SMOTE for class imbalance
validate_synthetic() Synthetic data quality validation

Quality Control & Clinical Reporting

Module Description
create_report_data() Assemble TRIPOD+AI report data from a benchmark result
generate_tripod_report() TRIPOD+AI compliant HTML/PDF reporting
os_bias_audit_report() Provenance/bias audit report for public-omics panels

Multi-Omics & Export

Module Description
stack_omics() Late integration of multi-omics data (stacked ensemble)
export_vetiver() Model export for deployment
export_bundle() Self-contained model bundle with metadata
export_onnx() ONNX export for cross-platform deployment

Citation

@article{stawiski2022omicselector,
  title={OmicSelector: automatic feature selection and deep learning
         modeling for omic experiments},
  author={Stawiski, Konrad and Kaszkowiak, Marcin and Mikulski, Damian and others},
  journal={bioRxiv},
  year={2022},
  doi={10.1101/2022.06.01.494299}
}

Publications Using OmicSelector

Selected peer-reviewed work built with OmicSelector:

  1. Stawiski K, Fortner RT, Pestarino L, et al. Validation of miRNA signatures for ovarian cancer earlier detection in the pre-diagnosis setting using machine learning approaches. Front Oncol. 2024;14:1389066. PMID: 38983926.

Authors

Supervised by: Prof. Wojciech Fendler, M.D., Ph.D.

Department of Biostatistics and Translational Medicine, Medical University of Lodz, Poland

Links

Reproducing the benchmark folds

Use os_make_benchmark_folds(y, group_id, k = 5L, seed = 101L) with the archived row order and seed to reproduce the benchmark partitions. It shuffles groups within outcome and assigns them round-robin. os_make_grouped_stratified_folds() retains its existing row-count balancing algorithm and produces different partitions.

Explore related projects