从短序列读取中推断异构体

### 从短序列读取中推断异构体 在基因研究中，从短序列读取推断异构体是一项重要任务。本文将详细介绍相关的理论基础、数据处理策略以及具体的推断算法。 #### 1. 二次规划公式化 - **基因与表达片段**：设 $G$ 为所有基因的集合，每个基因 $g$ 定义了一组表达片段 $S_g = \{s_1, s_2, \ldots, s_{|S_g|}\}$，这些表达片段根据它们在参考基因组中的位置排序。基因上的连接点是所有表达片段对 $(s_i, s_j)$，其中 $1 \leq i < j \leq |S_g|$，片段 $s_i$ 的长度为 $l_i$。 - **异构体表达水平**：基因 $g$ 的所有已知异构体集合记为 $F_g$，每个异构体 $f \in F_g$ 由表达片段的子集组成，异构体 $f$ 的表达水平（即每碱基的读取数）记为 $x_f$。所有基因的所有转录本长度之和，按其表达水平加权，为 $L_0 = C \cdot \sum_{g \in G} \sum_{s \in f, f \in F_g} l_s x_f$，其中常数 $C$ 定义了表达水平与对应异构体的转录本数量之间的线性关系，$C$ 可从数据中推断得出。 - **单端读取的分布**：设 $M$ 为映射到参考基因组的单端读取总数，$d_i$ 为落入表达片段 $s_i$ 的读取数。在均匀采样假设下，$d_i$ 是随机变量 $r_i$ 的观测值，$r_i$ 遵循二项分布 $B(M, p_i)$，其中 $p_i = C y_i l_i / L_0$，$y_i = \sum_{s_i \in f} x_f$。由于 $M$ 通常很大，$p_i$ 很小且 $M p_i$ 足够大，二项分布可近似为正态分布 $N(\mu_i, \sigma_i^2)$，其中 $\mu_i = M p_i$，$\sigma_i^2 = M p_i (1 - p_i) \approx M p_i = \mu_i$。 - **二次规划问题**：为了将读取映射到连接点，每个连接点 $(s_i, s_j)$ 被视为长度为 $2L_1 - 2$ 的片段，由 $s_i$ 的最后 $L_1 - 1$ 个碱基和 $s_j$ 的前 $L_1 - 1$ 个碱基组成。表达片段、异构体和单端读取之间的关系可以通过以下二次规划问题（QP）来描述： - 目标函数：$\min z = \sum_{s_i \in S \cup J} (\frac{\epsilon_i}{\sigma_i})^2$ - 约束条件： - $\sum_{s_i \in f} x_f l_i + \epsilon_i = d_i$，$s_i \in S \cup J$ - $x_f \geq 0$，$f \in F$ 其中 $\sigma_i$ 是正态分布 $N(\mu_i, \sigma_i^2)$ 的标准差，将从 $d_i$ 经验估计。该 QP 可通过简单变换证明为凸 QP，并可使用公共程序 QuadProg++ 在多项式时间内求解。由于 $\sigma_i$ 未知，用 $\sqrt{d_i}$ 近似代替 $\sigma_i$。 #### 2. 双端读取 双端读取由一对被间隙分隔的短（单端）读取组成。如果双端读取的跨度是一个遵循概率分布 $h(x)$ 的随机变量，有三种可能的生成双端读取的策略： - **策略 (a)**：双端读取的起始位置从所有表达异构体中均匀随机采样，然后根据分布 $h(x)$ 生成跨度。如果双端读取的结束位置超出异构体范围，则将其截断，使其结束位置位于异构体的末尾。 - **策略 (b)**：双端读取的中心位置从所有表达异构体中均匀随机采样，然后根据分布 $h(x)$ 生成跨度。如果起始（或结束）位置超出异构体范围，则将其截断，使其起始（或结束）位置位于异构体的起始（或结束）位置。 - **策略 (c)**：双端读取的结束位置从所有表达异构体中均匀随机采样，然后根据分布 $h(x)$ 生成跨度。如果双端读取的起始位置超出异构体范围，则将其截断，使其起始位置位于异构体的起始位置。 当应用这些策略生成一定数量的双端读取时，定理 1 给出了在第一个区间没有起始位置且在第三个区间没有结束位置的读取的概率的非平凡上界。 ```mermaid graph LR A[选择策略] -->|策略 (a)| B[随机采样起始位置] A -->|策略 (b)| C[随机采样中心位置] A -->|策略 (c)| D[随机采样结束位置] B --> E[生成跨度] C --> E D --> E E --> F{结束位置是否超出范围} F -->|是| G[截断读取] F -->|否| H[完成读取生成] G --> H ``` #### 3. 有效异构体 - **二进制向量表示**：对于具有表达片段 $S = \{s_1, s_2, \ldots, s_{|S|}\}$ 的基因，其异构体 $f$ 可以表示为长度为 $|S|$ 的二进制向量，当且仅当表达片段 $s_i$ 包含在 $f$ 中时，$f[i] = 1$。类似地，映射到表达片段子集 $S' \subseteq S$ 的单端或双端短读取可以表示为二进制向量 $r$，当且仅当 $s_i \in S'$ 时，$r[i] = 1$。 - **有效异构体的条件**： - **连接点信息**：连接点 $(s_i, s_