在转录本分析中,tophat2是近年来最常用的工具之一,可以将RNAseq的reads数据map到基因组上,找出外显子之间的结合位点。tophat2的下游软件是cufflinks,通过cufflinks我们就可以找到不同样本间的差异基因,进而进行后续的GO,KEGG分析等,所以说tophat2是转录本分析中的非常重要的一个中间工具。
tophat2的安装需要依赖三个软件:boost,samtools和Bowtie2。这三个软件安装完成后即可安装tophat2,命令如下:
cd tophat
./configure
make
make install
下面来为大家讲解写tophat2有参比对的用法:
tophat –p 20 –o output –G genome.gtf genome reads_1.fa reads_2.fa
-p 为安排运行tophat所需要的CPU线程数,根据服务器端的CPU总线程数来决定。
-o 为文件输出路径,一般建议每个样单独建立一个文件夹。
-G 后面跟着相应参考基因组的注释文件,在运行时会首先被tophat2调用bowtie2建立index,这个过程会占用一定的时间。
genome,样品物种的基因组index,可用bowtie2-build对物种的genome.fa文件进行建库,命令为:bowtie2-build genome.fa genome ,此命令建立index也需要花费较长的时间。
tophat既支持pair-end数据,也支持single-read数据。
此外,tophat的常用参数还有:
-r 成对的reads之间的平均inner距离,默认为50。
--library-type Tophat2处理的reads具有链特异性。一般Illumina数据的默认library-type为fr-unstranded。
输出的结果主要有以下几部分:
accepted_hits.bam reads排序的结果以bam格式生成文件,是后面cufflinks软件的输入文件。
junctions.bed deletions.bed insertions.bed 这三个分别是tophat处理的junctions,删除和插入的结果。
unmapped.bam 没有map上的序列
align_summary.txt 可以查看map上的reads所占的比例。
此类问题常见于第一次使用tophat时,部分tophat所调用的软件没有配置好。如上,将bowtie2-build在/usr/local/bin目录下建立软连接即可解决。
此问题的初心一般为注释文件.gtf和基因组的信息在染色体上格式的不匹配造成的,例如本例中的注释文件.gtf的1号染色体标识为chr1,而基因组.fa的1号染色体标识为1,导致两个文件相互无法识别。解决的办法是,将注释文件.gtf的chr1转换成1,操作代码为:sed 's/^chr//g' piggenome.gtf > piggenome1.gtf 。
关于tophat2程序后台运行的问题
可以使用nohup 命令 & 运行tophat2,这样可以后台同时挂多个tophat2程序,但是-p的值要根据服务器端的线程数进行设置,如总服务器CPU线程为56个,后台挂了6个tophat2程序,则单个-p值不应该超过9。可以用more nohup.out查看运行进度。
解决每次运行tophat2都要对注释文件建index的问题
每次用tophat2都要先进行注释文件建库会占用太多时间,所以我们可以提前建好index一劳永逸,使用命令:tophat –G XXX_genes.gtf --transcriptome-index=transcriptome_data/XXX genome(基因组index) 建立注释文件.gtf的index,对应的,使用的tophat命令改成:tophat –p 20 –o output --transcriptome-index=transcriptome_data/XXX genome reads_1.fa reads_2.fa