single-cell-rnaseq

Overview

本技能把单细胞 RNA-seq 数据从原始计数矩阵处理到可解释的细胞分群与细胞类型注释：自动识别 10x HDF5 / MTX / 稠密 CSV / h5ad 四类输入格式，做基于线粒体比例、基因数与 UMI 数的质控过滤，随后进行 CP10K + log1p 归一化、高变基因筛选、z-score 缩放、截断 SVD 主成分分析、kNN 邻接图构建与图聚类，最后输出 UMAP / t-SNE 二维嵌入、每簇 marker 基因表与可视化图表。

分析提供 双后端：--backend scanpy 调用官方 scanpy 实现，--backend numpy 使用内置的纯 NumPy/SciPy/scikit-learn 降级实现（自实现的 Louvain/UMAP 谱初始化与力优化），在没有 scanpy 的环境中仍能产出统计等价的结果。

When to Use This Skill

• 用户要做单细胞转录组分析（关键词：single-cell、scRNA-seq、单细胞、细胞图谱、细胞注释、UMAP、聚类）
• 用户有 10x Genomics 原始计数（matrix.h5 / filtered_feature_bc_matrix/）需要走标准 QC 与聚类流程
• 用户需要判断细胞质量（关键词：doublet、双细胞、线粒体比例、mitochondrial、QC、过滤阈值）
• 用户需要找每个细胞群的 marker 基因并据此命名细胞类型
• 用户环境里没有 scanpy，但仍需要可复现的降维与聚类结果
• 用户已有聚类结果，需要补充 t-SNE / marker 热图 / 分析报告

不适用：bulk RNA-seq（用 bulk-rnaseq）、空间转录组（CST/Slide-seq 需专门流程）、单细胞 ATAC、轨迹推断（trajectory inference）。

Quick Start

export PYTHONPATH="$PWD/packages"

# 完整流程（自动检测输入格式）
python skills/single-cell-rnaseq/scripts/run_pipeline.py \
    --input  /data/filtered_feature_bc_matrix \
    --output /results/scrna \
    --marker-file /ref/cell_marker_reference.csv

# 10x HDF5 + 强制使用降级后端（无 scanpy 环境）
python skills/single-cell-rnaseq/scripts/run_pipeline.py \
    --input /data/raw/filtered_feature_bc_matrix.h5 \
    --output /results/scrna \
    --backend numpy --n-neighbors 15 --resolution 1.0

# 只做质控，不做降维聚类
python skills/single-cell-rnaseq/scripts/quality_control.py \
    --input /data/filtered_feature_bc_matrix \
    --output /results/scrna_qc \
    --min-genes 200 --min-counts 500 --max-mt-pct 15

Workflow

1. 数据读取 — 自动嗅探输入类型：10x HDF5（read_10x_h5）、10x MTX 目录（read_10x_mtx）、稠密 CSV/TSV（read_dense_table）、h5ad（read_h5ad）；统一封装为轻量 SimpleAnnData 兼容容器，不强制依赖 anndata。
2. 质控指标 — 逐细胞计算 total_counts / n_genes_by_counts / pct_counts_mt（线粒体）/ pct_counts_hb（血红蛋白）及对数版本。
3. 质控过滤 — 按 min_genes、min_counts、max_mt_pct、可选 max_hb_pct 过滤细胞；按 min_cells_per_gene 过滤基因；输出过滤前后细胞数/基因数与去除率。
4. 归一化 — library size 归一化到 CP10K（--target-sum），随后 log1p 变换稳定方差。
5. 高变基因 — Seurat 风格分箱 z-score 法（20 个表达量分箱，对数期望值/标准差标准化，CV 校正），选出 --n-top-genes 个基因。
6. 缩放与 PCA — 逐基因 z-score 缩放并截断到 ±10；用截断 SVD（稀疏输入）/ 稠密 SVD 得到 --n-comps 个主成分，同时记录解释方差比例。
7. 邻接图 — 在 PCA 空间用 sklearn NearestNeighbors 找 k 近邻，用 Jaccard 相似度算边权，经 UMAP 式 prune 阈值剪枝后对称化并行归一化。
8. 聚类 — 优先 Leiden（leidenalg + igraph，RBConfiguration），缺失时回退到内置 Louvain（模块度增益的局部移动 + 多层图聚合）。--resolution 控制粒度。
9. 降维嵌入 — UMAP（umap-learn 可用时用官方实现；否则用内置谱初始化 + UMAP 真实力公式优化）与 t-SNE（sklearn TSNE）。
10. Marker 与注释 — Welch t 检验 + BH 校正对每簇排序 marker；用 marker 参考文件做基因重叠打分给出细胞类型标签，无参考文件时输出匿名标签。
11. 可视化与报告 — UMAP 分群图、UMAP 细胞类型图、簇大小条形图、marker 热图（PDF + PNG，300 DPI），并生成文字报告与 run_manifest.json。

Parameters

参数	默认	说明
--input	必填	输入矩阵：10x h5 / 10x mtx 目录 / csv / h5ad
--output	<results>/scrna	输出根目录，下设 7 个步骤子目录
--marker-file	无	marker 参考 CSV（含 class_label / subclass_label / marker 列）
--backend	auto	auto / scanpy / numpy（降级实现）
--min-genes	200	细胞至少检测到的基因数
--min-counts	100	细胞至少的 UMI 数
--max-mt-pct	20.0	线粒体基因比例上限（%）
--max-hb-pct	无	血红蛋白基因比例上限（%），默认不过滤
--min-cells-per-gene	3	基因最少在多少细胞中表达
--target-sum	1e4	归一化目标总计数（CP10K）
--n-top-genes	2000	参与 PCA 的高变基因数
--n-comps	50	PCA 主成分数
--n-neighbors	15	kNN 邻居数 k
--resolution	1.0	聚类分辨率（越大簇越细）
--min-dist	0.3	UMAP 最小嵌入距离
--n-markers	100	每簇保留的 marker 基因数
--seed	0	随机种子
--no-tsne	关	跳过 t-SNE（较慢）

Outputs

文件	说明
01_qc/obs_qc.csv	每细胞质控指标（counts / genes / mt% / hb%）
01_qc/qc_summary.json	过滤前后细胞数、基因数与所用阈值
02_normalized/var_hvg.csv	每基因的均值与高变标记
03_dimensionality/embeddings.csv	PCA / UMAP / t-SNE 坐标合并表
03_dimensionality/pca_variance.csv	各主成分解释方差比例
04_clustered/	聚类标签（写入 obs）
05_annotated/cluster_annotations.csv	每个 cluster 的细胞类型标签与 marker 重叠数
05_annotated/marker_genes.csv	每簇 marker 基因表（logFC / t 值 / p / padj）
05_annotated/cluster_summary.csv	每簇细胞数与 marker 数
06_figures/Figure{1..4}_*.{pdf,png}	UMAP 分群图、UMAP 细胞类型图、簇大小图、marker 热图
07_report/scrna_report.txt	完整分析报告
run_manifest.json	运行元数据（输入、后端、细胞/基因数、簇数）
qc_metrics.csv / qc_filtered_cells.csv / qc_summary.json	（quality_control.py 独立产出）

Validation

• 质控过滤率：细胞去除率通常在 5%–30%。超过 50% 说明阈值过严（常见于 --min-genes 偏高），低于 1% 说明阈值过松，应检查是否混入了大量低质量细胞。
• 线粒体比例：健康细胞通常 < 10%。若绝大多数细胞 > 30%，提示样本解离质量差或存在应激，此时上调 --max-mt-pct 只会引入噪声细胞。
• 簇大小：出现占比 < 1% 的小簇通常是双细胞或批次效应，可通过调高 --resolution 观察其是否稳定存在。
• UMAP 形态：簇应呈连续、互不重叠的团块。若所有细胞挤成一点，通常是 --n-top-genes 过小或 PCA 维度不足；若出现极端离群细胞，检查是否有未过滤的低质量细胞。
• marker 合理性：每簇的 top marker 应有明确生物学含义。若某簇的 marker 全是核糖体蛋白基因，多半是簇内缺少区分信号（分辨率过低）。
• 降维一致性：UMAP 与 t-SNE 上同一 cluster 的位置若完全不一致，说明邻接图不稳定，可增大 --n-neighbors 后重跑。
• 后端一致性：同一数据用 --backend scanpy 与 --backend numpy 各跑一次，聚类数与簇大小分布应大体一致（±1 个簇以内）。

References

• references/methodology.md — 各步骤方法学原理、公式推导与参数选择依据
• references/troubleshooting.md — 常见报错、结果异常与排错手册