本文我們來簡單介紹一下非??旖莺糜玫囊粋€RNAseq工具——Kallisto。Kallisto被我推薦的原因是其速度非???,在我的Mac Pro就可以運行使用,而且其結(jié)果也比較準,使用起來還十分簡單。
RNA-seq分析通常有以下幾種流程。
第一種是參考基因組,即先通過HISAR、STAR等軟件把序列比對到參考基因組然后再進行轉(zhuǎn)錄本鑒定及定量。根據(jù)有無GFF注釋可以分為兩種,如果沒有GFF注釋鑒定完之后再依據(jù)同源比對結(jié)果進行功能注釋。
第二種是今天要講的——參考轉(zhuǎn)錄組方法,直接將序列比對到轉(zhuǎn)錄組,然后進行轉(zhuǎn)錄本鑒定及定量。顯然,該方法的優(yōu)勢就是快捷,而缺點也很明顯,因為只和參考轉(zhuǎn)錄組進行非剪接比對所以無法鑒定出新的轉(zhuǎn)錄本或者是新的非編碼RNA包括lncRNA等。
第三種是無參考基因組,有時候我們做的物種比較小眾,所以還沒有參考基因組,所以只能先利用De Bruijin的方法對序列進行從頭拼接,然后再進行比對、定量,確定表達量。

因此,根據(jù)你的數(shù)據(jù)特點和你的需求可以選擇合適的方法。實際上,很多實驗室做RNA-seq可能暫時并不關(guān)注新的轉(zhuǎn)錄本,只想看一看不同條件下實驗組和對照組有哪些基因的表達量發(fā)生了變化,因此這時我們就可以選擇第二種方法,直接和轉(zhuǎn)錄組進行非剪接比對。今天我們就來講第二類方法中很優(yōu)秀的一個工具Kallisto。
Kallisto于2016年發(fā)表在Nature biotechnology,截至目前引用次數(shù)超過1300次。

Kallisto的安裝
#如果你的電腦是mac可以用以下的方式進行安裝
ruby -e "$(curl -fsSL https://raw.githubusercontent.com/Homebrew/install/master/install)"
brew install kallisto
#如果你已經(jīng)安裝了conda,可以用conda安裝:
conda install kallisto
#kallisto can also be installed on FreeBSD via the FreeBSD ports system using
pkg install kallisto
#**kallisto** binaries for Mac OS X, NetBSD, RHEL/CentOS and SmartOS can be installed on most POSIX platforms using pkgsrc:
pkgin install kallisto
安裝完成后,輸入kallisto:

Kallisto的使用
在正式開始之前我們需要準備以下數(shù)據(jù)
1、目標木中的參考轉(zhuǎn)錄組文件:cDNA文件
2、待分析的測序文件
本文我們以人的樣本為例下載相關(guān)的文件
準備工作
cDNA文件的下載:hg19(GRCh37)/hg38(GRCh38)


根據(jù)你需要的版本進行下載cDNA文件:
GRCh38:
ftp://ftp.ensembl.org/pub/release-96/fasta/homo_sapiens/cdna/
GRCh37:
ftp://ftp.ensembl.org/pub/grch37/current/fasta/homo_sapiens/cdna/
cDNA sequences for Ensembl or ab initio predicted genes,所以我們下載cdna.all.fa.gz的文件
這里,我下載的是GRCh37的版本。
第一步建立索引
這里要注意一下參數(shù)-i并不是指輸入文件,此處的i代表index,后面接的是你輸出的index名字。以后如果還是該物種,你可以直接使用本次建立的index,不用重復該步驟。
kallisto index ./Homo_sapiens.GRCh37.75.cdna.all.fa.gz -i Homo_sapiens.GRCh37.75.cdna.all.index
第二步轉(zhuǎn)錄本的鑒定及定量
#雙端測序
kallisto quant -i ./Homo_sapiens.GRCh37.75.cdna.all.index -o ./Result -t 4 -b 100 PATH/Sample_R1.fq.gz PATH/Sample_R2.fq.gz
#查看kallisto quant幫助
kallisto quant -h
Usage: kallisto quant [arguments] FASTQ-files
Required arguments:
-i, --index=STRING Filename for the kallisto index to be used for
quantification
-o, --output-dir=STRING Directory to write output to
Optional arguments:
--bias Perform sequence based bias correction
-b, --bootstrap-samples=INT Number of bootstrap samples (default: 0)
--seed=INT Seed for the bootstrap sampling (default: 42)
--plaintext Output plaintext instead of HDF5
--fusion Search for fusions for Pizzly
--single Quantify single-end reads
--fr-stranded Strand specific reads, first read forward
--rf-stranded Strand specific reads, first read reverse
-l, --fragment-length=DOUBLE Estimated average fragment length
-s, --sd=DOUBLE Estimated standard deviation of fragment length
(default: -l, -s values are estimated from paired
end data, but are required when using --single)
-t, --threads=INT Number of threads to use (default: 1)
--pseudobam Output pseudoalignments in SAM format to stdout
如果是單端測序還需要給-l參數(shù),后面跟估計的平均片段長度,-s參數(shù)后面跟估計的片段長度標準差。這兩個參數(shù)可以使用其他軟件如Agilent Bioanalyzer等確定。
#單端測序
kallisto quant -i index -o output --single -l length -s SD file.fq.gz
Kallisto的結(jié)果
然后就會生成三個文件:abundances.h5,abudances.tsv,run_info.json
abundance.h5
HDF5二進制格式的文件,包含了運行日志信息、表達豐度估計值、bootstrap估計和轉(zhuǎn)錄本長度信息。該文件可以直接用sleuth讀取處理,也可以使用kallisto h5dump命令將其轉(zhuǎn)變?yōu)榧兾谋镜膖sv格式文件
abundance.tsv
包含有表頭的純本文tsv格式文件,表頭是:target_id, length, eff_length, est_counts, tpm
run_info.json
一個json格式的日志文件
然后我們可以看各個轉(zhuǎn)錄本的TPM即其表達量。TPM具體的計算方式及其與RPKM、FPKM的差異可以看之前的日志RPM(CPM)/RPKM/FPKM/TPM
下一節(jié)我們將會講解如何用R包Sleuth對轉(zhuǎn)錄本進行差異表達分析等。