使用仿射间隙惩罚函数进行全局对齐

Question

我必须编写一个程序，使用仿射空位惩罚函数在两个序列之间进行全局对齐。动态算法（改进的 Needleman Wunsch）计算两个给定序列 s 和 t 的相似度（表示序列相似程度的最大分数）。它通过构建三个二维数组来考虑间隙，即序列中的连续空间块，这些空间比孤立空间更容易出现。数组可以不那么正式地描述为：

数组 C：保持以序列 s 字符结尾的块的最大分数与序列 t 的字符对齐；
数组 BS：保持最大分数对于以序列 t 字符结尾的块，与序列 s

中的空格对齐

数组 BT：保持以序列 s 字符结尾的块的最高分数，并与序列 t 中的空格对齐；
该算法具有以下递推关系：

C[i,j] = v(s[i],t[j]) + max{C[i-1][j-1], BS[i-1][j-1], BT [i-1][j-1]}
BS[i,j] = max{C[i][j-1]-(h+g), BS[i][j-1]-g, BT[i][j-1]-(h+ g)}
BT[i,j] = max{C[i-1][j]-(h+g), BS[i-1][j]-(h+g), BT[i-1][j] -g}

** v(s[i],t[i]) = 匹配值（当两个字符相同时）或不匹配（当字符不相同时） 相似度是每个数组的最后一个值中的最大值。 问题是当我运行程序时，它有一个奇怪的行为：
对于给定的序列对，如果我更改哪个是 t 或 s，我的程序会为同一对序列给出不同的值。那么，你能帮我找出为什么程序有这样的行为吗？你知道我做错了什么吗？关于代码，这里是：

    int main (void){

  int mat, mis, h, g,
      sim, i, j, m, n; 
  /* mat = match, mis = mismatch, h = open gap penalty, g = extend gap penalty */

  string s, t;
  s = malloc(1500);
  t = malloc(1500);

    scanf("%d %d %d %d", &mat, &mis, &h, &g);
    scanf("%s", s);
    scanf("%s", t);


    m = strlen(s);
    n = strlen(t);
    int C[m][n], BS[m][n], BT[m][n];  

    C[0][0] = 0;
    for(j = 1; j<= n; j++)
      C[0][j] = -32000;

    for(i = 1; i<= m; i++)
      C[i][0] = -32000;

    for(j = 1; j <= n; j++)
       BS[0][j] = -(h + g*j);


    for(i = 0; i <= m; i++)
      BS[i][0] = -32000;

      for(j = 0; j <= n; j++)
        BT[0][j] = -32000;

      for(i = 1; i <= m; i++)
        BT[i][0] = -(h + g*i);


for(i = 1; i <= m; i++){
    for(j = 1; j <= n; j++){
        C[i][j] = align(s[i-1],t[j-1],mat,mis) + max(C[i-1][j-1],BS[i-1][j-1],BT[i-1][j-1]);
        BS[i][j] = max((C[i][j-1]-(h+g)),(BS[i][j-1]-g),(BT[i][j-1])-(h+g));
        BT[i][j] = max((C[i-1][j]-(h+g)),(BS[i-1][j]-(h+g)),(BT[i-1][j]-g));
        }
    }
    printf("\n");
    printf("c[m][n]: %d bs[m][n]:%d bt[m][n]: %d\n", C[m][n], BS[m][n], BT[m][n]);
    sim = max(C[m][n], BS[m][n], BT[m][n]);
    printf("sim: %d\n", sim);
return 0;
}

Answer 1

好的，在尝试了许多 printfs 之后，我终于找到了问题，因为我不知道如何使用调试器。
我得到的第一个线索是 gcc 在我尝试阅读时告诉我的分段错误来自文件的序列（长度相对较大）。确实，分段错误可能有多种原因，但几乎所有我在程序中看到这个错误的时候，它的发生都是因为我试图访问一个实际上不存在于数组中的位置。
然后许多 printfs 在初始化步骤中显示了一些值，这些值与每个数组的第一行和第一列应该具有的值不同。我将分段错误与初始化步骤中的奇怪值相关联，并决定检查数组的声明步骤以及所有循环条件。就在那里！我只是忘记了任何数组的一个非常基本的特性：如果数组的大小为 n，则可以从零访问到 n-1。
为了解决我的新手错误，我在每个数组中为行和列添加了一个位置（因为位置 (0,0) 不与任何一对对齐字符相关联），因此每个数组都有大小 [ m][n]，其中m是序列s的（长度+1），n是（长度+1) 序列 t。此外，我更改了所有循环条件以访问直到 [m-1][n-1] 位置。现在程序运行良好。