前言

??今天來跟大家分享一款非常優(yōu)秀的桌面級(jí)的RNAseq定量軟件——kallisto。該軟可以同時(shí)用于bulk and single-cell RNA-Seq 以及一些常見的目標(biāo)捕獲的高通量測序數(shù)據(jù)的定量。與以往常用的STAR等依賴基因組序列的比對(duì)軟件不同，kallisto采用一種被稱作偽比對(duì)(pseudoalignment)的方式直接將測序片段直接比對(duì)到cDNA序列然后定量。之所以稱kallisto為桌面級(jí)軟件是因?yàn)閗allisto能夠使用個(gè)人筆記本在10分鐘內(nèi)完成一套含有30M reads 的雙端RNA-seq數(shù)據(jù)的定量分析，可謂是名副其實(shí)的桌面級(jí)分析軟件，而且速度也是相當(dāng)?shù)目?。這么好的軟件當(dāng)然要分享給大家一起使用了，下面我們就來看看如何使用該軟件做定量。

• bulk-RNAseq的定量
  ??首先來看一下如何使用kallisto來做常規(guī)RNAseq的定量，只需兩步即可，第一步是構(gòu)建比對(duì)參考基因組(該步驟只需做一次，后面可以直接使用)，第二步就是定量(該軟件無需單獨(dú)做一次比對(duì)，而是直接將reads比對(duì)后定量出結(jié)果，在步驟上一氣呵成)。下面來看看具體的代碼示例：

1. 構(gòu)建比對(duì)參考基因組index
   你可以自己準(zhǔn)備好參考轉(zhuǎn)錄組文件來構(gòu)建定量使用的index，或者你也可以省時(shí)省力去kallisto網(wǎng)站上下載已經(jīng)構(gòu)建好的index，目前準(zhǔn)備好的index有人、小鼠、大鼠等具體見網(wǎng)站kallisto-transcriptome-indices。

#構(gòu)建index
kallisto index -i gencode.vM24.transcripts gencode.vM24.transcripts.fa
# -i 指定生成的index文件的名字

2. bulk-RNAseq定量
   構(gòu)建好index后就可以來定量了，使用方法也是相當(dāng)?shù)暮唵?，不過需要注意單雙端測序的數(shù)據(jù)使用的命令不同，單端模式需要特定指定fragment的平均長度和標(biāo)準(zhǔn)差，而雙端模式軟件本身會(huì)自動(dòng)評(píng)估這兩個(gè)指標(biāo)所以無需提供，具體代碼如下：

#單端測序的定量
kallisto quant -i gencode.vM24.transcripts -t 5 --single -l 50 -s 1 -o . SRR9170683.fastq.gz
#雙端測序的定量
kallisto quant gencode.vM24.transcripts -t 5 -o . SRR11600559_R1.fastq.gz SRR11600559_R2.fastq.gz

比對(duì)完成后會(huì)得到三個(gè)文件，以h5結(jié)尾的文件是定量結(jié)果的HDF5格式文件，可以使用kallisto的子命令h5dump將其轉(zhuǎn)化為文本文件；以json為結(jié)尾的文件包含運(yùn)行結(jié)果的摘要信息；以tsv結(jié)尾的文件是重要的定量文件。目錄結(jié)構(gòu)如下：

quant
├── abundance.h5 
├── abundance.tsv  
└── run_info.json

如果你想保留中間的比對(duì)結(jié)果如bam文件，可以在上面的定量命令中添加兩個(gè)參數(shù)：--genomebam --gtf ${gtf_file}，這樣定量完成后會(huì)在結(jié)果目錄多生成一個(gè)bam文件。

• scRNAseq的定量
  kallisto支持很多單細(xì)胞測序技術(shù)如10X genomic(V1-V3)、CELSEQ(V1-V3)、DROPSEQ、INDROPS、SCRUBSEQ、SURECELL。你可以用如下的命令查看該軟件支持的所有測序技術(shù)：

kallisto bus -l
List of supported single-cell technologies

short name       description
----------       -----------
10xv1            10x version 1 chemistry
10xv2            10x version 2 chemistry
10xv3            10x version 3 chemistry
CELSeq           CEL-Seq
CELSeq2          CEL-Seq version 2
DropSeq          DropSeq
inDropsv1        inDrops version 1 chemistry
inDropsv2        inDrops version 2 chemistry
inDropsv3        inDrops version 3 chemistry
SCRBSeq          SCRB-Seq
SureCell         SureCell for ddSEQ

用kallisto軟件來定量單細(xì)胞數(shù)據(jù)使用的是bus命令，會(huì)生成BUS(Barcode,UMI,Set format)格式比對(duì)結(jié)果文件，然后借助bustools軟件來生成定量文件。具體分為三個(gè)步驟，第一步是創(chuàng)建比對(duì)用的參考index，第二步是生成BUS文件，第三步定量。具體過程如下：

1. 構(gòu)建index
   由于構(gòu)建index過程與bulk-RNAseq一致，這里就不再重復(fù)展示了，代碼同上方。

2. 生成BUS文件
   對(duì)于單細(xì)胞的數(shù)據(jù)定量前需要先把數(shù)據(jù)比對(duì)生成一個(gè)bus文件，具體代碼如下：

#生成bus文件
kallisto bus -i gencode.vM24.transcripts -o SRR9169420_out -x 10xv2 SRR9169420_1.fastq.gz SRR9169420_2.fastq.gz

bus命令結(jié)束后會(huì)在結(jié)果目錄生成四個(gè)文件，matrix.ec：表達(dá)量的矩陣等價(jià)類文件；output.bus：原始的BUS(Barcode,UMI,Set format)格式比對(duì)結(jié)果文件；run_info.json：比對(duì)結(jié)果的概要信息；transcripts.txt：定量生成的轉(zhuǎn)錄本文件；
目錄結(jié)構(gòu)如下：

SRR9169420_out
├── matrix.ec
├── output.bus
├── run_info.json
└── transcripts.txt

3. scRNAseq定量
   bustools在對(duì)bus文件做定量的時(shí)候?qū)嶋H上包含了三個(gè)小的步驟分別為correct、sort、count，首先是對(duì)barcode的矯正，然后是排序，最后才是定量。不過bustools可以接受標(biāo)準(zhǔn)輸入，所以可以將這三個(gè)小步驟可以用管道“|”連接起來，這樣就可以直接從bus文件到生成定量結(jié)果一氣呵成了，具體代碼如下：

#一步從bus到定量結(jié)果
bustools correct -w 10xv2_whitelist.txt -p SRR9169420_out/output.bus | bustools sort -T bus_out/tmp/ -t 2 -p - | bustools count -o bus_out/genecount -g transcripts_to_genes.txt -e bus_out/matrix.ec -t bus_out/transcripts.txt --genecounts -

需要提醒的是在這里需要提供兩個(gè)額外的文件：1、10xv2_whitelist.txt(10xv2試劑盒的barcode白名單文件，可以去10x官網(wǎng)下載，或者用bustools的whitelist命令從bus文件中生成一個(gè))；2、transcripts_to_genes.txt(轉(zhuǎn)錄本到基因的對(duì)應(yīng)關(guān)系文件)。定量完成后會(huì)生成一系列的結(jié)果文件，最重要的表達(dá)量文件是cells_x_genes.barcodes.txt、cells_x_genes.genes.txt、cells_x_genes.mtx(類似cellranger的結(jié)果)。

快捷方式

其實(shí)使用kallisto軟件來做定量分析有一個(gè)更為快捷方便的方式可選，那就是使用軟件kb_python來完成分析,該軟件包裝了kallisto和bustools兩個(gè)軟件，使用起來更為方便。就拿上面的scRNAseq定量來說，我們看看使用kb_python來分析有多方便，具體代碼如下：

#安裝kb軟件
pip3 install kb-python
#使用kb軟件來定量
kb count -i gencode.vM24.transcripts -g transcripts_to_genes.txt -x 10xv2 -o SRR9169420_out --filter bustools -t 5 SRR9169420_1.fastq.gz SRR9169420_2.fastq.gz

從上面的命令可以看出使用kb來做的定量可以省去中間的步驟，并且可以只給出試劑盒型號(hào)(這里是10xv2)而不提供barcode白名單文件。是不是省事省心，對(duì)于這么貼心的軟件應(yīng)該會(huì)心動(dòng)了吧！

最后

??kallisto軟件本身的特點(diǎn)就是對(duì)設(shè)備要求不高且速度快，如果只是做RNAseq的定量不妨試一試。這里主要介紹了定量方面的命令，對(duì)于其他的命令，大家可以自己去測試一下。emm，今天就分享到這里吧。