SURVIVOR 是一套用于模擬/評估 SV、合并和比較樣本內(nèi)及樣本間 SV 的工具集,其還可對SV進(jìn)行重新格式化及進(jìn)行統(tǒng)計(jì)總結(jié),其他詳細(xì)信息也可參考工具的wiki。
安裝
git clone https://github.com/fritzsedlazeck/SURVIVOR.git
cd SURVIVOR/Debug
make
快速使用
Program: SURVIVOR (Tools for Structural Variations in the VCF format)
Version: 1.0.7
Usage: SURVIVOR <command> [options]
Commands:
-- Simulation/ Evaluation
simSV 根據(jù)一個(gè)參考基因組模擬SV和SNP.
scanreads 從 mapped reads 中獲取錯(cuò)誤信息用于模擬.
simreads 模擬長 reads (Pacio or ONT).
eval 評估從模擬數(shù)據(jù)中鑒定的SV.
-- Comparison/filtering
merge 比較或合并VCF以生成一致或多個(gè)樣本的VCF.
filter 對VCF進(jìn)行特定SV大小或區(qū)域的過濾
stats 統(tǒng)計(jì)VCF中的多個(gè)信息
compMUMMer 使用MUMMer (Show-diff)找到的斷點(diǎn)來注釋VCF.
-- Conversion
bincov 將bins coverage vector 轉(zhuǎn)為bed文件來過濾低MQ區(qū)域的SV
vcftobed 將VCF格式轉(zhuǎn)為bed格式
bedtovcf 將bed格式轉(zhuǎn)為VCF格式
smaptovcf 將smap文件轉(zhuǎn)為VCF文件 (beta version)
bedpetovcf 將bedpe文件轉(zhuǎn)為VCF文件(beta version)
hapcuttovcf 使用原始提供給Hapcut2的SNP文件將Hapcut2最終結(jié)果轉(zhuǎn)為VCF文件
convertAssemblytics 將Assemblytics結(jié)果轉(zhuǎn)換為VCF文件
常用功能
這里介紹下SURVIVOR工具包中的幾個(gè)常用功能,其他功能及詳細(xì)參數(shù)介紹參見Wiki的Methods & Parmater部分。
模擬和評估SV
關(guān)于SV模擬的詳細(xì)參數(shù)介紹參見:Methods & Parmater。
首先,需要生成參數(shù)配置文件:
SURVIVOR simSV parameter_file
注:可以更改參數(shù)配置文件中的具體參數(shù),但不要更改每行的位置
然后,基于參考序列進(jìn)行SV模擬,參考序列需為fasta格式:
SURVIVOR simSV ref.fa parameter_file 0.1 0 simulated
該步驟使用的read模擬工具為Mason,也可以選擇其他模擬工具。
下一步,對模擬的SV進(jìn)行評估:
SURVIVOR eval caller.vcf simulated.bed 10 eval_res
該步驟會使用之前的模擬數(shù)據(jù)集對caller.vcf進(jìn)行評估,并允許斷點(diǎn)上下游有10bp的差異,評估結(jié)果保存在eval_res中。
最后會輸出一行評估的結(jié)果:
Overall: 20 11/0/0/0/9 0/0/0/0/0 0/0/0/0/0 1 0
這表明一共模擬了20個(gè)SV,然后的三組數(shù)據(jù)分別是真陽性結(jié)果(即模擬了并找到的),假陰性結(jié)果(模擬了但沒找到),假陽性結(jié)果(未模擬但找到了)。每組中的五個(gè)數(shù)字分別代表五種SV類型(DEL/DUP/INV/TRA/INS)。倒數(shù)第二個(gè)表示靈敏度,最后一個(gè)數(shù)字表示FDR概率??傊?,在這次模擬中,總共模擬了20個(gè)SV,從中再次發(fā)現(xiàn)11個(gè)缺失和9個(gè)插入,并且沒有丟失或報(bào)告其他的任何SV。
多樣本間SV的比較
為了降低一些SV鑒定工具的假陽性,同時(shí)保持較高的靈敏度,一般的建議是一個(gè)樣本使用多個(gè)鑒定工具。例如,可以在短reads數(shù)據(jù)上運(yùn)行Manta、GRIDSS和Lumpy。每個(gè)工具都會生成一個(gè)VCF文件作為結(jié)果。
這時(shí)SURVIVOR就得安排上了,它可比較這些VCF文件并生成一致性結(jié)果。注意,通常,SURVIVOR可以合并任何技術(shù)或其他工具生成的SV VCF文件。
首先, 將VCF文件名輸出到一個(gè)文件中:
ls *vcf > sample_files
然后, 生成合并的VCF結(jié)果:
SURVIVOR merge sample_files 1000 2 1 1 0 30 sample_merged.vcf
1000表示允許合并的SV間的距離最大為1000bp;
2表示僅輸出2個(gè)工具均鑒定出的SV;
1表示僅輸出2個(gè)工具鑒定出的同類型的SV;
1表示僅輸出2個(gè)工具鑒定出的同方向的SV;
30表示僅考慮長度在30bp以上的SV
識別SV鑒定錯(cuò)誤的區(qū)域
首先,從排好序的bam文件中提取比對質(zhì)量較差的reads:
samtools view -H our.sort.bam > lowMQ.sam
samtools view our.sort.bam | awk '$5<5 {print $0}' >> lowMQ.sam
samtools view -S -b -h lowMQ.sam > lowMQ.bam
現(xiàn)在就得到了MQ < 5 的reads的bam文件, 然后, 需要計(jì)算堿基覆蓋度:
samtools depth lowMQ.bam > lowMQ.cov
下一步, 將覆蓋度文件進(jìn)行聚類生成bed文件, 從而用于SV的過濾:
SURVIVOR bincov lowMQ.cov 10 2 > lowMQ.bed
該步驟允許將距離最大10bp的區(qū)域聚在一起,并且僅在覆蓋度大于2的情況下才考慮該區(qū)域。生成的bed文件可用于過濾SV時(shí)使用.