Triniy的原理及过程
Trinity是由 Broad Institute开发的,用于转录本的de novo拼接,主要由三个软件模块组成:Inchworm, Chrysalis and Butterfly,能处理大型的RNA数据。 摘自https://github.com/trinityrnaseq/trinityrnaseq/wiki:
Trinity partitions the sequence data into many individual de Bruijn graphs, each representing the transcriptional complexity at a given gene or locus, and then processes each graph independently to extract full-length splicing isoforms and to tease apart transcripts derived from paralogous genes. Briefly, the process works like so:
Inchworm assembles the RNA-seq data into the unique sequences of transcripts, often generating full-length transcripts for a dominant isoform, but then reports just the unique portions of alternatively spliced transcripts.(将RNA-seq reads拼接形成一个unique sequences,也就是contigs,这些contigs可能代表一个全长转录本,也可能为一个转录本中的一部分)
Chrysalis clusters the Inchworm contigs into clusters and constructs complete de Bruijn graphs for each cluster. Each cluster represents the full transcriptonal complexity for a given gene (or sets of genes that share sequences in common). Chrysalis then partitions the full read set among these disjoint graphs.(将contigs进行聚类形成de Bruijn graphs,然后将其区分形成一个个可能的转录本)
Butterfly then processes the individual graphs in parallel, tracing the paths that reads and pairs of reads take within the graph, ultimately reporting full-length transcripts for alternatively spliced isoforms, and teasing apart transcripts that corresponds to paralogous genes.
Trinity的安装
下载Trinity:https://github.com/trinityrnaseq/trinityrnaseq/releases
我选择的是2.4.0安装Trinity(先安装bowtie2)
2.1 下载bowties
http://bowtie-bio.sourceforge.net/bowtie2/index.shtml
然后选择版本后,安装并加到环境变量中
2.2 再切换到
trinityrnaseq-Trinity-v2.4.0
目录下,安装trinity及其相关部件(比如samtools、Trimmomatic等软件):make plugins
2.3 如有报错:
fatal error: zlib.h: No such file or directory
sudo apt-get install zlib1g-dev
测试Trinity是否正常安装
cd sample_data/test_Trinity_Assembly/ ./runMe.sh
Trinity的初步使用
查看Trinity相关参数,必须的参数为以下几个:
--seqType : 输入文件的格式,fq or fa
--max_memory : 设置trinity限制的最大使用内存
--left : left reads,多个文件以逗号隔开
--right : right reads,多个文件以逗号隔开
--samples_files : 如果文件过多,可以将信息写入samples文件中,如下:cond_A cond_A_rep1 A_rep1_left.fq A_rep1_right.fq cond_A cond_A_rep2 A_rep2_left.fq A_rep2_right.fq cond_B cond_B_rep1 B_rep1_left.fq B_rep1_right.fq cond_B cond_B_rep2 B_rep2_left.fq B_rep2_right.fq
--SS_lib_type : 如果是链特异性测序,则使用
--CPU : 设置CPU的使用数目(但是并不是trinity的所有程序都是能多线程使用最大的CPU的)还有些参数一般都是默认设定,如:
--min_contig_length : 拼接产生的contig的最小长度,默认200
--long_reads : fasta file containing error-corrected or circular consensus (CCS) pac bio reads
--genome_guided_bam : genome guided mode, provide path to coordinate-sorted bam file
--jaccard_clip : 减少拼接时产生融合转录子的可能。但是这个参数比较适合于真菌物种(由于其基因组比较稠密,转录子经常在UTR区域重叠),脊椎动物以及植物则不需要使用该参数
--no_normalize_reads : 不对reads进行silico normalization,但是默认是开启的(能有效减少trinity运行时间)
--output : 输出结果目录,默认是当前目录下的trinity_out_dir
文件下
--full_cleanup :only retain the Trinity fasta file, rename as ${output_dir}.Trinity.fasta如果要查看trinity的完整的全部参数:
Trinity --show_full_usage_info
从中我们可以看到Trinity其实还分别对Inchworm、Chrysalis以及Butterfly有各自的参数设定,以及一些特殊情况下的适用参数
对于大样本的mRNA的测序数据,现在最新的Trinity已经将对reads进行normalization设置为默认开启状态,所以如果需要关闭的话使用参数
--no_normalize_reads
,我觉得没必要关闭。。毕竟能节省计算资源嘛除了上述参数外,还有一个参数
--min_kmer_cov
min count for K-mers to be assembled by, Inchworm (default: 1),可以理解为k-mer 的coverage如果服务器配置不够,可以将
--min_kmer_cov
设置为2,这样做会去除unique occurring kmer,从而降低内存消耗,并且减少拼接过程中的噪音。但是,这会导致低表达量的转录本不再被拼接(在拼接过程中直接被舍弃),而且也会导致一些低表达量的转录本变得更加碎片化。因此看情况取舍吧此外
--min_glue
min number of reads needed to glue two inchworm contigs together. (default: 2),也就是在inchworm步骤中(在contig聚类成component的时候),两个contigs之间连接部位需要至少有多少个reads比对上。默认是2,但是可以人为设置以减少最后拼接出来的gene数目。对于以上两个参数,Trinity作者对此做了一定的解释,具体可见:https://github.com/trinityrnaseq/trinityrnaseq/issues/92
查看Trinity的FAQ以及别人提交的报错问题(几乎你的遇到的问题都能在里面找到别人提交并回答了的类似的问题)
https://github.com/trinityrnaseq/trinityrnaseq/wiki/Trinity-FAQ
https://github.com/trinityrnaseq/trinityrnaseq/issues
Trinity的拓展使用
通过上述Trinity的初步使用,并且基本按照默认参数运行后,会产生一个Trinity.fasta文件,这也是唯一我们需要后续分析使用的文件
Trinity官网不仅发布了Trinity对于无参转录组的组装使用教程,其还为组装后的接下来的转录本的评估,基因/转录本的定量,差异基因分析以及注释都做了流程化的搭建,对于每一步都有其建议以及相关代码,所有代码都是用perl写的,所以会perl的很容易就看懂哈
下面介绍几个配套流程的代码,有需要的话可以直接使用,免去自己写代码的过程哈
TrinityStats.pl
$TRINITY_HOME/util/TrinityStats.pl Trinity.fasta
了解这个代码前,需要理解一个概念:Trinity.fasta文件中转录本的id命名的规则,比如TRINITY_DN25847_c1_g2_i3,TRINITY_DN25847_c1表示Trinity read cluster(不同reads在拼接中聚类在一起的代号?),g2表示拼接出来可能的gene(2代表前面read cluster所所组装出的第2个可能的gene?),i3表示g2这个gene的第3个可能的转录本
因此这个代码会在基因水平和转录本水平上分别统计组装序列的一些统计信息,如:N50,序列个数、GC含量以及序列长度等
align_and_estimate_abundance.pl
$TRINITY_HOME/util/align_and_estimate_abundance.pl --transcripts Trinity.fasta --seqType fq --left reads_1.fq --right reads_2.fq --est_method RSEM --aln_method bowtie2 --trinity_mode --output_dir rsem_outdir
这个能在基因/转录本水平上进行表达丰度定量,主要调用bowtie/bowtie2比对到组装后的Trinity.fasta序列上,然后用RSEM进行表达丰度估算,最后得到在基因和转录本上的表达丰度定量值。当然比对软件以及定量软件都可以选择其他,比如kallisto,salmon等
run_DE_analysis.pl
$TRINITY_HOME/Analysis/DifferentialExpression/run_DE_analysis.pl \
--matrix counts.matrix --method edgeR --dispersion 0.1
这个主要用于差异分析,比如这个是调用edgeR包来做差异表达分析,trinity套件还提供其他方法,如DESeq2,voom,ROTS等
当然还有其他脚本代码,可以通过查看perl代码来学习其运行原理,然后再自己来写一遍符合自己分析情况的脚本哈
Trinity官网的教程也是非常棒的学习途径
https://github.com/trinityrnaseq/trinityrnaseq/wiki
本文出自于http://www.bioinfo-scrounger.com转载请注明出处