NGS的变异检测高效工具Sentieon使用体验

安装

以目前的20180806为例,下载并解压缩

wget https://s3.amazonaws.com/sentieon-release/software/sentieon-genomics-201808.06.tar.gz
tar xf sentieon-genomics-201808.06.tar.gz

之后需要将软件的安装位置加入到环境变量PATH中。

此外需要申请一个许可证文件,并且输出一个环境变量SENTIEON_LICENSE=license.lic

Sentieon读取文件的速度是 20-30Mb/每秒每核,因此为了获取极致性能,推荐使用SSD。

我的测试环境是:

  • Linux 3.10.0-862.14.4.el7.x86_64
  • Intel(R) Xeon(R) CPU E7-8890 v4 @ 2.20GHz
  • 普通硬盘

bwa-mem测试

下面的测试基于200M基因组,213011760 条 read,约100X测序,80个线程

一些测试结果:

  1. sentieon bwa建立的索引和原版的bwa索引可以共用
  2. 80个线程下,原版和sentieon版的运行用top看的时候都在8000%左右。

用time统计运行时间,结果如下:

# sentieon bwa比对
123772.31s user 3026.78s system 6601% cpu 32:00.84 total
# sentieon util排序
5775.55s user 659.32s system 319% cpu 33:35.59 total
# bwa比对
208744.12s user 2966.60s system 5942% cpu 59:22.51 total
# samtools sort排序
5457.31s user 377.73s system 84% cpu 1:55:41.84 total

最终输出的BAM,我用samtools idxstats比较结果时,两者输出一摸一样。

结论:是原版的bwa-mem的2倍以上速度。虽然现在有bwa-mem2了,但是据sentieon公司胡博士说,速度依旧比不上sentieon公司的优化版。同时,他们公司开发的配套排序和标记重复工具的速度也比目前开源软件快。

变异检测分析

单个样本

200M基因组,215050865条read,也是100X左右的深度,用了40个线程

标记重复: 49s + 139s

变异检测: 2402s

BSA项目

以我手头的几个BSA项目为例。

为参考序列(reference.fasta)建立索引文件

mkdir index && cd index
sentieon bwa index reference.fasta
samtools faidx reference.fasta
gatk CreateSequenceDictionary -R reference.fasta

然后创建输入文件

ls 0-raw-data | cut -d '_' -f 1 | uniq > samples.txt

编写批量运行的脚本

#!/bin/bash
# runing with bash align.sh
set -e
set -u
set -o pipefail

SAMPLES=samples.txt
REFERENCE=index/reference.fa
NUMBER_THREADS=80
SENTIEON_LICENSE=license.lic

export SENTIEON_LICENSE
export PATH=$PATH:/opt/biosoft/sentieon/sentieon-genomics-201808.06/bin

# read align
# \\t rather than \t
exec 0< $SAMPLES
while read sample;
do
(sentieon bwa mem -M -R "@RG\\tID:${sample}\\tSM:${sample}\\tPL:ILLUMINA" \
    -t $NUMBER_THREADS $REFERENCE 0-raw-data/${sample}_1.fq.gz 0-raw-data/${sample}_2.fq.gz  || echo -n 'error' ) \
    | sentieon util sort -r $REFERENCE -o 1-align/${sample}_sort.bam -t $NUMBER_THREADS --sam2bam -i - 

done

# Calculate data metrics
mkdir -p metrics
for BAM in 1-align/*_sort.bam;
do
    prefix=$(basename $BAM)
    # get metrics
    sentieon driver -t $NUMBER_THREADS -r $REFERENCE -i $BAM \
    --algo GCBias --summary metrics/${prefix}_GC_SUMMARY.txt metrics/${prefix}_GC_METRIC.txt \
    --algo MeanQualityByCycle metrics/${prefix}_MQ_METRIC.txt \
    --algo QualDistribution metrics/${prefix}_QD_METRIC.txt \
    --algo InsertSizeMetricAlgo metrics/${prefix}_IS_METRIC.txt \
    --algo AlignmentStat metrics/${prefix}_ALN_METRIC.txt
    # plotting
    sentieon plot QualDistribution -o metrics/${prefix}_QD_METRIC.pdf metrics/${prefix}_QD_METRIC.txt
    sentieon plot InsertSizeMetricAlgo -o metrics/${prefix}_IS_METRIC_PDF metrics/${prefix}_IS_METRIC.txt
done

# Remove duplicateds
for BAM in 1-align/*.bam;
do
    OUT=${BAM%_.*}
    sentieon driver -t $NUMBER_THREADS -i $BAM\
    --algo LocusCollector --fun score_info ${OUT}_SCORE.gz
    sentieon driver -t $NUMBER_THREADS -i $BAM \
    --algo Dedup --rmdup --score_info ${OUT}_SCORE.gz \
    --metrics ${OUT}_DEDUP_METRIC.txt ${OUT}_dup.bam
done

# Variant Calling
mkdir -p 3-variants
sentieon driver -t $NUMBER_THREADS -r $REFERENCE \
    $(for bam in 1-align/*_dup.bam; do echo "-i $bam"; done)\
    --algo Haplotyper 3-variants/final.vcf.gz

对于660M的两个亲本和两个子代重测序结果,换句话说就是要对4个样本进行变异建立。

时间记录:

  • 比对+排序平均耗时: 40min(2390s,2463s,2000s,1808s)
  • 标记重复平均耗时:6.4min (123 + 338, 63 + 193, 73 + 201, 74+ 272)
  • 变异检测耗时: 2h (7252s)

此外我还用96个线程测试了基因组大小200M的5个样本的变异检测,运行时间是71分钟(4309)

200个群体GBS分析

我还分别测试了200个GBS样本和157个GBS样本从BAM文件到最终的VCF文件的时间,分别用的是40个线程和56个线程。

运行的脚本代码如下

#!/bin/bash
# runing with bash align.sh
set -e
set -u
set -o pipefail

REFERENCE=index/reference.fa
NUMBER_THREADS=40
SENTIEON_LICENSE=license.lic

# Variant Calling
mkdir -p 3-variants
sentieon driver -t $NUMBER_THREADS -r $REFERENCE \
    $(for bam in 1-align/*_dup.bam; do echo "-i $bam"; done)\
    --algo Haplotyper 3-variants/final.vcf.gz

对于56线程的157个GBS样本,时间约为26小时(94141 s)

对于40线程的200个GBS样本,时间约为19小时(71528 s)

差不多1天就解决了我原本需要一周时间运行的任务。

总结

GATK里用到的模型使用C语言实现而非Java,会提高最终的运行效率,其实我不怎么惊讶。但是同样是通过C语言编写的bwa,Sentieon的开发者居然还能让bwa速度提高了2倍以上,速度还优于目前bwa-mem2,这就体现出他们团队的专业了。

前段时间和Sentieon公司的胡晋南博士交流后还了解到几件事情:Sentieon公司其实比broad研究所更熟悉GATK,因为broad研究所维护GATK软件的人流动性大,有段时间所有人都到谷歌公司去上班了。这也是为啥我用GATK那么难受的原因,很多命令行参数随着版本更新不知不觉就变了;GATK4.0为了追求效率,做了很多取舍,因此反而不如3.8的准确性高;GATK速度慢并且准的原因是因为有一步local assembly,因此只对SNP感兴趣的话,我个人认为可以不用GATK。还有GATK其实还有一些BUG,Sentieon公司为了保证一致性,在其中一个版本中保留了所有的bug,此外还开发另一个优化版本。

当然Sentieon是一个收费软件,对于公司或者一个大的研究机构而言,如果追求效率那么可以考虑一下。但就针对我而言,速度可能并不是我的刚需,手头有那么多项目,我可以切换个项目。不过Sentieon还支持广州超算那边的按核时计费模式,如果以后能够按样本付费的话,对于一些小的课题组会更加好一点吧。

©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 155,348评论 4 355
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 65,906评论 1 283
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 105,191评论 0 235
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 43,295评论 0 201
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 51,606评论 3 283
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 40,107评论 1 202
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 31,507评论 2 305
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 30,217评论 0 193
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 33,848评论 1 234
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 30,183评论 2 238
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 31,734评论 1 255
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 28,084评论 2 247
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 32,608评论 3 228
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 25,921评论 0 8
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 26,627评论 0 190
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 35,103评论 2 261
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 35,054评论 2 257

推荐阅读更多精彩内容