--- title: "Single-Cell RNA-seq Metabolic Flux Analysis" author: "Zaoqu Liu" date: "`r Sys.Date()`" output: rmarkdown::html_vignette: toc: true toc_depth: 3 vignette: > %\VignetteIndexEntry{Single-Cell RNA-seq Metabolic Flux Analysis} %\VignetteEngine{knitr::rmarkdown} %\VignetteEncoding{UTF-8} --- ```{r setup, include=FALSE} knitr::opts_chunk$set( collapse = TRUE, comment = "#>", fig.width = 7, fig.height = 5, warning = FALSE, message = FALSE ) ``` ## Overview This vignette demonstrates metabolic flux analysis for single-cell RNA-seq data using METAFLUX's community modeling approach. ## Biological Context ### Tumor Microenvironment Metabolism Single-cell metabolic analysis reveals: - **Metabolic heterogeneity** within cell populations - **Metabolic cooperation** between cell types - **Competition** for shared nutrients - **Metabolic crosstalk** via secreted metabolites ### Community Modeling Concept ``` ┌─────────────────────────────────────────────────────────────────┐ │ Tumor Microenvironment │ │ ┌──────────┐ ┌──────────┐ ┌──────────┐ │ │ │ Tumor │ ←→ │ T cell │ ←→ │ Macro- │ │ │ │ Cells │ │ │ │ phage │ │ │ └────┬─────┘ └────┬─────┘ └────┬─────┘ │ │ │ │ │ │ │ └────────────────┼────────────────┘ │ │ │ │ │ ┌────┴────┐ │ │ │ Shared │ │ │ │ Medium │ │ │ └─────────┘ │ └─────────────────────────────────────────────────────────────────┘ ``` ## Data Preparation ### Input Requirements | Component | Format | Description | |-----------|--------|-------------| | Seurat object | v4/v5 | Normalized counts in RNA assay | | Cell annotations | Metadata column | Cell type labels | | Medium profile | Data frame | Available nutrients | ### Load Example Data ```{r load-sc-data, eval=FALSE} library(METAFLUX) library(Seurat) # Load example single-cell data data("sc_test_example") data("human_blood") # Inspect Seurat object sc_test_example table(sc_test_example$Cell_type) ``` ## Step 1: Bootstrap Aggregation Single-cell data is sparse. Bootstrap aggregation reduces noise while preserving biological variation. ```{r bootstrap, eval=FALSE} # Calculate bootstrap-averaged expression avg_expr <- calculate_avg_exp( myseurat = sc_test_example, myident = "Cell_type", # Column with cell type labels n_bootstrap = 100, # Number of bootstrap iterations seed = 42 # For reproducibility ) # Check output dimensions dim(avg_expr) # Rows: genes # Columns: cell_type_1, cell_type_2, ..., cell_type_1, ... (repeated n_bootstrap times) ``` ### Bootstrap Algorithm ``` For each bootstrap iteration i: 1. Resample cells within each cell type (with replacement) 2. Compute average expression per cell type 3. Store result Output: Matrix of [genes × (n_celltypes × n_bootstrap)] ``` ## Step 2: Calculate MRAS ```{r sc-mras, eval=FALSE} # Calculate reaction activity scores mras <- calculate_reaction_score(avg_expr) # Verify dimensions dim(mras) # 13082 × (n_celltypes × n_bootstrap) ``` ## Step 3: Define Cell Type Fractions Cell type proportions weight the community model: ```{r fractions, eval=FALSE} # Calculate fractions from data cell_counts <- table(sc_test_example$Cell_type) fractions <- as.numeric(cell_counts / sum(cell_counts)) names(fractions) <- names(cell_counts) # Verify: must sum to 1 sum(fractions) # Should be 1 print(fractions) ``` ## Step 4: Community Flux Calculation ```{r sc-flux, eval=FALSE} # Number of cell types num_cell <- length(unique(sc_test_example$Cell_type)) # Compute community metabolic fluxes flux <- compute_sc_flux( num_cell = num_cell, fraction = fractions, fluxscore = mras, medium = human_blood, parallel = TRUE, # Enable parallel computing n_cores = 4, # Number of CPU cores use_cpp = TRUE # Use C++ acceleration ) # Check output dim(flux) head(rownames(flux), 20) ``` ### Output Structure The flux matrix contains: ``` celltype 1 HMR_0001 # Cell type 1 reactions celltype 1 HMR_0002 ... celltype 2 HMR_0001 # Cell type 2 reactions celltype 2 HMR_0002 ... external_medium HMR_XXXX # External exchange reactions ``` ## Step 5: Analyze Results ### Extract Cell Type-Specific Fluxes ```{r extract-flux, eval=FALSE} # Get reaction names reaction_names <- rownames(flux) # Extract fluxes for specific cell type get_celltype_flux <- function(flux_matrix, celltype_num) { pattern <- paste0("^celltype ", celltype_num, " ") idx <- grep(pattern, rownames(flux_matrix)) result <- flux_matrix[idx, , drop = FALSE] rownames(result) <- gsub(pattern, "", rownames(result)) return(result) } # Extract for each cell type tumor_flux <- get_celltype_flux(flux, 1) tcell_flux <- get_celltype_flux(flux, 2) ``` ### Compare Metabolic Phenotypes ```{r compare-flux, eval=FALSE} # Compare glucose uptake between cell types glucose_rxn <- "HMR_9034" tumor_glucose <- mean(tumor_flux[glucose_rxn, ]) tcell_glucose <- mean(tcell_flux[glucose_rxn, ]) cat(sprintf("Tumor glucose uptake: %.4f\n", tumor_glucose)) cat(sprintf("T cell glucose uptake: %.4f\n", tcell_glucose)) ``` ### External Medium Exchange ```{r external-exchange, eval=FALSE} # Extract external exchange reactions external_idx <- grep("^external_medium", rownames(flux)) external_flux <- flux[external_idx, , drop = FALSE] # Net nutrient consumption/production nutrient_balance <- rowMeans(external_flux) head(sort(nutrient_balance), 10) # Top uptakes tail(sort(nutrient_balance), 10) # Top secretions ``` ## Advanced Analysis ### Metabolic Competition Identify metabolites competed for by different cell types: ```{r competition, eval=FALSE} # Compare lactate exchange across cell types lactate_rxn <- "HMR_9135" lactate_by_celltype <- sapply(1:num_cell, function(ct) { ct_flux <- get_celltype_flux(flux, ct) mean(ct_flux[lactate_rxn, ]) }) names(lactate_by_celltype) <- names(cell_counts) print(lactate_by_celltype) ``` ### Metabolic Cooperation Identify metabolites exchanged between cell types: ```{r cooperation, eval=FALSE} # One cell type produces, another consumes producer <- which(lactate_by_celltype > 0) consumer <- which(lactate_by_celltype < 0) if (length(producer) > 0 && length(consumer) > 0) { cat(sprintf("%s produces lactate\n", names(lactate_by_celltype)[producer])) cat(sprintf("%s consumes lactate\n", names(lactate_by_celltype)[consumer])) } ``` ## Performance Tips ### Parallel Computing ```{r parallel-tips, eval=FALSE} # Check available cores parallel::detectCores() # Recommended: use physical cores - 1 n_cores <- max(1, parallel::detectCores(logical = FALSE) - 1) # Limit to 8 to prevent memory issues n_cores <- min(n_cores, 8) ``` ### Memory Management ```{r memory-tips, eval=FALSE} # For large datasets, reduce bootstrap iterations avg_expr <- calculate_avg_exp( myseurat = sc_test_example, myident = "Cell_type", n_bootstrap = 50, # Reduced from 100 seed = 42 ) # Or subsample cells first sc_subset <- subset(sc_test_example, downsample = 5000) ``` ## Biological Interpretation ### Warburg Effect ```{r warburg, eval=FALSE} # High glycolysis, low oxidative phosphorylation glycolysis_rxns <- c("HMR_4363", "HMR_4365", "HMR_4367") # Example oxphos_rxns <- c("HMR_6916", "HMR_6921") # Example tumor_flux <- get_celltype_flux(flux, 1) glycolysis_mean <- mean(abs(tumor_flux[glycolysis_rxns, ])) oxphos_mean <- mean(abs(tumor_flux[oxphos_rxns, ])) warburg_ratio <- glycolysis_mean / oxphos_mean cat(sprintf("Warburg ratio: %.2f\n", warburg_ratio)) ``` ### Immunometabolism ```{r immunometab, eval=FALSE} # T cell activation markers glutamine_rxn <- "HMR_9063" # Glutamine uptake fa_synthesis <- grep("fatty acid synthesis", rownames(tumor_flux), ignore.case = TRUE) tcell_flux <- get_celltype_flux(flux, 2) tcell_glutamine <- mean(abs(tcell_flux[glutamine_rxn, ])) tcell_fa <- mean(abs(tcell_flux[fa_synthesis, ])) cat(sprintf("T cell glutamine uptake: %.4f\n", tcell_glutamine)) cat(sprintf("T cell fatty acid synthesis: %.4f\n", tcell_fa)) ``` ## Output Files Save results for downstream analysis: ```{r save-output, eval=FALSE} # Save flux matrix write.csv(flux, "metabolic_flux_singlecell.csv") # Save cell type-specific fluxes for (ct in 1:num_cell) { ct_flux <- get_celltype_flux(flux, ct) filename <- sprintf("flux_celltype_%d.csv", ct) write.csv(ct_flux, filename) } ``` ## Next Steps 1. **Visualization**: See [Visualization Guide](visualization.html) 2. **Pathway analysis**: Aggregate fluxes by metabolic pathway 3. **Statistical testing**: Compare conditions 4. **Integration**: Combine with other omics data ## Author **Zaoqu Liu** - Package maintainer - GitHub: [Zaoqu-Liu](https://github.com/Zaoqu-Liu) - ORCID: [0000-0002-0452-742X](https://orcid.org/0000-0002-0452-742X) ## Session Information ```{r session, eval=FALSE} sessionInfo() ```