Helenr
写一个思路candidates = [ 'GGGAGCAGGCAAGGACTCTG', 'GCTCGGGCTTGTCCACAGGA', '...', # 被你看出来啦,这些其实人类基因的片段]
bg_db = [ 'CTGCTGACGGGTGACACCCA', 'AGGAACTGGTGCTTGATGGC', '...', # 这个更多,有十亿左右]因为你的数据其实是很有特点的,这里可以进行精简。因为所有的字符串都是20个字符长度,且都由ATCG四个字符组成。那么可以把它们变换为整数来进行比较。二进制表现形式如下A ==> 00T ==> 01C ==> 10G ==> 11因为一个字符串长度固定,每个字符可以由2个比特位表示,所以每个字符串可以表示为一个40位的整数。可以表示为32+8的形式,也可以直接使用64位整形。建议使用C语言来做。再来说说比较。因为要找到每一个candidate在bg_db中与之差异小于等于4的所有记录,所以只要两个整数一做^按位异或操作,结果中二进制中1不超过8个,且这不超过8个1最多只能分为4个组的才有可能是符合要求的(00^11=11,10^01=11)。把结果的40个比特位分作20个组,那么就是说最多只有4个组为b01 b10 b11这三个值,其余的全部为b00。那么比较算法就很好写了。可以对每个字节(四个组)获取其中有几个组是为三个非零值的,来简介获取整体的比较结果。因为每个字节只有256种可能的值,而符合条件的值只有3^4=81种,所以可以先将结果存储起来,然后进行获取。这里给出一个函数,来获取结果中有几个是非零组。/*****************下面table中值的生成******//**
int i;
for( i=0;i<256;++i){
int t =0;
t += (i&0x01 || i&0x02)?1:0;
t += (i&0x04 || i&0x08)?1:0;
t += (i&0x10 || i&0x20)?1:0;
t += (i&0x40 || i&0x80)?1:0;
printf("%d,",t);
if(i%10 ==9){putchar('\n');}
}
********************************************//
int table[] = {0,1,1,1,1,2,2,2,1,2,2,2,1,2,2,2,1,2,2,2,2,3,3,3,2,3,3,3,2,3,3,3,1,2,2,2,2,3,3,3,2,3,3,3,2,3,3,3,1,2,2,2,2,3,3,3,2,3,3,3,2,3,3,3,1,2,2,2,2,3,3,3,2,3,3,3,2,3,3,3,2,3,3,3,3,4,4,4,3,4,4,4,3,4,4,4,2,3,3,3,3,4,4,4,3,4,4,4,3,4,4,4,2,3,3,3,3,4,4,4,3,4,4,4,3,4,4,4,1,2,2,2,2,3,3,3,2,3,3,3,2,3,3,3,2,3,3,3,3,4,4,4,3,4,4,4,3,4,4,4,2,3,3,3,3,4,4,4,3,4,4,4,3,4,4,4,2,3,3,3,3,4,4,4,3,4,4,4,3,4,4,4,1,2,2,2,2,3,3,3,2,3,3,3,2,3,3,3,2,3,3,3,3,4,4,4,3,4,4,4,3,4,4,4,2,3,3,3,3,4,4,4,3,4,4,4,3,4,4,4,2,3,3,3,3,4,4,4,3,4,4,4,3,4,4,4};
int getCount(uint64_t cmpresult)
{
uint8_t* pb = &cmpresult; // 这里假设是小段模式,且之前比较结果是存在低40位
return table[pb[0]]+table[pb[1]]+table[pb[2]]+table[pb[3]]+table[pb[4]];
}