基礎医学生物学研究者のブログ

cellranger countによるfastqからカウントマトリックスの作成方法

シングルセルRNAシークエンスのリードをリファレンスシークエンスへマッピングし、リードカウントをやってくれるパイプラインがcellranger countである。ここでは、シークエンスから返ってきた、もしくはcellranger mkfastqなどで適切に作成したfastqファイルを、cellranger countで処理することでカウントマトリックスを作成する方法について記す。個人的に思うのが、10xGenomicsのシングルセルRNAシークエンスが普及したのは、この会社の宣伝やサポートに加えて、このcellranger countがすごく便利だったことも理由ではないかと思っている。

以下がcellrangerを行うための手順である。

目次

シングルセルRNAシークエンスでは十分なSWAP領域が必要である

SWAP領域の設定を256GBくらいにした方が良い。後のSeuratではそれでの足りなくなるケースがあるので、それ以上の512GBとかだと、なお良い。

過去の経験から、メモリの設定には要注意である。そもそもcellrangerはコンピューターのCPUが16スレッド以上、メモリが128GB以上を推奨している。しかし、個人で解析を行う場合はこのようなコンピューターを用意するのは簡単ではない。その一方で、その条件を満たしていないコンピューターでもcellrangerが実行できてしまう。

そしてcellrangerでは、リードのリファレンスシークエンスへのマッピングにSTARを利用している。このSTARというソフトは、大量のメモリを必要とする。最近ではGENCODEにアップロードされているリファレンスシークエンスのアノテーションがどんどん増えていくため、そのファイルサイズもどんどん増えており、それに従ってSTAR仕様時のメモリも大きくなってきていると個人的に感じている。少し前にはSTARのマニュアルには32GBのメモリが必要と書いてあったが、2019年くらいの時点で実際には60GBは必要だった。cellranger countで使用するGTFファイルは、元のGTFファイルから必要な情報を拾っていることから、元のGTFファイルよりも小さいものになっているが、それでもやはりメモリは必要である。そこで重重要なのは、Linux（今ではUbunntuかDebianだろうと思う。）でもWSL2（WIndows subsystem for Linux）でも、十分なSWAP領域を設定する必要があるということである。経験的に、STARを走らせているときにメモリが不足すると、そこでエラーがでて止まる。STAR単体でマッピングしているときは対処の使用もあるが、自分にはそのあたりをcellrangerで調整するスキルがないため、大きなSWAP領域を設定した方がよい。そして必要なメモリは、個人的にはcellrangerには256GBが必要と考えている。10xGenomicsの推奨通り、確かに128GBでも大丈夫だと思うが、メモリの使用量を見ていて、128GBではギリギリだった。そして余談だが、次にクラスタリングと発現解析に使用するSeuratでは、細胞数に応じて合計で400GBくらい（ということで、設定は合計で512GBくらいになると思う。）のメモリが必要になることを経験した。シングルセルRNAシークエンスはこのような解析になるので、大きな大学のようなサーバーやクラスターを組んでいるようなコンピューターが必要になるのだろう。ということで、SWAP領域は256GBくらいは確保したほうが良いと思う。

SWAP領域の設定手順

以下がSWAP領域の設定方法である。ここにも書いてある。以前の記事はDebian 9での設定方法であるが、Ubuntu 20.04でもそれは基本的に動く。

SWAPをONにするとき

以下がSWAPファイルの作成と有効化の手順である。

```sh
su # Rootに変更。
fallocate -l ２５６GB /swapfile # 使いたいサイズを設定する。
chmod 600 /swapfile #パーミッションを与える。
mkswap /swapfile #SWAP領域の作成
swapon /swapfile #SWAP領域を有効にする。
gedit /etc/fstab #テキストファイルを開く。最悪、viをつかう。

# テキストファイルに「/swapfile swap swap defaults 0 0」を追加してSave。

swapon --show #SWAP領域の状態を表示
free -h　#メモリとSWAPのサイズなどを表示
```

SWAPを消すとき

ディスク容量に対して設定したSWAPが大きすぎた場合など、OSを起動できなくなるときがあるので、それには要注意である。その場合はRecovery modeでOSを起動して、SWAPを消すなどの処置が必要である。それ以外にも、単に設定したSWAPが大きすぎた場合、逆に足りなかった場合などでは再設定が必要だろうと思う。そのときは、以下を使って設定したSWAP領域を消して、改めて設定する。

```sh
su
swapoff -v /swapfile
gedit /etc/fstab

# 「/swapfile swap swap defaults 0 0」を消してSave。

rm /swapfile
```

256GBのSWAPを設定したところ。

cellranger count

以下がcellranger countの基本的なコードである。–id=には任意の実験IDやサンプルIDを自分がわかりやすいように入れると良い。–fastqs=には入手したfastqファイルが入っているディレクトリを書く。–transcriptome=にはmkgtkとmkrefを使って作成したリファレンスシークエンスを入れる。–sample=、–lane=には、対象とするfastqのファイル名についているサンプル名とレーン番号を入れる。–expect-cells=にはライブラリ調整のまえのエマルジョン作成のときに使用した細胞数を入れるのが筋だろう。詳しくは10xGenomicsのページ（https://support.10xgenomics.com/single-cell-gene-expression/software/pipelines/latest/using/tutorial_ct?gclid=Cj0KCQiAgOefBhDgARIsAMhqXA5g0LJj5syFDDzbJCBRyOFVMoCjow1cGj9yoZ9wp4ULD-OJNcRLdZcaAqLYEALw_wcB#explorecount）を参照する。

```sh
cd /mnt/sda1/XXXX/YYYY/Data/ID/Analysis/count

cellranger count \
    --id=A1_1 \
    --fastqs=/mnt/sda1/XXXX/YYYY/Data/ID/Analysis/fastq \
    --transcriptome=/mnt/sda1/XXXX/YYYY/reference_seq_annotation/MsvM31 \
    --sample=A1 \
    --lanes=1 \
    --expect-cells=10000 \
    --localcores=24;

```

出力

出力は、なんかいっぱい出る。要はパイプラインを走らせているときに吐き出された中間生成物である。以下が個人的に重要と思うファイルもしくはディレクトリである。

/outs/filtered_feature_bc_matrix

最終的に必要なファイルは、生成されたディレクトリ内にある/outs/filtered_feature_bc_matrixというディレクトリに入っている3つのファイル（barcodes.tsv.gz, features.tsv.gz, matrix.mtx.gz）である。Seuratなどではこれを直接読み込むことでその後の解析を進めていく。

filtered_feature_bc_matrix.h5

このファイルも場合に依っては重要である。これも、Seurat、monocle3、Scanpyなどで読み込んで解析することができる。

possorted_gemome_bam.bam

これはRNA velocity解析に必要になるbamファイルである。

web_summary.html

これにはマッピング結果が書かれていて、そのマッピングが良かったのか、ノイズが多いのがが書かれている。しかし、これはほとんど参考にならない。なんやったら無視しても良い。その理由は、理想的とは言えない結果だったとして、それを解析しない、とかいう選択肢は無いためである。

cloupe.cloupe

これは10xGenomicsのLoupe Browserというソフトで開くことができるファイルである。というか、このLoupe Browserというソフト、超絶によく出来ていると思うのだが、10xGenomicsはこのソフトを論文に載せる結果に利用できるまでに開発しないのだろうか。そうしたらSeuratなんか使う必要ないのではないか？？どうなんだろう。

これらが無事に出力されれば、次はSeurat、monocle3、scanpyなどでクラスタリング、細胞のタイピング、DEGなんかを行う。それらが終われば結果によってはTrajectory解析やRNA velocity解析なんかを行う。ただし、解析していて思うのだが、こういったカッコいい解析を行うのは細胞のタイピングや発現解析でかなり良い結果がでた場合だと思う。