LCS最大公共子序列问题

来源:互联网 发布:如何备份网站数据库 编辑:程序博客网 时间:2024/05/16 19:38

在生物应用中,经常需要比较两个(或多个)不同生物体的DNA,

例如:某种生物的DNA可能为S1=ACCGGTCGAGTGCGCGGAAGCCGGCCGAA,

另一种生物的DNA可能为S2=GTCGTTCGGAATGCCGTTGCTCTGTAAA

我们比较两个DNA串的一个原因是希望确定它们的相似度,作为度量两种生物的近似程度指标


寻找第三个串S3,它所有碱基也都出现在S1和S2中,且三个串中的顺序都相同,但在S1和S2中不要求连续出现。

可以找到的S3越长,就可以认为S1和S2的相似度越高。在这个例子中最长的S3为GTCGTCGGAAGCCGGCCGAA


我们定义C[i, j]表示Xi和Yj的LCS长度。如果i = 0或j = 0,即一个序列长度为0,那么LCS的长度为0

根据LCS问题的最优子结构性质,可得如下公式:

C[i, j] = 0,若i = 0 或 j = 0

C[i, j] = C[i - 1, j - 1] + 1,若i,j > 0 且 Xi = Yj

C[i, j] = max(C[i, j - 1], C[i - 1, j]) ,若i, j > 0且Xi != Yj


代码如下:

package 动态规划;/** * Lcs即最长公共子序列问题(longest common subsequence problem) * @author wangdong20 * */public class Lcs {public static final int empty = 0;public static final int upLeft = 1;public static final int up = 2;public static final int left = 3;public static int[][][] lcsLength(String x, String y){int m = x.length();int n = y.length();int[][][] result = new int[2][m + 1][n + 1];  // result[0]表示子序列长度 result[1]表示LCS矩阵方向for(int i = 0; i < m + 1; i++){result[0][i][0] = 0;result[1][i][0] = empty;}for(int j = 0; j < n + 1; j++){result[0][0][j] = 0;result[1][0][j] = empty;}for(int i = 1; i <= m; i++){for(int j = 1; j <= n; j++){if(x.charAt(i - 1) == y.charAt(j - 1)){result[0][i][j] = result[0][i - 1][j - 1] + 1;result[1][i][j] = upLeft;}else if(result[0][i - 1][j] >= result[0][i][j - 1]){result[0][i][j] = result[0][i - 1][j];result[1][i][j] = up;}else{result[0][i][j] = result[0][i][j - 1];result[1][i][j] = left;}}}return result;}public static void printLcs(int[][][] b, String x, int i, int j){if(i == 0 || j == 0)return;if(b[1][i][j] == upLeft){printLcs(b, x, i - 1, j - 1);System.out.print(x.charAt(i - 1));}else if(b[1][i][j] == up){printLcs(b, x, i - 1, j);}else{printLcs(b, x, i, j - 1);}}/** * @param args */public static void main(String[] args) {// TODO 自动生成方法存根String s1 = "ACCGGTCGAGTGCGCGGAAGCCGGCCGAA";String s2 = "GTCGTTCGGAATGCCGTTGCTCTGTAAA";String s3 = "amputation";String s4 = "spanking";System.out.println("s1: " + s1);System.out.println("s2: " + s2);System.out.println("最长公共子序列: ");int result[][][] = lcsLength(s1, s2);printLcs(result, s1, s1.length(), s2.length());System.out.println("\ns3: " + s3);System.out.println("s4: " + s4);System.out.println("最长公共子序列: ");int result2[][][] = lcsLength(s3, s4);printLcs(result2, s3, s3.length(), s4.length());}}

实质上lcsLength(s3, s4)返回的是两个二维数组组成的三维数组


代码中result[0][i][j]保存的是图中显示的到字符串Xi, Yj目前的LCS长度

result[1][i][j]保存的是图中显示的字符串Xi, Yj的指引方向关系

得到这幅图我们就可以从中得出表b[m, n]

为了得出最后的LCS字符串,只需要从b[m, n]开始,按照箭头方向追踪下去即可。

当b[i, j]遇到upLeft左上时,意味着Xi = Yj是LCS的一个元素.

按照这种方法可以逆序依次构造出LCS的所有元素

public static void printLcs(int[][][] b, String x, int i, int j){if(i == 0 || j == 0)return;if(b[1][i][j] == upLeft){printLcs(b, x, i - 1, j - 1);System.out.print(x.charAt(i - 1));}else if(b[1][i][j] == up){printLcs(b, x, i - 1, j);}else{printLcs(b, x, i, j - 1);}}

最后运行结果: