2022年认证杯B题第三问思路

序列的遗传过程

序列同源性是DNA、RNA或蛋白质序列之间的生物学同源性，根据生命进化史中的共 同祖先定义。DNA、RNA或蛋白质之间的同源性通常是从它们的核苷酸或氨基酸序列 的相似性中推断出来的。显著的相似性有力地证明了两个序列与一个共同的祖先序列的进化变化相关。 考虑一个RNA序列的遗传过程，其中核苷酸碱基的突变是偶然发生的。为简单起 见，我们假设序列突变是由于单个碱基的变化（过渡或转换）、插入和删除而引起 的。所以我们可以用突变点的数量来测量两个序列的距离。紧密相连的多个碱基序列可以形成一个家族，它们被认为是同源的。 您的团队被要求开发一个合理的数学模型来完成以下问题。

3.如果一个家族中的多个碱基序列是从一个共同的祖先序列进化而来的，则设 计一种有效的算法来确定祖先序列，并绘制系谱树。

对于祖先序列，我们做出第一个假定，即祖先序列与其他的所有序列的差异突变点数目之和最少（当然也可取数量的平方和，立方和等等，这里简单考虑），这样便能保证祖先序列和其他所有序列相比，公共序列所含碱基数最多。

由第一问我们可以得到两两序列之间的碱基差异，即突变点的个数。突变点数目越多，两物种之间的差异性越大。同时由进化论的思想，突变是在亲代繁衍子代的过程中产生的，一旦突变点数目大于0，便存在繁衍过程。当突变点数目达到一定数量时，两物种之间形成生殖隔离，当数目过大时便无法保证物种间的亲子代关系。定义形成两物种之间的突变点的最大数量为两物种间的碱基差距，所以第二个假定便是两序列突变点序列小于碱基差距时，我们认为这两序列具有亲子代关系。

第三个假定：由第一条假定推广，第二代序列与其他所有序列对比，突变点数量之和介于祖先序列和第三代与其他所有序列对比后的突变点数量之和，之后的第N代子代序列以此类推。

第四个假定：当一条序列与所有序列对比后的突变点数量都大于碱基差距时，说明其既不存在亲代，也不存在子代，说明发生了断代现象。

做完上述假定后，我们便可以设计算法：

第一步先将问题简化，规定同时可检测的碱基序列为20条，同时限定只能检测出5代以内的序列亲子代关系。

第二步将碱基序列两两进行比较，并用矩阵进行储存。N条序列之间进行比对，可以得到一个N×N的矩阵，其中矩阵的第i行第j列，即A（i，j）表示第i条碱基序列和第j条碱基序列之间的不同碱基数。不难发现这个矩阵同时也是个实对称矩阵。且对角线上元素均为0（即某序列和本身做比对，突变点数为0），亲子关系就储存在这个矩阵中。

例如：随机输入一个9x9的矩阵如下

 

该矩阵表示共有9个碱基序列参与比对，且第一个与第二个序列之间突变点数为1，第一个与第三个序列之间突变点数为2，依次类推。

接着，可以任意选择判断亲缘关系时，突变点数的允许范围。

例如：保持例子的条件不变，输入碱基差距为2。它表示：当任意两个序列之间的突变点数小于等于2时，可以认为两者是相邻代

 

分别求出矩阵每一行的元素和，所得最小和所在行对应的序列即为祖先序列（可以有一个或多个祖先序列）。它表示：当一个或多个序列和所有序列比较后，所得突变点数的和最小时，该序列即为祖先序列。即按照生物遗传学基本规律，祖先序列与其他序列公共碱基数最多，突变点数最少。

例如：保持例子的条件不变。通过c++运行可得结果如下图所示。

             

确定完祖先序列后，由第二代开始依次寻找除祖先序列以外其他序列的代数。其基本思想为：在已知代数的序列所在的列中比较该序列与其他序列比较后的突变点数。判断条件仍然保持为“判断亲缘关系时，突变点数的允许范围”。

例如：保持例子的条件不变。可得c++运行结果如下图所示。

 

最后，根据已经判断的代系关系，应用“子代寻找父代”的基本思想，绘制出家谱树。

 

它的含义可以用matlab更清楚的表示出来：

 

附上代码：

#include<iostream>
using namespace std;
//除去int型预置的空间，并输出int 型数组的实际长度。
//大于等于0的原因是序列比对矩阵的行列元素均大于等于0，而未分配空间的元素为-858993460。
int number(int a[])
{int mylen=0;
	int i=0;
	for(;a[i]>=0;i++)
	{mylen=mylen+1;}
	return mylen;}


//除去数组里面重复的元素，并返回去除重复元素后的数组长度。
int remove(int a[],int n,int m)
{if (n==1){cout<<"第"<< a[0]<<"条碱基序列为第"<<m<<"代" <<endl;return n=1;}
else{	int  flag = 1,f=0;
	for (int i = 0; i<n; i++) {
		for (int k = 0; k<i; k++) {
			if (a[i] == a[k]) {
				flag = 0;
			}
		}
		if (flag) {
			cout <<"第"<< a[i]<<"条碱基序列为第"<<m<<"代" <<endl;f=f+1;
		}
		flag = 1;
	}n=f;return n;}}

//将两个数组的相同元素除去，并返回除去后数组的长度
int clear(int s[],int r[],int n){if(n==1){int k;r[0]=s[0];return k=0;}
else{
	int flag=1,f=0;
	for(int i=0;i<n;i++){
		for(int j=0;j<number(r);j++){
			if(s[i]==r[j]){flag=0;}}if(flag){s[f]=s[i];f=f+1;}
			flag=1;}return f;}}
//通过上一代，确定下一代
void findnext(int s[],int l,int fa[],int son[],int c[20][20],int d,int sum,int a)
   {int f=0;
    for(int j=0;j<l;j=j+1){
	for(int i=0;i<number(fa);i=i+1){
	if(c[s[j]-1][fa[i]-1]<=d)
            {son[f]=s[j];f=f+1;}
	                                              }
	                                              }
	         cout<<endl;}
//用于存放与祖先序列差距在5代或以上的剩余碱基序列。
void restplace(int r5[],int s[],int k){cout<<"剩余的碱基序列:";
	for(int j=0;j<k;j++){
		r5[j]=s[j];
		cout<<"第"<<r5[j]<<"条"<<" ";}cout<<"与祖先序列差距在5代或以上"<<endl;}
//用于检索与子代最接近的亲代，并绘制子代与亲代的对应关系，子代指向亲代。
void findclose(int fa[],int son[],int c[20][20]){int min;
if(number(fa)==0&&number(son)!=0)
       {cout<<"出现断代现象!!!无法检验!!!"<<endl;}
else{
	for(int i=0;i<number(son);i++){
		for(int j=0;j<number(fa);j++){
			if(j==0) {
			min=c[son[i]-1][fa[j]-1];	
			}
		else if(c[son[i]-1][fa[j]-1]<=min){
			min=c[son[i]-1][fa[j]-1];
			}}	//子代和最近的亲代之间，对应的突变点数量最小，并找出最小值
		for(int k=0;k<number(fa);k++){
		if(c[son[i]-1][fa[k]-1]==min){cout<<son[i]<<"→"<<fa[k]<<endl;}
		}
	}
}}

int main()
{cout<<"输入的突变点矩阵必须为正方形矩阵！"<<endl;
cout<<"A(i,j)表示第i个碱基与第j个碱基对比得到的突变点数"<<endl;
	int a,b,d,f=0,k,sum=0;cout<<"请输入行数"<<":";cin>>a;
cout<<"请输入列数"<<":";cin>>b;cout<<"请输入突变点矩阵:"<<endl;
int str[20],s[20],r1[20],r2[20],r3[20],r4[20],r5[20];
int c[20][20];for(int d=0;d<a;d++){for(int e=0;e<b;e++)
{cin>>c[d][e];}};cout<<"请输入形成两个物种的碱基差距:";cin>>d;
for(int d=0;d<a;d++){for(int e=0;e<b;e++)
{f=c[d][e]+f;str[d]=f;}f=0;};for(int i=0;i<a;i++)
{s[i]=i+1;};cout<<endl;
for (int i = 0; i < a; i++)
	{for (int j = 0;j<a-i-1; j++)
		{if (str[j]>str[j + 1])
			{int temp = str[j];str[j]=str[j+1];str[j+1] = temp;
int v=s[j];s[j]=s[j+1];s[j+1]=v;}}}
			for(int i=0; i<a;i++)
			{cout<<"第"<<s[i]<<"条碱基序列与其他碱基序列的突变位点之和为"<<str[i]<<endl;}cout<<endl;
			for(int i=0;str[i]==str[0];i++)
			{{r1[i]=s[i];cout<<"第"<<s[i]<<"条碱基序列为祖先序列"<<endl;}}
			sum=sum+remove(r1,number(r1),1);k=clear(s,r1,a);
			findnext(s,k,r1,r2,c,d,sum,a);
			sum=sum+remove(r2,number(r2),2);k=clear(s,r2,k);
			findnext(s,k,r2,r3,c,d,sum,a);
			sum=sum+remove(r3,number(r3),3);k=clear(s,r3,k);
			findnext(s,k,r3,r4,c,d,sum,a);
			sum=sum+remove(r4,number(r4),4);k=clear(s,r4,k);
			if (sum!=a){restplace(r5,s,k);}
			cout<<"家谱树对应的上下代关系为"<<endl;
            findclose(r1,r2,c);
				findclose(r2,r3,c);
				findclose(r3,r4,c);
				findclose (r4,r5,c);                 
				return 0;
}