# 论文要点提炼：《Enhancing Spatial Domain Detection in Spatial Transcriptomics with EnSDD》

# 1. 研究背景

• 空间转录组学（Spatial Transcriptomics, SRT）帮助理解器官功能和组织微环境，但空间域（spatial domains）的准确识别仍是挑战。
• 现有的空间域检测方法（如 BayesSpace、GraphST、STAGATE 等）存在局限，如忽略组织形态信息，或受数据复杂性影响效果不稳定。
• 传统方法难以选择最佳的检测方法，因此需要集成学习（ensemble learning）策略来整合多种方法的结果。

# 2. 研究目标

• 提出 EnSDD (Ensemble-learning for Spatial Domain Detection)，整合 8 种先进的空间域检测方法，并使用 动态加权机制 来优化每种方法的贡献。
• 该方法可 自动评估各方法的性能，无需依赖人工标注（ground truth）。
• EnSDD 还能够识别 特定空间域的差异表达基因（SVGs） 和 细胞类型空间分布，提供更深层次的组织异质性见解。

# 3. 方法创新

• 动态加权机制：EnSDD 通过优化模型自适应地分配各基方法的权重，并生成一个 共识相似度矩阵（consensus similarity matrix）。
• 自动调整聚类参数：采用 Louvain 算法的自适应分辨率调整策略，确保聚类结果的优化。
• 集成 8 种方法：包括 BayesSpace、DR-SC、GraphST、STAGATE、SpaGCN、stLearn、SiGra 和 spaVAE，结合空间坐标、基因表达数据和组织形态图像进行聚类分析。

# 4. 主要实验与结果

# (1) 在人类大脑 DLPFC 数据集上的表现

• EnSDD 在 12 个人类 DLPFC（前额叶皮层）切片上进行测试，结果显示：
  • 在大多数切片上表现最佳，准确划分了不同皮层层次。
  • 以 调整兰德指数（ARI）= 0.637 和 纯度（Purity）= 0.753 在部分切片中取得最高聚类准确度。
  • 与人工标注相比，部分区域的划分更细（如 Layer 3 进一步被细分），表明 EnSDD 具有较强的发现新结构的能力。

# (2) 在小鼠大脑数据集上的验证

• 采用 Allen Brain Atlas 作为参考标准，EnSDD 成功识别出 大脑皮层 6 层结构、海马 CA1 结构、齿状回（DG） 等重要区域。
• 相比其他方法更稳定，即使调整聚类数量，EnSDD 仍能保持空间域的连贯性，而其他方法如 BayesSpace、stLearn 在某些情况下生成不连续的区域。

# (3) 在癌症组织（乳腺癌、前列腺癌和卵巢癌）上的应用

• 乳腺癌：

  • EnSDD 识别出不同的 肿瘤区域（IDC、DCIS/LCIS） 及肿瘤边缘的细微异质性。
  • 通过差异基因分析，发现了 IDC_5 区域的关键基因（如 CSTA、FAM234B、HEBP1），揭示了潜在的肿瘤进展机制。
  • 肿瘤微环境分析 显示，IDC 5 区域的 巨噬细胞和树突状细胞 含量较高，可能与肿瘤免疫逃逸相关。

• 前列腺癌：

  • 成功识别 肿瘤区域 和 肿瘤边界（ecotone），相比其他方法，EnSDD 更精确地区分了肿瘤与周围组织的异质性。
  • 发现 KLK3、TP63、ATG9B 等基因在肿瘤边缘具有显著表达变化，可能影响前列腺癌的进展。

• 卵巢癌：

  • EnSDD 成功检测到 肿瘤 - 基质交界区（tumor-stroma interface），揭示了肿瘤微环境的复杂性。
  • 识别的关键基因 SFRP4、IGKN、THBS1 可能与肿瘤侵袭能力相关，提示 EnSDD 可用于肿瘤生物标志物的发现。

# 5. 贡献与未来展望

• EnSDD 比单一方法更稳定，可自动评估不同方法的贡献，适用于 不同类型的组织结构。
• 未来可扩展至 多切片分析（multi-slice analysis），探索组织在多个切片间的连续性。
• 结合外部参考数据 可能进一步提升 EnSDD 在未知数据集上的适用性。

# 总结

这篇论文提出了一种 基于集成学习的空间域检测方法 EnSDD，在多个数据集上表现优异，能更准确地识别 空间域结构、差异基因和细胞类型分布，并适用于 复杂组织结构及肿瘤微环境研究，为空间转录组学分析提供了新的工具。

STAMarker 和 SC3 都是利用集成学习思想的聚类方法，但它们应用的场景和侧重点有所不同：

• SC3 是一种用于单细胞 RNA 测序（scRNA-seq）数据聚类的方法。它通过整合多个聚类结果得到多个二元相似度矩阵，然后利用层次聚类生成最终的细胞分类结果。这种方法能够整合不同算法的信息，从而获得更稳健的聚类结果。

• STAMarker 则主要用于空间转录组学数据中的空间域检测。它计算不同聚类结果中任意两个空间点共同归入同一空间域的比例，形成一个连接矩阵，再利用层次聚类或 Louvain 算法来进行标签分配，从而识别出具有生物学意义的空间域。

这两种方法都是通过集成多种聚类结果来提高整体的准确性和稳健性，但 SC3 更侧重于单细胞数据的细胞分类，而 STAMarker 则专注于空间域的划分。

这是一个经典的带权 Frobenius 距离最小化问题，常出现在集成学习、矩阵逼近或多视图聚类中。你问的是下面这个优化目标：

# 📌 优化问题

其中：

• $S$ 是你要求解的集成矩阵
• $\omega_m$ 是每个 base 矩阵的权重（你暂时是常数或已知的情况下求偏导）
• $S^{(m)} \in \mathbb{R}^{n \times n}$ 是第 $m$ 个给定的相似度矩阵
• $\| \cdot \|_F$ 表示 Frobenius 范数

# 🎯 要求：

你现在是想对 $S$ 求偏导数，也就是：

# 🧮 计算步骤

我们知道：

所以每一项：

对 $S$ 求导的标准结果是：

# 🧾 于是整体偏导为：

# ✅ 最终结果：

# 💡 如果你要最小化它，只需让偏导为 0：

也就是说，最优解 $S$ 是所有矩阵的 加权平均！

需要我帮你把这推导写成 LaTeX 公式或加到你自己的代码推导里也可以继续说哈～

太好了！你这个问题非常经典，咱们现在就来一步步推导怎么用 Lagrange 乘子法（Lagrange multipliers）来求解权重 $\omega_m$ 的最优解。

# 🔢 你给的目标函数是：

带有约束：

其中 $S$ 是固定的（或者在交替优化里暂时固定），所以我们只针对 $\omega_m$ 来优化。

# ✏️ 第一步：写出 Lagrangian 拉格朗日函数

为了处理约束 $\sum \omega_m = 1$，引入拉格朗日乘子 $\mu$，构造：

# 🧮 第二步：对每个 $\omega_m$ 求偏导数

我们令每个偏导为 0：

整理一下：

再除以 $\lambda$ 并指数两边：

整理一下：

其中 $C = \exp\left( \frac{\mu - \lambda}{\lambda} \right)$ 是一个所有 $\omega_m$ 的共同常数。

# ✅ 第三步：用约束条件求出 $C$

由于所有 $\omega_m$ 要满足：

我们将上面形式代入得：

# 🎉 最终结果：

这其实就是一个 softmax 形式！每个 $\omega_m$ 是对 $S^{(m)}$ 的一个相似度评分（用负的距离衡量），用 softmax 来归一化。

# 📌 补充说明

• 如果 $\lambda \to 0$：权重趋向于 “winner takes all”（只选最小误差的那个）
• 如果 $\lambda \to \infty$：所有权重接近均匀平均
• 这个方法就是你之前代码中那个权重更新公式的数学基础

如果你想我也可以帮你推对 $S$ 的梯度，或者画图解释这个 softmax 权重的直觉含义 😄 想继续哪一块？

当然可以！我们就专注于这一步：如何对每个 $\omega_m$ 求偏导数。

你给的函数是：

我们现在只关注其中一个 $\omega_m$，要对它求偏导数。

# 🔍 拆开看，(\mathcal {L}) 对于单个 (\omega_m) 是：

选定一个 $m$，(\mathcal {L}$ 中关于 $\omega_m$ 的部分是：

这是三项：

1. 线性项：$\omega_m \|S - S^{(m)}\|_F^2$
2. 熵正则项：$\lambda \omega_m \log \omega_m$
3. 约束乘子项：$-\mu \omega_m$

# 🧮 逐项对 $\omega_m$ 求导

# 第一项：

这是因为 $\|S - S^{(m)}\|_F^2$ 是常数，对 $\omega_m$ 直接是乘法。

# 第二项（关键点）：

我们来慢慢推一下这个常见表达式的导数：

令：

用乘积法则：

所以：

# 第三项：

# 🔚 把三项加起来，得到偏导：

然后我们令偏导为 0，得到最优条件：

就是我们之前推导的式子啦！

# 🎯 总结一下你学到什么：

• 对 $\omega \log \omega$ 求导是 $\log \omega + 1$（这是熵类正则常见的导数）
• 在加权平均里用 softmax，其实就是从这个偏导数出发
• 拉格朗日乘子 $\mu$ 用于处理总和为 1 的约束

如果你想，我还可以做成一张「可视化推导卡片」，更方便你日后复习！要不要试试？

# 空间自相关（Spatial Auto-correlation, SA）的使用

在 EnSDD 研究中，空间自相关用于分析 基因表达在空间上的依赖性，即在组织切片中，基因的表达水平是否在空间上呈现聚集或分散模式。这有助于识别 空间可变基因（SVGs, Spatially Variable Genes），即在特定空间区域内具有特定表达模式的基因。

# 1. 什么是空间自相关？

空间自相关（SA）是指某一变量（这里是基因表达水平）在空间上的分布模式是否表现出 相似性或相关性。一般来说：

• 正空间自相关（Positive SA）：相邻位置的基因表达水平相似，呈现聚集模式。例如，在大脑皮层某个区域，如果某个基因的表达量高，其周围位置的表达量也较高。
• 负空间自相关（Negative SA）：相邻位置的基因表达水平相反，呈现分散模式。例如，某些基因在神经组织和周围支持细胞中表现出相反的表达模式。

EnSDD 采用 Local Getis and Ord’s Gi 统计量 来定量衡量某个基因在局部空间中的自相关性。

# 2. Local Getis and Ord’s $G_i$ 统计量

EnSDD 通过 Local Getis and Ord’s $G_i$ 统计量来评估基因表达的空间聚集性。这一统计量用于衡量每个基因在空间中的局部自相关性，数学表达如下：

其中：

• $x_j$：第 $j$ 个空间点的基因表达水平。
• $\mu$, $\sigma$：所有空间点（去除第 $i$ 个点）的基因表达的均值和标准差。
• $c_{ij}(k)$：空间权重矩阵，表示空间点 $i$ 与 $j$ 之间的距离权重（例如，最近的 k 个点的权重为 1，其余为 0）。
• $V_i = \sum_j c_{ij}$：第 $i$ 个点的总空间权重。
• $B_{1i} = \sum_j c_{ij}^2$：空间权重的平方和。
• $n$：所有空间点的总数。

解释：

• 如果 $G_i$ 值 大于 0，表示该基因的高表达区域在局部空间中聚集（正空间自相关）。
• 如果 $G_i$ 值 小于 0，表示该基因在该区域的表达呈现离散模式（负空间自相关）。

EnSDD 采用 $k = 6$ 作为最近邻窗口大小，计算基因在每个空间点上的 $G_i$ 统计量，并使用 Wilcoxon 检验 进一步确认哪些基因在特定空间域上表现出显著的空间可变性。

# 3. SA 在 EnSDD 中的应用

EnSDD 主要利用 Local Getis and Ord’s $G_i$ 统计量来：

1. 筛选空间可变基因（SVGs）：识别在不同空间域（Spatial Domain）中表达模式具有明显局部聚集或分散特征的基因。
2. 分析组织结构中的基因表达模式：
   • 例如，在 人类 DLPFC（前额叶皮层）数据 中，EnSDD 发现 Layer 3 被进一步细分，并通过 $G_i$ 统计量证实该区域内确实存在两个不同的基因表达聚集区域。
   • 在 乳腺癌数据 中，利用 $G_i$ 统计量找到了一些 与癌症进展相关的特异性基因（如 HEBP1、FAM234B、CSTA），它们在某些肿瘤区域内高度聚集。

# 4. 具体实现

# 步骤 1：计算局部空间自相关

对于每个基因：

1. 构建 空间权重矩阵，选取最近的 6 个相邻点 作为参考。
2. 计算 每个点的 $G_i$ 统计量，得到基因的局部空间聚集度信息。

# 步骤 2：Wilcoxon 检验

对比不同空间域之间的基因表达水平，检验在某个空间域中是否存在显著表达变化的基因。

# 步骤 3：可视化和分析

• 局部空间自相关（SA）热图：展示不同基因在不同空间点上的 $G_i$ 值，识别高表达聚集区域。
• 空间表达模式图：在组织切片上可视化基因的表达水平，观察其空间分布。

# 5. 实例解析

# 案例 1：人类前额叶皮层（DLPFC）

在 DLPFC 数据集中，EnSDD 发现：

• 在某些切片中，Layer 3 被自动分成两个子区域。
• 通过 $G_i$ 统计量计算，发现 129 个基因在这两个子区域的表达显著不同。
• 进一步分析发现这些基因 与神经功能和甲状腺癌相关，揭示了 Layer 3 内部可能存在未被人工标注的异质性。

# 案例 2：乳腺癌肿瘤区域

• EnSDD 识别了 肿瘤中心区域（IDC_5） 和 周围免疫浸润区域，并通过 $G_i$ 统计量找到了一组特异性基因：
  • 在 IDC_5 高表达的基因：HEBP1（与乳腺癌进展相关）
  • 在肿瘤边界高表达的基因：FAM234B（与肿瘤免疫浸润相关）

这些基因的空间表达模式揭示了 不同肿瘤区域的分子特征和潜在生物学机制，并可能成为新的肿瘤治疗靶点。

# 6. 总结

• 空间自相关（SA）是 EnSDD 识别空间可变基因（SVGs）的关键步骤。
• Local Getis and Ord’s $G_i$ 统计量用于衡量基因在局部空间的聚集或分散模式。
• 通过 SA 分析，EnSDD 发现了新的 组织异质性，并揭示了 肿瘤微环境 中的关键分子特征。
• 该方法在 人类大脑 DLPFC、小鼠大脑、乳腺癌、前列腺癌和卵巢癌 数据中均表现良好，能够有效提升空间转录组数据的生物学解释能力。

如果你对实现细节或代码感兴趣，可以让我帮你进一步解析如何在 Python 或 R 中实现 Local Getis and Ord’s $G_i$ 统计量计算！🚀

# RCTD（Robust Cell Type Decomposition）

在 EnSDD 研究中，RCTD（Robust Cell Type Decomposition） 主要用于 推测空间转录组数据（SRT）中的细胞类型组成，即利用 单细胞转录组测序（scRNA-seq）数据 作为参考，对空间转录组数据中的每个 spot 进行 细胞类型分解（deconvolution）。

# 1. 什么是 RCTD？

RCTD 是一种 基于权重优化的去卷积方法，用于从空间转录组数据中估计每个 spot（空间点）中 不同细胞类型的组成比例。

在空间转录组数据（SRT）中，每个 spot 包含多个细胞的混合信号，而单细胞 RNA 测序（scRNA-seq）数据可以提供 单细胞级别的基因表达信息。RCTD 通过匹配 SRT 和 scRNA-seq 数据，估算每个 spot 内不同细胞类型的比例。

关键思想：

• 输入数据：
  • 空间转录组数据（SRT）：基因表达矩阵，每个 spot 可能包含多个细胞的混合基因表达。
  • 单细胞 RNA-seq 数据（scRNA-seq）：提供每种细胞类型的基因表达模式。
• 输出数据：
  • 每个 spot 内部的 细胞类型组成比例，例如：

    1 2	Spot 1: 70% 肿瘤细胞, 20% 巨噬细胞, 10% 成纤维细胞 Spot 2: 50% 神经元, 30% 星形胶质细胞, 20% 少突胶质细胞

应用场景：

• 分析肿瘤微环境（TME），推测肿瘤区域内不同免疫细胞的分布情况。
• 探索大脑结构，研究不同皮层层次的细胞组成差异。
• 研究发育过程，分析不同组织区域的细胞类型变化。

# 2. RCTD 计算过程

# (1) 计算细胞类型的 “特征”

RCTD 通过 scRNA-seq 数据，获取每种细胞类型的特征基因表达模式。例如：

• 细胞类型 A：基因 1 高表达，基因 2 低表达。
• 细胞类型 B：基因 1 低表达，基因 2 高表达。

这些信息用于 构建细胞类型的特征矩阵。

# (2) 计算 SRT 数据中每个 spot 的基因表达模式

每个 spot 含有多个细胞，因此它的基因表达是 所有细胞类型的混合。RCTD 需要解出每个 spot 内 不同细胞类型的贡献比例。

# (3) 通过权重优化分解细胞类型

RCTD 假设 SRT 中每个 spot 的基因表达 $X$ 由 不同细胞类型的加权平均表达 组成：

其中：

• $X$ 是 SRT 数据的基因表达矩阵（待解）。
• $C$ 是 scRNA-seq 参考数据中每种细胞类型的基因表达模式（已知）。
• $W$ 是权重矩阵（待求解），表示每种细胞类型在每个 spot 的贡献比例。
• $\epsilon$ 是噪声项。

RCTD 采用 贝叶斯方法 + 线性优化 来求解 $W$，从而估计每个 spot 内部的细胞类型组成。

# (4) 归一化 & 计算细胞类型分布

为了确保所有细胞类型的比例总和为 1，RCTD 进行归一化处理，最终输出每个 spot 的细胞组成。

# 3. RCTD 在 EnSDD 研究中的应用

EnSDD 通过 RCTD 进行细胞类型推测，帮助分析 不同空间域的细胞组成，特别是在 肿瘤微环境分析 方面表现突出。例如：

# (1) 在乳腺癌数据中的应用

• 目标：分析肿瘤中心（IDC_5）与肿瘤边缘的细胞类型组成差异。
• 结果：
  • IDC_5 肿瘤中心：癌细胞比例高，免疫细胞（如巨噬细胞、树突状细胞）比例低。
  • 肿瘤边缘区域：免疫细胞（如巨噬细胞和树突状细胞）显著富集，可能与 肿瘤免疫逃逸机制 相关。

细胞比例结果（示例）：

# (2) 在前列腺癌数据中的应用

• 目标：研究肿瘤区域（Tumor）和肿瘤边界（Ecotone）的细胞类型差异。
• 结果：
  • 肿瘤区域：以 LE-KLK3 阳性上皮细胞 为主。
  • 肿瘤边界（Ecotone）：更高比例的 免疫细胞（MNP 和 B 细胞），可能影响肿瘤的进展和免疫治疗反应。

# (3) 在卵巢癌数据中的应用

• 目标：分析肿瘤、肿瘤 - 基质交界（Ecotone）和基质区域的细胞组成。
• 结果：
  • 肿瘤区域：EOC（上皮性卵巢癌细胞）占比高。
  • 肿瘤 - 基质交界：富含 上皮细胞、成纤维细胞（CAFs）、间皮细胞，可能在肿瘤侵袭过程中发挥作用。

# 4. RCTD 与其他方法的对比

相比 Seurat 和 Stereoscope，RCTD 对噪声的鲁棒性更强，适用于高复杂度的组织，如癌症组织和大脑皮层。

# 5. 结论

• RCTD 用于估计 SRT 数据中每个 spot 的细胞组成，帮助 EnSDD 解析不同空间域的细胞类型分布。
• 通过 scRNA-seq 参考数据，RCTD 能够揭示 肿瘤微环境、神经组织层次结构等。
• 在 乳腺癌、前列腺癌、卵巢癌 数据中，RCTD 发现了 肿瘤 - 免疫相互作用的关键特征，有助于癌症生物标志物发现和治疗策略优化。

# Purity（Pur，纯度）

Purity（纯度） 是一个用于评估聚类（Clustering）质量的指标，衡量聚类结果与真实类别标签（Ground Truth）的一致性。它的取值范围为 $$\begin{align}0, 1\end{align}$$，值越高，表示聚类结果越接近真实分类。

# 1. Purity 公式

其中：

• $N$：样本总数（数据点总数）。
• $L$：预测聚类的数量（即 EnSDD 识别出的空间域数）。
• $T$：真实类别（Ground Truth）数量（如人工标注的组织结构）。
• $u_l$：第 $l$ 个预测聚类（即 EnSDD 识别出的某个空间域）。
• $v_t$：第 $t$ 个真实类别（如已知组织切片的层级）。
• $|u_l \cap v_t|$：预测聚类 $u_l$ 内属于真实类别 $v_t$ 的样本数量。
• $\max_t |u_l \cap v_t|$：对于每个预测的聚类 $u_l$，找出其中最多的真实类别样本数量。

解释：

1. 对于每个预测的聚类 $u_l$，找出其中占比最多的真实类别 $v_t$。
2. 求所有预测聚类的占比最多真实类别样本数的总和，并除以总样本数 $N$，得到纯度值。

# 2. Purity 计算示例

# (1) 示例数据

假设我们有 10 个数据点，它们的真实类别（Ground Truth）和聚类算法（如 EnSDD）给出的预测聚类如下：

# (2) 计算每个聚类的最大类别

• 聚类 1（u₁）：

  • 真实类别 A：2 个
  • 最大类别：A（2）

• 聚类 2（u₂）：

  • 真实类别 A：1 个
  • 真实类别 B：2 个
  • 最大类别：B（2）

• 聚类 3（u₃）：

  • 真实类别 B：1 个
  • 真实类别 C：4 个
  • 最大类别：C（4）

# (3) 计算 Purity

# 3. Purity 在 EnSDD 论文中的应用

在 EnSDD 研究中，Purity 主要用于 评估空间域划分（Spatial Domain Detection）的准确性，即检查 EnSDD 生成的空间域与人工标注的组织结构（如大脑层次结构或肿瘤组织分区）的匹配程度。

# 实验结果示例

• 在 人类大脑 DLPFC 数据集 中：
  • EnSDD 纯度（Purity）= 0.704，优于 STAGATE（0.678）和 spaVAE（0.672）。
• 在 乳腺癌数据 中：
  • EnSDD 纯度（Purity）= 0.668，高于其他方法。

这些结果表明，EnSDD 生成的空间域划分与真实组织结构的匹配度较高，能够更好地识别不同的组织区域。

# 4. Purity 与 Adjusted Rand Index（ARI）的区别

总结：

• Purity 直观易懂，适合快速评估聚类质量，但它 不考虑错误分类的影响，如果一个大类被错误拆分成多个小类，Purity 仍然可能较高。
• ARI 更加严格，它不仅考虑正确分类，还会对错误分类进行惩罚，因此对于聚类方法的评估更为精确。

在 EnSDD 论文中，Purity 和 ARI 一起使用，以确保空间域检测方法的可靠性。

# 5. 结论

• Purity 衡量聚类质量，数值越高，说明 EnSDD 预测的空间域更接近真实组织结构。
• 与 ARI 结合使用时，Purity 提供直观的匹配度衡量，而 ARI 可检测错误分类的影响。
• 在 EnSDD 研究中，Purity 证明了该方法比 STAGATE、spaVAE 等方法具有更好的空间域检测能力。

如果你想了解 代码如何计算 Purity 或者 如何优化空间域划分算法，可以告诉我，我可以进一步解析！

当然可以，我们来逐行解读这个 apply_BayesSpace 函数。它的功能是：在一个 Seurat 对象上使用 BayesSpace 包进行空间转录组聚类分析。你这段代码流程很完整，也有些细节挺有意思，我会分块讲解，穿插解释函数作用和每一步背后的意图。

# 🧠 函数定义和参数说明：

1
2

apply_BayesSpace <- function(Seurat.data, platform = "Visium", n.PCs = 15, n.HVGs = 2000,
                             n.setting = 7, nrep = 50000){

• Seurat.data ：传入一个处理过的 Seurat 对象，假设已经做完 SCTransform。
• platform ：平台类型，如 "Visium" 或 "ST" （决定空间结构等处理方式）。
• n.PCs ：PCA 主成分数量，默认为 15。
• n.HVGs ：选择前多少个高变基因（Highly Variable Genes）。
• n.setting ：聚类数量，即 q 值，默认是 7 类。
• nrep ：MCMC 采样的迭代次数，BayesSpace 聚类时使用，默认 50000。

# 📦 步骤一：从 Seurat 对象中提取数据

1
2

exp_data <- list(counts = as.matrix(Seurat.data@assays$SCT@counts),
                 logcounts = as.matrix(Seurat.data@assays$SCT@data))

• 提取标准化前后的表达矩阵：
  • counts : 原始计数
  • logcounts : SCT 标准化后的 log-scale 表达值

1

df <- Seurat.data@meta.data[, c("x", "y", "pixel_x", "pixel_y")]

• 从元数据中提取空间坐标。

1
2
3
4

temp_colnames <- colnames(df)
temp_colnames <- sub("x", "row", temp_colnames)
temp_colnames <- sub("y", "col", temp_colnames)
colnames(df) <- temp_colnames

• BayesSpace 需要用 row , col 命名空间位置，所以这里做了列名替换。

# 📦 步骤二：构造 SingleCellExperiment 对象

1

sce <- SingleCellExperiment::SingleCellExperiment(assay = exp_data, colData = df)

• 创建一个符合 Bioconductor 标准的数据容器，包含表达矩阵和元数据。
• BayesSpace 是以 SingleCellExperiment 为输入格式。

# 📈 步骤三：寻找高变基因并 PCA

1
2

dec <- scran::modelGeneVar(sce)
top <- scran::getTopHVGs(dec, n = n.HVGs)

• 用 scran 估计每个基因的变异程度，选择前 n.HVGs 个变异最大的基因（通常对降维和聚类效果好）。

1
2

set.seed(102)
sce <- scater::runPCA(sce, subset_row=top)

• 在 top HVGs 上做 PCA 降维。

# 🧼 步骤四：预处理空间信息

1

sce <- BayesSpace::spatialPreprocess(sce, platform=platform, skip.PCA=TRUE)

• 执行 BayesSpace 所需的空间预处理：
  • skip.PCA=TRUE 表示我们已经做过 PCA，不重复。
  • 添加空间网格、邻接矩阵等结构信息。

# 🎚️ 步骤五：设置 gamma 超参数（平台特异）

1
2
3
4
5
6

if(platform == "Visium"){
  gamma <- 3
}
if(platform == "ST"){
  gamma <- 2
}

• gamma 是 BayesSpace 聚类模型的平滑参数，控制聚类边界的 “平滑程度”。

# 🔁 步骤六：执行 BayesSpace 聚类

1
2
3

set.seed(150)
sce <- BayesSpace::spatialCluster(sce = sce, q = n.setting, d = n.PCs,
                                  platform = platform, nrep = nrep, gamma = gamma)

• 开始空间聚类：
  • q : 聚类的簇数
  • d : 使用的 PCA 数量（如前 15 个）
  • nrep : MCMC 迭代次数（越多越稳定但越慢）

# 📊 步骤七：提取聚类结果

1
2

res.cluster <- sce$spatial.cluster
names(res.cluster) <- rownames(sce@colData)

• 聚类结果是一个向量（每个 spot 的聚类标签）。
• 给结果加上对应 spot 的名字。

# 🔁（注释掉的代码）自动寻找最佳聚类数

1
2
3

# sce <- BayesSpace::qTune(...) 
# opt_k <- ...
# sce <- BayesSpace::spatialCluster(q = opt_k, ...)

• 这段注释代码原本用于自动寻找最优 q 值（聚类数），通过比较对数似然。但作者这里选择手动设定 n.setting 。

# 🔚 返回结果

1

return(res.cluster)

• 返回一个聚类结果向量，用户可以根据它进行可视化、下游分析等。

# ✅ 总结一句话：

这个函数将 Seurat 对象转换为 SingleCellExperiment ，挑选高变基因做 PCA，然后用 BayesSpace 进行空间聚类分析，并返回每个 spot 的聚类标签。

你可以配合 BayesSpace::clusterPlot(sce) 进行可视化，也可以把这个结果加到 Seurat 对象里，继续做差异分析。

如果你还想了解聚类原理或者对某个包函数不理解，可以继续问我！

当然可以！下面我来逐行解读这段 apply_SpaGCN 函数的 R 代码。这段代码的功能是：调用 Python 实现的 SpaGCN 空间聚类算法，并返回聚类结果。它通过 reticulate 接口从 R 中运行 Python 脚本。

# 🧠 函数定义及参数说明

1
2
3
4
5
6
7
8
9
10

apply_SpaGCN <- function(Seurat.data,
                         SpaGCN_beta = 55,
                         SpaGCN_alpha = 1,
                         SpaGCN_p = 0.5,
                         SpaGCN_l = 0.5,
                         num.cluster = 7,
                         res.setting = NULL,
                         SpaGCN_lr = 0.05,
                         SpaGCN_epoches = 200,
                         platform = "Visium"){

• Seurat.data : Seurat 对象（空间数据），包含图像路径信息
• SpaGCN_beta , SpaGCN_alpha , SpaGCN_p , SpaGCN_l : 是 SpaGCN 聚类算法的超参数，控制空间权重、邻域强度等
• num.cluster : 聚类簇数量
• res.setting : 聚类分辨率设置（用户自定义）
• SpaGCN_lr : 学习率
• SpaGCN_epoches : 训练迭代次数
• platform : 平台，如 "Visium"，决定空间坐标解读方式

# 📦 1. 提取图像路径

1
2
3
4

temp <- Seurat.data@images
counts_path <- temp$counts_path
meta_path <- temp$loc_path
img_path <- temp$img_path

假设 Seurat.data@images 存放了一些自定义字段（比如用 $counts_path 来表示 count 矩阵路径、meta 信息和图像路径）：

• counts_path : 表达矩阵（spot × gene）的文本路径
• meta_path : spot 的空间位置信息
• img_path : H&E 图像文件路径

⚠️ 这些字段可能是你自己前期写入 SeuratObject@images 的，不是默认字段。

# 🐍 2. 加载 Python 脚本

1

reticulate::source_python(system.file("python", "apply_spaGCN.py", package = "EnSDD"))

• source_python(...) ：通过 R 的 reticulate 包加载 Python 脚本
• system.file(...) ：定位 R 包 EnSDD 中的 apply_spaGCN.py 脚本路径
  • 假设你安装了一个 R 包叫 EnSDD ，它里面封装了运行 SpaGCN 所需的 Python 脚本

📝 被注释的 EnSDD_pyexec(...) 也可能是某个自定义封装，但你实际使用的是 reticulate 。

# 🚀 3. 调用 Python 的 run_spaGCN 函数

1
2
3
4
5
6
7
8
9

run_spaGCN(counts_path, meta_path, img_path,
           platform = platform,
           res_setting_by_user = res.setting,
           b = as.integer(SpaGCN_beta),
           a = SpaGCN_alpha,
           p = SpaGCN_p,
           n_clusters = as.integer(num.cluster),
           lr = SpaGCN_lr,
           max_epochs = as.integer(SpaGCN_epoches))

• 这个函数在 apply_spaGCN.py 中定义，接收参数并执行 SpaGCN 空间聚类任务
• SpaGCN 是一个 GCN（图卷积网络）框架，结合空间坐标和图像信息进行聚类
• 参数解析：
  • b = beta ：控制空间光滑性
  • a = alpha 、 p ：控制边权的计算
  • n_clusters ：聚类个数
  • res_setting_by_user : 分辨率
  • lr : 学习率
  • max_epochs : 迭代次数

✅ 此处训练过程全在 Python 中执行，输出为文件结果。

# 📁 4. 加载结果并处理

1
2
3
4

res.path <- paste0(dirname(img_path), "/spaGCN_res", "/spaGCN_label.txt")
res_temp <- read.table(res.path, sep = "\t")
res <- res_temp$V2
names(res) <- res_temp$V1

• Python 会在图片所在文件夹下创建 spaGCN_res/spaGCN_label.txt
• 文件中记录的是： spot名称 \t 聚类标签
• 把 V2 作为聚类标签， V1 为 spot 名，转换为一个 R 命名向量

# 🔢 5. SpaGCN 聚类结果中可能含有 "0" 类 → 调整为从 1 开始

1
2
3

if(sum(res == 0) != 0){
  res <- res + 1
}

• 如果存在 0 类，就全部 +1，使得聚类编号从 1 开始（防止误差）

# ✅ 6. 根据 Seurat 的列顺序排列聚类结果，并清理中间文件夹

1
2
3

res <- res[colnames(Seurat.data)]
unlink(paste0(dirname(img_path), "/spaGCN_res"), recursive=TRUE)
return(res)

• 确保返回的聚类向量顺序和 Seurat 中列顺序一致
• 删除 SpaGCN 临时输出文件夹，避免污染
• 返回：spot-level 聚类标签

# 🎉 总结一句话：

apply_SpaGCN() 是一个 R 函数，用于封装和调用 Python 的 SpaGCN 空间聚类过程，输入是 Seurat 对象，输出是每个 spot 的聚类标签。

如果你想看 apply_spaGCN.py 文件里 run_spaGCN() 是如何实现的，也可以把那个 Python 文件贴给我，我可以继续帮你逐行解读～

或者你想把这个结果接回 Seurat 进行可视化的话，也可以问我哦。

当然可以！这段代码是一个 STAGATE 方法的调用和结果处理模块，通常嵌套在一个更大的函数体内，负责：

运行 STAGATE 空间聚类，记录运行时间，处理聚类结果并保存。

我会按逻辑结构详细解释每个部分👇：

# ✅ 1. 判断是否执行 STAGATE 聚类

1

if(STAGATE){

• 如果传入的 STAGATE 参数为 TRUE ，就会执行 STAGATE 聚类分析。

# 🧪 2. 调用 apply_STAGATE() 函数并计时

1
2
3
4
5

cat("Run STAGATE...\n")
pred.STAGATE <- Sys.time()
temp.STAGATE <- apply_STAGATE(Seurat.data, number.setting.by.user, STAGATE_alpha)
end.STAGATE <- Sys.time()
time.STAGATE <- difftime(end.STAGATE, pred.STAGATE, units = "mins")

• 打印提示语： Run STAGATE...
• 记录开始时间 pred.STAGATE
• 调用 STAGATE 聚类函数 apply_STAGATE() ，传入参数：
  • Seurat.data ：Seurat 格式空间转录组对象
  • number.setting.by.user ：用户设定的聚类簇数
  • STAGATE_alpha ：控制聚类图平滑程度的参数（STAGATE 超参数）
• 记录结束时间 end.STAGATE
• 计算运行时间（分钟） time.STAGATE

# 📝 3. 保存 STAGATE 聚类结果到列表对象中

1
2
3
4

cat("Run time for STAGATE: ", time.STAGATE, "min","\n")
time.methods <- c(time.methods, time.STAGATE)
res.clustering$STAGATE <- temp.STAGATE
k <- k + 1

• 输出 STAGATE 的运行时间
• 把运行时间加入时间列表 time.methods
• 把聚类标签结果 temp.STAGATE 保存进聚类结果列表 res.clustering$STAGATE
• k 是用来统计或记录方法数量的变量，这里自增（可能用于 ensemble 阶段）

# 🧼 4. 检查是否有包含 NA 的聚类方法，清理无效结果

1
2
3
4
5
6
7
8

names(time.methods) <- Methods.used
index.list <- c()
for (i in 1:length(res.clustering)) {
  na.index <- sum(is.na(res.clustering[[i]]))
  if(na.index != 0){
    index.list <- c(index.list, i)
  }
}

• 将时间向量命名（与方法对应）
• 检查每种聚类方法的结果中是否存在 NA 值
• 如果某方法含 NA，就将其索引加入 index.list

# 🧽 5. 删除含 NA 的聚类方法

1
2
3
4
5
6

if(length(index.list) != 0){
  methods.names.tmp <- names(res.clustering)[index.list]
  for (i in methods.names.tmp) {
    res.clustering[[i]] <- NULL
  }
}

• 若有 NA，获取含 NA 方法的名称
• 从结果列表中将其移除（设置为 NULL）

# 📐 6. 聚类向量转为矩阵（后续可用于 ensemble）

1

res.clustering.matrix <- lapply(res.clustering, clust2Mat)

• clust2Mat ：将每种方法的聚类标签向量转换为 二值指示矩阵（每行一个样本，每列一个簇，1 表示属于某簇）

# 📦 7. 打包最终结果为列表对象

1
2
3
4
5

res <- list(
  clustering.vec = res.clustering,
  clustering.mat = res.clustering.matrix,
  time = time.methods
)

结构说明：

• clustering.vec : 不同方法对应的聚类标签向量（含 STAGATE）
• clustering.mat : 聚类标签的二值表示矩阵（用于后续集成）
• time : 每种方法的运行时间

# 💾 8. 可选保存结果到 .RData

1
2
3
4
5

if(saving_results)
{
  cat("Saving individual clustering analysis results...\n")
  save(res, file = "res.clustering.RData")
}

• 如果 saving_results = TRUE ，则保存结果到 res.clustering.RData 文件

# 🔚 9. 返回最终聚类结果对象

1

return(res)

# 🔁 总结一下这段代码做了什么？

如果你还想看看 apply_STAGATE() 的具体内容（即 STAGATE 是如何被调用、处理数据的），可以发我代码，我也可以逐行讲解！

或者你也可以告诉我，你是在哪个项目 / 代码库里看到的这些，比如是 EnSDD、SpaGCN 或其他空间分析工具～我可以更有针对性地解释模块设计逻辑。

clust2Mat: Function to compute adjacency matrix of cluster graph given a vector of cluster memberships
Description
Given a list of cluster memberships, this function computes the adjacency matrix for the associated cluster graph. The adjacency matrix is a matrix whose rows and columns are the nodes of the cluster graph and whose entries are 0's or 1's. A 1 entry indicates that the corresponding nodes are connected, and a 0 indocates that they are not.
Usage
clust2Mat(memb)
Arguments
memb
A numeric vector, with each entry representing a node, the entry's value being the number of the cluster to which that node belongs.
Value
Details
Given a vector of cluster membership numbers, where the number of entries is the number of nodes n, the function computes an nxn ``adjacency matrix'' for the corresponding cluster graph. The cluster graph is the graph in which two nodes are connected by an edge if and only if they are members of the same cluster. The adjacency matrix for the graph has rows and columns representing the nodes, in the same order as the input vector. The (i,j) entry is 1 if and only if node i and node j are in the same cluster. Otherwise, the entry is 0. By convention, diagonal entries are 0.

当然可以！我们一起来逐行解读这段代码，它的作用是：

基于多个聚类数量设置运行聚类分析，并可视化每种聚类结果（例如用于 EnSDD 的聚类 ensemble 处理）。

它主要涉及：

• 加载 Python 脚本（可视化）
• 读取多个不同聚类数下的聚类结果
• 做 ensemble 聚类（融合多个方法 / 模型结果）
• 用 Louvain 算法在相似度矩阵上聚类
• 存储结果并生成聚类可视化图像

# 🔧 准备阶段

1
2
3

## setting the number of clustering equal to the annotation provided by SEDR
library(reticulate)
reticulate::source_python(system.file("python", "Visualization_HBC_clustering.py", package = "EnSDD"))

• reticulate : 用于 R 调用 Python
• source_python(...) : 加载 Python 文件（用于空间聚类可视化 visualization_clustering() ）

1
2

label_setting <- c(10, 16, 17, 20, 25)
sankey_df <- list()

• label_setting : 聚类数的多个设置（不同聚类簇个数）
• sankey_df : 用来保存每个设定下的聚类标签，可能用于后续的 Sankey plot 可视化流向图

# 🔁 遍历多个聚类数设定

1
2

for (i in 1:length(label_setting)) {
  message(paste0("Setting the number of clusters: ", label_setting[i]))

逐个设置聚类簇数，准备执行后续流程。

# 📦 加载已有的聚类结果 .RData

1
2

res.tmp <- load(paste0("../data/res_individual_", label_setting[i], "_all.RData"))
assign("res", get(res.tmp))

• 加载保存好的结果文件
• 结果中保存了 clustering.mat （聚类标签矩阵）和 clustering.vec （单一聚类标签向量）

# 🧠 聚类结果融合（Ensemble）

1

res.ensemble <- solve_ensemble(res$clustering.mat, lambda = NULL, prob.quantile = 0.5, niter = 100, epsilon = 1e-5)

使用 solve_ensemble() 函数对聚类矩阵做 集成聚类融合：

• clustering.mat : 多模型、多参数下的聚类标签（行：spot，列：不同聚类方案）
• 输出 res.ensemble$H : 相似度矩阵（用于后续 Louvain 聚类）

# 📊 Louvain 聚类（基于集成的相似度网络）

1
2

sim_net_louvain <- res.ensemble$H
res.ensemble.clustering.louvain <- Louvain_clustering_binary(Seurat.data, similarity_mat = sim_net_louvain, setting_k = label_setting[i], max.res = 8)

• 使用 Louvain_clustering_binary() 对相似度矩阵 H 进行社区检测（聚类）
• setting_k : 目标聚类数
• max.res : 最大分辨率（控制 Louvain 聚类是否更细或更粗）

# 📝 聚类结果格式化并保存

1
2
3
4

meta_data <- Seurat.data@meta.data
res.vec <- res$clustering.vec
res.vec.df <- as.data.frame(res.vec)
res.vec.df$EnSDD <- res.ensemble.clustering.louvain[rownames(meta_data)]

• 构造一个包含两个聚类标签的 data frame：
  • 原始标签： res.vec
  • 融合聚类标签（EnSDD 集成结果）： res.ensemble.clustering.louvain

# 🗂 结果写入文件并加入总表

1
2
3
4
5

domain_name <- paste0('Domain_', label_setting[i])
sankey_df[[domain_name]] <- res.ensemble.clustering.louvain[rownames(meta_data)]

meta_data_all <- cbind(meta_data, res.vec.df)
write.table(meta_data_all, file = paste0("../data/clustering.res.", label_setting[i], ".txt"), sep = "\t")

• 把每个聚类设定的集成标签加入 sankey_df ，方便后续绘制 Sankey 流图
• 将聚类结果保存为 .txt 文件，以便后续可视化或分析

# 🖼 调用 Python 进行聚类可视化

1
2
3
4
5

meta_data_path <- paste0("../data/clustering.res.", label_setting[i], ".txt")
visualization_clustering(counts_path, meta_data_path, img_path,
                         cluster_num = as.integer(label_setting[i]),
                         spot_diameter_fullres = 150,
                         spot_size = 1.5)

使用之前加载的 Visualization_HBC_clustering.py 中的 visualization_clustering() ：

• 输入：
  • 表达数据（ counts_path ）
  • 聚类结果路径
  • 图像路径
• 输出：空间聚类的图像（spot 被聚类标签着色）

# 🧹 清理（注释掉）

1

# unlike('res.clustering.txt')

这行可能是个笔误或作者写错了，可能本意是 unlink() 删除文件。

# ✅ 总结一句话：

该代码是为了测试多个聚类簇数设置下的空间聚类效果，结合 EnSDD 的聚类融合方法和 Python 的可视化函数，输出空间聚类图和标签，为后续分析（如 Sankey plot）做准备。

如果你有任何一个函数想继续详细分析（比如 solve_ensemble() 或 visualization_clustering() 的内部逻辑），可以发我，我继续给你拆解！😎

好呀！我们来详细解读一下这个 solve_ensemble 函数 —— 它实现的是一个聚类集成（Ensemble Clustering）方法，用于融合多个聚类结果，输出一个共识相似度矩阵 H ，并自动学习每个聚类结果的权重 w 。

# 🧠 目标简述

1

solve_ensemble(Results.clustering, ...)

• 输入：多个聚类方法生成的样本 - 样本相似度矩阵（比如你前面问过的 clust2Mat 转换出来的）
• 输出：
  • H ：最终融合后的相似度矩阵（作为新的聚类输入）
  • w ：每个方法的权重（表示该方法对最终共识的贡献程度）

该方法基于一个优化目标，使用迭代优化来寻找最合适的加权平均。

# 📥 输入参数说明

# 🔧 算法结构详解

# 1. 参数初始化

1
2
3

num.methods <- length(Results.clustering)
num.spots <- nrow(Results.clustering[[1]])
num.cell.type <- ncol(Results.clustering[[1]])

拿到方法数量、样本数、细胞类型数（实际没用）。

# 2. 初始化权重和共识矩阵

1
2

w <- rep(1 / num.methods, num.methods)  # 初始权重平均
H <- Reduce("+", Map("*", Results.clustering, w))  # 加权平均得到初始共识矩阵

📌 共识矩阵 H 是多个相似度矩阵的加权平均。

# 3. 计算初始 loss 和 lambda（如未指定）

1
2
3
4

temp2 <- sapply(Results.clustering, L2_norm, Y = H)
if(is.null(lambda)){
    lambda <- quantile(temp2, probs = prob.quantile)
}

• temp2 是每个方法和当前共识矩阵之间的 L2 距离（差异）。
• 如果没有指定 lambda ，就根据这些距离的分布自适应选一个中位数。

# 4. 主循环：迭代优化权重和共识矩阵

1
2
3

while (k <= niter) {
  ...
}

每一轮迭代包括 3 步：

# 🌀 步骤 1：更新权重

1
2

temp2 <- sapply(Results.clustering, L2_norm, Y = H)
w <- exp(-temp2/lambda)/sum(exp(-temp2/lambda))

• 越接近 H 的聚类方法，其权重 w 越大。
• 形式上是 Softmax，受 lambda 控制平滑程度。

# 🌀 步骤 2：更新共识矩阵

1

H <- Reduce("+", Map("*", Results.clustering, w))

• 用新的权重对所有相似度矩阵加权平均。

# 🌀 步骤 3：计算损失函数

1
2
3

loss_main <- sum(sapply(Results.clustering, L2_norm, Y = H) * w)
loss_entropy <- sum(w * log(w))  # 熵惩罚项
loss_all <- loss_main + lambda * loss_entropy

• 主损失 loss_main ：共识矩阵和原始方法之间的加权差异
• 正则项 loss_entropy ：鼓励权重分布均匀，避免过拟合
• 总损失是两者的加权和

# 5. 收敛判断与跳出循环

1
2
3

diff_iter <- abs(loss_all - loss_all_temp)
if(diff_iter < epsilon | k >= niter)
  break

# 6. 返回结果

1

return(list(H = H, w = w))

# ✅ 最终输出解释

1
2
3
4

list(
  H = 聚类共识相似度矩阵,
  w = 各方法权重向量
)

你可以把 H 当作一个新的样本之间的相似度矩阵，再送入例如 Louvain 聚类、谱聚类等进一步聚类分析中。

# 📌 总结一张图（逻辑结构）

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16

多种聚类结果（A, B, C...） => [相似度矩阵1, 2, 3...]

初始化：
    权重 w = 平均
    共识矩阵 H = w 加权求和

迭代：
    1. 计算每个方法与 H 的距离（L2）
    2. Softmax 更新权重 w
    3. 用 w 更新共识矩阵 H
    4. 计算损失函数：L2 + 熵
    5. 检查收敛条件（loss 差值 < ε）

输出：
    - 最终共识矩阵 H
    - 每个方法的权重 w

如果你也想看 L2_norm 这个函数我也可以帮你查或分析下～

这段 R 代码定义了一个用于空间转录组数据预处理的函数 data_process 。它使用了 Seurat 包来创建 Seurat 对象并执行标准的数据预处理流程（如标准化、变量基因筛选、PCA）。下面是对这段代码逐行的详细解析：

# 🔧 函数定义

1

data_process <- function(counts_path, loc_path, img_path, n_HVG = 2000, n_PCA = 20)

定义了一个函数 data_process ，接受以下参数：

• counts_path ：表达矩阵的路径（行是基因，列是空间点）
• loc_path ：空间位置信息的路径（spot 坐标）
• img_path ：图片路径，暂时没被实际读取
• n_HVG ：选择多少个高变异基因（默认 2000）
• n_PCA ：使用多少主成分进行 PCA（默认 20）

# 📂 读取表达矩阵

1
2

counts = read.table(counts_path, sep = '\t')
counts_ST <- as.matrix(counts)

读取表达矩阵文件，并转为矩阵形式。假设行为基因、列为空间点（spot）。

# 📍 读取坐标信息

1

coordinates <- data.frame(read.table(loc_path, sep = '\t'))

读取空间坐标信息，通常包含每个 spot 的 x , y 坐标。

# 🧬 对齐 counts 和坐标

1

colnames(counts_ST) <- rownames(coordinates)

将 counts_ST 的列名设置为空间坐标的 rownames，确保表达矩阵的列与坐标信息对齐。

# 🔬 创建 Seurat 对象

1
2
3
4

suppressWarnings({
  Seurat.ST <- Seurat::CreateSeuratObject(counts = counts_ST, assay = "Spatial",
                                          meta.data = coordinates)
})

用表达矩阵和空间坐标创建一个 Seurat 对象，使用 Spatial 作为 assay 名称，并把坐标作为元数据添加进去。

# 📊 标准化处理

1

Seurat.ST <- Seurat::SCTransform(Seurat.ST, assay = "Spatial", verbose = FALSE)

使用 SCTransform 对数据进行标准化，这是一种较新的标准化方法，比传统的 NormalizeData 更鲁棒。

# 📈 高变异基因选择

1

Seurat.ST <- Seurat::FindVariableFeatures(Seurat.ST, nfeatures = n_HVG, verbose = FALSE)

从标准化后的数据中选择 n_HVG 个高变异基因（默认 2000 个）。

# 📐 PCA 降维

1
2

Seurat.ST <- Seurat::ScaleData(Seurat.ST, verbose = FALSE)
Seurat.ST <- Seurat::RunPCA(Seurat.ST, npcs = n_PCA, verbose = FALSE)

先对数据标准化（center & scale），再进行 PCA 降维，保留 n_PCA 个主成分。

# 🖼️ 处理图片信息（未使用）

1
2

data.list <- list(counts_path = counts_path, loc_path = loc_path, img_path = img_path)
Seurat.ST@images <- data.list

将文件路径信息存储到 Seurat.ST@images 槽中。注意这里并没有真正读入图像，只是把路径作为一个列表存进去。如果后续不使用这个结构，可能没太大意义。

# 🔁 返回结果

1

return(Seurat.ST)

返回预处理后的 Seurat 对象。

# ✅ 总结

这是一个用于空间转录组数据预处理的函数，完成了以下几步：

1. 读取表达矩阵与坐标
2. 创建 Seurat 对象
3. SCTransform 标准化
4. 选择高变异基因
5. PCA 降维
6. 保存原始路径信息到 Seurat 对象中

如你后续要加上图像处理功能（例如和 histology 图像配准），可以使用 Read10X_Image() 或 ReadImage() 之类的方法来真正地载入图像。

如果你有具体的输入数据结构（比如坐标文件长什么样），我也可以帮你判断这段代码是否能跑通，或者怎么修改以适应不同的数据来源。你想继续深入看看这部分吗？

这段 R 代码定义了一个名为 clust2Mat 的函数，用于把聚类标签（cluster membership）向量转换成一个 共聚类关系矩阵。下面我们来逐步分析：

# 🔧 函数定义

1

clust2Mat <- function(memb){

• memb 是一个向量，表示每个样本的聚类标签。
  • 可以是一个具名向量（比如 names(memb) 是样本名）。
  • 每个元素的值是对应样本所属的聚类编号。

# 📏 获取样本数量

1

N <- length(memb)

• 获取一共多少个样本（即多少个聚类成员）。

# 📐 构造共聚类矩阵

1

mat_cn <- as.numeric(outer(memb, memb, FUN = "==")) - outer(1:N, 1:N, "==")

这是核心逻辑。逐步拆解如下：

# ✅ outer(memb, memb, FUN = "==")

• 构造一个 N x N 的矩阵，表示每对样本是否属于同一个聚类。
• 例如，如果 memb = c("A" = 1, "B" = 1, "C" = 2) ，则结果是：

1
2
3
4

     A     B     C
A  TRUE  TRUE  FALSE
B  TRUE  TRUE  FALSE
C  FALSE FALSE  TRUE

# ✅ as.numeric(...)

• 把逻辑值 TRUE/FALSE 转为 1/0 。

# ✅ outer(1:N, 1:N, "==")

• 构造一个对角矩阵，主对角线上是 TRUE ，其他是 FALSE ，表示 i == j ，也就是样本与自身是否相等。

# ✅ 两者相减：

1

as.numeric(outer(memb, memb, FUN = "==")) - outer(1:N, 1:N, "==")

这样可以将对角线上的 “自我比较” 从 1 减为 0 ，避免误判为共聚类（因为我们通常不关心 “自己和自己” 是否在一类中）。

# ✅ 举个例子：

如果：

1

memb <- c("A" = 1, "B" = 1, "C" = 2)

则：

1
2
3
4
5

outer(memb, memb, "==") =
     A     B     C
A  TRUE  TRUE  FALSE
B  TRUE  TRUE  FALSE
C  FALSE FALSE  TRUE