# 第 09 章 降维聚类与细胞类型注释

> 本章配套脚本：`code/04_annotation.R`（含 cell2location 调用说明）；本章图：`figures/fig09_celltype.png`

第 08 章结束，我们有了一个干净、归一化、带聚类的 Seurat 对象。但聚类编号（0、1、2……）本身没有生物学含义——它是“表达相似的 spot 群落”。本章要回答的问题只有一个：**每个 cluster 里到底是什么细胞？** 这一步叫细胞类型注释（cell type annotation）。它是承上启下的一环：没有注释，第 10 章的生态位无法命名，第 11 章的通讯分析没有“说话的主体”。本章给你三套策略，从最轻量到最严谨，并教你遇到问题怎么排查。

## 本章目标

1. 理解为什么 Visium 的 spot 必须注释：每个 spot 是 1–10 个细胞的混合信号，不是单细胞。
2. 掌握策略 A：用 marker 基因（标志基因）直接注释，会看 marker 的空间分布。
3. 理解策略 B：单细胞参考去卷积（deconvolution）的原理，会用 cell2location 估计每个 spot 的细胞类型丰度。
4. 了解策略 C：与公共单细胞图谱整合注释。
5. 能正确解读注释结果图（fig09），并处理“marker 不特异”“与病理不符”等常见问题。

## 9.1 为什么 spot 需要注释：混合信号的真相

先想清楚一个底层事实：**Visium 的 spot 不是单细胞**。每个 spot 直径约 55 μm，覆盖约 1–10 个细胞（第 02 章）。因此一个 spot 的计数矩阵是多种细胞表达的叠加——比如一个 spot 里同时有 3 个肿瘤细胞和 2 个巨噬细胞，它的表达谱就是“肿瘤 + 巨噬”的混合体。

这带来两个推论：

1. **聚类 ≠ 细胞类型**。单细胞数据里一个 cluster 通常≈一种细胞；Visium 里一个 cluster 可能是“某种细胞占主导”的区域，也可能是几种细胞的稳定组合。
2. **注释的对象是“spot 的主要成分”或“细胞类型比例”**，而不是“这个 spot 就是某某细胞”。所以注释策略天然分两派：直接派（这个 spot 更像谁）和比例派（这个 spot 里各占多少），本章的 A/B/C 三套策略正是这两派的体现。

理解了这一点，你就不会在看到“巨噬细胞 spot”时误以为那个 spot 里只有巨噬细胞。

## 9.2 策略 A：marker 基因直接注释（最轻量，先做它）

marker 基因（标志基因）是某类细胞特征性高表达、其他细胞低表达的基因。乳腺癌 Visium 的常见 cell type 与 marker 如下（与 00-TECH 速查第 7 节完全一致，全书统一）：

| 细胞类型 | marker 基因 | 备注 |
| --- | --- | --- |
| 肿瘤细胞（上皮） | EPCAM / KRT8 / KRT19 | 上皮细胞标志，肿瘤组织里通常是癌细胞 |
| 癌症相关成纤维细胞 CAF | COL1A1 / COL3A1 / DCN / PDGFRA | 成纤维细胞活化后高表达胶原与 PDGFRA |
| CD8 T 细胞 | CD3D / CD8A / GZMB | CD3D 是所有 T 细胞的通用标志，CD8A/GZMB 定位到杀伤性 T |
| 巨噬细胞 | CD68 / CD163 / C1QA | 肿瘤相关巨噬细胞（TAM）的标志 |
| B 细胞 | MS4A1 / CD79A | 淋巴细胞标志 |
| 内皮细胞 | PECAM1 / VWF | 血管内皮标志 |

直接注释的操作流程：**把每个 marker 基因画到组织切片上（SpatialFeaturePlot）→ 肉眼观察哪些区域高表达 → 给 cluster 或 spot 打分 → 命名**。先看空间分布：

```r
# 把 6 类细胞的代表 marker 画到组织切片上，一次看 6 张
SpatialFeaturePlot(obj,
  features = c("EPCAM", "COL1A1", "CD3D", "CD68", "MS4A1", "PECAM1"),
  ncol = 3, pt.size.factor = 1.2)
```

**运行结果解读**：你会看到 6 张切片图。理想情况下：EPCAM 高表达区域成片出现在肿瘤主体；COL1A1 环绕在肿瘤周围（基质）；CD3D/CD68 呈散点状散布（免疫浸润）；PECAM1 呈细条状网络（血管）。**如果某个 marker 全片低表达，先别慌**，可能是数据问题也可能是真实的生物学（见 9.6 节）。

肉眼观察之后，把观察变成可重复的“打分”：用 `AddModuleScore` 把同一细胞类型的多个 marker 合成一个模块分数，再按分数给每个 spot 判断主要细胞类型。

```r
# 定义 6 类细胞的 marker 列表（与上方表格一致）
marker_list <- list(
  Tumor   = c("EPCAM", "KRT8", "KRT19"),
  CAF     = c("COL1A1", "COL3A1", "DCN", "PDGFRA"),
  CD8T    = c("CD3D", "CD8A", "GZMB"),
  Macro   = c("CD68", "CD163", "C1QA"),
  Bcell   = c("MS4A1", "CD79A"),
  Endo    = c("PECAM1", "VWF")
)

# 每个 spot 计算 6 个模块分数（列名自动为 score1~score6，与 marker_list 顺序对应）
obj <- AddModuleScore(obj, features = marker_list, name = "score")

# 每个 spot 取分数最高的细胞类型作为"主要类型"
score_cols <- paste0("score", 1:length(marker_list))
obj$celltype <- names(marker_list)[apply(obj[[]][score_cols], 1, which.max)]

# 看看各类型 spot 数量
table(obj$celltype)

# 主要细胞类型画到切片上
SpatialDimPlot(obj, group.by = "celltype", label = TRUE)
```

**运行结果解读**：`table()` 给出每种主要类型的 spot 数——如果“CD8T”占了 60%，八成有问题（T 细胞不可能占乳腺癌切片的一半），需要检查打分逻辑或回到 marker 分布图。`SpatialDimPlot` 把“每个 spot 的主要细胞类型”映射回切片，你会看到肿瘤区域一片 Tumor 色、基质区域一片 CAF 色——这就是策略 A 的最终产物。

**策略 A 的适用场景**：快速出图、初步判断、以及作为策略 B/C 的交叉验证。它的局限也很明显：只回答“谁占主导”，答不了“比例是多少”，而且对混合严重的 spot 容易误判。

## 9.3 策略 B：单细胞参考去卷积（更定量）

去卷积（deconvolution）解决策略 A 回答不了的问题：**每个 spot 里各种细胞类型各占多少**。思路是：先有一个参考——高质量的同一组织类型单细胞图谱（已知每个细胞是什么类型）——然后反推：把每个 spot 的混合表达谱拆解成各细胞类型的比例。

**cell2location 原理**（本教程主推）：它分两步走，都是贝叶斯模型。

- **第一步**：用参考单细胞数据训练一个负二项回归模型，学到每种细胞类型的“平均表达谱”与基因离散度参数——相当于给每种细胞建一份“基因指纹”。
- **第二步**：把第一步学到的指纹固定住，用来解释 Visium 每个 spot 的总 UMI：模型假设 spot 的计数 = 若干种细胞类型的指纹按“每种类型有多少个细胞”加权求和，然后用变分推断（variational inference）反推出每个 spot 每种细胞类型的**绝对细胞数**（丰度）。

输出的核心矩阵是 `spot × 细胞类型` 的丰度表，存于 `adata.obsm["q05_cell_abundance_w_sf"]`。下面是核心代码（Python，cell2location ≥ 0.8；需要先准备好参考单细胞 h5ad 与 Visium h5ad，见第 07 章数据清单）：

```python
import cell2location
from cell2location.models import RegressionModel, Cell2location
import scanpy as sc

# 输入准备：参考单细胞图谱（含细胞类型标签列 "cell_type"）与 Visium 数据
adata_ref  = sc.read_h5ad("data/ref_breast_sc.h5ad")    # 参考单细胞
adata_vis  = sc.read_h5ad("data/visium_breast.h5ad")    # Visium（scanpy 版）

# 第一步：用参考单细胞训练"指纹"模型
ref_model = RegressionModel(adata_ref)
ref_model.train(max_epochs=250, use_gpu=False)
adata_ref = ref_model.export_posterior(adata_ref)

# 第二步：固定指纹，估计每个 spot 的细胞类型丰度
# N_cells_per_location：每个 spot 平均细胞数，Visium 取 10（覆盖 1-10 个细胞）
mod = Cell2location(adata_vis, adata_ref, N_cells_per_location=10)
mod.train(max_epochs=30000, use_gpu=False)

# 导出丰度矩阵并写回 adata
adata_vis = mod.export_posterior(adata_vis)
abundance = adata_vis.obsm["q05_cell_abundance_w_sf"]   # spot × 细胞类型
abundance.head()
```

**运行结果解读**：`abundance` 是数据框，行是 spot，列是细胞类型，数值是“该 spot 里大约有几个这种细胞”。把丰度最高的类型作为 spot 的主要细胞类型，可以画出和策略 A 类似的切片图；把丰度按生态位汇总，就是 9.5 节堆叠条形图的原料。注意两点：训练要跑很久（教程里 `max_epochs=30000` 是保守设置，实际常用早停，机器建议 ≥ 16 GB 内存，见第 06 章）；`use_gpu=False` 表示纯 CPU 也能跑，只是慢。

**一句话认识另外两个去卷积工具**：RCTD（R 包 spacexr）用受限二次规划把参考表达谱拟合到 spot 上，速度快、文档友好；SPOTlight 先用非负矩阵分解（NMF）从参考数据提炼“细胞类型特征基”，再用非负最小二乘估计比例。三者结论通常大体一致；cell2location 的卖点是能估计**绝对细胞数**（而非只给相对比例），更适合我们后面“按生态位比较组成”的需求。

## 9.4 策略 C：与公共单细胞图谱整合注释

如果你手头没有自己的单细胞参考，可以用公共图谱。乳腺癌方向有多个公开的单细胞图谱（如 Broad 的 single cell portal、CellxGene 上的乳腺癌图谱，见第 07 章数据来源）。整合注释的通用做法是 Seurat 的标签转移（label transfer）：把 Visium 和参考图谱放进同一个低维空间，让参考的细胞类型标签“传导”到空间 spot 上。

```r
# 伪代码：Seurat 标签转移流程（参考图谱 adata_ref 需先转成 Seurat 对象）
# anchors <- FindTransferAnchors(reference = ref, query = obj, dims = 1:30)
# pred <- TransferData(anchorset = anchors, refdata = ref$cell_type, dims = 1:30)
# obj$celltype_integrated <- pred$predicted.id
```

> 上面是伪代码：真实运行需要先下载并整理参考图谱（第 07 章的下载清单）。本书以策略 A/B 为主，策略 C 的价值在于**没有自备参考时的备选**，以及多策略交叉验证。

**三套策略怎么选？** 实践建议：先用策略 A 快速出图建立直觉 → 有条件就用策略 B 拿到定量比例（本课题主线）→ 用策略 C 或文献结果交叉验证，三者一致时结论才敢写进论文。

## 9.5 结果解读：注释结果图

图 9-1 是注释结果示例（模拟数据演示）：

![图 9-1 细胞类型注释结果（模拟数据演示）。左：每个 spot 的主要细胞类型映射回组织切片；右：按空间生态位汇总的平均细胞组成堆叠条形图。](figures/fig09_celltype.png){width=85%}

**怎么读这张图**：左图看“空间”——理想结果里，Tumor 色块成片位于切片中央，CAF 色块包裹在外周，免疫细胞（CD8T、Macro、Bcell）散布其间，Endo 呈细线穿插——这与乳腺癌“肿瘤巢 + 间质 + 免疫浸润”的组织学结构相符。右图看“组成”——每个生态位（柱）里各细胞类型的占比（堆叠色块），一眼就能比较“免疫浸润区”和“肿瘤核心”的细胞构成差异。**注意：本图为模拟数据演示，仅用于展示“真实分析会长什么样”。** 右图的生态位划分来自第 10 章，这里先眼熟它——下一章我们就要亲手做出来。

**统计提示**：如果后续要对生态位间的细胞组成做差异检验（比如“侵袭前沿的巨噬细胞比例显著更高吗”），建议用置换检验（permutation test）或广义线性模型，并对多重比较做 BH-FDR 校正（见第 00-SPEC 第 8 节），不要直接用 `t.test` 反复比。

## 9.6 常见问题排查

**问题 1：marker 不特异怎么办？**

这是新手注释翻车的第一大原因。比如 `CD68` 在巨噬细胞和树突状细胞里都表达；`EPCAM` 在发生上皮间质转化（EMT）的肿瘤细胞里会下调，导致“肿瘤区域测不到 EPCAM”。对策：

- 用**多个 marker 联合打分**（策略 A 的 `AddModuleScore` 已经这么做了），单个 marker 高不算数；
- 换更特异的 marker（如巨噬细胞加看 `CD163`、`MRC1`；T 细胞加看 `CD4`/`CD8A` 细分）；
- 用去卷积结果交叉验证：如果 cell2location 说某个区域巨噬细胞丰度很高，而 marker 打分却低，多半是 marker 选择问题而不是数据问题；
- 查阅文献：以“乳腺癌 + 空间转录组 + marker”为关键词检索（见第 07 章文献检索记录），用同行验证过的 marker 列表。

**问题 2：注释结果与病理（H&E）不符怎么办？**

比如病理切片上明明是大片肿瘤，注释却标成基质。可能的原因按可能性排序：

1. **EMT 或异质性**：部分肿瘤细胞下调上皮 marker，看起来像基质（生物学真相，不是 bug）；
2. **混合信号稀释**：肿瘤细胞与基质细胞共处一个 spot，marker 信号被稀释到阈值以下（Visium 分辨率的固有限制）；
3. **参考不匹配**：策略 B/C 的参考图谱来自不同亚型或不同平台，批次效应污染了注释。

排查顺序建议：先回到原始表达（`SpatialFeaturePlot` 直接看 marker），再换一套 marker 或换策略，最后才考虑“是不是数据读错了”。**切记：注释是计算推断，最终解释权在病理与实验证据**，教程的注释只能作为假说来源。

## 9.7 动手练习

1. 用策略 A 的 `AddModuleScore` 流程注释你的对象，把 `table(obj$celltype)` 的结果与 9.2 节的示例对照，解释各类 spot 数量的合理性。
2. 把 cell2location 的 `abundance` 结果按“丰度最高的类型”画成切片图，与策略 A 的结果叠加比较：哪些 spot 两种策略一致？哪些不一致？试着解释不一致的原因。
3. （进阶）下载一个公共乳腺癌单细胞图谱，用 9.4 节的标签转移流程注释，与策略 A/B 三方交叉验证。

## 本章小结

本章回答了“spot 里是什么细胞”：先讲清 spot 是 1–10 个细胞的混合信号这个底层事实，然后给出三套注释策略——marker 直接注释（快、直观）、cell2location 去卷积（定量、给出比例）、公共图谱整合（无参考时的备选与验证）。同时教了两种常见问题的排查思路。注释完成后，每个 spot 都有了“主要细胞类型”或“丰度向量”，这正是下一章识别空间生态位的输入——因为生态位的本质，就是“细胞组成有特色、空间位置连成片”的区域。

## 练习与思考

1. **概念题**：为什么说“Visium 聚类得到的 cluster ≠ 细胞类型”？混合信号对注释策略的选择有什么影响？
2. **概念题**：cell2location 与策略 A 的根本区别是什么？（提示：一个回答“谁占主导”，一个回答“各占多少”）
3. **动手题**：对 `EPCAM`、`KRT8`、`KRT19` 三个 marker 分别画 `SpatialFeaturePlot`，观察它们的高表达区域是否完全重合。如果不重合，说明了什么？
4. **动手题**：把 cell2location 的丰度结果与策略 A 的主要类型结果做一个一致性表格（如列联表），统计一致率，并讨论不一致 spot 的特征。
5. **思考题**：如果你的数据里 B 细胞 marker（MS4A1/CD79A）全片几乎不表达，你会先怀疑数据质量，还是先怀疑生物学？请列出你验证判断的具体步骤。
