kallisto使ってみた！

Psudoalignmentというある意味姑息な手法を使って高速にRNA-seq data解析を行うソフトkallisto使ってみた。

インストールその他は、チュートリアルGetting Started及びKallistoInstallが詳しい！

リファレンスゲノムファイル（cDNA）をEnsembleのFTPサイトで落とし、macのterminalで

$ Kallisto index -i myMouseIndex Mus_musculus.GRCm38.cdna.all.fa

としてインデックスファイルをつくる。

その後kallistoのマニュアルを参考に、kallistoの定量のコマンドを打ち込めば良い！

今回pair endedのシーケンスだったため、XXX_SYYY_L001_R1_001.fastq, XXX_SYYY_L001_R2_001.fastqという２種のファイルがあったのだけれど、それぞれR1, R2というディレクトリを作り移動！

$ mkdir R1, R2

$ mv *_R1_* R1
$ mv *_R2_* R2

それから実際のコマンドをうちこむ。ファイル名は共通部分をbasenameコマンドで、変数として切り取り利用すれば良い！

変数の定義の際に=の前後でスペースを入れないこと
変数をテキストとして使う場合に${変数%%.test}を使うところに気づかずちょっと手間取ったが、

以下で、カウンティングファイルとともにBAMファイルも作ってくれる！（前もってsamtoolsのインストールが必要である）。

$ for file in R1/*.fastq

do

echo $file
basef=`basename $file _L001_R1_001.fastq`
R1file=${basef%%.text}_L001_R1_001.fastq
R2file=${basef%%.text}_L001_R2_001.fastq
echo $R1file
echo $R2file
kallisto quant -i myMouseIndex -o ${basef%%.text} --pseudobam R1/$R1file R2/$R2file | samtools view -bS - > ${basef%%.text}.bam

done

で各サンプルごとのフォルダにカウントファイルほか、BAMファイルもつくってくれる！

次に

for folder in *_S*
do
echo $folder
cp ${folder%%.txt}/abundance.tsv ${folder%%.txt}_abundance.tsv
done

で各フォルダーに保存してあるカウントファイルにサンプル名をつけたファイルができる。

kallistoの文献はこちら
kallistoの紹介記事Kallisto, a new ultra fast RNA-seq quantitation method