最大子序列、最长递增子序列、最长公共子串、最长公共子序列、字符串编辑距离

最大子序列

最大子序列是要找出由数组成的一维数组中和最大的连续子序列。比如{5,-3,4,2}的最大子序列就是 {5,-3,4,2}，它的和是8,达到最大；而 {5,-6,4,2}的最大子序列是{4,2}，它的和是6。你已经看出来了，找最大子序列的方法很简单，只要前i项的和还没有小于0那么子序列就一直向后扩展，否则丢弃之前的子序列开始新的子序列，同时我们要记下各个子序列的和，最后找到和最大的子序列。

代码如下：

 

#include<iostream>
#include<vector>
using namespace std;
int maxSubSum(const vector<int> & arr,int &begin,int &end){
    int maxSum=0;
    int currSum=0;
    int newbegin=0;
    for(int i=0;i<arr.size();++i){
        currSum+=arr[i];
        if(currSum>maxSum){
            maxSum=currSum;
            begin=newbegin;
            end=i;
        }
        if(currSum<0){
            currSum=0;
            newbegin=i+1;
        }
    }
    return maxSum;
}

int main(){
    int len;
    cout<<"Input array length"<<endl;
    cin>>len;
    cout<<"Input an integer vector"<<endl;
    vector<int> arr;
    int a;
    for(int i=0;i<len;++i){
        cin>>a;
        arr.push_back(a);
    }
    int begin,end;
    cout<<maxSubSum(arr,begin,end)<<endl;
    for(int i=begin;i<=end;++i)
        cout<<arr[i]<<" ";
    cout<<endl;
    return 0;
}

最长递增子序列

 比如arr={1,5,8,2,3,4}的最长递增子序列是1,2,3,4

动态规划的思想，考虑{arr[0],...,arr[i]}的最长递增子序列时需要找到所有比arr[i]小的arr[j]，且j<i，结果应该是所有{arr[0],...,arr[j]}的最长递增子序列中最长的那一个再加1。即我们需要一个辅助的数组aid_arr，aid_arr[i]的值是{arr[0],...,arr[i]}的最长递增子序列的长度，aid_arr[0]=1。

 

#include<iostream>
#include<stack>
#include<vector>

using namespace std;

int main(){
	const int len=14;
	int arr[len]={1,9,3,8,11,4,5,6,4,1,9,7,1,7};
	vector<int> vec(&arr[0],&arr[len]);

	vector<int> monoseqlen(len,1);
	vector<int> preindex(len,-1);
	int maxmonoseqlen=-1;
	int maxmonoindex=-1;

	for(int i=1;i<len;++i){
		int curr=vec[i];
		for(int j=0;j<i;++j){
			if(vec[j]<vec[i]){
				int msl=monoseqlen[j]+1;
				if(msl>monoseqlen[i]){
					monoseqlen[i]=msl;
					preindex[i]=j;
				}
			}
		}
	}

	for(int i=0;i<len;++i){
		if(monoseqlen[i]>maxmonoseqlen){
			maxmonoseqlen=monoseqlen[i];
			maxmonoindex=i;
		}
	}

	stack<int> st;
	while(maxmonoindex>=0){
		st.push(vec[maxmonoindex]);
		maxmonoindex=preindex[maxmonoindex];
	}

	vector<int> rect;
	while(!st.empty()){
		rect.push_back(st.top());
		st.pop();
	}	

	vector<int>::iterator itr=rect.begin();
	while(itr!=rect.end()){
		cout<<*itr++<<"\t";
	}
	cout<<endl;
	return 0;
}

最长公共子串（LCS）

找两个字符串的最长公共子串，这个子串要求在原字符串中是连续的。其实这又是一个序贯决策问题，可以用动态规划来求解。我们采用一个二维矩阵来记录中间的结果。这个二维矩阵怎么构造呢？直接举个例子吧："bab"和"caba"(当然我们现在一眼就可以看出来最长公共子串是"ba"或"ab")

　　 b　　a　　b

c　　0　　0　　0

a　　0　　1　　0

b　　1　　0　　1

a　　0　　1　　0

我们看矩阵的斜对角线最长的那个就能找出最长公共子串。

不过在二维矩阵上找最长的由1组成的斜对角线也是件麻烦费时的事，下面改进：当要在矩阵是填1时让它等于其左上角元素加1。

　　 b　　a　　b

c　　0　　0　　0

a　　0　　1　　0

b　　1　　0　　2

a　　0　　2　　0

这样矩阵中的最大元素就是 最长公共子串的长度。

在构造这个二维矩阵的过程中由于得出矩阵的某一行后其上一行就没用了，所以实际上在程序中可以用一维数组来代替这个矩阵。

代码如下：

 

#include<iostream>
#include<cstring>
#include<vector>
using namespace std;
//str1为横向，str2这纵向
const string LCS(const string& str1,const string& str2){
	int xlen=str1.size();		//横向长度
	vector<int> tmp(xlen);		//保存矩阵的上一行
	vector<int> arr(tmp);		//当前行
	int ylen=str2.size();		//纵向长度
	int maxele=0;				//矩阵元素中的最大值
	int pos=0;					//矩阵元素最大值出现在第几列
	for(int i=0;i<ylen;i++){
		string s=str2.substr(i,1);
		arr.assign(xlen,0);		//数组清0
		for(int j=0;j<xlen;j++){
			if(str1.compare(j,1,s)==0){
				if(j==0)
					arr[j]=1;
				else
					arr[j]=tmp[j-1]+1;
				if(arr[j]>maxele){
					maxele=arr[j];
					pos=j;
				}
			}		
		}
//		{
//			vector<int>::iterator iter=arr.begin();
//			while(iter!=arr.end())
//				cout<<*iter++;
//			cout<<endl;
//		}
		tmp.assign(arr.begin(),arr.end());
	}
	string res=str1.substr(pos-maxele+1,maxele);
	return res;
}
int main(){
	string str1("21232523311324");
	string str2("312123223445");
	string lcs=LCS(str1,str2);
	cout<<lcs<<endl;
	return 0;
}

 

最长公共子序列

最长公共子序列与最长公共子串的区别在于最长公共子序列不要求在原字符串中是连续的，比如ADE和ABCDE的最长公共子序列是ADE。

我们用动态规划的方法来思考这个问题如是求解。首先要找到状态转移方程：

符号约定，C1是S1的最右侧字符，C2是S2的最右侧字符，S1‘是从S1中去除C1的部分，S2'是从S2中去除C2的部分。

LCS(S1,S2)等于下列3项的最大者：

（1）LCS（S1，S2’）

（2）LCS（S1’，S2）

（3）LCS（S1’，S2’）--如果C1不等于C2； LCS（S1'，S2'）+C1--如果C1等于C2；

边界终止条件：如果S1和S2都是空串，则结果也是空串。

下面我们同样要构建一个矩阵来存储动态规划过程中子问题的解。这个矩阵中的每个数字代表了该行和该列之前的LCS的长度。与上面刚刚分析出的状态转移议程相对应，矩阵中每个格子里的数字应该这么填，它等于以下3项的最大值：

（1）上面一个格子里的数字

（2）左边一个格子里的数字

（3）左上角那个格子里的数字（如果 C1不等于C2）； 左上角那个格子里的数字+1（ 如果C1等于C2）

举个例子：

 　　　　  G　　C　　T　　A

　　 0　　0　　0　　0　　0

G　 0　　1　　1　　1　　1

B　 0　　1　　1　　1　　1

T　 0　　1　　1　　2　　2

A    0　　1　　1　　2　　3

填写最后一个数字时，它应该是下面三个的最大者：

（1）上边的数字2

（2）左边的数字2

（3）左上角的数字2+1=3,因为此时C1==C2

所以最终结果是3。

在填写过程中我们还是记录下当前单元格的数字来自于哪个单元格，以方便最后我们回溯找出最长公共子串。有时候左上、左、上三者中有多个同时达到最大，那么任取其中之一，但是在整个过程中你必须遵循固定的优先标准。在我的代码中优先级别是左上>左>上。

下图给出了回溯法找出LCS的过程：

奉上代码：

 

#include<iostream>
#include<cstring>
#include<stack>
#include<utility>
#define LEFTUP 0
#define LEFT 1
#define UP 2
using namespace std;
int Max(int a,int b,int c,int *max){			//找最大者时a的优先级别最高，c的最低.最大值保存在*max中
	int res=0;			//res记录来自于哪个单元格
	*max=a;
	if(b>*max){
		*max=b;
		res=1;
	}
	if(c>*max){
		*max=c;
		res=2;
	}
	return res;
}
//调用此函数时请注意把较长的字符串赋给str1,这主要是为了在回溯最长子序列时节省时间。如果没有把较长的字符串赋给str1不影响程序的正确执行。
string LCS(const string &str1,const string &str2){
	int xlen=str1.size();				//横向长度
	int ylen=str2.size();				//纵向长度
	if(xlen==0||ylen==0)				//str1和str2中只要有一个为空，则返回空
		return "";
	pair<int,int> arr[ylen+1][xlen+1];	//构造pair二维数组，first记录数据，second记录来源
	for(int i=0;i<=xlen;i++)			//首行清0
		arr[0][i].first=0;
	for(int j=0;j<=ylen;j++)			//首列清0
		arr[j][0].first=0;
	for(int i=1;i<=ylen;i++){
		char s=str2.at(i-1);
		for(int j=1;j<=xlen;j++){
			int leftup=arr[i-1][j-1].first;
			int left=arr[i][j-1].first;
			int up=arr[i-1][j].first;
			if(str1.at(j-1)==s)			//C1==C2
				leftup++;
			int max;
			arr[i][j].second=Max(leftup,left,up,&arr[i][j].first);
//			cout<<arr[i][j].first<<arr[i][j].second<<" ";
		}
//		cout<<endl;
	}		/*矩阵构造完毕*/
	//回溯找出最长公共子序列
	stack<int> st;
	int i=ylen,j=xlen;
	while(i>=0&&j>=0){
		if(arr[i][j].second==LEFTUP){
			if(arr[i][j].first==arr[i-1][j-1].first+1)
				st.push(i);
			--i;
			--j;
		}
		else if(arr[i][j].second==LEFT){
			--j;
		}
		else if(arr[i][j].second==UP){
			--i;
		}
	}
	string res="";
	while(!st.empty()){
		int index=st.top()-1;
		res.append(str2.substr(index,1));
		st.pop();
	}
	return res;
}
int main(){
	string str1="GCCCTAGCG";
	string str2="GCGCAATG";
	string lcs=LCS(str1,str2);
	cout<<lcs<<endl;
	return 0;
}

下面给一个Java版本

public static <E> List<E> longestCommonSubsequence(E[] s1,E[] s2){
		int[][] num=new int[s1.length+1][s2.length+1];
		for(int i=1;i<s1.length+1;i++){
			for(int j=1;j<s2.length+1;j++){
				if(s1[i-1].equals(s2[j-1])){
					num[i][j]=1+num[i-1][j-1];
				}
				else{
					num[i][j]=Math.max(num[i-1][j],num[i][j-1]);
				}
			}
		}
		System.out.println("lenght of LCS= "+num[s1.length][s2.length]);
		int s1position=s1.length,s2position=s2.length;
		List<E> result=new LinkedList<E>();
		while(s1position>0 && s2position>0){
			if(s1[s1position-1].equals(s2[s2position-1])){
				result.add(s1[s1position-1]);
				s1position--;
				s2position--;
			}
			else if(num[s1position][s2position-1]>=num[s1position-1][s2position]){
				s2position--;
			}
			else{
				s1position--;
			}
		}
		Collections.reverse(result);
		return result;
	}

std::endl是一个特殊值，称为操纵符（manipulator），将它写入输出流时具有输出换行的效果，并刷新与设备相关联的缓冲区(buffer)。通过刷新缓冲区用户可立即看到写入到流中的输出。

我在调试以上代码时在45行（cout<<endl;）处设置断点，结果发现“43行（cout<<arr[i][j].first<<arr[i][j].second<<" ";） 没有执行”，这就是因为43行末尾没有加endl，所以用户没有立即看到输出流中的数据。

python代码 

def lcs(arr, brr):
    l1 = len(arr)
    l2 = len(brr)
    if l1 == 0 or l2 == 0:
        return 0

    tmp_list_1 = [0] * (l1 + 1)
    tmp_list_2 = [0] * (l1 + 1)

    for i in xrange(l2):
        for j in xrange(l1):
            ele2 = brr[i]
            ele1 = arr[j]
            if ele2 == ele1:
                tmp_list_2[j + 1] = tmp_list_1[j] + 1
            else:
                tmp_list_2[j + 1] = max(tmp_list_2[j], tmp_list_1[j])
        tmp_list_1 = tmp_list_2
        tmp_list_2 = [0] * (l1 + 1)

    return tmp_list_1[l1]


arr = list('GCCCTAGCG')
brr = list('GCGCAATG')
print lcs(arr, brr)

　　

字符串编辑距离

 要想把字符串S1变成S2，可以经过若干次下列原子操作：

1.删除一个字符

2.增加一个字符

3.更改一个字符

字符串S1和S2的编辑距离定义为从S1变成S2所需要原子操作的最少次数。

解法跟上面的最长公共子序列十分相似，都是动态规划，把一个问题转换为若干个规模更小的子问题，并且都借助于一个二维矩阵来实现计算。

约定：字符串S去掉最后一个字符T后为S'，T1和T2分别是S1和S2的最后一个字符。

则dist(S1,S2)是下列4个值的最小者：

1.dist(S1',S2')--当T1==T2

2.1+dist(S1',S2)--当T1!=T2，并且删除S1的最后一个字符T1

3.1+dist(S1,S2')--当T1!=T2，并且在S1后面增加一个字符T2

4.1+dist(S1',S2')--当T1!=T2，并且把S1的最的一个字符T1改成T2

 

把问题转换为二维矩阵：

arr[i][j]表示S1.sub(0,i)和S2.sub(0,j)的编辑距离，则

arr[i][j]=min{1+arr[i][j-1], 1+arr[i-1][j], 1+arr[i-1][j-1](当S1[i]!=S2[j]), arr[i-1][j-1](当S1[i]==S2[j])}

边界情况：arr[0][j]=j, arr[i][0]=i

代码就不写了，跟最长公共子序列很相似。

 

计算两个字条串的相似度除了Edit Distance，还有一种方法是计算Jaro Distance。具体怎么算读者可以搜一下。