Seurat Integration Guide
Zaoqu Liu
2026-01-24
Source:vignettes/seurat-integration.Rmd
seurat-integration.RmdIntroduction
SCEVAN seamlessly integrates with Seurat, the most widely used R toolkit for single-cell analysis. This guide demonstrates how to:
- Extract count matrices from Seurat objects
- Add SCEVAN results back to Seurat
- Visualize CNA information in Seurat plots
SCEVAN supports both Seurat v4 and v5 data structures.
Extracting Data from Seurat
Using getCountMtxFromSeurat()
SCEVAN provides a helper function that works with both Seurat v4 and v5:
# Load your Seurat object
seurat_obj <- readRDS("path/to/seurat_object.rds")
# Extract raw counts (v4/v5 compatible)
count_mtx <- getCountMtxFromSeurat(
seurat_obj,
assay = "RNA", # Default assay
layer = "counts" # For Seurat v5
)
# Check dimensions
dim(count_mtx)Manual Extraction
For more control, you can extract manually:
# Seurat v4
count_mtx <- as.matrix(GetAssayData(seurat_obj, slot = "counts"))
# Seurat v5
count_mtx <- as.matrix(GetAssayData(seurat_obj, layer = "counts"))Running SCEVAN
Standard Analysis
# Run SCEVAN pipeline
results <- pipelineCNA(
count_mtx,
sample = "MySample",
par_cores = 4,
organism = "human",
SUBCLONES = TRUE
)
# View results
head(results)Using Seurat Clusters as Prior
You can use Seurat’s clustering to identify potential normal cells:
# Get cluster assignments
clusters <- Idents(seurat_obj)
# Identify clusters with normal cell markers
# (e.g., immune clusters, fibroblast clusters)
normal_clusters <- c("T_cells", "Macrophages", "Fibroblasts")
potential_normals <- names(clusters)[clusters %in% normal_clusters]
# Run SCEVAN with known normals
results <- pipelineCNA(
count_mtx,
sample = "MySample_guided",
norm_cell = potential_normals,
SUBCLONES = TRUE
)Adding Results to Seurat
Add Cell Classifications
# Ensure cell names match
results_matched <- results[colnames(seurat_obj), , drop = FALSE]
# Add all SCEVAN columns
seurat_obj <- AddMetaData(seurat_obj, metadata = results_matched)
# Check added columns
head(seurat_obj@meta.data)Add CNA Scores
# Load CNA matrix
load("./output/MySample_CNAmtx.RData")
load("./output/MySample_count_mtx_annot.RData")
# Calculate global CNA score per cell
global_cna <- colMeans(abs(CNA_mtx_relat))
seurat_obj$CNA_score <- global_cna[colnames(seurat_obj)]
# Calculate chromosome-specific scores
chr7_cna <- colMeans(CNA_mtx_relat[count_mtx_annot$seqnames == 7, ])
seurat_obj$chr7_cna <- chr7_cna[colnames(seurat_obj)]
chr10_cna <- colMeans(CNA_mtx_relat[count_mtx_annot$seqnames == 10, ])
seurat_obj$chr10_cna <- chr10_cna[colnames(seurat_obj)]Visualization in Seurat
Feature Plot with CNA Scores
# Global CNA burden
FeaturePlot(
seurat_obj,
features = "CNA_score",
cols = c("lightgray", "red")
) +
ggtitle("Global CNA Score")
# Chromosome-specific
FeaturePlot(
seurat_obj,
features = c("chr7_cna", "chr10_cna"),
cols = c("blue", "white", "red"),
min.cutoff = -0.3,
max.cutoff = 0.3
)Advanced Integration
Region-Specific CNA Visualization
# Define genomic regions of interest
regions <- list(
EGFR = c(chr = 7, start = 55019017, end = 55211628),
CDKN2A = c(chr = 9, start = 21967751, end = 21995301),
PTEN = c(chr = 10, start = 87863113, end = 87971930)
)
# Calculate CNA for each region
for(gene in names(regions)) {
r <- regions[[gene]]
region_cna <- colMeans(CNA_mtx_relat[
count_mtx_annot$seqnames == r["chr"] &
count_mtx_annot$start >= r["start"] &
count_mtx_annot$end <= r["end"], , drop = FALSE
])
seurat_obj[[paste0(gene, "_cna")]] <- region_cna[colnames(seurat_obj)]
}
# Visualize
FeaturePlot(seurat_obj, features = c("EGFR_cna", "PTEN_cna"))Combined Visualization
# Create combined plot with classification and clustering
p1 <- DimPlot(seurat_obj, group.by = "class") +
ggtitle("SCEVAN Classification") + NoLegend()
p2 <- DimPlot(seurat_obj, group.by = "seurat_clusters") +
ggtitle("Seurat Clusters") + NoLegend()
p3 <- FeaturePlot(seurat_obj, features = "CNA_score") +
ggtitle("CNA Score")
p4 <- DimPlot(seurat_obj, group.by = "subclone", na.value = "lightgray") +
ggtitle("Subclones")
# Combine
(p1 | p2) / (p3 | p4)Downstream Analysis
Differential Expression by Subclone
# Set identity to subclone
Idents(seurat_obj) <- "subclone"
# Find markers for each subclone
subclone_markers <- FindAllMarkers(
seurat_obj,
only.pos = TRUE,
min.pct = 0.25,
logfc.threshold = 0.25
)
# Top markers per subclone
top_markers <- subclone_markers %>%
group_by(cluster) %>%
top_n(10, avg_log2FC)
print(top_markers)Pathway Analysis per Subclone
# Use marker genes for pathway analysis
# Example with enrichR or clusterProfiler
# Get genes from subclone 1
s1_genes <- subclone_markers %>%
filter(cluster == 1, p_val_adj < 0.05) %>%
pull(gene)
# Run enrichment (requires clusterProfiler)
# enrichGO(s1_genes, OrgDb = org.Hs.eg.db, keyType = "SYMBOL")Creating New Seurat Object
From Scratch with SCEVAN Results
# Create new Seurat object with SCEVAN metadata
seurat_new <- CreateSeuratObject(
counts = count_mtx,
meta.data = results,
project = "SCEVAN_analysis"
)
# Standard Seurat workflow
seurat_new <- NormalizeData(seurat_new)
seurat_new <- FindVariableFeatures(seurat_new)
seurat_new <- ScaleData(seurat_new)
seurat_new <- RunPCA(seurat_new)
seurat_new <- RunUMAP(seurat_new, dims = 1:30)
# Visualize with SCEVAN results
DimPlot(seurat_new, group.by = "class")Best Practices
Memory Management
# For large datasets, process in chunks
# or use sparse matrices
# Check memory usage
format(object.size(seurat_obj), units = "GB")
# Remove intermediate objects
rm(CNA_mtx_relat)
gc()Troubleshooting
Session Info
sessionInfo()
#> R version 4.4.0 (2024-04-24)
#> Platform: aarch64-apple-darwin20
#> Running under: macOS 15.6.1
#>
#> Matrix products: default
#> BLAS: /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/lib/libRblas.0.dylib
#> LAPACK: /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/lib/libRlapack.dylib; LAPACK version 3.12.0
#>
#> locale:
#> [1] C
#>
#> time zone: Asia/Shanghai
#> tzcode source: internal
#>
#> attached base packages:
#> [1] stats graphics grDevices utils datasets methods base
#>
#> loaded via a namespace (and not attached):
#> [1] digest_0.6.39 desc_1.4.3 R6_2.6.1 fastmap_1.2.0
#> [5] xfun_0.56 cachem_1.1.0 knitr_1.51 htmltools_0.5.9
#> [9] rmarkdown_2.30 lifecycle_1.0.5 cli_3.6.5 sass_0.4.10
#> [13] pkgdown_2.1.3 textshaping_1.0.4 jquerylib_0.1.4 systemfonts_1.3.1
#> [17] compiler_4.4.0 tools_4.4.0 ragg_1.5.0 bslib_0.9.0
#> [21] evaluate_1.0.5 yaml_2.3.12 otel_0.2.0 jsonlite_2.0.0
#> [25] rlang_1.1.7 fs_1.6.6 htmlwidgets_1.6.4