大家好,今天来介绍hisat2和bowtie2(hisat2建立索引需要多久)的问题,以下是渲大师小编对此问题的归纳和整理,感兴趣的来一起看看吧!
DNA/RNA序列比对软件整理
文章仅是记录自己的学习使用,有错误请指出,我立刻改正
在对比对工具进行比较时,通常将其分为DNA比对工具(DNA-seq)和RNA比对工具(RNA-seq)。它们的区别在于是否会考虑跨外显子的比对,即:是否会将没有比对上的reads劈开,对劈开后的两部分再次比对)。
随着现在各种seq测序的出现,我们已经不能简单的根据是比对DNA还是RNA来判断。比对工具的选择主要依据reads的比对是否需跨外显子。(PRO-seq/GRO-seq,它们虽然在建库时捕获的RNA,但是它们的比对并不需要考虑跨外显子。)
常用工具:
DNA-seq:BWA;bowtie&bowtie2
RNA-seq:STAR;HISAT2;Tophat&Tophat2
BWA主要应用二代测序后的大量短小片段与参考基因组之间的定位比对。需要先嫌燃对参考序列建建立索引,BWA也是基于 BWT和 FM-Index 理论来对参考基因组做索引。根据测序方法的不同,有单末端序列(Single-end,SE)比对和双末端序列(Pair-end,PE)比对。
bowtie出现在测序行业还不成熟的时候,序列长度普遍在50bp以下,bowtie的只满足长度在50bp以下的reads的比对。官方称其可以把短的DNA序列(35bp)快速的比对到人类基因组上。
Bowtie2 是一款经典的短读长序列( 50-100 bp,最多可到1000 bp ) 比对软件,节约内存且灵活与成熟的短序列比对软件,比较适合下一代测序技术。支持单端测序(unpaired) 和双端测序的比对。支持全局比对(end-to-end align ) 和 局部比对( local align )。其通常使用全文分索引(FM-index)以及Burrows-Wheeler 变换(BWT)索引基因组使得比对非常快速且内存高效,但是这种方法不适合于找到较长的、带缺口的序列比对
结论:bowtie和bowtie2,是两个不同类型的比对工具,bowtie2并非是bowtie的升级。尺有所长寸有所短,bowtie适合长度在50b长度以内的reads比对,而bowtie2适合50-100b,甚至更长的reads比对。但是这两个都属DNA-seq比对工具
RNA-Seq测序的特性,天然的会有一部分数据延伸到内含子区,这部分跨越外显子和内含子的reads就称为『junction reads』,所以RNA-Seq比对软件需要针对此进行优化。
( junction:转录组reads比对不同于基因组reads比对(如ChIP-seq、WES等)的地方在于,比对的reads可能来源于2个被内含子隔开的外显子区域,导致reads一端比对在第一个外显子的后面部分,另一端比对在第二个外显子的前面部分,芹衡虚即跨剪切位点,从而形成exon-exon junction (剪接点)。这些reads又称为junction reads,对转录本的拼接、鉴定和差异分析具有重要的意义。)
(soft-clip事件: 即reads末端存在低质量碱基或接头导致比对不上的, STAR会自动尝试截去未比对部分,只保留比对上的部分。)
STAR是ENCODE皇家御用的RNA-seq比对工具,ENCODE计划(ENCyclopedia Of DNA Elements)又称人类基因组DNA元件百科全书计划,是2003年在人类基因组计划完成之后紧接着的又一个大型国际科研项目。
Tophat2的原作者们也不知道是出于什么考虑,不再更新Tophat2,转而开发了一个新的比对工具HISAT2,更是推荐人们使用HISAT2,声称其速度更快,内存占用率更小,准确率更高。
此外,HISAT2不仅支持RNA-seq的比对还支持DNA-seq比对,唯一需要做的就是加上一个参数–no-spliced-alignment。但是就目前拦散来看,大部分人都是使用HISAT2做RNA-seq,没人使用它做DNA-seq
Tophat/Tophat2工具本身不能进行比对,它是通过调用bowtie/bowtie2进行比对的。划重点,bowtie2不是bowtie的升级版,但是Tophat2是Tophat2的升级版。因此Tophat只可以调用bowtie,而Tophat2不仅可以调用bowtie2(默认)还可以更改设置调用bowtie。
Tophat/Tophat2调用bowtie/bowtie2后,会首先使用bowtie/bowtie2对序列进行比对,对于那些没有比对上的,会考虑其跨外显子的可能性,将reads劈开重新比对。
全长转录组(Full-length transcriptome)是基于PacBio和Nanopore三代测序平台,无需打断拼接,直接获得包含5’UTR、3’UTR、polyA尾的mRNA全长序列及完整结构信息,从而准确分析有参考基因组物种可变剪接及融合基因等结构信息,克服无参考基因组物种转录本拼接较短、信息不完整的难题。同时还可以借助二代测序数据,进行转录本特异性表达分析,获得更加全面的注释信息。
传统的使用比较多的长读长比对软件是GMAP,05年发表公布,最开始是用来比对低通量的est序列的,后来也有进一步升级为GSNAP支持高通量的二代测序。PacBio测序技术出现后,常用于Iso-seq转录本的鉴定,目前仍是相关研究引用量最高的比对软件,该软件也一直在持续更新升级。其可以将转录本序列与参考基因组序列比对,输出gff文件,比对速度稍慢。
Minimap2是生信大牛李恒18年用C语言开发的可以用于三代数据(subreads、iso-seq)比对的长序列比对软件,与传统的三代比对工具GMAP相比,其速度有非常显著的提升,当然同时消耗的内存也比较大。使用方法也比较简单,近几年引用次数增长的也很迅速,所以大家可以试试用minimap2进行Iso-seq的比对。
HISAT2 建立索引警告和比对时报错解决方案
tags: HISAT2 RNA-seq
HISAT2 发表的文章中强调了它的速度很快,我就测试了一下这个工具。
HISAT2 建立索引:
然而没多久就看到这样的警告:
只是警告,并没有报错。
HISAT2 建参考索引很慢,等 HISAT2 建完索引(建索引花了 16 个小时),然后用 HISAT2 比对 RNA-seq 数据测试。
几分钟之后报错了:
github 上有人在 HISAT2 项目中报告过这个错误,虽然没有最终讨论出解决办法,但是都觉得跟建索引不完整有关,或许与建索引时候的警告有关。我查了一些资料,综合 biostar 和 SEQanswer 中的讨论,建立索引时遇到的警告是由于参考序列中存在大段的 n 碱基导致的,例如其中一条 fasta 中 n (我遇到的是小写的 n)太多。解决办法也很简单,过滤掉参考序列中长度小于 50bp 的 contig 和序列中连续 n 碱基超过 40bp 的contig 。然后重新建索引,就没有任何警告了,但是会损失部分参考序列。这样处理之后建的索引再睁唤庆用 HISAT2 比对 RNA-seq 数据,就没有问题了。
HISAT2 比对速度确实很快,一个样本转录组数据比对 2G 的核糖体 RNA 参考基因组约 25 分钟,bowtie2 需要 170 分钟(线程数都是 4)。 bowtie2 的 local 模式比对出 21% 的 rRNA 污染链宴,而 HISAT2 比对 2% 的 rRNA 污染,差异也挺大的……
但是,HISAT2 的线程数不能设置太高,用户手册中建议对人类参考基因组进行比对时,线程数设置在 1-8 之间。我用文献悉握 Transcript-level expression analysis of RNA-seq
experiments with HISAT, StringTie and Ballgown 中的数据测试也发现当线程数超过 12 时整个比对步骤消耗的时间反而会增加。
bowtie2可以更改错配碱基数么
1、Sam文件各标签含义(tophat/hisat2)
•NH:i:: N=1时 为unique。常用于tophat/hisat2产生的sam文件unique read筛选。
•CC:Z: 当为‘=’为map到同一条基因上,一般在map基因组时由于内含子存在而容易出现,他只代表两种不同的方式,计数时应记为1。此处一般为其他基因的名字。CP:i 和HI:i标签为map到第i条基因及起始位置。则游
•YT:Z::代表的含义与bowtie产生的sam也不同。具体还未知!其他标签AS,XN,XM,XO,XG,NM,MD等如下图可以看出都相同。
•对于tophat/hisat2比对产生的sam文件我们可以直接筛选NH标签。
grep‘NH:i:1’ out.sam >unique.sam
2、Sam文件各标签含义(bowtie2)
•AS:i:Alignmentscore.可以为负的,在local下可以为正的。 只有当Align≥1 time才出现
•XS:i:Alignmentscorefor second-best alignment. 当Align>1 time出现
•YS:i:Alignmentscorefor opposite mate in the paired-end alignment. 当该read是双末端测序中的另一条时出现
•XN:i:Thenumber of ambiguous bases in the reference covering this alignment.(推测是指不知道错配发生在哪个位置,推测是针对于和缺失,待查证)
•XM:i:错配碱基的数目
•XO:i:Thenumberof gap opens(针对于比对中的和缺失)
•XG:i:Thenumberof gap extensions(针对于比对中的和缺失延伸数目)
•NM:i:Theeditdistance。(edits:/缺失/替换数目)
•YF:Z:该reads被过滤掉的原因。可能为LN(错配数太多,待查证)、NS(read中包含N或者.)、SC(match bonus低于设定的阈值)、QC(failing quality control,待证)
•YT:Z:值为UU表示不是pair中一部分、CP是pair且可以完美匹配、DP是pair但不能很好的匹配、UP是pair但是无法比对到参考序列上。
•MD:Z:比对上的错配碱基的字符串表示。
由于bowtie2产生的sam文件并没有NH标签帆耐,所以提取uniqueread可能比较麻烦。首先提取“AS”标签表示能比对上的read(>=1 time),然后利用grep反正则表达式过态盯春滤掉XS标签得到我们需要的unique read。
grep “AS:” aligned.sam grep –v“XS:” >unique_alignments.sam
对于双端测序用bowtie2比对筛选unique concordant pair时则需要在上一步的基础上增加如下命令:
grep ‘YT:Z:CP’ unique.sam>pair-end_unique.sam
3、Bwa获取unique
samtoolsviewbwa.bam grep “XT:A:U”
如何批量去除基因组中内含子
批量去除基因组中内含子的方法如下:
1、使用基因注释工具:许多基因注释工具如Glimmer,GeneMark,AUGUSTUS和MAKER等可以用来预测基因和内含子的位置,可以使用这些工具来标记内含子并将其删除。
2、基于比对的方首迟法:使用比对算法,如BLAT,Bowtie2,BWA和HISAT2等,将基因组序列与已注释的基因组序列核汪比对,然后通过比对结果识别和删除内含子序列。
3、基于RNA-Seq数据的方法:通过RNA-Seq数据的拼接转录本来预测内含子和外显子的位置,然后根据这些位置去除内含子序列。
4、基于机器学习的方法:使用机器学习算法来训练模型,预测基因和内含子的位置,然后根据这些位置改芹仔去除内含子序列。
conda环境报错: Solving environment: failed with initial frozen solve
刷了机,安装好conda环境后conda又疯了。
conda install -y samtools hisat2 bowtie2 bwa htseq bowtie
报错开中冲带始
先查询当前conda版本:conda -V
更新conda到最新版本:conda update -n base conda
再查一下conda版本:conda -V
然后更新一下所有:conda update –all
还是一样的报错
但 conda install -c conda-forge opencv=4.2.0可以被顺利执行并完成,
查询发现是清华镜像的问题( https://www.cnblogs.com/YlnChen/p/12710893.html ),立马换中科大镜像, 没删除清华镜像,把中科大的放在最前面:
先少量安装尝试–最好一次放三个,一次性放判仿太多,中间卡住就gg了
conda install -y samtools hisat2 bowtie2
成功!卖芦
或者你可以用pip来安装 https://pypi.org/ ,如: