2016-03-03 64 views
0

我有一個爆炸輸出,並想獲得倒數最佳匹配,即,尋找最好通過確保查詢是一個參考的最佳匹配,反之亦然基於值從第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[$1]="0";b[$1]="";c[$2]="0";d[$2]=""; 
     if (e[$1,$2]==0) 
     e[$1,$2]=$12; 
     else { 
     score=e[$1,$2]+$12; 
     e[$1,$2]=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 
+0

這個條件不滿足任何行'$ 3> 40 && $ 11>中2E-04' – karakfa

回答

1

你不會從單一的爆炸運行的輸出得到RBH。 收集來自目標數據庫的最佳得分命中數,然後對查詢數據庫進行反擊。 如果來自第二次爆炸的最佳匹配與第一次爆炸的查詢序列匹配,則查詢和目標序列是「相互最佳匹配」

只有一次爆炸輸出可以做的最多的是從目標數據庫爆炸回查詢分貝。

希望有幫助。

也較小E(xpect)值是更好的,所以你會想column_11 < 2E-04

+0

或者使用更高bitscore其是更好的表現數(象可能有任何浮點舍入問題)。 – peterjc

0

只有一個輸出文件的最佳匹配:

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