使用动态规划算法进行 DNA 序列比对
4.93/5 (16投票s)
使用动态规划算法比较两个序列
引言
序列比对是通过寻找在所有序列中以相同顺序排列的一系列字符来比较两个（成对比对）或更多（多序列比对）序列的过程。可以通过将两个序列写在一页上两行来对齐它们。相同或相似的字符放在同一列中，而不同的字符可以作为不匹配放在同一列中，或者与另一个序列中的间隙（-）对齐。以这种方式对齐的序列被称为相似。序列比对对于发现生物序列中的功能、结构和进化信息很有用。考虑以下DNA序列GACGGATTAG和GATCGGAATAG。请注意，当我们把它们一个放在另一个上面时
GA-CGGATTAG
GATCGGAATAG
唯一的区别在上述序列中用颜色标记。观察到间隙（-）被引入第一个序列中，以使相等的碱基完美对齐。本文的目的是介绍一种有效的算法，该算法接收两个序列并确定它们之间的最佳比对。比对的总得分取决于比对的每一列。如果该列有两个相同的字符，它将收到+1的值（匹配）。不同的字符将为该列提供-1的值（不匹配）。最后，一列中的间隙会将其值降低到-2（间隙惩罚）。最佳比对将是总得分最高的那个。上述比对将给出总分：9 × 1 + 1 × (-1) + 1 × (-2) = 6。
这些参数的匹配、不匹配和间隙惩罚可以根据序列的选择或实验结果调整为不同的值。
计算两个序列之间相似性的一种方法是生成所有可能的比对并选择最佳比对。但是，两个序列之间的比对数量是指数级的，这将导致算法变慢，因此，动态规划被用作一种产生更快比对算法的技术。动态规划试图通过使用已经计算出的相同问题较小实例的解来解决问题的实例。给定两个序列Seq1和Seq2，动态规划不是确定序列之间的整体相似性，而是通过确定两个序列的任意前缀之间的所有相似性来构建解决方案。该算法从较短的前缀开始，并使用先前计算的结果来解决较大前缀的问题。
设 M = Seq1 的大小，N= Seq2 的大小，计算被安排到一个 (N+1) × (M+1) 的数组中，其中条目 (j,i) 包含 Seq2[1.....j] 和 Seq1[1.....i] 之间的相似性。该算法通过查看前三个条目来计算条目 (j,i) 的值
(j-1,i-1) Diagonal Cell to entry (j,i)
(j-1,i) Above Cell to entry (j,i)
(j,i-1) Left Cell to entry (j,i)
j-1,i
j,i
j-1,i-1
j,i-1
条目 (j,i) 的值可以通过以下公式计算
 |
公式 (1.1) |
其中 p(j,i)= +1 如果 Seq2[j]=Seq1[i] (匹配得分)，p(j,i)= -1 如果 Seq2[j]!=Seq1[i]。
比对得分的最大值位于单元格 (N-1,M-1) 中，算法将从该单元格回溯到第一个条目单元格 (1,1) 以产生结果比对。如果单元格 (j,i) 的值是使用对角单元格的值计算的，则比对将包含 Seq2[j] 和 Seq1[i]。如果该值是使用上面的单元格计算的，则比对将包含 Seq2[j] 和 Seq1[i] 中的一个间隙 ('-')。如果该值是使用左侧单元格计算的，则比对将包含 Seq1[i] 和 Seq2[j] 中的一个间隙 ('-')。结果比对将通过遍历单元格 (N-1,M-1) 回溯到单元格 (1,1) 的初始条目而完全产生。
Using the Code
我的代码有两个类，第一个名为DynamicProgramming.cs，第二个名为Cell.cs。我将在接下来的几行中讨论DynamicProgramming.cs类的细节，因为它描述了我的文章的主要思想。
第一个类包含三个方法，描述了动态规划算法的步骤。第一个方法名为Intialization_Step，此方法准备矩阵 a[i,j]，该矩阵保存两个序列的任意前缀之间的相似性。该算法从较短的前缀开始，并使用先前计算的结果来解决较大前缀的问题。
public static Cell[,] Intialization_Step
(string Seq1, string Seq2,int Sim,int NonSimilar,int Gap)
{
int M = Seq1.Length;//Length+1//-AAA
int N = Seq2.Length;//Length+1//-AAA
Cell[,] Matrix = new Cell[N, M];
//Intialize the first Row With Gap Penalty Equal To i*Gap
for (int i = 0; i < Matrix.GetLength(1); i++)
{
Matrix[0, i] = new Cell(0, i, i*Gap);
}
//Intialize the first Column With Gap Penalty Equal To i*Gap
for (int i = 0; i < Matrix.GetLength(0); i++)
{
Matrix[i, 0] = new Cell(i, 0, i*Gap);
}
// Fill Matrix with each cell has a value result from method Get_Max
for (int j = 1; j < Matrix.GetLength(0); j++)
{
for (int i = 1; i < Matrix.GetLength(1); i++)
{
Matrix[j, i] = Get_Max(i, j, Seq1, Seq2, Matrix,Sim,NonSimilar,Gap);
}
}
return Matrix;
}
第二个方法名为Get_Max，通过公式 1.1 计算单元格 (j,i) 的值。
public static Cell Get_Max(int i, int j, string Seq1,
string Seq2, Cell[,] Matrix,int Similar,int NonSimilar,int GapPenality)
{
Cell Temp = new Cell();
int Sim;
int Gap = GapPenality;
if (Seq1[i] == Seq2[j])
Sim = Similar;
else
Sim = NonSimilar;
int M1, M2, M3;
M1 = Matrix[j - 1, i - 1].CellScore + Sim;
M2 = Matrix[j, i - 1].CellScore + Gap;
M3 = Matrix[j - 1, i].CellScore + Gap;
int max = M1 >= M2 ? M1 : M2;
int Mmax = M3 >= max ? M3 : max;
if (Mmax == M1)
{ Temp = new Cell(j, i, M1, Matrix[j - 1, i - 1],
Cell.PrevcellType.Diagonal); }
else
{
if (Mmax == M2)
{ Temp = new Cell(j, i, M2, Matrix[j, i - 1], Cell.PrevcellType.Left); }
else
{
if (Mmax == M3)
{ Temp = new Cell(j, i, M3, Matrix[j - 1, i],
Cell.PrevcellType.Above); }
}
}
return Temp;
}
第三个方法名为Traceback_Step。此方法将通过遍历单元格矩阵 (N-1,M-1) 回溯到单元格矩阵 (1,1) 的初始条目来产生比对。
public static void Traceback_Step(Cell[,] Matrix,
string Sq1, string Sq2, List<char> Seq1, List<char> Seq2)
{
//List<char> Seq1 = new List<char>();
//List<char> Seq2 = new List<char>();
Cell CurrentCell = Matrix[Sq2.Length - 1, Sq1.Length - 1];
while (CurrentCell.CellPointer != null)
{
if (CurrentCell.Type == Cell.PrevcellType.Diagonal)
{
Seq1.Add(Sq1[CurrentCell.CellColumn]);
Seq2.Add(Sq2[CurrentCell.CellRow]);
}
if (CurrentCell.Type == Cell.PrevcellType.Left)
{
Seq1.Add(Sq1[CurrentCell.CellColumn]);
Seq2.Add('-');
}
if (CurrentCell.Type == Cell.PrevcellType.Above)
{
Seq1.Add('-');
Seq2.Add(Sq2[CurrentCell.CellRow]);
}
CurrentCell = CurrentCell.CellPointer;
}
}
我的代码中的第二个类名为Cell.cs。此类操作矩阵的单元格。每个单元格都有
- 一个由行索引和列索引指示的位置
- 一个由比对得分表示的值
- 一个指向用于计算当前单元格得分的前一个单元格的指针 [注意：指针值“对角线、上方和左侧”]
参考
- “计算分子生物学导论” 作者：JOÃO SETUBAL 和 JOÃO MEIDANIS