之前寫過一篇如何使用blast+套件進行本地blast庫的創(chuàng)建及比對袱饭,今天跟大家聊聊比對結(jié)果的輸出格式敲才。
比對命令
blastn -db test -query test.fa -outfmt 0 -out test.o0
通過outfmt參數(shù)指定輸出格式晦墙,官方提供的輸出格式是19種,以下是具體的介紹歧斟。
第 1 種:outfmt 0
軟件默認的輸出格式纯丸,與網(wǎng)頁版展示的比對結(jié)果類似。包括比對統(tǒng)計結(jié)果静袖、比對序列等信息觉鼻,示例如下。
Query= test1
Length=50
Score
Sequences producing significant alignments: (Bits)
23 dna:primary_assembly primary_assembly:ARS-UCD1.2:23:1:524986... 93.5
4 dna:primary_assembly primary_assembly:ARS-UCD1.2:4:1:12000060... 75.0
X dna:primary_assembly primary_assembly:ARS-UCD1.2:X:1:13900914... 54.7
> 23 dna:primary_assembly primary_assembly:ARS-UCD1.2:23:1:52498615:1
REF
Length=52498615
Score = 93.5 bits (50), Expect = 5e-18
Identities = 50/50 (100%), Gaps = 0/50 (0%)
Strand=Plus/Plus
Query 1 GAAACAGACTCACAGACGTAGAAAACAGACATGTGGCTGCCAAGGTGGTG 50
||||||||||||||||||||||||||||||||||||||||||||||||||
Sbjct 24604334 GAAACAGACTCACAGACGTAGAAAACAGACATGTGGCTGCCAAGGTGGTG 24604383
從上到下每條序列結(jié)果依次展示队橙,首先是全部比對結(jié)果的信息坠陈,只顯示序列信息以及比對得分score,接下來是每部分的詳細信息捐康,包括evalue仇矾、gap、正負鏈以及比對序列等信息解总。
第 2 種:outfmt 1
Query= test1
Length=50
Score E
Sequences producing significant alignments: (Bits) Value
23 dna:primary_assembly primary_assembly:ARS-UCD1.2:23:1:524986... 93.5 5e-18
4 dna:primary_assembly primary_assembly:ARS-UCD1.2:4:1:12000060... 75.0 2e-12
X dna:primary_assembly primary_assembly:ARS-UCD1.2:X:1:13900914... 54.7 2e-06
Query_1 1 GAAACAGACTCACAGACGTAGAAAACAGACATGTGGCTGCCAAGGTGGTG 50
22 24604334 .................................................. 24604383
3 4383526 .............................T............. 4383568
29 108133828 ............................. 108133856
這種格式只保留了序列比對信息贮匕,給出輸入序列,庫中一致的用.表示花枫,不一致的會顯示不一致的堿基刻盐,可讀性更強。
第 3 種: outfmt 2
與第二種一致劳翰,只是序列展示有區(qū)別敦锌。
Query= test1
Length=50
Score E
Sequences producing significant alignments: (Bits) Value
23 dna:primary_assembly primary_assembly:ARS-UCD1.2:23:1:524986... 93.5 5e-18
4 dna:primary_assembly primary_assembly:ARS-UCD1.2:4:1:12000060... 75.0 2e-12
X dna:primary_assembly primary_assembly:ARS-UCD1.2:X:1:13900914... 54.7 2e-06
Query_1 1 GAAACAGACTCACAGACGTAGAAAACAGACATGTGGCTGCCAAGGTGGTG 50
22 24604334 GAAACAGACTCACAGACGTAGAAAACAGACATGTGGCTGCCAAGGTGGTG 24604383
3 4383526 AAACAGACTCACAGACGTAGAAAACAGACTTGTGGCTGCCAAG 4383568
29 108133828 GAAACAGACTCACAGACGTAGAAAACAGA 108133856
所有比對的序列都保留堿基序列格式。
第 4-5 種:outfmt 3/4
兩種輸出格式與第二佳簸、三種格式完全一致供屉。
第 6 種:outfmt 5
這種格式為xml格式的文件,很多種編程語言都可以直接解析溺蕉。
<?xml version="1.0"?>
<!DOCTYPE BlastOutput PUBLIC "-//NCBI//NCBI BlastOutput/EN" "http://www.ncbi.nlm.nih.gov/dtd/NCBI
<BlastOutput>
<BlastOutput_program>blastn</BlastOutput_program>
<BlastOutput_version>BLASTN 2.4.0+</BlastOutput_version>
<BlastOutput_reference>Zheng Zhang, Scott Schwartz, Lukas Wagner, and Webb Miller (2000), "
<BlastOutput_db>test</BlastOutput_db>
<BlastOutput_query-ID>Query_1</BlastOutput_query-ID>
<BlastOutput_query-def>test1</BlastOutput_query-def>
<BlastOutput_query-len>50</BlastOutput_query-len>
<BlastOutput_param>
<Parameters>
<Parameters_expect>10</Parameters_expect>
<Parameters_sc-match>1</Parameters_sc-match>
<Parameters_sc-mismatch>-2</Parameters_sc-mismatch>
<Parameters_gap-open>0</Parameters_gap-open>
<Parameters_gap-extend>0</Parameters_gap-extend>
<Parameters_filter>L;m;</Parameters_filter>
</Parameters>
</BlastOutput_param>
<BlastOutput_iterations>
<Iteration>
<Iteration_iter-num>1</Iteration_iter-num>
<Iteration_query-ID>Query_1</Iteration_query-ID>
<Iteration_query-def>test1</Iteration_query-def>
<Iteration_query-len>50</Iteration_query-len>
<Iteration_hits>
<Hit>
<Hit_num>1</Hit_num>
每一個輸入序列及比對結(jié)果都會以嵌套格式進行整理伶丐。
第 7 種:outfmt 6
test1 23 100.000 50 0 0 1 50 24604334 24604383 4.81e-18 93.5
test1 4 97.674 43 1 0 2 44 4383526 4383568 1.74e-12 75.0
test1 X 100.000 29 0 0 1 29 108133828 108133856 2.27e-06 54.7
test2 6 100.000 50 0 0 1 50 85451248 85451297 4.81e-18 93.5
test3 1 100.000 50 0 0 1 50 56773372 56773323 4.81e-18 93.5
test3 1 100.000 50 0 0 1 50 56777851 56777900 4.81e-18 93.5
test3 1 100.000 50 0 0 1 50 56790245 56790294 4.81e-18 93.5
這種格式會輸出比較多的統(tǒng)計信息,包括比對位置疯特、得分哗魂、evalue等,可讀性很強漓雅,一般都會選擇輸出這種格式录别。
第 8 種:outfmt 7
這種格式基本與第7種一致朽色,只不過每條序列加了說明信息,便于查看結(jié)果组题,相當于加了一個表頭葫男。
# BLASTN 2.4.0+
# Query: test1
# Database: test
# Fields: query id, subject id, % identity, alignment length, mismatches, gap opens, q. start, q. end, s. start, s. end, evalue, bit score
# 3 hits found
test1 23 100.000 50 0 0 1 50 24604334 24604383 4.81e-18 93.5
test1 4 97.674 43 1 0 2 44 4383526 4383568 1.74e-12 75.0
test1 X 100.000 29 0 0 1 29 108133828 108133856 2.27e-06 54.7
如果只是一條序列比對,建議使用這種格式崔列。
第 9-10 種:outfmt 8/9
兩種格式的編碼方式是ASN.1
data align {
{
type partial,
dim 2,
score {
{
id str "score",
value int 50
},
{
id str "blast_score",
value int 50
},
{
id str "e_value",
value real { 480860855605031, 10, -32 }
},
{
id str "bit_score",
value real { 934527768506114, 10, -13 }
},
{
id str "num_ident",
value int 50
},
{
id str "hsp_percent_coverage",
value real { 1, 10, 2 }
}
},
格式與xml格式類似梢褐,都是采用嵌套的方式進行展示,不過這種格式我沒用過赵讯,有興趣可以嘗試一下盈咳。
第 11 種:outfmt 10
test1,23,100.000,50,0,0,1,50,24604334,24604383,4.81e-18,93.5
test1,4,97.674,43,1,0,2,44,4383526,4383568,1.74e-12,75.0
test1,X,100.000,29,0,0,1,29,108133828,108133856,2.27e-06,54.7
test2,6,100.000,50,0,0,1,50,85451248,85451297,4.81e-18,93.5
test3,1,100.000,50,0,0,1,50,56773372,56773323,4.81e-18,93.5
這種格式的數(shù)據(jù)與第七、八種基本一致边翼,只不過分隔符由tab鍵換成了逗號鱼响。
第 12 種:outfmt 11
這種格式也是一種ASN.1格式的文件,只不過是用blast編碼對數(shù)據(jù)進行了處理组底。
seq {
id {
local str "Query_3"
},
descr {
user {
type str "CFastaReader",
data {
{
label str "DefLine",
data str ">test3"
}
}
},
title "test3"
},
inst {
repr raw,
mol na,
length 50,
seq-data ncbi2na '7F34C88518E691FDB5D1C4AF90'H
}
},
序列信息使用一種叫NCBI2na
的方法進行重新編碼丈积,具體可參照https://ncbi.github.io/cxx-toolkit/pages/ch_datamod
。
其他格式 outfmt 12-18
多數(shù)使用json或者xml格式的文件進行編碼债鸡,具體包括以下幾個江滨。
12 = Seqalign (JSON),
13 = Multiple-file BLAST JSON,
14 = Multiple-file BLAST XML2,
15 = Single-file BLAST JSON,
16 = Single-file BLAST XML2,
17 = Sequence Alignment/Map (SAM),
18 = Organism Report
上述格式中,如果輸入文件中存在多條序列娘锁,標注Multiple的會單獨輸出每條序列的比對結(jié)果。其中17輸出sam格式的文件饺鹃,18一般的比對時不能使用莫秆。
# outfmt 17
@HD VN:1.2 SO:coordinate GO:reference
@SQ SN:Query_1 LN:50
@SQ SN:Query_2 LN:50
@SQ SN:Query_3 LN:50
@SQ SN:Query_4 LN:50
@SQ SN:Query_5 LN:50
@PG ID:0 VN:2.4.0+ CL:blastn -db test -query test.fa -outfmt 17 -out test.o17 PN:blastn
BL_ORD_ID:22 0 Query_1 1 255 24604333H50M27894232H * 0 0 AS:i:50 EV:f:4.80861e-18 NM:i:0 PI:f:100.00 BS:f:93.4528
BL_ORD_ID:29 0 Query_1 1 255 108133827H29M30875288H * 0 0 AS:i:29 EV:f:2.26911e-06 NM:i:0 PI:f:100.00 BS:f:54.6731
BL_ORD_ID:3 0 Query_1 2 255 4383525H43M115617033H * 0 0 AS:i:40 EV:f:1.74176e-12 NM:i:0 PI:f:97.67 BS:f:74.9863
BL_ORD_ID:5 0 Query_2 1 255 85451247H50M32355043H * 0 0 AS:i:50 EV:f:4.80861e-18 NM:i:0 PI:f:100.00 BS:f:93.4528
BL_ORD_ID:0 16 Query_3 1 255 101760738H50M56773322H * 0 0 AS:i:50 EV:f:4.80861e-18 NM:i:0 PI:f:100.00 BS:f:93.4528
特殊使用
上述格式中outfmt為6、7悔详、10三種格式的文件在輸出時還可以指定輸出特定數(shù)據(jù)镊屎,具體代碼如下。
blastn -db test -query test.fa \
-outfmt "6 qseqid sseqid pident evalue bitscore qseq sseq" -out test.o6.target
以上代碼輸出了輸入及數(shù)據(jù)庫中的序列名稱茄螃、evalue以及兩條序列的堿基序列等信息缝驳。
test1 23 100.000 4.81e-18 93.5 GAAACAGACTCACAGACGTAGAAAACAGACATGTGGCTGCCAAGGTGGTG GAAACAGACTCACAGACGTAGAAAACAGACATGTGGCTGCCAAGGTGGTG
test1 4 97.674 1.74e-12 75.0 AAACAGACTCACAGACGTAGAAAACAGACATGTGGCTGCCAAG AAACAGACTCACAGACGTAGAAAACAGACTTGTGGCTGCCAAG
test1 X 100.000 2.27e-06 54.7 GAAACAGACTCACAGACGTAGAAAACAGA GAAACAGACTCACAGACGTAGAAAACAGA
test2 6 100.000 4.81e-18 93.5 CCACCACAGGGGTTTGAGTAAGAGGAGGGATGTTTTGTGGGAGGCTGTTA CCACCACAGGGGTTTGAGTAAGAGGAGGGATGTTTTGTGGGAGGCTGTTA
test3 1 100.000 4.81e-18 93.5 CTTTATCATAGAGACCACGATGCGGCACTTTCGTCCTCACTACAGGTTGC CTTTATCATAGAGACCACGATGCGGCACTTTCGTCCTCACTACAGGTTGC
test3 1 100.000 4.81e-18 93.5 CTTTATCATAGAGACCACGATGCGGCACTTTCGTCCTCACTACAGGTTGC CTTTATCATAGAGACCACGATGCGGCACTTTCGTCCTCACTACAGGTTGC
test3 1 100.000 4.81e-18 93.5 CTTTATCATAGAGACCACGATGCGGCACTTTCGTCCTCACTACAGGTTGC CTTTATCATAGAGACCACGATGCGGCACTTTCGTCCTCACTACAGGTTGC
指定時可以設定包括qseqid在內(nèi)的共40種統(tǒng)計結(jié)果及相關信息,具體見參考文獻[1]归苍。