【算法详解】多切片空间分辨转录组数据部分对齐算法PASTE2，PASTE算法升级版

最新推荐文章于 2024-08-05 14:31:47 发布

yaoyao255

最新推荐文章于 2024-08-05 14:31:47 发布

阅读量621

点赞数 4

分类专栏：算法详解文章标签：机器学习人工智能算法论文阅读论文笔记

本文链接：https://blog.csdn.net/weixin_60734652/article/details/133068417

版权

算法详解专栏收录该内容

17 篇文章 0 订阅

订阅专栏

0 参考文献

[1] PASTE论文：Alignment and integration of spatial transcriptomics data

[2] PASTE2论文：PASTE2: Partial Alignment of Multi-slice Spatially Resolved Transcriptomics Data.

1 PASTE回顾

1.1 算法概述

在之前的博客中，我已经介绍过了PASTE算法。

PASTE（probabilistic alignment of Spacial Transcriptomics experiments）是一种通过同时利用数据中提供的基因表达水平信息和位置信息，将相邻组织切片通过所获得的数据进行对齐或整合的算法。PASTE提供了两种模式，即数据对齐（alignment）和数据整合（integration）。

“对齐”模式是在两个切片的点与点之间找到具有高度相似性的点（这里的“点”由技术手段决定）进行对齐，即寻找从一个切片中的点到另一个切片中的点的合理概率映射，从而能够从2D切片数据恢复出3D组织结构。

“整合”模式是将多个切片整合为单一的“中心切片”，该中心切片与独立切片在基因表达水平和空间关系上具有高度的相似性。其目的是克服由于不同的测序覆盖、组织解剖或捕获区上的组织放置而导致的单个切片的可变性。

PASTE2基于PASTE的“对齐”模式进行改进。

1.2 数据的符号定义

这部分在介绍PASTE时已经介绍过了，PASTE2与PASTE仅有细微差别，但为了方便，这里再写一次对数据的符号定义。

ST实验得到的结果数据是一对矩阵(X, Z)，其中 $X=[x_{ij}]$ 是一个n*p的矩阵，p为基因数量，n为捕获区的点数， $x_{ij}$ 代表第i个点处第j个基因的表达水平（数量），向量 $x_{i\cdot }$ 为第i个点的基因表达水平。矩阵Z是一个2*n矩阵，其中第j列的二维向量 $z_{\cdot j}$ 代表第j个点在捕获区的二维坐标，即空间位置。

由于在进行ST实验时组织切片放置方向比较随机，因此直接使用Z矩阵作为空间位置并不高效，因此使用相对位置矩阵D， $D=[d_{ij}]$ 是一个n*n的矩阵， $d_{ij}=||z_{\cdot i}-z_{\cdot j}||$ 为第i个点与第j个点的相对位置。

此外，对每一个点i，设置一个严格大于零的权重 $g_i$ 用于表示该点相对于其他点的重要程度，并且有限制 $\sum_{i} g_{i}=1$ 。如果没有先验知识指导权重的设置，那么就令每个点的权重相同。

设置代价函数c: $\mathbb{R}_+^p\times \mathbb{R}_+^p\rightarrow \mathbb{R}_+$ ，该函数用于衡量两个点处基因表达水平的相似程度，例如使用KL散度，应该越相似函数值小。

那么现在的两个切片数据就变为(X, D, g)和(X', D', g')，其中两个切片分别具有n和n'个点，点数不一定相同是因为只考虑捕获区接触到组织切片的点。

“对齐”模式要寻找两个切片之间的概率映射 $\pi =[\pi_{ij}]\in \mathbb{R}^{n\times n'}$ ，其中 $\pi_{ij}$ 代表第一个切片中的第i个点与第二个切片中第j个点对齐的概率。

1.3 PASTE“对齐”模式实现

由于PASTE2仅对PASTE的“对齐”模式进行改进，因此这里只回顾PASTE“对齐”模式的算法实现。

PASTE通过最小化下面这个函数来找到这个映射

其中F函数包含两项，第一项用于衡量两个点之间基因表达水平的相似程度，第二项用于衡量两对点（每一对分别属于一个切片）在切片内距离的相似程度，均为相似程度越高，该项的值越小。具体内容可以参考PASTE算法的介绍。

第一项被称为Wasserstein距离，第二项被称为Gromov-Wasserstein距离，函数F被称为Fused Gromov-Wasserstein (FGW) 最优传输目标。

所求的该映射要符合如下约束

这意味着概率映射矩阵 $\pi$ 每一行进行求和得到的向量应该等于g，每一列求和得到的向量转置后等于g’。这使得两个切片中的所有点都必须参与对齐。

1.4 PASTE“对齐”模式存在的问题

PASTE在进行切片对齐时假设两个切片可以在整体的2D切片内（所有点）进行对齐。即两个切片整体在生理上、技术上都具有高度的相似性。然而这通常并不是一个合理的假设，因为组织解剖和阵列（捕获区）放置在一致性上存在技术困难，并且相邻切片之间组织形态可能也难免存在差异。有时可能相邻切片之间沿z轴（垂直于切片的轴）只有部分组织重叠。如下图所示。

PASTE还假设两个切片在细胞类型上具有相似和数量和细胞比例。然而有时相邻切片会存在细胞类型的差异，有的细胞类型可能只存在于某一切片中。

在上述两种情况下，每个切片中都存在部分不应该参与对齐的点，而PASTE并没有考虑到这一点，因此会导致这些点进行了较为随机的、错误的对齐。

2 PASTE2算法

2.1 PASTE2的改进工作

针对如上所述的PASTE算法的问题，PASTE2进行了如下改进：

1）实现了一种“部分对齐”算法，使得两个切片只有一定比例的点需要进行对齐。

2）提出了一种模型选择步骤，用于估计两切片的重叠比例，用于指导进行切片对齐的点的比例系数的设置。

3）使用了组织学图像进行辅助（这在PASTE中是没有的）。

4）提出了广义普鲁克分析（generalized Procrustes analysis）使用部分对齐切片实现3D结构重建。

2.2 “部分对齐”算法

“部分对齐”实际上是一个比较简单的实现，PASTE2使用与PASTE相同的目标函数，即最小化下面函数

PASTE2的不同在于使用了如下不同的约束条件

这里给出一个超参数 $s\in [0,1]$ ，用于表示从g到g’的传输质量，上面的等式约束意味着概率映射矩阵所有元素之和为s，用于保证只有以s为比列的概率质量进行传输。同时取消了原约束条件中的等式约束，改为不等式约束。s被认为是表示两切片之间重叠占比的参数，当 $g_i=1/n$ 时，将有占比大约为s的细胞参与了对齐。

2.3 切片重叠占比的估计

我们已经知道了PASTE2是依靠超参数s来实现部分对齐的，那么如何来估计这个超参数应该设置的值呢？

作者发现，当高估了s的值时，大部分参与对齐的点将相邻；当低估了s的值时，在区域内将会有随机分布的未参与对齐的点。

因此文中定义一个edge inconsistency score，这个值越大，参与对齐的点之间约不相邻，我们希望这个分数较小。

具体而言，对于一个图 $G=(V,E)$ ，其中 $V$ 是节点的集合， $E$ 是边的集合。这个图有一组标签 $L=[l(i)]$ ，即每个节点具有一个标签 $l(i)$ ， $l(i)\in \{1,2\}$ 。定义 $E'$ 为 $E$ 的一个子集，其中包含所有两端节点标签不同的变（即一端为1，另一端为2）。定义edge inconsistency score如下

$H(G,L)=|E'|/|E|$

这代表着 $E'$ 中边的数量在总边数中的占比。

给定两个切片 $(X,D),(X',D')$ ，逐步改变s，运行PASTE2,。对于每一个s，在两个切片参与对齐的点的凸包中计算edge inconsistency score。在s大于s*时，这个分数保持较低水平，在s小于s*时，这个分数会增加，并在s=s*/2时达到最大值。因此分别计算出使得两个切片edge inconsistency score达到最大值的s分别为 $s'_1,s'_2$ ，则估计的s值为 $\widehat{s}=2min\{s'_1,s'_2\}$ 。

（可以理解为在凸包内包含的未参与对齐的点越少，那么 $E'$ 所包含的边越少，所计算出的分数越小）

2.4 概率映射矩阵的优化

PASTE2使用迭代条件梯度算法对概率映射矩阵进行优化。

首先将目标函数F写成矩阵形式如下

其中 $C=[c_{ij}]\in \mathbb{R}^{n\times n'}$ ， $c_{ij}=c(x_{i\cdot},x'_{j\cdot})$ 。

$L(D,D')=[L_{i,j,k,l}]\in\mathbb{R}^{n\times n'\times n\times n'},L_{i,j,k,l}=(D_{ik}-D'_{jl})^2$ 。

$\bigotimes$ 是张量-矩阵乘法算子， $L\bigotimes \pi$ 是一个 $n\times n'$ 的矩阵，其中第i行第j列的元素为 $\sum_{k,l}L_{i,j,k,l}\cdot \pi_{k,l}$ 。

$\left \langle \cdot ,\cdot \right \rangle$ 是Frobenius点积，即相同大小的矩阵对位元素相乘再求和。

先随机初始化概率映射 $\pi^{(0)}$ ，随后每个循环迭代执行以下三个步骤：

步骤1：

通过求解下面的线性规划问题找到梯度下降方向

其中目标函数的梯度为

这对于确定的概率映射矩阵而言是一个常量。

步骤2：

确定了下降方向 $\widetilde{\pi}^{(k)}$ 后，我们需要求一个步长 $\gamma ^{(k)}\in [0,1]$ ，使得新的概率映射矩阵处目标函数达到最小值

求解时，我们定义 $E^{(k)}=\widetilde{\pi}^{(k)}-\pi^{(k)}$ 和函数 $\Phi :[0,1]\rightarrow \mathbb{R}$

将 $\Phi (\gamma )$ 按照下面步骤进行重写

其中常数a、b、c如下

这样最小化 $\Phi (\gamma )$ 的任务就变成了在[0,1]区间上求一个一元二次函数的最小值。这样即可解出步长。

步骤3：

根据步骤1和步骤2计算出的梯度下降方向和步长，即可通过下面式子更新概率映射矩阵

2.5 组织学图像的应用

每个切片数据通常伴随有H&E图像，同时使用转录组数据和H&E图像可以提升对齐效果。

图像数据矩阵为 $H\in\mathbb{N}^{n\times 3}$ ，这里n为点数，每一行对应该点处的RGB图层数值。事实上，每个点处会有很多像素，因此每一行的向量是对应点中的所有像素RGB向量的平均值。

这样两个切片的数据就变成了 $(X,D,H),(X',D',H')$ 。

此时需要将目标函数更改为以下形式

其中 $C_{gene}$ 就是之前的 $C$ 矩阵， $C_{image}$ 是用于衡量来自两个切片中的点在H&E图像上的相似关系的矩阵，其中第i行第j列的元素为 $||h_{i\cdot}-h'_{j\cdot}||_2$ 。由于两矩阵元素值会存在数值规模的差异，因此改变 $C_{image}$ 使其最大元素与 $C_{gene}$ 最大元素相等。

2.6 3D结构重建

给出一组同组织、连续的切片数据 $(X^{(1)},Z^{(1)}),...,(X^{(t)},Z^{(t)})$ ，可以重建组织的3D结构。

对于 $k=1,2,...,t-1$ ，我们将第k+1个切片的坐标投影到第k个切片的坐标上。这种投影方式由一个旋转矩阵 $R\in \mathbb{R}^{2\times 2}$ 和一个平移向量 $t\in\mathbb{R}^2$ 来定义。

给定两个切片的坐标矩阵 $Z\in\mathbb{R}^{2\times n},W\in\mathbb{R}^{2\times n'}$ 和一个概率映射矩阵 $\pi$ （按照上文方法求得），我们希望为 $W$ 找到 $R$ 和 $t$ 使得下面式子达到最小值

求解时，我们先假设无需平移（t=0）并且R是固定的，那么我们用Q对t求偏导使其等于0如下

其中

这样我们有

然后我们使用这个t替换式子中的t得到

由于 $Zp/s$ 和 $Wq/s$ 不依赖于R，因此如果我们将 $Z_{\cdot i}$ 和 $W_{\cdot j}$ 分别替换为 $Z_{\cdot i}-Zp/s$ 和 $W_{\cdot j}-Wq/s$ ，那么我们就得到的给定一个t时Q的最小化形式，只需要求一个最优的R即可。

将Q重新写成矩阵形式如下

其中α是一个与R无关的常数。令 $U\Sigma V^T$ 为矩阵 $W\pi^TZ^T$ 的SVD分解，则有

由于V、R、U是标准正交矩阵， $\Sigma$ 是对角元素大于0的对角矩阵，因此当 $V^TRU=I$ 时Q有最小值，因此可以求得 $R=VU^T$ 。

逐步完成每一层的投影变换即可重建组织的3D结构。

yaoyao255

关注

4
点赞
踩
2

收藏

觉得还不错? 一键收藏
2
评论
【算法详解】多切片空间分辨转录组数据部分对齐算法PASTE2，PASTE算法升级版

在之前的博客中，我已经介绍过了PASTE算法。PASTE（probabilistic alignment of Spacial Transcriptomics experiments）是一种通过同时利用数据中提供的基因表达水平信息和位置信息，将相邻组织切片通过所获得的数据进行对齐或整合的算法。PASTE提供了两种模式，即数据对齐（alignment）和数据整合（integration）。
复制链接

扫一扫

专栏目录