scRNA解析に有用 - SingleCellExperiment クラス

SingleCellExperimentクラス

single cell RNA (scRNA-seq) 解析のためのS4 Classesです。
様々なscRNA-seq解析のパッケージで採用されており、scRNA-seq解析に有用なデータを保持しながら解析を進めることができます。

Bioconductor
Bioconductor - SingleCellExperiment

チュートリアル
An introduction to the SingleCellExperiment class

scRNAseqパッケージ

オープンソース分化って素晴らしいですね。
NCBIのSRA や、EBIのArrayExpressからも入手可能ですが、
scRNAパッケージは、代表的な3つのscRNAのデータセットをRのオブジェクトとしてパッケージ化してくれています。 以下がその内容です。

変数名 内容
fluidigm FLudigm C1で取得した異なるdepthの神経細胞65 細胞分(Pollen et al. 2014)(SRA: SRP041736)
th2 Tヘルパー細胞 96 細胞分(Mahata et al. 2014)
allen マウス視覚皮質379 細胞分 (Tasic et al. 2016)

チュートリアル: Overview of the scRNAseq dataset collection

例えば、allenのデータを読み込みたい際は、 ライブラリを読み込んで

> library(scRNAseq)
> data(allen)

のように実行するだけです。 データを見てみましょう。

> allen
class: SummarizedExperiment
dim: 20908 379
metadata(2): SuppInfo which_qc
assays(4): tophat_counts cufflinks_fpkm rsem_counts rsem_tpm
rownames(20908): 0610007P14Rik 0610009B22Rik ... Zzef1 Zzz3
rowData names(0):
colnames(379): SRR2140028 SRR2140022 ... SRR2139341 SRR2139336
colData names(22): NREADS NALIGNED ... Animal.ID passes_qc_checks_s

このオブジェクトは、SummarizedExperimentというクラスのオブジェクトになっていることがわかります。

オブジェクトの変換

SingleCellExperimentクラスのオブジェクトに変換したい際は、as関数を使用します。

> sce <- as(allen, "SingleCellExperiment")
> sce
class: SingleCellExperiment
dim: 20908 379
metadata(2): SuppInfo which_qc
assays(4): tophat_counts cufflinks_fpkm rsem_counts rsem_tpm
rownames(20908): 0610007P14Rik 0610009B22Rik ... Zzef1 Zzz3
rowData names(0):
colnames(379): SRR2140028 SRR2140022 ... SRR2139341 SRR2139336
colData names(22): NREADS NALIGNED ... Animal.ID passes_qc_checks_s
reducedDimNames(0):
spikeNames(0):

SingleCellExperimentに変換されました!

データをざっくり見てみる

サンプルの細胞ごとのメタデータを見たい場合、colData関数を使用します。 先っちょだけみたいので、[0:5, 0:5]で左上の5 × 5行列だけ取ってきます。

> colData(sce)[0:5,0:5]
DataFrame with 5 rows and 5 columns
              NREADS  NALIGNED    RALIGN TOTAL_DUP     PRIMER
           <numeric> <numeric> <numeric> <numeric>  <numeric>
SRR2140028  13743900  13011100   94.6681     51.11  0.0148148
SRR2140022  14078700  13521900   96.0454   55.9157 0.00853083
SRR2140055   5842930   5135980   87.9008   59.1126  0.0356112
SRR2140083  16784400  15585800   92.8587   55.3076  0.0209695
SRR2139991  11558600  10864300   93.9929   50.2258   0.016408

リード数、アライメントされたリード数等、サンプルに関する細胞ごとのデータがデータフレーム形式で格納されています。

> colnames(colData(sce))
 [1] "NREADS"                       "NALIGNED"
 [3] "RALIGN"                       "TOTAL_DUP"
 [5] "PRIMER"                       "PCT_RIBOSOMAL_BASES"
 [7] "PCT_CODING_BASES"             "PCT_UTR_BASES"
 [9] "PCT_INTRONIC_BASES"           "PCT_INTERGENIC_BASES"
[11] "PCT_MRNA_BASES"               "MEDIAN_CV_COVERAGE"
[13] "MEDIAN_5PRIME_BIAS"           "MEDIAN_3PRIME_BIAS"
[15] "MEDIAN_5PRIME_TO_3PRIME_BIAS" "driver_1_s"
[17] "dissection_s"                 "Core.Type"
[19] "Primary.Type"                 "Secondary.Type"
[21] "Animal.ID"                    "passes_qc_checks_s"

可視化してみる

試しにコーディングRNAの塩基数の割合(PCT_CODING_BASES)と、メッセンジャーRNAの塩基数の割合(PCT_MRNA_BASES)を2軸にとってプロットしてみましょう。

> plot(x=colData(sce)$PCT_CODING_BASES, 
       y=colData(sce)$PCT_MRNA_BASES,
       xlab="PCT_CODING_BASES",
       ylab="PCT_MRNA_BASES")

f:id:kimoppy126:20180927085233p:plain

正に相関してそうですね。

ポアソンの相関係数を出してみましょう。

> cor(x=colData(sce)$PCT_CODING_BASES, y=colData(sce)$PCT_MRNA_BASES)
[1] 0.9679943

確かに相関が見られますね!

おわりに

簡単な例でしたがSingleCellExperiment クラスの使い方でした。
scaterやscPipe、scranといったscRNA-seq解析用のパッケージではこのクラスのオブジェクトをインプットとしてデータ解析を行います。
時間があれば各パッケージの解析例も紹介しますねー。