LinearPartition

LinearPartition: linear-time approximation of RNA folding partition function and base-pairing probabilities

Year: 2020
Authors: He Zhang, Liang Zhang, David H. Mathews and Liang Huang
Journal Name: Bioinformatics

Motivation

传统分割方法的复杂度与序列长度呈三次方关系。

Research Objective

设计一种线性时间的算法来近似分割函数和碱基对概率

Dataset

定义 $\text{x}=x_1,...,x_n$ 为长度为 n 的 RNA 序列， $\mathcal{Y}(\text{x})$ 代表 $\mathbf{x}$ 所有可能的二级结构。

Method

分割函数为
$$Q(\text{x})=\sum _{\text{y}\in \mathcal{Y}(\text{x})}{e}^{-\frac{\Delta G^{\circ }(y)}{RT}}$$
其中， $\Delta G^{\circ }(y)$ 为结构 $y$ 的 Gibb 自由能变化， $R$ 是气体常量， $T$ 是热动力学温度。
$\Delta G^{\circ }(y)$ 使用逆 Nussinov–Jacobson 能量模型计算。在位置 $j$ 的未配对碱基的自由能变化表示为 $\delta (\text{x},j)$ ，碱基对 $(i,j)$ 的自由能变化表示为 $\xi (\text{x},i,j)$ 。
因此， $\Delta G^{\circ }(y)$ 表示为
$$\Delta G^{\circ }(y)=\sum _{j\in unpaired(y)}\delta (\text{x},j)+\sum _{(i,j)\in pairs(y)}\xi (\text{x},i,j)$$
其中， $unpaired(y)$ 表示 $y$ 中未配对的碱基集合， $pairs(y)$ 表示 $y$ 中的碱基对集合。
分割函数可表示为
$$Q(\text{x})=\sum _{\text{y}\in \mathcal{Y}(\text{x})}(\prod _{j\in unpaired(y)}{e}^{-\frac{\delta (\text{x},j)}{RT}}\prod _{(i,j)\in pairs(y)}{e}^{-\frac{\xi (\text{x},i,j)}{RT}})$$
作者定义 $[i,j]$ 为子序列 $x_i,...,x_j$ （ $[j,j-1]$ 为空）。
$[i,j]$ 对应的分割函数 $Q_{i,j}={\sum }_{y\in \mathcal{Y}(x_i,...,x_j)}{e}^{-\frac{\Delta G^{\circ }(y)}{RT}}$
线性切割算法的伪代码如下所示

1. function LinearPartition($\text{x}$, b)
2. n = len( $\text{x}$ )
3. 令所有的 $Q_{j,j-1}=1$
4. for j in (1, …, n) do
5.   for each $[i,j-1]$ in $Q$ do
6.     $Q_{i,j}+=Q_{i,j-1}\cdot e^{-\frac{\delta (\text{x},j)}{RT}}$
7.     if $x_{i-1}{x}_j$ in {AU, UA, CG, GC, GU, UG} then
8.       for each $[k,i-2]$ in $Q$ do
9.         $Q_{k,j}+=Q_{k,i-2}\cdot Q_{i,j-1}\cdot e^{-\frac{\xi (\text{x},i-1,j)}{RT}}$
10.     BeamPrune($Q$, j, b)
11. return $Q$

在 第 5-6 行，通过添加碱基 $x_j$ 将 $[i,j-1]$ 扩展为 $[i,j]$ 。在 8-9 行，如果 $x_{i-1}$ 和 $x_j$ 为碱基对，就结合 $[i,j-1]$ 和左侧的 $[k,i-2]$ 并更新 $[k,j]$ ，如下图所示。

1-9 行算法的复杂度为 $O(n^3)$ 。于是，作者在每次循环 $j$ 时加入了 BeamPrune($Q$, j, b) ，伪代码如下所示

1. function BeamPrune($Q$, j, b)
2. for [i, j] in Q do
3.   candidates[i] = $Q_{1,i-1}\cdot Q_{i,j}$
4. candidates = SelectTopB(candidates, b) # 选出前 b 个最大的值
5. for $[i,j]$ in $Q$ do
6.   if $i$ not in candidates then
7.     delete $Q_{i,j}$

这样就使得复杂度降低为 $O(n{b}^2)$ ， $b$ 通常取 100 。
核苷酸 $i$ 与 核苷酸 $j$ 配对的概率为
$$\begin{gathered} p_{i,j}=\sum _{y\in \mathcal{Y}(\text{x}),(i,j)\in pairs(y)}p(y) \\ p(y)=\frac{e^{-\frac{\Delta G^{\circ }(y)}{RT}}}{Q(\text{x})} \end{gathered}$$

Limitations

1. 本方法无法保证近似质量，因为 b 为固定值。但如果以阈值作为筛选条件的话又不能保证运行时间。
2. LinearPartition 在随机序列，特别是存在很多冲突的配对选择时表现不好。