RLE 正規化

DESeq/DESeq2 の正規化は size factor と呼ばれる係数を計算し、この係数を利用してカウントデータを正規化している。この size factor は、edgeR に実装されている TMM 正規化法で計算される normalization factor と異なるものである。

m 個のライブラリーが存在するとき、ライブラリー j の size factor は次のように求める。ここで、kij はライブラリー j における遺伝子 i のカウントデータを表す。

\[ s_{j} = median \frac{k_{ij}}{\left(\prod_{v=1}^{m}k_{iv}\right)^{\frac{1}{m}}} \]

このライブラリーサイズを利用して、データを正規化するとき、次の手順により正規化後のカウントデータを qij を計算する。

\[q_{i\rho} = \frac{k_{ij}}{s_{j}} \]

実際に R で実行した際のスクリプトは次のようになる。サンプルデータは A 群と B 群の 2 群間を想定し、biological replicate が 2 とする。また、遺伝子数を 5 個とする。

x <- data.frame(
        A1 = c(20, 40, 30, 120, 80),
        A2 = c(25, 50, 30, 100, 90),
        B1 = c(25, 20, 30, 420, 90),
        B2 = c(20, 25, 30, 390, 80)
)

sizefactors <- apply(x / apply(x, 1, function(y){(prod(y))^(1/4)}), 2, median)
sizefactors
##       A1       A2       B1       B2 
## 0.942809 1.060660 1.060660 0.942809 

library(DESeq)
d <- newCountDataSet(x, condition=c(1,1,2,2))
d <- estimateSizeFactors(d)
d@phenoData@data$sizeFactor
## [1] 0.942809 1.060660 1.060660 0.942809

normalization factors と size factors

発現変動遺伝子同定を行うパッケージとして edgeR と DESeq がよく知られている。両者では検定法が異なるほか、正規化法も異なっている。edgeR は TMM 法による正規化を行う。normalization factors とよばれる係数が計算され、正規化はこの係数とライブラリーサイズ(総リード数)の積を利用して正規化を行う。一方、DESeq/DESeq2 は size factors とよばれる係数が計算され、上述の方法により正規化を行う。

normalization factors と size factors はともに正規化に関わる係数であるが、両者は異なる概念である。normalization factors はライブラリーサイズに依存するが、size factors は依存しない。従って、両者を直接に比較できない。もし、比較する場合は、normalization factors とライブラリーサイズの積を size factors と比較すべきである(参照:grokbase)。

References

  • Anders S, Huber W. Differential expression analysis for sequence count data. Genome Biol. 2010, 11(10):R106. PubMed Abstract
  • Robinson MD, Oshlack A. A scaling normalization method for differential expression analysis of RNA-seq data. Genome Biol. 2010, 11(3):R25. PubMed Abstract