三代測(cè)序因?yàn)殚L(zhǎng)讀長(zhǎng)的特點(diǎn),在基因組組裝方面有明顯的優(yōu)勢(shì)。今天給大家介紹一款基因組組裝軟件,miniasm。該軟件可用于Pacbio和ONT兩個(gè)平臺(tái)的數(shù)據(jù)組裝,分析速度很快。
軟件安裝
- minimap2安裝
miniasm軟件組裝時(shí),需要從比對(duì)結(jié)果中獲取reads之間的相對(duì)位置以及overlap情況。相對(duì)于二代測(cè)序,三代reads較長(zhǎng),bwa等軟件效率太低,這里使用minimap2進(jìn)行比對(duì),速度很快。
$ git clone https://github.com/lh3/minimap2 && (cd minimap2 && make)
- miniasm安裝
組裝軟件在github可以獲取源碼,直接安裝即可
$ git clone https://github.com/lh3/miniasm && (cd miniasm && make)
基因組組裝
數(shù)據(jù)過(guò)濾
三代測(cè)序的原始數(shù)據(jù)錯(cuò)誤率比較高,需要先進(jìn)性數(shù)據(jù)的矯正,不同平臺(tái)有不同的數(shù)據(jù)格式,對(duì)應(yīng)不同的軟件,可以自行搜索進(jìn)行過(guò)濾。
數(shù)據(jù)比對(duì)
使用minimap2軟件對(duì)fq數(shù)據(jù)比對(duì)。因?yàn)橐@得reads之間的overlap關(guān)系,比對(duì)時(shí)只需要進(jìn)行fq自身數(shù)據(jù)間的比對(duì)即可。
$ minimap2 -x ava-pb -t8 pb-reads.fq pb-reads.fq > pb-reads.minimap2.out
其中-x表明數(shù)據(jù)類型及組裝方法,最常用的包含4中,其中map-pb/map-ont表示Pacbio和ONT平臺(tái)數(shù)據(jù)與參考基因組進(jìn)行比對(duì),ava-pb/ava-ont表示Pacibo和ONT平臺(tái)fq數(shù)據(jù)自身數(shù)據(jù)間的比對(duì);-t為線程數(shù)
$ head -5 pb-reads.minimap2.out
read1 204 35 204 - read4 419 144 313 119 169 0 tp:A:S cm:i:1s1:i:119 dv:f:0.0042 rl:i:181
read2 204 35 204 + read5 414 64 232 119 169 0 tp:A:S cm:i:1s1:i:119 dv:f:0.0042 rl:i:181
read3 204 35 204 + read6 231 36 204 119 169 0 tp:A:S cm:i:1s1:i:119 dv:f:0.0021 rl:i:181
以上是minimap2軟件的輸出結(jié)果。前12列為比對(duì)數(shù)據(jù)統(tǒng)計(jì)結(jié)果,各列依次是query名、query長(zhǎng)度、query起始位置、query終止位置、鏈方向、target名、target長(zhǎng)度、target起始位置、target終止位置、比對(duì)時(shí)match堿基數(shù)、比對(duì)堿基數(shù)(包括gap等)、比對(duì)質(zhì)量。后面的各列與bam文件中信息類似,每個(gè)tag對(duì)應(yīng)不同的統(tǒng)計(jì)結(jié)果,具體可查詢官方文檔。
基因組組裝
使用miniasm主程序?qū)?shù)據(jù)進(jìn)行組裝,要求輸入fastq以及minimap2比對(duì)結(jié)果,示例如下。
$ miniasm -1 -2 -m 50 -s 1000 -c 5 -f ../../X5.HB5.extract.fq X5.extract.minimap2.sam > X5.all.extractsam.miniasm.out.new
上述參數(shù)中,-1和-2表示數(shù)據(jù)在組裝時(shí)跳過(guò)對(duì)reads的篩選,直接使用(正常情況下默認(rèn)沒(méi)有這兩個(gè)參數(shù),也沒(méi)有必要特意加上),-m為最小的匹配長(zhǎng)度,-s為最小的span長(zhǎng)度,-c為測(cè)序深度。
對(duì)于基因組組裝而言,最重要的幾個(gè)參數(shù)基本就是這些,需要給定overlap長(zhǎng)度,span長(zhǎng)度以及測(cè)序深度,有的軟件還需要添加輸出最小的contigs長(zhǎng)度等。在組裝之前,應(yīng)該對(duì)自己的數(shù)據(jù)有一個(gè)基本的認(rèn)識(shí),比如測(cè)序數(shù)據(jù)中reads的長(zhǎng)度分布范圍、測(cè)序深度等信息,方便組裝不理想時(shí)調(diào)整參數(shù)。
S utg000002c ACTTAGACCTACCGTTCACCTAATGACT...CTTAGACCTACCGTTCACTAATGACTTG LN:i:125
L utg000002c + utg000002c + 0M
L utg000002c - utg000002c - 0M
a utg000002c 0 read1:64-361 + 124
a utg000002c 124 read2:498-743 - 1
S utg000003c TAAGGCTAATGGCACTCAGTAAC...TAACAATAAGGCTAATAG LN:i:224
L utg000003c + utg000003c + 0M
L utg000003c - utg000003c - 0M
a utg000003c 0 read3:52-294 + 106
a utg000003c 106 read4:341-600 + 118
S utg000004l ACGCCAGTTGCTATGGAGCCATCC...CTCAACGGATAAAAGGGTACTCCACA LN:i:302
a utg000004l 0 read5:738-982 - 12
a utg000004l 12 read6:492-762 + 12
a utg000004l 24 read7:929-1194 + 15
a utg000004l 39 read8:76-337 + 69
a utg000004l 108 read9:796-989 + 194
x utg000002c 125 2
x utg000003c 224 2
x utg000004l 302 5 0 0 read12:738-982 - read15:796-989 -
其中S表示組裝后的序列,包括名字、序列以及序列長(zhǎng)度,a表示最佳組裝路徑,L表示overlap信息,x為各組裝序列的統(tǒng)計(jì)結(jié)果
參考文獻(xiàn)
[1] https://www.ncbi.nlm.nih.gov/pmc/articles/PMC4937194/ PMID: 27153593