本文我們來(lái)簡(jiǎn)單介紹一下非逞芬疲快捷好用的一個(gè)RNAseq工具——Kallisto。Kallisto被我推薦的原因是其速度非常快玫膀，在我的Mac Pro就可以運(yùn)行使用，而且其結(jié)果也比較準(zhǔn)爹脾，使用起來(lái)還十分簡(jiǎn)單帖旨。
RNA-seq分析通常有以下幾種流程。
第一種是參考基因組灵妨，即先通過(guò)HISAR解阅、STAR等軟件把序列比對(duì)到參考基因組然后再進(jìn)行轉(zhuǎn)錄本鑒定及定量。根據(jù)有無(wú)GFF注釋可以分為兩種泌霍，如果沒(méi)有GFF注釋鑒定完之后再依據(jù)同源比對(duì)結(jié)果進(jìn)行功能注釋货抄。
第二種是今天要講的——參考轉(zhuǎn)錄組方法，直接將序列比對(duì)到轉(zhuǎn)錄組朱转，然后進(jìn)行轉(zhuǎn)錄本鑒定及定量蟹地。顯然，該方法的優(yōu)勢(shì)就是快捷藤为，而缺點(diǎn)也很明顯怪与，因?yàn)橹缓蛥⒖嫁D(zhuǎn)錄組進(jìn)行非剪接比對(duì)所以無(wú)法鑒定出新的轉(zhuǎn)錄本或者是新的非編碼RNA包括lncRNA等。
第三種是無(wú)參考基因組缅疟，有時(shí)候我們做的物種比較小眾分别，所以還沒(méi)有參考基因組遍愿，所以只能先利用De Bruijin的方法對(duì)序列進(jìn)行從頭拼接，然后再進(jìn)行比對(duì)耘斩、定量沼填，確定表達(dá)量。
因此括授，根據(jù)你的數(shù)據(jù)特點(diǎn)和你的需求可以選擇合適的方法倾哺。實(shí)際上，很多實(shí)驗(yàn)室做RNA-seq可能暫時(shí)并不關(guān)注新的轉(zhuǎn)錄本刽脖，只想看一看不同條件下實(shí)驗(yàn)組和對(duì)照組有哪些基因的表達(dá)量發(fā)生了變化羞海，因此這時(shí)我們就可以選擇第二種方法，直接和轉(zhuǎn)錄組進(jìn)行非剪接比對(duì)曲管。今天我們就來(lái)講第二類方法中很優(yōu)秀的一個(gè)工具Kallisto却邓。
Kallisto于2016年發(fā)表在Nature biotechnology，截至目前引用次數(shù)超過(guò)1300次院水。
Kallisto的安裝
#如果你的電腦是mac可以用以下的方式進(jìn)行安裝
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的使用
在正式開(kāi)始之前我們需要準(zhǔn)備以下數(shù)據(jù)
1檬某、目標(biāo)木中的參考轉(zhuǎn)錄組文件：cDNA文件
2撬腾、待分析的測(cè)序文件
本文我們以人的樣本為例下載相關(guān)的文件
準(zhǔn)備工作
cDNA文件的下載：hg19(GRCh37)/hg38(GRCh38)
根據(jù)你需要的版本進(jìn)行下載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漏隐，不用重復(fù)該步驟喧半。
kallisto index ./Homo_sapiens.GRCh37.75.cdna.all.fa.gz -i Homo_sapiens.GRCh37.75.cdna.all.index
第二步轉(zhuǎn)錄本的鑒定及定量
#雙端測(cè)序
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
如果是單端測(cè)序還需要給-l參數(shù)，后面跟估計(jì)的平均片段長(zhǎng)度青责，-s參數(shù)后面跟估計(jì)的片段長(zhǎng)度標(biāo)準(zhǔn)差挺据。這兩個(gè)參數(shù)可以使用其他軟件如Agilent Bioanalyzer等確定。
#單端測(cè)序
kallisto quant -i index -o output --single -l length -s SD file.fq.gz
Kallisto的結(jié)果
然后就會(huì)生成三個(gè)文件：abundances.h5,abudances.tsv,run_info.json
abundance.h5
HDF5二進(jìn)制格式的文件脖隶，包含了運(yùn)行日志信息扁耐、表達(dá)豐度估計(jì)值、bootstrap估計(jì)和轉(zhuǎn)錄本長(zhǎng)度信息浩村。該文件可以直接用sleuth讀取處理做葵，也可以使用kallisto h5dump命令將其轉(zhuǎn)變?yōu)榧兾谋镜膖sv格式文件
abundance.tsv
包含有表頭的純本文tsv格式文件，表頭是：target_id, length, eff_length, est_counts, tpm
run_info.json
一個(gè)json格式的日志文件
然后我們可以看各個(gè)轉(zhuǎn)錄本的TPM即其表達(dá)量心墅。TPM具體的計(jì)算方式及其與RPKM酿矢、FPKM的差異可以看之前的日志RPM(CPM)/RPKM/FPKM/TPM
下一節(jié)我們將會(huì)講解如何用R包Sleuth對(duì)轉(zhuǎn)錄本進(jìn)行差異表達(dá)分析等。
