# 第 10 章 空间生态位识别

> 本章配套脚本：`code/05_spatial_domain.R`；本章图：`figures/fig10_niche.png`、`figures/fig11_nhood.png`

前面两章我们完成了“每个 spot 是什么”的拼图：第 08 章清洗数据，第 09 章注释细胞类型。但从“spot 有细胞类型”到“组织里有几个生态位”，还差关键一步——**把散落的 spot 组织成有生物学意义的空间区域**。这就是本章的主题：空间生态位识别（spatial niche identification）。它是本课题的核心章节：没有生态位，就没有“分区通讯”，也就没有后面的“通讯重编程”故事。本章的三步法会带你从分割、验证到命名，完整走一遍。

## 本章目标

1. 理解生态位识别“三步法”的整体框架：空间域分割 → 邻域富集验证 → 命名与特征化。
2. 学会用 BayesSpace 做空间域分割，理解其马尔可夫随机场原理与聚类数 q 的选择。
3. 学会用 Squidpy 的邻域富集（nhood enrichment）验证细胞类型的空间共定位。
4. 掌握生态位命名规则与特征化方法（marker 基因 + Hallmark 通路富集）。
5. 能正确解读生态位结果图，并判断“几个生态位合适”。

## 10.1 三步法总览

生态位识别不是一个算法能完成的，而是一条验证链。本书采用三步法：

1. **空间域分割（spatial domain segmentation）**：用考虑空间位置信息的聚类算法，把相邻且表达相似的 spot 归为同一“空间域”。这一步产生“候选生态位”。
2. **邻域富集验证（neighborhood enrichment validation）**：用统计检验确认候选生态位的细胞组成在空间上确实“扎堆”——即某些细胞类型倾向于互为邻居，而不是随机散布。这一步验证分割结果不是统计假象。
3. **生态位命名与特征化（niche annotation & characterization）**：根据每个生态位的细胞组成与特征通路，结合组织学先验（比如肿瘤核心、侵袭前沿），给它们命名，并用差异表达与通路富集描述“这个生态位是什么、擅长干什么”。

三步缺一不可：第一步给地图，第二步验地图，第三步给地图写图例。下面依次实操。

## 10.2 第一步：空间域分割（BayesSpace 实操）

### BayesSpace 原理

第 08 章的 Seurat 聚类只用了表达信息（PCA 空间），完全无视 spot 的物理位置，所以它可能把组织两侧“表达相似但位置不相邻”的 spot 聚到一起——这对识别空间区域是缺陷。BayesSpace 的做法是引入**空间先验**：它假设空间上相邻的 spot 更可能属于同一个空间域。具体地，BayesSpace 在聚类模型里加入一个**马尔可夫随机场（Markov random field, MRF）**，用一个“平滑”参数约束相邻 spot 的聚类标签趋于一致，再通过贝叶斯推断同时估计每个 spot 的域归属和基因表达模式。**聚类数 q 是唯一需要你拍板的关键参数**——它直接决定组织被切成几块。别担心，10.2 节和 10.6 节会教你怎么选。

### 实操代码

BayesSpace 的输入是 SingleCellExperiment（SCE）对象，所以第一步从 Seurat 转换（本教程以 q = 7 为例，参数与 00-TECH 速查第 3 节一致）：

```r
library(BayesSpace)
library(SingleCellExperiment)

# 从 Seurat 对象转成 SCE（BayesSpace 的输入格式）
sce <- as.SingleCellExperiment(obj)

# 空间预处理：对数标准化 + 提取前 15 个主成分（供聚类使用）
sce <- spatialPreprocess(sce, platform = "Visium", n.PCs = 15)

# 选 q：用 qTune 对多个 q 值各跑一轮，比较结果（数据量大时较慢，可先用 4-8）
sce <- qTune(sce, qs = seq(4, 8), platform = "Visium", d = 15)
# qTune 会为每个 q 存一组聚类标签（colData 里的 q4~q8 列），
# 可用 clusterPlot(sce, label = "q4") 依次查看每个 q 的分割效果

# 正式空间聚类（此处 q 取自 qTune 的折中值，实践中按 10.6 节标准定）
sce <- spatialCluster(sce, q = 7, platform = "Visium", d = 15)

# 可视化：空间域映射回组织切片
clusterPlot(sce)

# 把空间聚类标签传回 Seurat 对象，便于后续统一操作
obj$BayesSpace_cluster <- sce$spatial.cluster
```

**运行结果解读**：`clusterPlot` 画出一张组织切片，不同颜色代表不同空间域。与第 08 章的 `SpatialDimPlot` 对比，你会看到两个明显区别：一是分块更“整”——颜色区域边界平滑、成片分布，因为 BayesSpace 强制相邻 spot 同域；二是分块往往对应组织学的解剖结构。`qTune` 的输出是每个 q 的一组分域结果，肉眼扫一遍就能发现：q 太小会把明显不同的区域硬并在一起，q 太大则出现碎斑。**注意：spatialCluster 基于 MCMC 采样，每次运行结果有细微随机性，正式分析建议设随机种子（`set.seed(123)`）保证可复现。**

### SpaGCN 与 STAGATE 一句话对比

BayesSpace 只用表达 + 坐标；另外两个流行工具思路不同，按场景选：

- **SpaGCN**：把空间坐标、表达和（可选的）组织学图像颜色融合成图卷积网络的输入，再聚类。**如果手头有高质量 H&E 图像，SpaGCN 能利用染色信息，通常加分**——但参数更多、更吃调参。
- **STAGATE**：基于图注意力自编码器（graph attention autoencoder），把 spot 邻居关系编码进嵌入，再聚类。对组织边界敏感，但同样依赖超参数。

**本书选 BayesSpace 的理由**：原理直观（新手容易讲清楚）、参数少（主要就 q 和 d）、结果稳定，是教学与入门的最优解。真实课题若追求极致分域，可以三家都跑再取一致区。

## 10.3 第二步：邻域富集验证（Squidpy 实操）

分割之后要回答一个尖锐的问题：**这些空间域是真的“细胞扎堆”的结果，还是聚类算法的巧合？** 邻域富集分析（neighborhood enrichment）就是干这个的。它的逻辑是：把每个 spot 的邻居（比如六边形网格上相邻的 6 个 spot）统计出来，计算“细胞类型 A 的邻居里出现类型 B”的频率，再与随机置换（permutation）下的期望频率比较，得到 z 分数。z 分数显著为正 = A 与 B 倾向共定位（富集）；显著为负 = 倾向互斥。

```python
import squidpy as sq

# 第一步：构建空间邻居图（Visium 是六边形网格，coord_type="generic"，取 6 邻域）
sq.gr.spatial_neighbors(adata, coord_type="generic", n_neighs=6)

# 第二步：邻域富集检验（基于置换检验，cluster_key 是细胞类型标签列）
sq.gr.nhood_enrichment(adata, cluster_key="celltype")

# 第三步：画热图——z 分数，红=共定位富集，蓝=互斥
sq.pl.nhood_enrichment(adata, cluster_key="celltype")
```

**运行结果解读**：热图的行列都是细胞类型，格子的 z 分数表示两类细胞的共定位程度（图 10-2 就是它的成品，见 10.5 节）。典型模式：肿瘤细胞 × 肿瘤细胞对角线上深红——肿瘤细胞确实彼此紧挨成片；肿瘤细胞 × CD8 T 细胞显著为正——免疫细胞浸润到肿瘤巢；肿瘤 × B 细胞接近 0 或为负——它们通常不直接相邻。**这套“谁和谁挨着”的证据，正是第 11 章通讯分析的前置证据：配体-受体互作要有意义，前提是两类细胞真的空间相邻。**

**co-occurrence 补充**：除了“是不是邻居”，有时你还想知道“随着距离增加，共现如何衰减”。co-occurrence 分析给出随距离变化的共现概率曲线，能区分“紧密接触”与“松散共处”：

```python
# 共现分析：随距离变化的共现概率（补充 nhood_enrichment 的视角）
sq.gr.co_occurrence(adata, cluster_key="celltype")
sq.pl.co_occurrence(adata, cluster_key="celltype")
```

**运行结果解读**：输出一组曲线，横轴是距离（spot 数），纵轴是共现概率。如果肿瘤与巨噬细胞的曲线在近距离处明显高于随机线，说明二者紧密相邻——这个信息在后面的通讯重编程故事里会反复用到。

## 10.4 第三步：生态位命名与特征化

### 命名规则

分割出的空间域（BayesSpace cluster）没有自带名字，需要我们把“细胞组成 + 位置 + 先验知识”综合起来命名。本书采用的教学命名规则如下（与 00-TECH 速查第 4 节一致；**注意：名称仅为教学约定，真实课题应基于数据 + 先验命名并给出依据**）：

| 生态位 | 典型组成 | 特征通路 / marker |
| --- | --- | --- |
| 肿瘤核心（tumor core） | 肿瘤细胞为主、血管内皮 | EPCAM/KRT、VEGF 高 |
| 侵袭前沿（invasive front） | 肿瘤 + CAF 交界、EMT 特征 | VIM/ZEB1/TGF-β 高 |
| 免疫浸润区（immune-infiltrated region） | T/B/巨噬细胞富集 | CD3D/CD8A/MS4A1 高 |
| 基质区（stromal region） | CAF、ECM 富集 | COL1A1/DCN/PDGFRA 高 |

**命名的实操流程**：先按第 09 章的细胞组成给每个空间域画堆叠条形图（每个域 = 一根柱子，看各细胞类型占比），再对照上表“对号入座”。比如一个域的巨噬细胞和 T 细胞占比最高，就命名为免疫浸润区；一个域肿瘤细胞占比 90% 且位于切片中央，就是肿瘤核心。命名要**可证伪**：给出你依据的 marker 证据，而不是凭感觉。

### 特征化：marker 基因与通路富集

命名之后，用两个分析把每个生态位的“个性”量化出来：一是差异表达基因（每个生态位相对其他生态位高表达的基因），二是通路富集（这些基因富集在哪些生物学通路）。

```r
# 把生态位作为分组变量，找每个生态位的特征基因（BH-FDR 校正，FindAllMarkers 内置）
Idents(obj) <- obj$niche   # niche 列 = 命名后的生态位标签
niche_markers <- FindAllMarkers(obj, only.pos = TRUE,
                                min.pct = 0.25, logfc.threshold = 0.5)
# 查看侵袭前沿的 top 特征基因
head(niche_markers[niche_markers$cluster == "侵袭前沿", ], 10)
```

```r
# Hallmark 通路富集（ORA，过表达分析）
# 先获取 Hallmark 基因集（msigdbr 包：BiocManager::install("msigdbr")）
library(msigdbr); library(clusterProfiler)
hallmark <- msigdbr(species = "Homo sapiens", category = "H") %>%
  dplyr::select(gs_name, gene_symbol)

# 对"侵袭前沿"的特征基因做富集（其他生态位同理，循环即可）
res <- enricher(gene = niche_markers$gene[niche_markers$cluster == "侵袭前沿"],
                TERM2GENE = hallmark, pvalueCutoff = 0.05)
head(res@result[, c("ID", "p.adjust", "GeneRatio")])
```

**运行结果解读**：`FindAllMarkers` 输出每个生态位的特征基因表；`enricher` 输出富集到的 Hallmark 通路与校正后 p 值。教学示例里你常会看到：侵袭前沿富集到 EMT（上皮间质转化）、TGF-β 信号；免疫浸润区富集到干扰素 γ 响应、炎症反应；肿瘤核心富集到血管生成（VEGF 通路）。**注意：通路富集结论要注明统计方法（这里是超几何检验 + BH-FDR 校正），并且“富集到 EMT 通路”≠“细胞正在发生 EMT”，前者是计算证据，后者需要实验验证**（第 00-SPEC 第 8 节）。

## 10.5 结果解读：生态位图与邻域富集热图

图 10-1 是生态位识别的最终结果（模拟数据演示）：

![图 10-1 空间生态位识别最终结果（模拟数据演示）。左：组织切片按四个生态位着色（肿瘤核心、侵袭前沿、免疫浸润区、基质区）；右：各生态位占比图例与 spot 数。](figures/fig10_niche.png){width=85%}

**怎么读这张图**：左图是成品的“生态位地图”——四种颜色对应四个生态位，边界清晰、成片分布。检查它是否合理，看三点：① 肿瘤核心是否位于切片中央的大块区域；② 侵袭前沿是否夹在肿瘤核心与基质区之间（“前沿”的字面意思）；③ 免疫浸润区是否零散分布但成簇。如果某生态位碎成几十个孤立小点，说明分割过细或数据质量有问题，回 10.2 节调 q。右图给出各生态位的 spot 数，方便评估样本量——spot 数太少的生态位（比如 < 50）后面做分区通讯时统计功效不足。

图 10-2 是邻域富集验证结果（模拟数据演示）：

![图 10-2 细胞类型邻域共现富集热图（模拟数据演示，Squidpy nhood_enrichment 风格）。行列为细胞类型，颜色为 z 分数：红色表示两类细胞显著共定位（互为邻居的频率高于随机），蓝色表示显著互斥。](figures/fig11_nhood.png){width=85%}

**怎么读这张图**：对角线格子（A×A）深红说明该类细胞“自己挨自己”，是成片分布的标志；非对角线格子红色说明两类细胞倾向相邻。把图 10-2 与图 10-1 对照：图 10-1 里侵袭前沿的组成是“肿瘤 + CAF 交界”，那么图 10-2 里肿瘤细胞 × CAF 就应该显著富集——两张图互相印证，生态位的划分才站得住脚。**注意：本图为模拟数据演示，仅用于展示“真实分析会长什么样”。**

## 10.6 生态位数量的生物学判断

q 选多少？统计上可以看模型拟合，但更重要的判断标准是**生物学可解释性**。参考原则：

- **4–7 个较常见**：太少（< 4）会把生物学上不同的区域并到一起，后面“重编程”比较无从谈起；太多（> 10）则出现无法解释的碎片。
- **宁少勿多**：每个生态位都要能回答“它是什么、为什么存在”。一个说不清的生态位是分析的负担而不是收获。
- **三个硬性检查**：① 每个生态位有独特、可解释的细胞组成（堆叠条形图差异明显）；② 空间上连成片（不是碎点）；③ 能在 H&E 切片上找到对应的组织学结构（肿瘤巢、纤维间质、淋巴浸润灶）。三条都过，这个 q 才敢用。

实践技巧：把 q 从 4 跑到 8，画出每个 q 的生态位地图和组成条形图，放在一起看“哪个 q 让每个域都说得通”——这就是 10.7 节练习要做的事。

## 10.7 动手练习：跑 BayesSpace，比较 q = 4 / 6 / 8

```r
library(patchwork)                  # 用 + 拼接多张 ggplot 图
set.seed(123)                       # 保证 MCMC 可复现
sce <- as.SingleCellExperiment(obj)
sce <- spatialPreprocess(sce, platform = "Visium", n.PCs = 15)

# 对三个 q 分别聚类，结果存进不同列
sce$niche_q4 <- spatialCluster(sce, q = 4, platform = "Visium",
                               d = 15)$spatial.cluster
sce$niche_q6 <- spatialCluster(sce, q = 6, platform = "Visium",
                               d = 15)$spatial.cluster
sce$niche_q8 <- spatialCluster(sce, q = 8, platform = "Visium",
                               d = 15)$spatial.cluster

# 三个结果并排画在组织切片上
clusterPlot(sce, label = "niche_q4") + clusterPlot(sce, label = "niche_q6") +
  clusterPlot(sce, label = "niche_q8")
```

**运行结果解读**：三张切片图并排，你能直观看到 q 增大时生态位如何“分裂”：q = 4 时基质区可能只有一大片；q = 8 时同一个区域被切成两三块。接着按 10.4 节的流程给每套结果命名、画细胞组成条形图，回答：**q = 6 相比 q = 4 多分出来的域，有独特的细胞组成吗？q = 8 相比 q = 6 多出来的域，能解释吗？** 把回答写进你的分析笔记——这就是你为最终 q 值准备的“决策记录”，论文方法部分可以直接引用。

## 本章小结

本章完成了课题的核心一步：用三步法把组织切成有生物学意义的生态位。第一步用 BayesSpace（马尔可夫随机场空间先验）做空间域分割并选定 q；第二步用 Squidpy 邻域富集验证细胞共定位不是巧合；第三步用 marker 与 Hallmark 通路给每个生态位命名和画像。至此你手里有了“生态位 × 细胞类型”的完整空间地图——下一章开始，我们将在这张地图上做细胞通讯分析，回答本课题真正的问题：不同生态位的通讯网络有什么不同，谁在“重编程”谁。

## 练习与思考

1. **概念题**：BayesSpace 与第 08 章 Seurat 聚类（FindClusters）的本质区别是什么？（提示：从“是否使用空间位置”和“如何保证相邻 spot 同域”两个角度回答）
2. **概念题**：为什么说“邻域富集验证”是生态位识别中不可省略的一步？如果跳过它，分割结果可能有什么问题？
3. **动手题**：按 10.7 节跑 q = 4/6/8 三套结果，用细胞组成堆叠条形图比较，说明你最终选择哪个 q 以及理由。
4. **动手题**：用 `SpatialFeaturePlot` 把侵袭前沿的两个特征基因（如 VIM、ZEB1）画到切片上，验证它们是否真的在“侵袭前沿”生态位高表达。
5. **思考题**：你的数据里免疫浸润区被分成了两个生态位，一个 T 细胞多、一个巨噬细胞多。你会把它们合并成一个“免疫浸润区”，还是保留为两个？请给出你的判断标准（提示：联系第 11–13 章“分区通讯”对生态位数量的需求）。
