東京大学 新領域創成科学研究科 メディカル情報生命専攻 2024年1月実施 問題12

标签：

• Tokyo-University
• Sequence-Alignment

Author

zephyr, 祭音Myyura

Description

ナノポアシークエンサーは DNA 配列を電流値の変化を通して読み取ることができる。 ある DNA 配列を読み取ったところ、時間 $t = 1, 2, \ldots , n$ でそれぞれ読み取られた電流値列は $\{A_t\} \ (0 \leq A_t \leq 1)$ であった。 この DNA 配列には基準配列が含まれていることが分かっており、基準配列を単独で読み取ったときの電流値列は時間 $t = 1, 2, \ldots, m$ に対して $\{B_t\} \ (0 \leq B_t \leq 1)$ であった（基準波形）。 このとき、$A_t$ の中で基準配列に相当する時間区間を見つけたい。

DNA 配列を読み取る速度にはばらつきがあり急な変動があるため、基準配列を探す場合には時間軸方向へのずれを許したい。 このため、$1 \leq p_1 \leq p_2 \leq \cdots \leq p_m \leq n$ (ただし $1 \leq t < m$ に対して $p_{t+1} - p_t \leq W$, $W$ は正の整数) が存在して $B_t$ と $A_{p_t}$ が対応していると仮定する。 読み取りノイズ等が無ければ $B_t = A_{p_t}$ が成立して欲しいが、実際には読み取りノイズや電流値読み取りタイミングのずれを考慮して、相違度 $\sum_{t=1}^m (B_t - A_{p_t})^2$ が最小値をとるような対応関係 $p_1, p_2, \ldots, p_m$ が基準波形とナノポアシークエンサー出力の一部分との正しい対応関係を示していると仮定する。

(1) $m = 1$ とする。相違度が最小となる $p_1$ を１つ求めるアルゴリズムを示せ。

(2) $m = 2$ とする。相違度が最小となる $p_1$, $p_2$ の組を１つ求めるアルゴリズムを示せ。

(3) $W$ はどのような物理条件に対応していると考えられるか。１行で説明せよ。

(4) $D_{i,j}$ は $B_i$ と $A_j$ が対応しているときの、基準配列 $\{B_1, B_2, \ldots, B_i\}$ に対する $\{A_1, A_2, \ldots, A_j\}$ の相違度の最小値である。$D_{x,y} \ (x < i, y \leq j)$ を用いて $D_{i,j}$ を表せ。

(5) 基準配列 $\{B_t\}$ に対する $\{A_t\}$ の相違度の最小値を計算するアルゴリズムを示し、その時間計算量を $n$, $m$, $W$ を用いて示せ。

题目描述

纳米孔测序仪以电流变化读取 DNA。某待测序列在时刻 $t=1,\ldots,n$ 的电流为

$$\{A_t\},\qquad0\le A_t\le1,$$

已知其中包含一段参考序列；单独读取参考序列时，在 $t=1,\ldots,m$ 得到参考波形

$$\{B_t\},\qquad0\le B_t\le1.$$

由于读取速度波动，需要允许时间轴错位。假设存在对应位置

$$1\le p_1\le p_2\le\cdots\le p_m\le n,$$

且对相邻参考时刻 $1\le t<m$

$$p_{t+1}-p_t\le W,$$

其中 $W$ 为正整数。把参考波形与测序输出部分的差异度定义为

$$\sum_{t=1}^{m}(B_t-A_{p_t})^2,$$

并以使其最小的 $p_1,\ldots,p_m$ 作为正确对应。结合上图回答：

1. 当 $m=1$ 时，给出求一组最优 $p_1$ 的算法。
2. 当 $m=2$ 时，给出求一组最优 $(p_1,p_2)$ 的算法。
3. 用一行说明参数 $W$ 对应何种物理限制。
4. 定义 $D_{i,j}$：在 $B_i$ 与 $A_j$ 对应的条件下，参考前缀 $\{B_1,\ldots,B_i\}$ 与输出前缀 $\{A_1,\ldots,A_j\}$ 的最小差异度。用满足 $x<i,y\le j$ 的先前状态 $D_{x,y}$ 表示 $D_{i,j}$，并体现单调对应与相邻位置差不超过 $W$ 的约束。
5. 设计计算完整参考波形 $\{B_t\}$ 对 $\{A_t\}$ 最小差异度的算法，并用 $n,m,W$ 表示时间复杂度。

Kai

1. Finding $p_1$ for $m = 1$

When $m = 1$, the problem reduces to finding a single time point $p_1$ such that the dissimilarity $\sum_{t=1}^1 (B_t - A_{p_t})^2$ is minimized. The algorithm can be described as follows:

1. Initialize $p_1$ to 1.
2. For each $i = 1, 2, \dots, n$, compute the dissimilarity $d_i = (B_1 - A_i)^2$.
3. Find the index $i$ that minimizes $d_i$.

The resulting $p_1$ will be the index that minimizes the dissimilarity.

Algorithm

Input: A sequence {A_t}, reference value B_1
Output: Index p_1 that minimizes dissimilarity

1. Initialize min_dissimilarity = infinity
2. Initialize p_1 = 1
3. For i = 1 to n do
     d_i = (B_1 - A_i)^2
     If d_i < min_dissimilarity then
         min_dissimilarity = d_i
         p_1 = i
4. Return p_1

2. Finding $p_1, p_2$ for $m = 2$

When $m = 2$, the problem involves finding two indices $p_1$ and $p_2$ such that $1 \leq p_1 \leq p_2 \leq n$ and $p_{t+1} - p_t \leq W$. The goal is to minimize the dissimilarity $\sum_{t=1}^2 (B_t - A_{p_t})^2$. The algorithm can be described as follows:

1. Initialize minimum dissimilarity to infinity and $p_1, p_2$ to 1.
2. For each $i = 1, 2, \dots, n$, do:
   • For each $j = i, i+1, \dots, \min(i + W, n)$, do:
     • Compute the dissimilarity $d_{i,j} = (B_1 - A_i)^2 + (B_2 - A_j)^2$.
     • If $d_{i,j} < \text{min\_dissimilarity}$, update $p_1 = i$, $p_2 = j$, and $\text{min\_dissimilarity} = d_{i,j}$.

Algorithm

Input: A sequence {A_t}, reference values {B_1, B_2}, window W
Output: Indices p_1, p_2 that minimize dissimilarity

1. Initialize min_dissimilarity = infinity
2. Initialize p_1, p_2 = 1
3. For i = 1 to n do
     For j = i to min(i + W, n) do
         d_ij = (B_1 - A_i)^2 + (B_2 - A_j)^2
         If d_ij < min_dissimilarity then
             min_dissimilarity = d_ij
             p_1 = i
             p_2 = j
4. Return p_1, p_2

3. Physical Meaning of $W$

The parameter $W$ bounds local time stretching: at most $W$ output sampling intervals may elapse between two adjacent reference samples.

4. Recursive Formula for $D_{i,j}$

The minimum dissimilarity $D_{i,j}$ can be computed using the following recurrence relation:

$$\begin{aligned} D_{i,j} &=(B_i-A_j)^2+ \min_{\max(1,j-W)\le k\le j}D_{i-1,k} \qquad(i\ge2). \end{aligned}$$

5. Algorithm for Minimum Dissimilarity for $m > 2$

To calculate the minimum dissimilarity between the sequence $\{A_t\}$ and the reference signal $\{B_t\}$, we use a dynamic programming approach.

1. Initialize a 2D array $D$ of size $m \times n$ with infinity, where $D_{i,j}$ represents the minimum dissimilarity for the first $i$ elements of the reference sequence mapped to the first $j$ elements of the output sequence.
2. Set $D_{1,j} = (B_1 - A_j)^2$ for all $j = 1, 2, \dots, n$.
3. For each $i = 2, 3, \dots, m$, do:
   • For each $j = 1, 2, \dots, n$, do:
     • Set $D_{i,j} = (B_i-A_j)^2+\min_{x = \max(1, j - W)}^{j}D_{i-1,x}$.
4. Return $\min_{1\le j\le n}D_{m,j}$.

Time Complexity

The time complexity of this algorithm is $O(m \cdot n \cdot W)$, where:

• $m$ is the length of the reference signal,
• $n$ is the length of the sequence $\{A_t\}$,
• $W$ is the maximum allowed shift between adjacent points.

Knowledge

难点思路

本题的主要难点在于第 4 和第 5 问，需要正确地推导出动态规划的递推关系，并设计出高效的算法。关键是要理解如何处理时间轴偏移的约束条件，以及如何在保证正确性的同时优化算法的时间复杂度。

解题技巧

对于序列对齐问题，动态规划是一种有效的解决方案。特别是在处理时间序列中的波动和对齐问题时，可以灵活调整对齐窗口 $W$ 来达到最优匹配。

重点词汇

• Dissimilarity: 相异度
• Reference Signal: 参考信号
• Window: 窗口

参考资料

1. Dynamic Time Warping Algorithm: Introduction to the dynamic time warping algorithm, Wikipedia.
2. Time Series Analysis and Its Applications, Third Edition, Robert H. Shumway and David S. Stoffer.