metabolomics-de

name	metabolomics-de
description	Metabolomics differential analysis using univariate tests (t-test, FDR), multivariate methods (PCA, PLS-DA, OPLS-DA, sPLS-DA), Random Forest, and ROC analysis for biomarker discovery.
version	0.1.0
author	OmicsClaw
license	MIT
tags	["metabolomics","differential","PLS-DA","volcano","biomarker","ROC"]
metadata	{"omicsclaw":{"domain":"metabolomics","emoji":"📈","trigger_keywords":["metabolomics differential","PLS-DA","volcano plot","biomarker","OPLS-DA"],"allowed_extra_flags":["--group-a-prefix","--group-b-prefix"],"legacy_aliases":["met-diff"],"saves_h5ad":false}}

📈 Metabolomics Differential Analysis

Univariate and multivariate statistical analysis for identifying differentially abundant metabolites and biomarker discovery.

Core Capabilities

1. Univariate testing: t-test, Wilcoxon, ANOVA with BH FDR correction
2. Fold change: Log2FC calculation with volcano plot visualization
3. PCA: Unsupervised exploration and outlier detection
4. PLS-DA / sPLS-DA: Supervised classification with VIP scores and sparse feature selection
5. OPLS-DA: Orthogonal PLS-DA with S-plot for predictive features
6. Random Forest: Non-linear feature importance ranking
7. ROC analysis: Diagnostic performance evaluation of candidate biomarkers

CLI Reference

python omicsclaw.py run met-diff --demo
python omicsclaw.py run met-diff --input <feature_table.csv> --output <dir>

Algorithm / Methodology

Univariate Analysis (R)

library(tidyverse)

data <- read.csv('normalized_data.csv', row.names = 1)
groups <- factor(read.csv('sample_info.csv')$group)

# T-test for each feature
ttest_results <- apply(data, 2, function(x) {
    test <- t.test(x ~ groups)
    c(pvalue = test$p.value,
      fc = mean(x[groups == levels(groups)[2]]) - mean(x[groups == levels(groups)[1]])

ttest_results  as.data.frametttest_results
ttest_resultsfdr  p.adjustttest_resultspvalue method  

sig_features  ttest_resultsttest_resultsfdr   

Volcano Plot (R)

library(ggplot2)

results$significant <- results$fdr < 0.05 & abs(results$log2fc) > 1

ggplot(results, aes(x = log2fc, y = -log10(pvalue), color = significant)) +
    geom_point(alpha = 0.6) +
    scale_color_manual(values = c('gray', 'red')) +
    geom_hline(yintercept = -log10(0.05), linetype = 'dashed') +
    geom_vline(xintercept = c  linetype   
    labsx   y   
    theme_bw

PCA (R)

library(pcaMethods)

pca_result <- pca(data, nPcs = 5, method = 'ppca')
scores <- as.data.frame(scores(pca_result))
scores$group <- groups

ggplot(scores, aes(x = PC1, y = PC2, color = group)) +
    geom_point(size = 3) +
    stat_ellipse(level = 0.95) +
    labs(x = paste0('PC1 (', round(pca_result@R2[1] * 100, 1), '%)'
         y  paste0 pca_resultR2     
    theme_bw

PLS-DA with VIP Scores (mixOmics)

library(mixOmics)

plsda_result <- plsda(as.matrix(data), groups, ncomp = 3)

# Cross-validation
perf_plsda <- perf(plsda_result, validation = 'Mfold', folds = 5, nrepeat = 50)
ncomp_opt <- perf_plsda$choice.ncomp['BER', 'centroids.dist']

# Final model + VIP scores
final_plsda <- plsda(as.matrix(data), groups, ncomp = ncomp_opt)
plotIndiv(final_plsda, group = groups, ellipse = TRUE, legend = TRUE)

vip <- vip(final_plsda)
top_vip <- sort(vip ncomp_opt decreasing  

sPLS-DA (Sparse — Feature Selection)

tune_splsda <- tune.splsda(as.matrix(data), groups, ncomp = 3,
                            validation = 'Mfold', folds = 5, nrepeat = 50,
                            test.keepX = c(5, 10, 20, 50, 100))
optimal_keepX <- tune_splsda$choice.keepX

splsda_result <- splsda(as.matrix(data), groups, ncomp = ncomp_opt, keepX = optimal_keepX)
selected_features <- selectVar(splsda_result, comp = 1)$name

OPLS-DA (ropls)

library(ropls)

oplsda <- opls(data, groups, predI = 1, orthoI = NA)
plot(oplsda, typeVc = 'x-score')  # Scores plot
plot(oplsda, typeVc = 'x-loading')  # S-plot

vip_scores <- getVipVn(oplsda)
top_vip <- sort(vip_scores, decreasing = TRUE)[1:20]

Random Forest Feature Importance

library(randomForest)

rf_model <- randomForest(x = data, y = groups, importance = TRUE, ntree = 500)
importance <- importance(rf_model)
top_features <- rownames(importance)[order(importance[, 'MeanDecreaseAccuracy'], decreasing = TRUE)[1:20]]
varImpPlot(rf_model, n.var = 20)

ROC Analysis (pROC)

library(pROC)

top_feature <- 'feature_123'
roc_result <- roc(groups, data[, top_feature])
plot(roc_result, main = paste('AUC =', round(auc(roc_result), 3)))

Heatmap

library(pheatmap)

top_features <- rownames(sig_features)[1:50]
annotation_row <- data.frame(Group = groups)
rownames(annotation_row) <- rownames(data)

pheatmap(t(data[, top_features]), annotation_col = annotation_row,
         scale = 'row', clustering_method = 'ward.D2')

Parameters

Parameter	Default	Description
--method	ttest	Statistical test method
--fdr-cutoff	0.05	FDR significance threshold
--fc-cutoff	1.0	Log2FC threshold
--multivariate	pca	pca, plsda, oplsda, splsda

Why This Exists

• Without it: Complex metabolomic variance is hard to parse; minor batch effects can mimic true biological signal
• With it: Powerful multivariate (PLS-DA) and univariate testing models separate true biomarkers from noise
• Why OmicsClaw: Automates strict R-based multidimensional workflows with single Python calls

Workflow

1. Calculate: Log transformations and scaling.
2. Execute: Run PCA/OPLS-DA and extract variable importance (VIP) scores.
3. Assess: Execute parallel univariate FDR-corrected models.
4. Generate: Output candidate biomarker lists.
5. Report: Synthesize VIP plots, score plots, and hierarchical heatmaps.

Example Queries

• "Run PLS-DA and find significant metabolites"
• "Perform differential analysis on this normalized metabolomics matrix"

Output Structure

output_directory/
├── report.md
├── result.json
├── statistics.csv
├── figures/
│   ├── plsda_scores.png
│   ├── vip_plot.png
│   └── volcano_plot.png
├── tables/
│   └── differential_results.csv
└── reproducibility/
    ├── commands.sh
    ├── requirements.txt
    └── checksums.sha256

Safety

• Local-first: Strict offline processing without external upload.
• Disclaimer: Requires OmicsClaw reporting structures and disclaimers.
• Audit trail: Hyperparameters and operational flow states are logged fully.

Integration with Orchestrator

Trigger conditions:

• Automatically invoked dynamically based on tool metadata and user intent matching.

Chaining partners:

• met-normalize — Upstream data scaling
• met-annotate — Downstream identification of biomarkers

Version Compatibility

Reference examples tested with: mixOmics 6.24+, ropls 1.32+

Dependencies

Required: numpy, pandas, scipy, scikit-learn Optional: mixOmics (R), ropls (R), randomForest (R), pROC (R), pheatmap (R)

Citations

• mixOmics — Rohart et al., PLoS Computational Biology 2017
• ropls — Thévenot et al., Journal of Proteome Research 2015
• MetaboAnalyst — Pang et al., Nucleic Acids Research 2021

Related Skills

• met-normalize — Data normalization before analysis
• met-pathway — Pathway enrichment of significant metabolites
• xcms-preprocess — Feature extraction upstream