相似基因序列


/*众所周知,人类基因可以看作一个碱基对序列,它包含了4种核苷酸,简记为A,C,G,T。
让我们观察这样一段基因序列 “ACCAGGTT”,这段序列共由8个核苷酸构成,其中第1位和第4位是碱基“A”
第2位和第3位是碱基“C”,第5位和第6位是碱基“G”,第7位和第8位是碱基“T”。
Tom构造了这样一个0,1矩阵:
1, 0, 0, 1, 0, 0, 0, 0
0, 1, 1, 0, 0, 0, 0, 0
0, 1, 1, 0, 0, 0, 0, 0
1, 0, 0, 1, 0, 0, 0, 0
0, 0, 0, 0, 1, 1, 0, 0
0, 0, 0, 0, 1, 1, 0, 0
0, 0, 0, 0, 0, 0, 1, 1
0, 0, 0, 0, 0, 0, 1, 1
如果第i位的碱基与第j位的碱基一样,那么0,1矩阵的i行j列为1,否则为0。
如果基因序列X与基因序列Y等长且具有相同的0,1矩阵,Tom就会认为X与Y是相似的基因序列。 
现在的问题是:给你两段长度为N的基因序列,请你帮助Tom判断它们是否相似。
Input
可以有多组测试数据,每组数据第1行输入一个正整数N(1≤N≤1000)
第2行和第3行分别输入两段长度为N的基因序列(只由A,C,G,T四种字符构成)
Output
每组数据输出仅一行,如果相似则输出 YES,否则输出 NO。
解
看规律:对称的矩阵这个规律不是主要的。而且此稀疏矩阵矩阵的利用率很差,
所以不要用矩阵,在观察时发现一些结论,先举个例子:
AACCTTGG与
CCAATTGG
是相似的 因为A C T G 个数相等而且 第二个基因序列可以看作 第一个基因序列 中A与C做了全换,T与G做了全换---所以:第一个字符串中只有个数相等的才能完全替换;
AAATTCCCG
GGGTTAAAC
AAAGGTTTC······
也是相似的,因为对应的位置做了全换
上述的四个字母的个数是相同的,我们可以根据字母的在序列的位置判断出 谁与谁对应 
AAACCTTG
CCCTATAG
不是相似的 虽然第一个字符串中的CT 与第二个中的TA 个数相等 但对应的位置不同。
AAACCTTAGT
CCCTTGGCA
如果我们还是只看匹配的下标,那么会发现总是在接近末尾才能判断是不是相似序列
所以我们上来就看4个字母的数量是否一一匹配 数量匹配了 则进行全匹配
那么效率会高一些
用"A","C","G","T"与两个基因序列匹配得到的下标
假设用二维数组aaa[8][N]存第一个字符串匹配的下标存在aaa[0--3],
第二个aaa[4--7])满足aaa[0--3]与aaa[4--7]下标值某两个一一对应共四组;
假设:
AACCTTTG与
CCTTAAAG
分别用"A","C","G","T"替换
则对于第一个基因序列的各个下标以数组下标为准:
A  0 1      C  2 3    T 4 5 6  G   7   
第二个:
A  4 5 6    C  0 1    T 2 3   G   7   
我们该怎样找到以上四个的一一对应关系呢? 
如果我以第一个字符串的字母为准 
则第一个字符串中的A匹配了两个 C 2个 T 3个  G 1个
第二个字符串中 A 3  C 2  T 2 G 7  
将以上的各项匹配的下标与个数合并 得到下表:
AACCTTTG:
---------------A  0 1   2     个
---------------C  2 3   2     个
---------------T  4 5 6 3     个
---------------G  7 1         个
CCTTAAAG
---------------A  4 5 6 3     个
---------------C  0 1   2     个
---------------T  2 3   2     个
---------------G  7 1         个
然后看上表中 字母 与 "个"之间的数字  如果我们忽略字母(AGCT)
那么中间的数字如果我们得到了 并对上下两个字符串找到一一对应关系 这道题就解决了
第一个中的A 与第二个中的 C 相对应   ,第一个中的G与第二个中的 G相对应   ,第一个中的C 与第二个中的 T相对应   ,第一个中的T 与第二个中的 A 相对应   
这两个得到相似
大体流程:
输入两个基因序列;
对两个序列用  A C  T G 匹配
将两个基因序列找到对应关系
得到对应之后 再全部检验;
*/
# include <stdio.h>
# define  N 10001    //如果想录入x个字节那么就把N的数值改成x+1
int gainchar(char *a,int min,int max);//返回字符串长度
void BF(char a[],char b[],int *c);//a为主串,b为被检验的串 c保存匹配的下标
int main(){
	int i,j,LEN,*p[4]={NULL,NULL,NULL,NULL};
	int  pipei[8][N]={0,0,0,0,0,0,0,0};//初始化
	char a[2][N];
    char SS[][3]={"一","二"},moban[][2]={"A","C","G","T"};
	do{
		 do{
			 LEN=N+1;
			 printf("请输入基因片断的长度:(4--%d):\n",N-1);
			 scanf("%d",&LEN);
			 while(getchar()!='\n');
		 }while(LEN<4||LEN>N);
		for(i=0;i<2;i++)
		{
			printf("输入第%s个基因大写序列:长度(%d--%d):\n",SS[i],LEN,LEN);
			gainchar(a[i],LEN,LEN+1);
			for(j=0;j<LEN;j++)
				if(a[i][j]!=moban[0][0]&&a[i][j]!=moban[1][0]&&a[i][j]!=moban[2][0]&&a[i][j]!=moban[3][0])
				{
					printf("输入的基因字母有误!\n请重新");
					i--;
					break;
				}
				
		}
		printf("基因序列:\n%s\n----------------------------\n%s\n",a[0],a[1]);
		for(i=0;i<2;i++)
			for(j=0;j<4;j++)   //将ACGT的匹配的下标存储在pipei数组里有两组字符串 共八个数组
			  BF(a[i],moban[j],pipei[i*4+j]);
		for(i=0;i<4;i++)
			for(j=4;j<8;j++)   //将p[i]指向 第二个基因序列中与moban[i]匹配的下标相同的pipei[j]的地址
			 if(pipei[i][1]==pipei[j][1]&&pipei[i][0]==pipei[j][0])
				p[i]=&pipei[j][0];
		if(p[0]&&p[1]&&p[2]&&p[3])   //如果第一个和第二个字符串是一一对应则继续检查
		{
			for(i=0;i<4;i++)
				for(j=0;j<pipei[i][0];j++)
					if(pipei[i][j]!=p[i][j])
					{
						printf("NO!\n");
						break;
					}
			if(i==4)  //如果检查完了没问题  则为YES
				printf("YES!\n");
		}
		else
			printf("NO!\n");//如果不是一一对应则NO
	p[0]=p[1]=p[2]=p[3]=NULL;//再初始化
	}while(LEN);
	return 0;
}
int gainchar(char *a,int min,int max)//对*a输入范围[min,max-1]
{
	int c,k;
	do{
		c=-1;k=0;
		fgets(a,max,stdin);
		while(a[++c]);
		c=a[c-1]=='\n'&&c<max?c-1:c;
		if(c>=max-1)
			while(getchar()!='\n')
				k++;
			else
				a[c]='\0';
			if(k||c&&c<min)//如果用户只输入'\n'则不提示输入错误,否则提示错误
				printf("您输入字符串长度为:%d字节,请重新输入!\n注:只录入(%d--%d)字节:\n",c+k,min,max-1);
	}while(k||c<min);
    return c;
}
void BF(char a[],char b[],int *c)//BF算法 一般o(strlen(A)+strlen(B))  最坏o(strlen(A)*strlen(B))
{                     
	int i=0,j=0,k=0; 
    do{
		if (b[j]&&a[i++]==b[j]) // 继续比较后继字
			++j;
		else
		{
			if(!b[j])//如果b字符串已经结束
				c[++k]=i-j;
			else
			    i-=j;  //否则i回溯
			j=0;      //无论b字符串是否结束因为(b[j]&&a[i++]==b[j])已经不满足 所以j要回到初位置
		}
	}while(a[i-1]);
		c[0]=k;
}

---------------------------------------------------多日之后的改进

# include <stdio.h>  
# define  N 10001    //如果想录入x个字节那么就把N的数值改成x+1  
int gainchar(char *A,int min,int max);//返回字符串长度  
void BF(char a[],char b[],int *c);//a为主串,b为被检验的串 c保存匹配的下标  
int main(){  
    int i,j,LEN,*p[4]={NULL},pipei[8][N]={0};//初始化  
    char a[2][N];  
    char SS[][3]={"一","二"},moban[][2]={"A","C","G","T"};  
    do{  
         do{  
             LEN=N+1;  
             printf("请输入基因片断的长度:(4--%d):\n",N-1);  
             scanf("%d",&LEN);  
             while(getchar()!='\n');  
         }while(LEN<4||LEN>N);  
        for(i=0;i<2;i++)  
        {  
            printf("输入第%s个基因大写序列:长度(%d--%d):\n",SS[i],LEN,LEN);  
            gainchar(a[i],LEN,LEN);  
            for(j=0;j<LEN;j++)  
                if(a[i][j]!=moban[0][0]&&a[i][j]!=moban[1][0]&&a[i][j]!=moban[2][0]&&a[i][j]!=moban[3][0])  
                {  
                    printf("输入的基因字母有误!\n请重新");  
                    i--;  
                    break;  
                }  
        }    
        for(i=0;i<2;i++)  
            for(j=0;j<4;j++)   //将ACGT的匹配的下标存储在pipei数组里有两组字符串 共八个数组  
              BF(a[i],moban[j],pipei[i*4+j]);  
        for(i=0;i<4;i++)  
            for(j=4;j<8;j++)   //将p[i]指向 第二个基因序列中与moban[i]匹配的下标相同的pipei[j]的地址  
             if(pipei[i][1]==pipei[j][1]&&pipei[i][0]==pipei[j][0])  
                p[i]=&pipei[j][0];  
        if(p[0]&&p[1]&&p[2]&&p[3])   //如果第一个和第二个字符串是一一对应则继续检查  
        {  
            for(i=0;i<4;i++)  
                for(j=0;j<pipei[i][0];j++)  
                    if(pipei[i][j]!=p[i][j])  
                    {  
                        printf("No!\n");  
                        break;  
                    }  
            if(i==4)  //如果检查完了没问题  则为YES  
                printf("Yes!\n");  
        }  
        else  
            printf("No!\n");//如果不是一一对应则NO  
    p[0]=p[1]=p[2]=p[3]=NULL;//再初始化  
    }while(LEN);  
    return 0;  
}  
int gainchar(char *A,int min,int max)//长度在[min,max]  <闭区间>  之间时 函数结束 返回字符串A的长度        
{        
    int B,C;      
 do{        
        A[max]=B=C=0;    
        while((A[B++]=getchar())!='\n'&&B<max);    
        if(A[B-1]!='\n')while(getchar()!='\n'&&++C);        
        else A[--B]=0;     
    if(C||B&&B<min)    
       printf("您录入的字符串长度:%d字节\n只录入[%d,%d]个字节!\n",B+C,min,max);        
    }while(C||B<min);        
    return B;      
}   
void BF(char a[],char b[],int *c)//a为主串 b为次串 c[0]存储匹配的个数 c[1]及后面存储匹配的下标   
{                         
    int i=0,j=0,k=0;     
    do{    
        if (b[j]&&a[i++]==b[j])++j;    
        else    
        {    
            b[j]?(i-=j):(c[k++]=i-j);     
            j=0;                
        }    
    }while(a[i-1]);    
    c[0]=k;    
}  



评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值