算法笔记学习000——Smith-Waterman算法寻找两个字符串中匹配度最高的子串
来源:互联网 发布:商标域名注册管理局 编辑:程序博客网 时间:2024/04/30 08:19
<span style="font-family: Arial, Helvetica, sans-serif; background-color: rgb(255, 255, 255);"></span>
<span style="font-family: Arial, Helvetica, sans-serif; background-color: rgb(255, 255, 255);"> </span>
<span style="font-family: Arial, Helvetica, sans-serif; background-color: rgb(255, 255, 255);">//第一次在CSDN上写笔记,是想着身为一个程序员,有些时候有些知识是转瞬即逝的,现在理解了说不定会很快忘记,所以记录一下加深印象,也有助于以后复习~</span>
<span style="font-family: Arial, Helvetica, sans-serif; background-color: rgb(255, 255, 255);">//本文参考了著名经典论文<em>Identification of Common Molecular Subsequences,</em>【其实根本就是照着人家论文学的嘛喂!</span>
<span style="font-family: Arial, Helvetica, sans-serif; background-color: rgb(255, 255, 255);">Smith-Waterman算法与1981年提出,出发点是信息生物学,意欲寻找两个碱基序列中匹配度最高的两个片段,以提取出有价值的生物学信息。</span>
Smith-Waterman算法的基本思路是,对于给定的长度为m的串s1(a1,a2,a3,...,am)和长度为n的串s2(b1,b2,...,bn),构造一个(m+1)*(n+1)的得分矩阵H。H中的元组H(i,j)的含义是,匹配子串(a1,a2,..,ai)和(b1,b2,...bj)的最高匹配度。
下面我们来讨论匹配度的定义。
我们知道,当匹配进行到(ai,bj)元组时,我们有三种选择:
1. 在bj元组前插入一个空格,用这个空格去匹配ai,匹配后的形式如下,此种匹配方式得分为-4/3
...,ai,...
...,--,bj...
2.在ai元组前插入一个空格,用这个空格去匹配bj,匹配后的形式如下,此种匹配方式得分为-4/3
...,--,ai,...
...,bj,...
3.直接匹配(ai,bj),如果ai与bj相同,则得分为1,否则得分为-1/3。
那么,匹配度H(i,j),就应该是在保证匹配度大于0的前提下这三种匹配方式中得分最高者。于是我们有了一个递归的定义:
H(i,j) = max{H(i-1,j-1) + MATCH(或者DISMATCH), H(i-1,j) - 4, H(i, j-1) - 4,0}
其中0<=i<=m, 0<=j<=n.
矩阵的第一行和第一列都为0,即H(0,j)和H(i,0)都等于零,相当于尚未开始匹配的匹配度。
有了匹配度矩阵H,我们的目标就是找到H中得分最大的元组H(i,j),并用回溯的方法找到一个轨迹,按照上述的评分规则沿轨迹从头到尾走一遍,就可以得到那个最高匹配度H(i,j),也就是找到s1,s2的匹配过程,自然也就找到了两个匹配度最高的子序列。
下面举例子来说明,这个例子也来源于Smith-Waterman的论文原文。我们假设需要匹配的两个序列分别为s1=AAUGCCAUUGACGG,S2=ACAGCCUCGCUUAG。
首先,计算匹配度矩阵H:
找到矩阵中得分最大(3.3)的元组H(10,8),开始回溯的过程。
其实回溯的思路很简单,就是检查位于该元组上方,左方,和左上方的元组,看它的得分是等于上-4/3,还是左-4/3,还是左上+1,还是左上-1/3。简而言之,就是看看这个元组是“从谁那儿走过来的”。
回溯终止的临界条件是,某个元组的得分为0,这意味着我们尚未找到匹配这两个串的子串头。
整个回溯过程结束后,找到的子串如下:
下面是用java语言写的源代码,大家可以到自己机器上试着跑一下~
import java.io.BufferedReader;import java.io.IOException;import java.io.InputStreamReader;import java.util.ArrayList;import java.util.Stack;public class SWSq { private int[][] H; private int[][] isEmpty; private static int SPACE ; //空格匹配的得分 private static int MATCH ; //两个字母相同的得分 private static int DISMACH; //两个字母不同的得分 private int maxIndexM, maxIndexN; private Stack<Character> stk1, stk2; public String subSq1, subSq2; //相似度最高的两个子串 public SWSq(){ stk1 = new Stack<Character>(); stk2 = new Stack<Character>(); SPACE = -4; MATCH = 3; DISMACH = -1; } private int max(int a, int b, int c){ int maxN; if(a >= b) maxN = a; else maxN = b; if(maxN < c) maxN = c; if(maxN < 0) maxN = 0; return maxN; } private void calculateMatrix(String s1, String s2, int m, int n){//计算得分矩阵 if(m == 0) H[m][n] = 0; else if(n == 0) H[m][n] = 0; else{ if(isEmpty[m - 1][n - 1] == 1) calculateMatrix(s1, s2, m-1, n-1); if(isEmpty[m][n - 1] == 1) calculateMatrix(s1, s2, m, n-1); if(isEmpty[m - 1][n] == 1) calculateMatrix(s1, s2, m-1, n); if(s1.charAt(m-1) == s2.charAt(n-1)) H[m][n] = max(H[m - 1][n - 1] + MATCH, H[m][n - 1] + SPACE, H[m - 1][n] + SPACE); else H[m][n] = max(H[m - 1][n - 1] + DISMACH, H[m][n - 1] + SPACE, H[m - 1][n] + SPACE); } isEmpty[m][n] = 0; } private void findMaxIndex(int[][] H, int m, int n){//找到得分矩阵H中得分最高的元组的下标 int curM, curN, i, j, max; curM = 0; curN = 0; max = H[0][0]; for(i = 0; i < m; i++) for(j = 0; j < n; j++) if(H[i][j] > max){ max = H[i][j]; curM = i; curN = j; } maxIndexM = curM; maxIndexN = curN; } private void traceBack(String s1, String s2, int m, int n){//回溯 寻找最相似子序列 if(H[m][n] == 0) return; if(H[m][n] == H[m-1][n] + SPACE) { stk1.add(s1.charAt(m-1)); stk2.add('-'); traceBack(s1, s2, m - 1, n); } else if(H[m][n] == H[m][n-1] + SPACE) { stk1.add('-'); stk2.add(s2.charAt(n-1)); traceBack(s1, s2, m, n - 1); } else { stk1.push(s1.charAt(m - 1)); stk2.push(s2.charAt(n-1)); traceBack(s1, s2, m - 1, n - 1); } } public void find(String s1, String s2){ //initMatrix(s1.length(), s2.length()); int i, j; H = new int[s1.length() + 1][s2.length() + 1]; isEmpty = new int[s1.length() + 1][s2.length() + 1]; for(i = 0; i<=s1.length(); i++) for(j = 0; j<=s2.length(); j++) isEmpty[i][j] = 1; calculateMatrix(s1, s2, s1.length(), s2.length()); findMaxIndex(H, H.length, H[0].length); traceBack(s1, s2, maxIndexM, maxIndexN); ArrayList<Character> arr1 = new ArrayList<Character>(); ArrayList<Character> arr2 = new ArrayList<Character>(); while(!stk1.empty()) arr1.add(stk1.pop()); subSq1 = arr1.toString(); while(!stk2.empty()) arr2.add(stk2.pop()); subSq2 = arr2.toString(); } public static void main(String[] args) throws IOException { SWSq x = new SWSq(); String s1,s2; BufferedReader br = new BufferedReader(new InputStreamReader(System.in)); System.out.println("Please enter one String:"); s1 = br.readLine(); System.out.println("Enter another String:"); s2 = br.readLine(); x.find(s1, s2); System.out.println("The subsequences with greatest similarity are " + x.subSq1 + " and " + x.subSq2); }}
- 算法笔记学习000——Smith-Waterman算法寻找两个字符串中匹配度最高的子串
- 简单的Smith Waterman算法实现
- Needleman-Wunsch 算法和Smith-Waterman算法
- Smith-waterman算法 openmp+mpi实现
- Needleman-wunsch 和 Smith-Waterman 比对算法
- [算法学习笔记]寻找子字符串第一次出现位置
- 寻找最长回文子串Manacher算法学习笔记
- 动态编程之序列比对:Needleman-Wunsch 算法和Smith-Waterman算法
- 两个字符串模式匹配的算法
- 朴素匹配算法-子字符串的查找
- 数据结构——算法之(028)( 寻找其中的一个子字符串个数)
- 数据结构——算法之(028)( 寻找其中的一个子字符串个数)
- 匹配两个字符串的最大子串
- 笔试——字符串算法题——寻找最大回文子串
- 关于寻找两个字符串中最长子序列的问题
- 字符串匹配-KMP算法学习笔记
- 算法题-两个字符串的最大公共子串
- 获取两个字符串所有公共的子串算法
- Linux下安装telnet服务
- HDOJ题目2870 Largest Submatrix(动态规划)
- 在蓝鸥的日子 2014.10.19
- Multiply String
- android 学习笔记3——WebView的使用
- 算法笔记学习000——Smith-Waterman算法寻找两个字符串中匹配度最高的子串
- 小说明
- 《数据结构》第一章知识导图
- 直击用户大脑——用户研究新方法(眼动与脑电数据分析)
- delete后仍然可以调用问题
- ajax简介
- 如何在swift中自定义基本类型Bool
- Cocos2d-x中精灵图片的渲染顺序
- hdu-oj 3351 Seinfeld