查找查询的相互最佳匹配

Question

我有一个 blast 输出并且想要获得相互最佳命中，即通过确保查询是参考的最佳命中来找到最佳命中，反之亦然基于第 3 列和第 11 列的值（对于阈值第 3 列 >40 和第 11 列 < 2e-04。

MnCG00500.1 Ma.15G248500.1  27.78   180 60  5   9   188 1   110 6e-04   41.2
MnCG00510.1 Ma.15G003800.1  87.88   33  4   0   6   38  3   35  2e-13   60.1
MnCG00510.1 Ma.10G208900.1  84.85   33  5   0   6   38  3   35  2e-12   57.0
MnCG00510.1 Ma.17G173700.1  84.85   33  5   0   6   38  3   35  9e-12   55.5
MnCG00280.1 Ma.11G114700.1  97.17   106 3   0   22  127 22  127 8e-65   210
MnCG00280.1 Ma.05G074900.1  98.18   55  1   0   420 474 11  65  3e-29   111
MnCG00280.1 Ma.20G242300.1  80.36   56  11  0   419 474 95  150 3e-22   95.1
MnCG00890.1 Ma.05G094500.1  89.55   67  7   0   321 387 4   70  1e-34   125
MnCG00890.1 Ma.01G201500.1  91.07   56  5   0   332 387 1   56  3e-28   107

我尝试单独使用基于第 12 列的类似衬垫，并尝试修改以适合我的条件

 awk '{
       a[]="0";b[]="";c[]="0";d[]="";
       if (e[,]==0) 
         e[,]=; 
       else {
         score=e[,]+; 
         e[,]=score
       }
     }
     END{
       for (i in a) 
         for (j in e) {
           split(j,f,SUBSEP); 
           if (f[1]==i && e[j]>a[i]) {
             a[i]=e[j];b[i]=f[2]
           }
         }; 
         for (i in c) for (j in e) {
           split(j,f,SUBSEP); 
           if (f[2]==i && e[j]>c[i]) {
            c[i]=e[j];d[i]=f[1]
         }
       }; 
       for (i in b) 
         if (b[i] in d && d[b[i]]==i) 
           print i"\t"b[i]"\t"a[i]"\t"c[b[i]]
     }' result.blast

但是不起作用。

最后我需要类似的东西

MnCG00500.1, Ma.15G248500.1, MnCG00500.1, no_match  
MnCG00510.1, Ma.15G003800.1, match  
MnCG00510.1, Ma.10G208900.1, match  
MnCG00510.1, Ma.17G173700.1, match  
MnCG00280.1, Ma.11G114700.1, match  
MnCG00280.1, Ma.05G074900.1, match  
MnCG00280.1, Ma.20G242300.1, match  
MnCG00890.1, Ma.05G094500.1,MnCG00890.1, no match   
MnCG00890.1, Ma.01G201500.1, match

Answer 1

您不会从单次爆炸的输出中获得 RBH 运行。收集来自目标数据库的最佳得分命中，然后返回查询数据库。如果第二次爆炸的最佳命中与第一次爆炸的查询序列匹配，则查询和目标序列为 "reciprocal best hits"

你最多只能用一个 blast 输出来从目标数据库中获取 id 以对查询数据库进行 blast 返回。

希望对您有所帮助。

E(xpect) 值越小越好，因此您需要 column_11 < 2e-04

Answer 2

只有输出文件中的最佳匹配：

blastx -query myFile.fasta -db myDB -evalue 0.00001 -num_alignments 1 -outfmt "7 qseqid sseqid"  | uniq > myResult.txt

查找查询的相互最佳匹配

find reciprocal best hits for a query

awk

grep

sed

biopython

blast