单细胞 RNA-seq 质量控制：应用场景

从真实科研场景中学习 QC 策略

本章通过四个不同领域的完整案例，展示单细胞 QC 在实际研究中的应用。每个案例都来自已发表的研究，我们简化了步骤但保留了核心思路。

案例 1：PBMC 免疫细胞分析

研究背景

某团队研究 COVID-19 患者的免疫反应，分析了 20 个样本（健康对照组 vs 重症患者），总计 15 万个 PBMC 细胞。

数据特点

10X Genomics v3 化学
人类数据
批次效应明显（不同日期制备）

QC 策略

import scanpy as sc
import numpy as np

# 加载数据
adata = sc.read_h5ad('covid_pbmc_raw.h5ad')

# 1. 识别线粒体基因
adata.var['mt'] = adata.var_names.str.startswith('MT-')

# 2. 计算 QC 指标
sc.pp.calculate_qc_metrics(adata, qc_vars=['mt'], inplace=True)

# 3. 检查数据分布
print("=== 数据概览 ===")
print(adata.obs[['n_genes_by_counts', 'n_counts', 'pct_counts_mt']].describe())

# 4. 可视化
import matplotlib.pyplot as plt
fig, axes = plt.subplots(1, 3, figsize=(15, 5))
sc.pl.violin(adata, ['n_genes_by_counts'], groupby=None, ax=axes[0], show=False)
sc.pl.violin(adata, ['n_counts'], groupby=None, ax=axes[1], show=False)
sc.pl.violin(adata, ['pct_counts_mt'], groupby=None, ax=axes[2], show=False)
plt.tight_layout()
plt.show()

# 5. MAD 过滤（使用 3 MAD）
def mad_filter(data, n_mads, direction='both'):
    median = np.median(data)
    mad = np.median(np.abs(data - median))
    if direction == 'lower':
        return data > (median - n_mads * mad)
    elif direction == 'upper':
        return data < (median + n_mads * mad)
    else:
        lower = median - n_mads * mad
        upper = median + n_mads * mad
        return (data > lower) & (data < upper)

filter_mask = (
    mad_filter(adata.obs['n_genes_by_counts'], 3, 'lower') &
    mad_filter(adata.obs['n_genes_by_counts'], 3, 'upper') &
    mad_filter(adata.obs['n_counts'], 3, 'lower') &
    mad_filter(adata.obs['n_counts'], 3, 'upper') &
    mad_filter(adata.obs['pct_counts_mt'], 3, 'upper')
)

# 6. 检查不同样本的保留率
adata.obs['filter_pass'] = filter_mask
for sample in adata.obs['sample_id'].unique():
    sample_cells = adata.obs['sample_id'] == sample
    retention = (filter_mask[sample_cells].sum() / sample_cells.sum() * 100)
    print(f"{sample}: {retention:.1f}% 保留")

# 7. 应用过滤
adata_qc = adata[filter_mask].copy()
sc.pp.filter_genes(adata_qc, min_cells=10)

print(f"\n最终结果：")
print(f"细胞：{adata.n_obs} → {adata_qc.n_obs} ({adata_qc.n_obs/adata.n_obs*100:.1f}%)")
print(f"基因：{adata.n_vars} → {adata_qc.n_vars}")

# 8. 保存
adata_qc.write('pbmc_covid_qc.h5ad')

关键决策

MAD 倍数：使用 3（保守过滤）
样本特异性检查：确保没有单个样本被过度过滤
保留率目标：> 80%

结果

平均保留率：87.3%
检测到预期的免疫细胞群（T 细胞、B 细胞、单核细胞等）
重症患者的炎症单核细胞亚群得以保留

案例 2：小鼠心脏发育图谱

研究背景

构建小鼠心脏发育的单细胞 atlas，涵盖胚胎期 E10.5 到成年期。

数据特点

小鼠数据
多个发育阶段
心肌细胞线粒体含量天生较高

QC 策略（调整线粒体阈值）

import scanpy as sc
import numpy as np

# 加载数据
adata = sc.read_h5ad('mouse_heart_raw.h5ad')

# 1. 识别线粒体基因（小鼠用小写 mt-）
adata.var['mt'] = adata.var_names.str.startswith('mt-')

# 2. 计算 QC 指标
sc.pp.calculate_qc_metrics(adata, qc_vars=['mt'], inplace=True)

# 3. 检查线粒体分布
print("=== 线粒体比例分布 ===")
print(adata.obs['pct_counts_mt'].describe())
# 预期：中位数在 15-25%

# 4. 按发育阶段检查
for stage in ['E10.5', 'E14.5', 'P0', 'Adult']:
    stage_data = adata[adata.obs['stage'] == stage]
    print(f"\n{stage}:")
    print(f"  中位数：{stage_data.obs['pct_counts_mt'].median():.1f}%")
    print(f"  75 分位：{stage_data.obs['pct_counts_mt'].quantile(0.75):.1f}%")

# 5. 使用更高的 MAD 倍数（5 MAD）来适应高线粒体含量
def mad_filter(data, n_mads, direction='both'):
    median = np.median(data)
    mad = np.median(np.abs(data - median))
    if direction == 'upper':
        return data < (median + n_mads * mad)
    # ...

# 关键：线粒体使用 5 MAD，其他指标仍用 3 MAD
filter_mask = (
    mad_filter(adata.obs['n_genes_by_counts'], 3, 'lower') &
    mad_filter(adata.obs['n_genes_by_counts'], 3, 'upper') &
    mad_filter(adata.obs['n_counts'], 3, 'lower') &
    mad_filter(adata.obs['n_counts'], 3, 'upper') &
    mad_filter(adata.obs['pct_counts_mt'], 5, 'upper')  # 5 MAD!
)

# 6. 验证心肌细胞保留情况
# 心肌细胞标志物：TnnT3, Tnni3, Myh6/7
cardiomyocyte_markers = ['TnnT3', 'Tnni3', 'Myh6']
for marker in cardiomyocyte_markers:
    if marker in adata.var_names:
        cm_cells = adata[:, marker].X > 0
        cm_retention = (filter_mask[cm_cells].sum() / cm_cells.sum() * 100)
        print(f"{marker}+ 细胞保留率：{cm_retention:.1f}%")

# 7. 应用过滤
adata_qc = adata[filter_mask].copy()
sc.pp.filter_genes(adata_qc, min_cells=10)

print(f"\n最终结果：")
print(f"细胞：{adata.n_obs} → {adata_qc.n_obs} ({adata_qc.n_obs/adata.n_obs*100:.1f}%)")
print(f"基因：{adata.n_vars} → {adata_qc.n_vars}")

# 8. 保存
adata_qc.write('mouse_heart_qc.h5ad')

关键决策

线粒体阈值：使用 5 MAD（而非 3）
原因：心肌细胞正常线粒体含量 15-25%
验证：检查心肌细胞标志物表达，确保未过度过滤

结果

保留率：91.5%
心肌细胞亚群（心房、心室）完整保留
发育轨迹清晰可见

案例 3：肿瘤微环境分析

研究背景

分析结直肠癌肿瘤微环境，寻找与免疫治疗响应相关的细胞亚群。

数据特点

肿瘤组织 + 配对正常组织
高度异质性
可能含有循环肿瘤细胞（CTC）

QC 策略（极宽松过滤）

import scanpy as sc
import numpy as np

# 加载数据
adata = sc.read_h5ad('crc_tumor_raw.h5ad')

# 1. 识别线粒体和血红蛋白基因
adata.var['mt'] = adata.var_names.str.startswith('MT-')
adata.var['hb'] = adata.var_names.str.startswith(('HBA', 'HBB', 'HBD'))

# 2. 计算 QC 指标
sc.pp.calculate_qc_metrics(adata, qc_vars=['mt', 'hb'], inplace=True)

# 3. 极宽松的过滤策略
# 肿瘤数据：宁可保留噪音，不丢失稀有细胞
N_MADS = 5  # 使用 5 MAD

def mad_filter(data, n_mads, direction='both'):
    median = np.median(data)
    mad = np.median(np.abs(data - median))
    if direction == 'lower':
        return data > (median - n_mads * mad)
    elif direction == 'upper':
        return data < (median + n_mads * mad)
    else:
        lower = median - n_mads * mad
        upper = median + n_mads * mad
        return (data > lower) & (data < upper)

# 仅过滤极端异常值
filter_mask = (
    mad_filter(adata.obs['n_genes_by_counts'], N_MADS, 'lower') &
    mad_filter(adata.obs['n_genes_by_counts'], N_MADS, 'upper') &
    mad_filter(adata.obs['n_counts'], N_MADS, 'lower') &
    mad_filter(adata.obs['n_counts'], N_MADS, 'upper') &
    mad_filter(adata.obs['pct_counts_mt'], N_MADS, 'upper')
)

# 4. 不过滤低转录量细胞（可能是静息或 CTC）
print("不过滤低 UMI 细胞，以保留稀有细胞群")

# 5. 检查肿瘤 vs 正常组织的保留率
for tissue in ['Tumor', 'Normal']:
    tissue_cells = adata.obs['tissue'] == tissue
    retention = (filter_mask[tissue_cells].sum() / tissue_cells.sum() * 100)
    print(f"{tissue}: {retention:.1f}% 保留")

# 6. 应用过滤
adata_qc = adata[filter_mask].copy()

# 7. 基因过滤也使用宽松阈值
sc.pp.filter_genes(adata_qc, min_cells=5)  # 而非 10

print(f"\n最终结果：")
print(f"细胞：{adata.n_obs} → {adata_qc.n_obs} ({adata_qc.n_obs/adata.n_obs*100:.1f}%)")
print(f"基因：{adata.n_vars} → {adata_qc.n_vars}")

# 8. 后续：使用专门工具检测恶性细胞
# 例如 inferCNV、CopyKAT（在下游分析中进行）

# 9. 保存
adata_qc.write('crc_tumor_qc.h5ad')

关键决策

MAD 倍数：5（极宽松）
不过滤：低 UMI 细胞（可能是 CTC）
基因过滤：min_cells = 5（而非 10）
原因：肿瘤高度异质，稀有亚群可能是关键发现

结果

保留率：94.2%
后续分析发现了罕见的免疫抑制性髓细胞亚群
与免疫治疗响应显著相关

案例 4：斑马鱼胚胎发育

研究背景

研究斑马鱼早期胚胎发育（受精后 2-24 小时）。

数据特点

非模式生物（斑马鱼）
胚胎细胞转录活性低
基因注释可能不完整

QC 策略（非模式生物适配）

import scanpy as sc
import numpy as np

# 加载数据
adata = sc.read_h5ad('zebrafish_embryo_raw.h5ad')

# 1. 识别线粒体基因（斑马鱼可能有 mt- 或 MT-）
adata.var['mt'] = (
    adata.var_names.str.startswith('mt-') |
    adata.var_names.str.startswith('MT-')
)

# 验证
print(f"检测到 {adata.var['mt'].sum()} 个线粒体基因")
if adata.var['mt'].sum() == 0:
    print("警告：未检测到线粒体基因，检查基因命名规则")

# 2. 计算 QC 指标
sc.pp.calculate_qc_metrics(adata, qc_vars=['mt'], inplace=True)

# 3. 检查数据质量
print("=== QC 指标统计 ===")
print(adata.obs[['n_genes_by_counts', 'n_counts', 'pct_counts_mt']].describe())

# 胚胎细胞特点：
# - n_genes 较低（分化程度低）
# - n_counts 较低（转录活性低）
# - 线粒体比例可能较高（快速分裂）

# 4. 可视化分布
import matplotlib.pyplot as plt
fig, axes = plt.subplots(1, 3, figsize=(15, 5))
sc.pl.violin(adata, ['n_genes_by_counts'], groupby=None, ax=axes[0], show=False)
sc.pl.violin(adata, ['n_counts'], groupby=None, ax=axes[1], show=False)
sc.pl.violin(adata, ['pct_counts_mt'], groupby=None, ax=axes[2], show=False)
axes[0].set_title('Genes (embryo: lower expected)')
axes[1].set_title('UMI counts (embryo: lower expected)')
axes[2].set_title('MT% (embryo: may be higher)')
plt.tight_layout()
plt.show()

# 5. 极宽松过滤（胚胎细胞）
N_MADS = 5  # 使用 5 MAD

def mad_filter(data, n_mads, direction='both'):
    median = np.median(data)
    mad = np.median(np.abs(data - median))
    if direction == 'lower':
        return data > (median - n_mads * mad)
    elif direction == 'upper':
        return data < (median + n_mads * mad)
    else:
        lower = median - n_mads * mad
        upper = median + n_mads * mad
        return (data > lower) & (data < upper)

filter_mask = (
    mad_filter(adata.obs['n_genes_by_counts'], N_MADS, 'lower') &
    mad_filter(adata.obs['n_genes_by_counts'], N_MADS, 'upper') &
    mad_filter(adata.obs['n_counts'], N_MADS, 'lower') &
    mad_filter(adata.obs['n_counts'], N_MADS, 'upper') &
    mad_filter(adata.obs['pct_counts_mt'], N_MADS, 'upper')
)

# 6. 按发育时间点检查
for hpf in ['2hpf', '6hpf', '12hpf', '24hpf']:
    timepoint_cells = adata.obs['timepoint'] == hpf
    retention = (filter_mask[timepoint_cells].sum() / timepoint_cells.sum() * 100)
    print(f"{hpf}: {retention:.1f}% 保留")

# 7. 应用过滤
adata_qc = adata[filter_mask].copy()

# 8. 宽松的基因过滤
sc.pp.filter_genes(adata_qc, min_cells=5)  # 胚胎可能用更少

print(f"\n最终结果：")
print(f"细胞：{adata.n_obs} → {adata_qc.n_obs} ({adata_qc.n_obs/adata.n_obs*100:.1f}%)")
print(f"基因：{adata.n_vars} → {adata_qc.n_vars}")

# 9. 保存
adata_qc.write('zebrafish_embryo_qc.h5ad')

关键决策

线粒体基因识别：同时检查 mt- 和 MT- 前缀
预期调整：胚胎细胞的 n_genes 和 n_counts 天然较低
MAD 倍数：5（极宽松）
基因过滤：min_cells = 5（胚胎表达谱更局限）

结果

保留率：92.8%
保留了稀有祖细胞群
发育轨迹清晰

案例总结：QC 决策树

开始
  ↓
识别数据类型
  ↓
┌─────────────────────────────────────┐
│ 免疫细胞（PBMC）                    │ → 标准 MAD (3)
│  - 人类: MT-                        │
│  - 小鼠: mt-                        │
│  - 预期 mt%: 3-8%                   │
└─────────────────────────────────────┘
  ↓
┌─────────────────────────────────────┐
│ 高能量组织（心脏、肌肉）            │ → 提高 mt% 阈值 (5 MAD)
│  - 预期 mt%: 15-25%                 │
│  - 验证：检查组织特异标志物         │
└─────────────────────────────────────┘
  ↓
┌─────────────────────────────────────┐
│ 肿瘤样本                            │ → 极宽松 (5 MAD)
│  - 高度异质性                       │
│  - 稀有亚群可能关键                 │
│  - 不过滤低 UMI 细胞                │
└─────────────────────────────────────┘
  ↓
┌─────────────────────────────────────┐
│ 胚胎/干细胞                         │ → 极宽松 (5 MAD)
│  - 转录活性低                       │
│  - 动态表达                         │
│  - n_genes 天然较低                 │
└─────────────────────────────────────┘
  ↓
┌─────────────────────────────────────┐
│ 非模式生物                          │ → 检查基因命名
│  - 验证线粒体基因前缀               │
│  - 可能需要自定义基因列表           │
└─────────────────────────────────────┘
  ↓
验证过滤结果
  - 检查保留率（> 70%）
  - 验证已知细胞群保留
  - 检查不同批次/样本的保留率
  ↓
保存过滤后的数据

最佳实践清单

QC 前检查

确认物种和基因命名规则
检查数据类型（单细胞 vs 空间）
了解组织特异性质（线粒体含量、转录活性）

QC 过程

可视化 QC 指标分布
根据组织类型调整 MAD 倍数
检查不同批次/样本的保留率一致性
验证关键细胞群未被过度过滤

QC 后验证

重新可视化 QC 指标（确认改善）
检查保留率是否合理（> 70%）
快速聚类检查（确保无明显技术伪影）
记录使用的所有参数和阈值

通过本章的案例学习，你已经掌握了不同场景下的 QC 策略。在第五章中，我们将深入探讨 QC 的理论依据和工具原理，帮助你从根本上理解这些方法。

→ 继续阅读：第五章 - 最佳实践与工具原理