sailfish:不需要比对的转录本定量软件

欢迎关注”生信修炼手册”!

传统的转录组定量分析都是基于比对的结果去定量，根据比对上基因/转录本的reads的数目来确定其表达量，是能够想到的最直接的定量方式。基于这样的策略，有很多经典的工具被开发来执行这样的任务，虽然软件的运行速度和硬件消耗不断被优化，但是整个比对定量这条流程的运行时间还是比较长的。

科研人员又开发出一个新的定量思路，叫做Alignment-free RNA quantification, 也就是不需要要比对，直接根据测序reads的特征来定量。

核心思想是将转录本划分成kmer, 以kmer出现的次数作为转录本丰度的衡量证据。当然实际处理时，会非常的复杂。因为同一个kmer会对应多个转录本，一个转录本会出现相同的kmer等情况，kmer的次数如何合理分配就很关键。

由于不需要比对，相比比对后再定量的思路，其运行速度提高了非常多。目前基于这个思路的定量软件有以下几种

sailfish
salmon
kallisto

本文主要介绍sailfish这个软件, 官网如下

http://www.cs.cmu.edu/~ckingsf/software/sailfish/index.html

sailfish的运行需要两步，第一步对转录本建索引，第二步定量。

1. 对转录本的fasta文件建立索引

将转录本划分成kmer, 记录转录本和kmer之间的对应关系，用法如下

代码语言：javascript复制

sailfish index 
-t hg19_refGeneMrna.fa  
-k 31 
-p 10 
-o hg19_saifish_db

-t参数指定参考基因组转录本的fasta格式的文件，-k参数指定kmer的长度，-p参数指定线程数，-o参数指定输出目录，在该目录下，会生成如下文件

代码语言：javascript复制

├── hash.bin
├── header.json
├── logs
├── quasi_index.log
├── rsd.bin
├── sa.bin
├── txpInfo.bin
└── versionInfo.json

2. 定量

索引构建好后，就可以定量了。sailfish支持单端/双端测序数据，用法如下

代码语言：javascript复制

sailfish quant 
-i hg19_saifish_db 
-g hg19.gtf 
-p 10 
-l IU 
-1 R1.fq 
-2 R2.fq 
-o out_dir

-i参数指定转录本索引所在的目录，-g指定转录本和基因之间的对应关系，可以是GTF格式的文件，也可以是t分隔的两列，第一列为转录本ID, 第二列为基因ID, 如果指定了这个文件，会给出基因和转录本两个水平的定量结果，否则只给出转录本水平的定量结果；-1和-2参数用于指定双端测序的reads, 单端测序的reads用-r参数指定，注意sailfish不支持读取压缩文件，-p参数指定线程数，-o参数指定输出结果的目录。

-l参数指定文库类型library type, sailfish采用了自定义的字母来表示，常见的文库类型如下

具体设置可以参考下图

在输出目录会生成以下两个文件

quant.sf
quant.genes.sf

quant.sf 是转录本水平的定量结果，内容如下

Name	Length	EffectiveLength	TPM	NumReads
NR_024540	1769	1556.3	12.3822	505.16
NR_110022	3187	2974.3	0.413261	32.2216
NR_033380	2915	2702.3	1.81983	128.915

quant.genes.sf是基因水平的定量结果，内容如下

Name	Length	EffectiveLength	TPM	NumReads
PLEKHG4B	12614	12401.3	0.368822	119.901
COL6A6	9563	9350.3	0.371259	91
ZNF660-ZNF197	1524.1	1480.69	1.57437	299.411

sailfish运行一个样本的定量只需要几分钟，传统的定量思路，hisat stringTie 也需要半个小时左右，速度上快了20倍以上。

·end·

—如果喜欢，快分享给你的朋友们吧—

0 人点赞