Scanpy实战:从数据预处理到细胞类型注释的完整单细胞分析流程

Scanpy实战:从数据预处理到细胞类型注释的完整单细胞分析流程 1. 单细胞分析入门为什么选择Scanpy单细胞RNA测序技术已经成为生命科学研究的重要工具它能让我们在单个细胞水平上观察基因表达模式。但面对海量的单细胞数据如何高效处理和分析成了关键问题。Scanpy作为Python生态中的单细胞分析利器凭借其出色的性能和易用性已经成为许多研究人员的首选工具。我第一次接触Scanpy是在处理一批10X Genomics数据时当时被它流畅的分析流程和丰富的可视化功能所吸引。相比R语言中的SeuratScanpy在Python环境下运行更加高效特别是处理大规模数据集时优势明显。而且Python生态中有大量优秀的机器学习库可以无缝衔接Scanpy的分析结果。Scanpy的核心优势在于完整的分析流程从原始数据读取到最终细胞注释提供一站式解决方案出色的性能基于AnnData数据结构内存占用低计算效率高丰富的可视化内置多种专业级绘图函数轻松生成发表级图表强大的扩展性支持harmonypy等整合工具方便多样本分析2. 数据准备与预处理2.1 10X数据读取与整理处理10X Genomics数据是单细胞分析的常见起点。Scanpy提供了专门的函数来读取这种格式的数据。我通常会先检查数据目录结构确保barcodes.tsv.gz、features.tsv.gz和matrix.mtx.gz三个文件齐全。import scanpy as sc # 读取单个10X样本 adata sc.read_10x_mtx( path/to/sample, # 样本目录路径 var_namesgene_symbols, # 使用基因符号作为变量名 cacheTrue # 缓存数据加速后续读取 )对于多样本项目我习惯先单独处理每个样本再进行合并。这样可以保留样本来源信息方便后续批次校正samples [sample1, sample2, sample3] # 样本列表 adatas [] for sample in samples: adata sc.read_10x_mtx(fpath/{sample}, var_namesgene_symbols) adata.obs[sample] sample # 添加样本标签 adatas.append(adata) # 合并所有样本 combined adatas[0].concatenate(adatas[1:], index_unique_)2.2 数据质量控制数据质量直接影响后续分析结果。我通常会从三个维度进行质控细胞层面过滤每个细胞检测到的基因数过低可能是空液滴过高可能是双细胞每个细胞的总UMI数反映测序深度基因层面过滤在至少一定数量细胞中表达的基因去除噪音基因线粒体基因比例高比例通常表示细胞状态不佳或破裂# 计算质控指标 adata.var[mt] adata.var_names.str.startswith(MT-) # 标记线粒体基因 sc.pp.calculate_qc_metrics(adata, qc_vars[mt], inplaceTrue) # 绘制质控图 sc.pl.violin(adata, [n_genes_by_counts, total_counts, pct_counts_mt], jitter0.4, multi_panelTrue) # 应用过滤 sc.pp.filter_cells(adata, min_genes200) # 至少检测到200个基因 sc.pp.filter_genes(adata, min_cells3) # 至少在3个细胞中表达 adata adata[adata.obs.pct_counts_mt 20, :] # 线粒体基因比例20%3. 数据标准化与特征选择3.1 数据标准化单细胞数据存在技术噪音和测序深度差异标准化是必不可少的步骤。Scanpy提供了几种标准化方法我最常用的是CPMCounts Per Million标准化# 保存原始计数 adata.layers[counts] adata.X.copy() # 标准化到中位数1e4 sc.pp.normalize_total(adata, target_sum1e4) # 对数转换 sc.pp.log1p(adata) # 保存标准化数据 adata.raw adata3.2 高变基因筛选不是所有基因都携带有用的生物学信息。筛选高变基因可以降低噪音提高分析效率# 识别高变基因 sc.pp.highly_variable_genes( adata, min_mean0.0125, max_mean3, min_disp0.5, batch_keysample # 考虑批次效应 ) # 可视化高变基因 sc.pl.highly_variable_genes(adata) # 可选仅保留高变基因 # adata adata[:, adata.var.highly_variable]4. 数据整合与批次校正4.1 使用harmonypy整合多样本数据当处理来自不同批次或实验条件的样本时批次效应会影响分析结果。harmonypy是处理这类问题的利器import scanpy.external as sce # 运行harmonypy整合 sce.pp.harmony_integrate(adata, sample) # 按样本批次整合 # 检查整合效果 sc.pp.neighbors(adata, n_pcs15, use_repX_pca_harmony) sc.tl.umap(adata) sc.pl.umap(adata, color[sample], wspace0.5)4.2 细胞周期效应校正细胞周期阶段会影响基因表达模式可能导致虚假聚类。我们可以通过评分和回归来校正这种影响# 定义细胞周期基因集 s_genes [...] # S期基因 g2m_genes [...] # G2/M期基因 # 计算细胞周期评分 sc.tl.score_genes_cell_cycle(adata, s_geness_genes, g2m_genesg2m_genes) # 回归去除细胞周期效应 sc.pp.regress_out(adata, [S_score, G2M_score]) # 可视化校正效果 sc.pl.umap(adata, colorphase)5. 降维、聚类与细胞注释5.1 降维与可视化主成分分析PCA是单细胞分析的标准降维方法# 缩放数据 sc.pp.scale(adata, max_value10) # PCA降维 sc.tl.pca(adata) # 选择主成分数量 sc.pl.pca_variance_ratio(adata, n_pcs50, logTrue) # UMAP可视化 sc.pp.neighbors(adata, n_neighbors15, n_pcs20) sc.tl.umap(adata) sc.pl.umap(adata, color[sample, n_genes])5.2 细胞聚类Leiden算法是目前单细胞分析中最常用的聚类方法# 尝试不同分辨率 for res in [0.2, 0.5, 0.8, 1.0]: sc.tl.leiden(adata, resolutionres, key_addedfleiden_{res}) # 可视化不同分辨率下的聚类结果 sc.pl.umap(adata, color[leiden_0.5, leiden_1.0], legend_locon data)5.3 双细胞检测与去除双细胞两个细胞被捕获在一个液滴中会影响分析结果。scrublet是检测双细胞的常用工具# 运行scrublet检测双细胞 sc.pp.scrublet(adata, batch_keysample) # 可视化双细胞评分 sc.pl.umap(adata, color[doublet_score, predicted_doublet]) # 去除预测的双细胞 adata adata[adata.obs.predicted_doublet False, :]5.4 细胞类型注释基于已知标记基因进行细胞类型注释是单细胞分析的关键步骤# 定义标记基因字典 marker_genes { T细胞: [CD3D, CD3E], B细胞: [CD79A, MS4A1], 髓系细胞: [CD14, FCGR3A], # 其他细胞类型标记... } # 可视化标记基因表达 sc.pl.dotplot(adata, marker_genes, groupbyleiden_1.0) # 差异表达分析识别cluster特征基因 sc.tl.rank_genes_groups(adata, leiden_1.0, methodwilcoxon) sc.pl.rank_genes_groups_dotplot(adata, n_genes5) # 手动注释细胞类型 celltype_mapping { 0: T细胞, 1: B细胞, # 其他cluster映射... } adata.obs[celltype] adata.obs[leiden_1.0].map(celltype_mapping)6. 高级分析与结果导出6.1 亚群分析有时我们需要对特定细胞群体进行更精细的分析# 提取T细胞亚群 t_cells adata[adata.obs.celltype T细胞].copy() # 重新聚类分析 sc.pp.pca(t_cells) sc.pp.neighbors(t_cells) sc.tl.leiden(t_cells, resolution0.3) sc.tl.umap(t_cells) # T细胞亚群标记基因 tcell_markers { CD4 T细胞: [CD4, IL7R], CD8 T细胞: [CD8A, CD8B], # 其他T细胞亚型... } sc.pl.dotplot(t_cells, tcell_markers, groupbyleiden)6.2 结果保存与可视化完成分析后保存结果非常重要# 保存完整分析结果 adata.write(final_analysis.h5ad) # 导出关键结果表格 import pandas as pd pd.DataFrame(adata.obs).to_csv(cell_metadata.csv) pd.DataFrame(adata.var).to_csv(gene_metadata.csv) # 生成发表级图表 with plt.rc_context({figure.figsize: (8, 6)}): sc.pl.umap(adata, colorcelltype, legend_locon data, frameonFalse, save_celltypes.pdf) sc.pl.dotplot(adata, marker_genes, groupbycelltype, save_marker_expression.pdf)在实际项目中我通常会建立一个完整的分析脚本包含从原始数据处理到最终结果可视化的所有步骤。这样不仅方便复现也便于与其他研究人员分享。Scanpy的另一个优势是它与Jupyter Notebook的完美集成可以实时查看分析结果和图表大大提高了分析效率。