Tooluniverse Metabolomics Analysis

Metabolomics Analysis

Comprehensive analysis of metabolomics data from metabolite identification through quantification, statistical analysis, pathway interpretation, and integration with other omics layers.

Domain Reasoning

Metabolomics quantification depends critically on normalization. Total ion current (TIC) normalization corrects for sample-loading variation and works well for global abundance changes; internal standard normalization is more accurate for targeted analysis where specific metabolite concentrations matter. Missing values in a peak table may reflect signal below the detection limit — not true absence — and should be imputed or handled explicitly rather than treated as zero. Failing to account for batch effects across instrument runs is a frequent source of spurious differential metabolites.

LOOK UP DON'T GUESS

• Metabolite identities: use Metabolite_search and Metabolite_get_info to confirm names, CIDs, and HMDB IDs; never assume identity from m/z alone.
• Pathway memberships: query KEGG, MetaCyc, or Reactome tools; do not list pathways from memory.
• Disease associations: retrieve from CTD via Metabolite_get_diseases; do not infer clinical relevance without database evidence.
• CV thresholds and QC criteria: apply the values defined in this workflow (CV < 30%, blank ratio > 3x); do not override with guesses.

---

When to Use This Skill

Triggers:

• User has metabolomics data (LC-MS, GC-MS, NMR)
• Questions about metabolite abundance or concentrations
• Differential metabolite analysis requests
• Metabolic pathway analysis
• Multi-omics integration with metabolomics
• Metabolic biomarker discovery
• Flux balance analysis or metabolic modeling
• Metabolite-enzyme correlation

Example Questions: 1. "Analyze this LC-MS metabolomics data for differential metabolites" 2. "Which metabolic pathways are dysregulated between conditions?" 3. "Identify metabolite biomarkers for disease classification" 4. "Correlate metabolite levels with enzyme expression" 5. "Perform pathway enrichment for differential metabolites" 6. "Integrate metabolomics with transcriptomics data"

---

Core Capabilities

---

Workflow Overview

Input: Metabolomics Data (Peak Table or Spectra)
    |
    v
Phase 1: Data Import & Metabolite Identification
    |-- Load peak table or process raw spectra
    |-- Match features to HMDB, KEGG (accurate mass +/- 5 ppm)
    |-- Confidence scoring (Level 1-4)
    |
    v
Phase 2: Quality Control & Filtering
    |-- CV in QC samples (<30%)
    |-- Blank subtraction (sample/blank > 3)
    |-- Remove features with >50% missing
    |
    v
Phase 3: Normalization
    |-- Sample-wise: TIC, PQN, or internal standards
    |-- Transformation: log2, Pareto, or auto-scaling
    |-- Batch effect correction (if multi-batch)
    |
    v
Phase 4: Exploratory Analysis
    |-- PCA for sample clustering
    |-- PLS-DA for supervised separation
    |-- Outlier detection
    |
    v
Phase 5: Differential Analysis
    |-- t-test / ANOVA / Wilcoxon
    |-- Fold change + FDR correction
    |-- Volcano plots, heatmaps
    |
    v
Phase 6: Pathway Analysis
    |-- Metabolite set enrichment (MSEA)
    |-- KEGG/Reactome pathway mapping
    |-- Pathway topology (hub/bottleneck metabolites)
    |
    v
Phase 7: Multi-Omics Integration
    |-- Metabolite-enzyme Spearman correlation
    |-- Pathway-level concordance scoring
    |-- Metabolic flux inference
    |
    v
Phase 8: Generate Report
    |-- Summary statistics, differential metabolites
    |-- Pathway diagrams, biomarker panel

---

Phase Summaries

Phase 1: Data Import & Identification

Load peak tables (CSV/TSV) or process raw spectra (mzML). Match features to HMDB by accurate mass (+/- 5 ppm). Assign confidence levels: L1 (standard match), L2 (MS/MS), L3 (mass only), L4 (unknown).

Phase 2: Quality Control

Assess CV in QC samples (reject >30%), compute blank ratios (keep >3x blank), filter features with >50% missing values. Check internal standard recovery (95-105% acceptable).

Phase 3: Normalization

Three methods available: TIC (simple, assumes similar total abundance), PQN (robust to large changes, recommended), Internal Standard (most accurate with spiked standards). Follow with log2 transform or Pareto scaling.

Phase 4: Exploratory Analysis

PCA reveals sample grouping and batch effects. PLS-DA provides supervised separation (report R2 and Q2 for model quality). Flag and investigate outliers.

Phase 5: Differential Analysis

Welch's t-test (two groups) or ANOVA (multiple groups) with Benjamini-Hochberg FDR correction. Significance thresholds: adj. p < 0.05 and |log2FC| > 1.0.

Phase 6: Pathway Analysis

Map differential metabolites to KEGG compound IDs. Perform MSEA for pathway enrichment. Consider topology: metabolites at pathway hubs (high degree/betweenness centrality) have greater impact.

Phase 7: Multi-Omics Integration

Correlate metabolite levels with enzyme expression (Spearman). Expected: substrate-enzyme negative correlation (consumption), product-enzyme positive correlation (production). Score pathway dysregulation using combined metabolite + gene evidence.

Phase 8: Report

See report_template.md for full example output.

---

Integration with ToolUniverse

---

Quantified Minimums

---

Limitations

• Identification: Many features remain unidentified (Level 4)
• Coverage: Cannot detect all metabolites (depends on method)
• Quantification: Relative abundance (not absolute without standards)
• Isomers: Difficult to distinguish structural isomers
• Ion suppression: Matrix effects can affect quantification
• Dynamic range: Limited compared to targeted methods

---

References

Methods:

• MetaboAnalyst: https://doi.org/10.1093/nar/gkab382
• XCMS: https://doi.org/10.1021/ac051437y
• MSEA: https://doi.org/10.1186/1471-2105-11-395

Databases:

• HMDB: https://hmdb.ca
• KEGG Compound: https://www.genome.jp/kegg/compound/
• Reactome: https://reactome.org

---

Reference Files

• code_examples.md - Python code for all phases (data loading, QC, normalization, statistics, pathway analysis)
• report_template.md - Full example report (LC-MS disease vs control)

Metabolomics Analysis Code Examples

Phase 1: Data Import & Metabolite Identification

Data Loading

def load_metabolomics_data(file_path, file_type='peak_table'):
    """
    Load metabolomics data.
    file_type: 'peak_table' (CSV/TSV), 'mzml' (raw LC-MS), 'nmr'
    """
    import pandas as pd

    if file_type == 'peak_table':
        data = pd.read_csv(file_path, index_col=0)
        # Rows = samples, Columns = metabolites
        return data
    elif file_type == 'mzml':
        # Process raw MS data (requires pymzml)
        pass

Metabolite Identification

def identify_metabolites(feature_data, mass_list, rt_list=None):
    """Match features to metabolite databases (HMDB, KEGG, PubChem)."""
    from tooluniverse import ToolUniverse
    tu = ToolUniverse()

    identified_metabolites = []
    for i, mass in enumerate(mass_list):
        # Mass-based annotation via Metabolomics Workbench. NOTE: `mass` is the
        # observed m/z and the result `data` is a TSV string (parse it); set
        # `adduct` to your ionization mode (M+H positive, M-H negative).
        hmdb_result = tu.run_one_function({
            "name": "MetabolomicsWorkbench_search_by_mz",
            "arguments": {"mz_value": mass, "adduct": "M+H", "tolerance": 0.005}
        })
        if hmdb_result and hmdb_result.get('status') == 'success':
            rows = _parse_tsv(hmdb_result['data'])  # split TSV into list[dict]
            if rows:
                identified_metabolites.append({
                    'feature_id': i,
                    'metabolite_name': rows[0].get('Name'),
                    'formula': rows[0].get('Formula'),
                    'matched_mz': rows[0].get('Matched m/z'),
                    'confidence': calculate_confidence(rows[0])
                })
    return identified_metabolites

Confidence scoring:

Level 1: Confirmed with authentic standard (MS + RT match)
Level 2: Probable structure (accurate mass + MS/MS)
Level 3: Tentative match (accurate mass only)
Level 4: Unknown metabolite

Phase 2: Quality Control & Filtering

def metabolomics_qc(data, sample_metadata):
    """
    QC metrics: CV in QC samples (<30%), blank ratios (>3x),
    missing values (<50%), total ion current per sample.
    """
    qc_samples = sample_metadata['sample_type'] == 'QC'
    qc_data = data[qc_samples]
    cv_per_metabolite = qc_data.std() / qc_data.mean()

    blank_samples = sample_metadata['sample_type'] == 'Blank'
    blank_data = data[blank_samples]
    blank_means = blank_data.mean()
    sample_means = data[~blank_samples & ~qc_samples].mean()
    blank_ratio = sample_means / blank_means

    keep_metabolites = blank_ratio > 3
    missing_per_metabolite = (data == 0).sum() / data.shape[0]
    keep_metabolites &= (missing_per_metabolite < 0.5)

    return data.loc[:, keep_metabolites]

Phase 3: Normalization

Total Ion Current (TIC)

def normalize_tic(data):
    """Assumes total metabolite abundance is similar across samples."""
    tic = data.sum(axis=1)
    median_tic = tic.median()
    norm_factors = median_tic / tic
    return data.multiply(norm_factors, axis=0)

Probabilistic Quotient Normalization (PQN)

def normalize_pqn(data, reference_sample=None):
    """More robust than TIC to large metabolite changes."""
    import numpy as np

    if reference_sample is None:
        reference = data.median(axis=0)
    else:
        reference = data.loc[reference_sample]

    quotients = data.div(reference, axis=1)
    norm_factors = quotients.median(axis=1)
    return data.div(norm_factors, axis=0)

Internal Standard Normalization

def normalize_internal_standard(data, is_metabolite):
    """Most accurate if added before sample processing."""
    is_abundance = data[is_metabolite]
    norm_factors = is_abundance.median() / is_abundance
    normalized = data.multiply(norm_factors, axis=0)
    return normalized.drop(columns=[is_metabolite])

Transformation

def transform_data(data, method='log'):
    """Methods: 'log' (log2), 'pareto' (mean-center/sqrt(std)), 'auto' (z-score)."""
    import numpy as np

    if method == 'log':
        return np.log2(data + 1)
    elif method == 'pareto':
        mean, std = data.mean(axis=0), data.std(axis=0)
        return (data - mean) / np.sqrt(std)
    elif method == 'auto':
        mean, std = data.mean(axis=0), data.std(axis=0)
        return (data - mean) / std

Phase 4: Exploratory Analysis

def perform_pca_metabolomics(data, sample_groups):
    """PCA for sample clustering and outlier detection."""
    from sklearn.decomposition import PCA
    import matplotlib.pyplot as plt

    pca = PCA(n_components=2)
    pca_result = pca.fit_transform(data)

    plt.figure(figsize=(8, 6))
    for group in sample_groups.unique():
        mask = sample_groups == group
        plt.scatter(pca_result[mask, 0], pca_result[mask, 1], label=group)
    plt.xlabel(f'PC1 ({pca.explained_variance_ratio_[0]:.1%})')
    plt.ylabel(f'PC2 ({pca.explained_variance_ratio_[1]:.1%})')
    plt.legend()
    plt.title('PCA - Metabolomics Data')


def plsda_analysis(X, y, n_components=2):
    """PLS-DA for supervised separation (better than PCA for classification)."""
    from sklearn.cross_decomposition import PLSRegression
    from sklearn.preprocessing import LabelEncoder

    le = LabelEncoder()
    y_encoded = le.fit_transform(y)
    pls = PLSRegression(n_components=n_components)
    X_pls = pls.fit_transform(X, y_encoded)[0]
    return X_pls

Phase 5: Differential Metabolite Analysis

def differential_metabolites(data, group1_samples, group2_samples):
    """Identify differential metabolites with FDR correction."""
    from scipy import stats
    from statsmodels.stats.multitest import multipletests
    import numpy as np
    import pandas as pd

    results = []
    for metabolite in data.columns:
        group1 = data.loc[group1_samples, metabolite]
        group2 = data.loc[group2_samples, metabolite]
        fold_change = group2.mean() / group1.mean()
        t_stat, p_value = stats.ttest_ind(group1, group2, equal_var=False)
        results.append({
            'metabolite': metabolite,
            'fold_change': fold_change,
            'log2FC': np.log2(fold_change),
            'mean_group1': group1.mean(),
            'mean_group2': group2.mean(),
            'p_value': p_value,
            't_statistic': t_stat
        })

    results_df = pd.DataFrame(results)
    results_df['adj_p_value'] = multipletests(results_df['p_value'], method='fdr_bh')[1]
    results_df['significant'] = (
        (results_df['adj_p_value'] < 0.05) & (np.abs(results_df['log2FC']) > 1.0)
    )
    return results_df

Phase 6: Metabolic Pathway Analysis

def pathway_enrichment_metabolites(metabolite_list, organism='human'):
    """Perform pathway enrichment using KEGG metabolic pathways."""
    from tooluniverse import ToolUniverse
    tu = ToolUniverse()

    # Metabolite-set pathway enrichment takes metabolite NAMES directly (it
    # resolves them to KEGG internally) -- no separate ID lookup needed.
    enrichment = tu.run_one_function({
        "name": "MetaboAnalyst_pathway_enrichment",
        "arguments": {"metabolites": metabolite_list, "organism": organism}
    })
    return enrichment

Phase 7: Multi-Omics Integration

def correlate_metabolite_enzyme(metabolite_data, enzyme_expression):
    """Correlate metabolite levels with enzyme expression (Spearman)."""
    from scipy.stats import spearmanr

    correlations = {}
    for metabolite in metabolite_data.columns:
        enzymes = find_metabolite_enzymes(metabolite)
        for enzyme in enzymes:
            if enzyme in enzyme_expression.index:
                r, p = spearmanr(
                    metabolite_data[metabolite],
                    enzyme_expression.loc[enzyme]
                )
                correlations[f'{metabolite}_{enzyme}'] = {
                    'r': r, 'p': p,
                    'relationship': 'product' if r > 0 else 'substrate'
                }
    return correlations

Metabolomics Analysis Report Template

Example Report: LC-MS Disease vs Control

# Metabolomics Analysis Report

## Dataset Summary
- **Platform**: LC-MS/MS (Orbitrap)
- **Method**: Untargeted metabolomics
- **Samples**: 40 (20 disease, 20 control)
- **Metabolites Identified**: 324 (Level 1/2 confidence)
- **Metabolites Quantified**: 298 (after QC)

## Quality Control
- **CV in QC samples**: 18% median (acceptable: <30%)
- **Blank ratios**: All metabolites > 3x blank signal
- **Missing values**: 8% average per metabolite
- **Internal standard**: Recovery 95-105% across samples

## Normalization
- **Method**: Probabilistic Quotient Normalization (PQN)
- **Transformation**: log2
- **Batch correction**: Not required (single batch)

## Exploratory Analysis
- **PCA**: Clear separation between groups (PC1: 28%, PC2: 18%)
- **PLS-DA**: Excellent discrimination (R2=0.89, Q2=0.75)
- **Outliers**: 1 sample removed (technical failure)

## Differential Metabolites
- **Significant metabolites**: 87 (adj. p < 0.05, |log2FC| > 1)
  - Increased: 52 metabolites
  - Decreased: 35 metabolites

### Top Increased Metabolites
1. **Lactate** (log2FC=3.2, p=1e-12) - Glycolysis
2. **Glutamine** (log2FC=2.8, p=1e-10) - Amino acid metabolism
3. **Palmitate** (log2FC=2.5, p=1e-9) - Fatty acid synthesis

### Top Decreased Metabolites
1. **Citrate** (log2FC=-2.9, p=1e-11) - TCA cycle
2. **ATP** (log2FC=-2.3, p=1e-9) - Energy metabolism
3. **NAD+** (log2FC=-2.1, p=1e-8) - Redox balance

## Pathway Enrichment
### Top Dysregulated Pathways
1. **Glycolysis/Gluconeogenesis** (p=1e-15)
   - 12 metabolites: glucose, pyruvate, lactate, etc.
   - Direction: Increased flux to lactate (Warburg effect)
2. **TCA Cycle** (p=1e-12)
   - 8 metabolites: citrate, succinate, malate, etc.
   - Direction: Decreased activity
3. **Glutaminolysis** (p=1e-10)
   - 6 metabolites: glutamine, glutamate, a-KG, etc.
   - Direction: Increased glutamine consumption

## Multi-Omics Integration
### Metabolite-Enzyme Correlations
- **LDHA (lactate dehydrogenase)**
  - Expression: 3.5-fold increased (RNA + protein)
  - Lactate: 3.2-fold increased
  - Correlation: r=0.85 (p<0.001) - Concordant upregulation
- **IDH1 (isocitrate dehydrogenase)**
  - Expression: 2.1-fold decreased
  - Citrate: 2.9-fold decreased
  - Correlation: r=0.78 (p<0.001) - TCA cycle suppression

### Metabolic Phenotype
Integration with RNA-seq and proteomics reveals:
- **Warburg effect**: Shift from oxidative to glycolytic metabolism
- **Glutamine addiction**: Increased glutaminolysis for anaplerosis
- **Redox imbalance**: Decreased NAD+/NADH ratio, oxidative stress

## Biomarker Discovery
### Top 10 Metabolites for Classification
Random Forest model (10-fold CV):
- **AUC**: 0.96 +/- 0.03
- **Accuracy**: 92%

**Biomarker Panel**:
1. Lactate
2. Glutamine
3. Citrate
4. ATP
5. Palmitate
6. Pyruvate
7. Succinate
8. NAD+
9. a-ketoglutarate
10. Glucose-6-phosphate

## Biological Interpretation
Metabolomics reveals fundamental metabolic reprogramming in disease state:

1. **Glycolytic switch**: Increased glycolysis with lactate accumulation despite oxygen
   availability (Warburg effect), driven by LDHA upregulation.
2. **TCA cycle suppression**: Decreased citrate and TCA intermediates, consistent with
   IDH1 downregulation. Shunts carbon to biosynthesis.
3. **Glutamine dependence**: Elevated glutamine consumption provides alternative carbon
   source for anaplerosis and NADPH for biosynthesis.
4. **Biosynthetic activation**: Increased palmitate indicates active fatty acid synthesis,
   supporting membrane production for proliferation.
5. **Energy stress**: Despite active glycolysis, ATP levels are decreased, suggesting high
   energy demand outpacing production.

## Clinical Relevance
- **Therapeutic targets**: LDHA inhibitors, glutaminase inhibitors
- **Biomarkers**: Lactate/citrate ratio as metabolic activity marker
- **Drug response**: Metabolic phenotype may predict sensitivity to metabolic inhibitors