フローサイトメトリー(FCM)という実験がある。これは、細胞表面に抗原があるのを抗体反応で標識して、標識された細胞に光を当てて抗体反応が起きていたか、つまりその抗原が存在するかを識別する実験である。
実験機器が許す限り、抗原標識は可能で、たいてい2〜3の抗体、多くて何十できるとかできないとか。
正常な細胞でも白血病などの異常な細胞でも、ある種の抗原パターンがあって、それらの組み合わせで理論的には細胞を分離できる。
flowCoreパッケージでFCMデータを扱ってみる。
# Required once per R installation source("http://bioconductor.org/biocLite.R") biocLite("flowCore")
library(flowCore) file.name <- system.file("extdata","0877408774.B08", package="flowCore") x <- read.FCS(file.name, transformation=FALSE) # 定量データの部分だけ別に取っておく。時系列も抜いておく。 dd <- x@exprs[, grep("-", colnames(x@exprs))] time0 <- x@exprs[, "Time"]
サンプルファイルを読み込んでみると、7つの抗体で標識された50 sec ほどのデータであることが分かる。
x
flowFrame object '0877408774.B08' with 10000 cells and 8 observables: name desc range minRange maxRange $P1 FSC-H FSC-H 1024 0 1023 $P2 SSC-H SSC-H 1024 0 1023 $P3 FL1-H 1024 0 1023 $P4 FL2-H 1024 0 1023 $P5 FL3-H 1024 0 1023 $P6 FL1-A <NA> 1024 0 1023 $P7 FL4-H 1024 0 1023 $P8 Time Time (51.20 sec.) 1024 0 1023 147 keywords are stored in the 'description' slot
実験データは1細胞当たり、7つの抗体それぞれの光量の実数データを持つ。
head(x@exprs)
FSC-H SSC-H FL1-H FL2-H FL3-H FL1-A FL4-H Time [1,] 382 77 618 0 225 55 286 1 [2,] 628 280 245 431 259 0 371 1 [3,] 1023 735 699 448 215 143 638 1 [4,] 373 128 202 354 94 0 149 1 [5,] 1023 1023 618 742 408 61 866 1 [6,] 489 292 179 374 154 0 363 1
とりあえずかっこよく散布図を描く。なんとなくクラスターっぽいのがある。
library(GGally) colnames(dd) <- gsub("-", "_", colnames(dd)) # - が入っているとなぜかエラー ggpairs(dd, upper = list(continuous = "density", combo = "box"), lower = list(continuous = "points", combo = "dot"), params=c(cex=0.3, col=memb), )
ある抗原があるかないかは、gating という閾値を決めてその値を超えるものはさらに次のgatingを満たすものは以下略というようにして絞っていく。
ここでmultiple experienceを読み込んでgatingしてみる。
# loading 5 experiment together into an object fs fs <- read.flowSet(path=system.file("extdata", "compdata", "data", package="flowCore"), name.keyword="SAMPLE ID", phenoData = list(name="SAMPLE ID", Filename="$FIL")) pData(phenoData(fs)) # meta data for 5 samples summary(fs) x1 <- fs[[1]] # same structure as x above
Gatingを行う。
plot(x1@exprs[, c("FL1-H", "FL2-H")], cex=0.5) # どこに垂直・水平の2直線を引きたい? How do you draw vertical and holizontal lines? abline(v=c(0, 12), h=c(0, 12), lty=2, col=2, lwd=3) # 枠内に入った割合は? Fractio within the gated rectangle? rectGate <- rectangleGate(filterId="Fluorescence Region", "FL1-H"=c(0, 12), "FL2-H"=c(0, 12)) result <- filter(x1, rectGate) summary(result)
Fluorescence Region+: 9811 of 10000 events (98.11%)
98%が右下にgatingされている。この区切り方で本当にいいかと思ったりしたら、クラスタリングしてみる。4つくらいに分かれるかなと思ったら
dist.1.2 <- dist(x1@exprs[, c("FL1-H", "FL2-H")]) cl <- hclust(dist.1.2) plot(cl) memb <- cutree(cl, k = 4) # クラスタリングの枝葉切り plot(x1@exprs[, c("FL1-H", "FL2-H")], col=memb, pch=16, cex=1)
この条件でgatingすると、陽性の98%の細胞クラスターと陰性の2%の細胞クラスターが抽出される。
# extracting a subset with the gate fs2 <- Subset(fs, rectGate)
fs[[1]]
flowFrame object 'NA' with 10000 cells and 7 observables: name desc range minRange maxRange $P1 FSC-H FSC-Height 1024 0 1023 $P2 SSC-H SSC-Height 1024 0 1023 $P3 FL1-H <NA> 1024 1 10000 $P4 FL2-H <NA> 1024 1 10000 $P5 FL3-H <NA> 1024 1 10000 $P6 FL1-A <NA> 1024 0 1023 $P7 FL4-H <NA> 1024 1 10000 141 keywords are stored in the 'description' slot
fs2[[1]]
flowFrame object 'NA' with 9811 cells and 7 observables: name desc range minRange maxRange $P1 FSC-H FSC-Height 1024 0 1023 $P2 SSC-H SSC-Height 1024 0 1023 $P3 FL1-H <NA> 1024 1 10000 $P4 FL2-H <NA> 1024 1 10000 $P5 FL3-H <NA> 1024 1 10000 $P6 FL1-A <NA> 1024 0 1023 $P7 FL4-H <NA> 1024 1 10000 141 keywords are stored in the 'description' slot
x1split <- split(x1, rectGate) show(x1split) # split into two subsets with the gate
$`Fluorescence Region+` flowFrame object 'NA (Fluorescence Region+)' with 9811 cells and 7 observables: name desc range minRange maxRange $P1 FSC-H FSC-Height 1024 0 1023 $P2 SSC-H SSC-Height 1024 0 1023 $P3 FL1-H <NA> 1024 1 10000 $P4 FL2-H <NA> 1024 1 10000 $P5 FL3-H <NA> 1024 1 10000 $P6 FL1-A <NA> 1024 0 1023 $P7 FL4-H <NA> 1024 1 10000 141 keywords are stored in the 'description' slot $`Fluorescence Region-` flowFrame object 'NA (Fluorescence Region-)' with 189 cells and 7 observables: name desc range minRange maxRange $P1 FSC-H FSC-Height 1024 0 1023 $P2 SSC-H SSC-Height 1024 0 1023 $P3 FL1-H <NA> 1024 1 10000 $P4 FL2-H <NA> 1024 1 10000 $P5 FL3-H <NA> 1024 1 10000 $P6 FL1-A <NA> 1024 0 1023 $P7 FL4-H <NA> 1024 1 10000 141 keywords are stored in the 'description' slot
kmeans法でクラスターするのもいいかも。
k0 <- kmeans(x1@exprs, 4) memb <- k0$cluster # ggpairs の色付けがうまくいかない plot(as.data.frame(x1@exprs), col=memb, cex=0.5)
外れ値解析を混ぜてもいいかもしれない。
library(Rlof) lof0 <- lof(dd, 20) # 95% tile より外れ値スコアが大きい細胞は外れ値と考える。 plot(as.data.frame(dd), col=(lof0 > quantile(lof0, 0.95))+1, cex=0.2)