原创: hxj7
前言: 本文介绍的是采用仿射罚分模型的Smith-Waterman算法。
关于全局联配,局部联配以及仿射罚分模型的介绍可参见前文:
序列比对(一)全局比对Needleman-Wunsch算法
序列比对(二)Needleman-Wunsch算法之仿射罚分
序列比对(三)局部联配Smith-Waterman算法
其中,《序列比对(二)Needleman-Wunsch算法之仿射罚分》一文中对**X(i, j)以及Y(i, j)**的计算需要补充说明一下(以X(i, j)为例,Y(i, j)与之类似):
在上文中X(i, j)的计算方式如下:
图片引自https://fengcong97.cn/sequece-alignment/
其前提是假设一个缺失和一个插入不会连在一起。否则,X(i, j)的计算就会变成下面这样:
局部联配线性罚分模型与仿射罚分模型结果的区别
(alignSW代表局部联配线性罚分模型,alignSW2代表局部联配仿射罚分模型):
C代码
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#define MAXSEQ 1000
#define GAP_CHAR '-'
// 对空位的罚分是仿射的
struct Unit {
int W1; // 是否往上回溯一格
int W2; // 是否往左上回溯一格
int W3; // 是否往左回溯一格
float X;
float Y;
float M;
float O; // 得分矩阵第(i, j)这个单元的分值,即序列s(1,...,i)与序列r(1,...,j)比对的最高得分
};
typedef struct Unit *pUnit;
void strUpper(char *s);
float max2(float a, float b);
float max3(float a, float b, float c);
float getFScore(char a, char b);
void printAlign(pUnit** a, const int i, const int j, char* s, char* r, char* saln, char* raln, int n);
void align(char *s, char *r);
int main() {
char s[MAXSEQ];
char r[MAXSEQ];
printf("The 1st seq: ");
scanf("%s", s);
printf("The 2nd seq: ");
scanf("%s", r);
align(s, r);
return 0;
}
void strUpper(char *s) {
while (*s != '\0') {
if (*s >= 'a' && *s <= 'z') {
*s -= 32;
}
s++;
}
}
float max2(float a, float b) {
return a > b ? a : b;
}
float max3(float a, float b, float c) {
float f = a > b ? a : b;
return f > c ? f : c;
}
// 替换矩阵:match分值为5,mismatch分值为-4
// 数组下标是两个字符的ascii码减去65之后的和
float FMatrix[] = {
5, 0