iDEGES正規化

iDEGES正規化法とは

RPKM正規化法の章でもお話ししましたが、発現変動遺伝子(Diffrentially Expressed Gene: DEG)は正規化時には邪魔になります。そこで、正規化する前にDEGを取り除いてしまうDEG elimination strategy(DEGES:でげす)法が考案されました(Sun et al., BMC Bioinfomatics, 2013)。

2群間比較・3反復のデータの場合

1列目から遺伝子名, サンプルA_rep1, サンプルA_rep2, サンプルA_rep3, サンプルB_rep1, サンプルB_rep2, サンプルB_rep3のカウント数が並んでいるデータ(count.txt)を想定しています。

in_f <- "count.txt"
out_f1 <- "non_vs_snf_iDEGESedgeR_DE.txt"
out_f2 <- "non_vs_snf_iDEGESedgeR_DE_FDR.png"

param_G1 <- 3
param_G2 <- 3
param_DEmethod <- "edger"
param_FDR <- 0.001
param_fig <- c(400, 380)

library(TCC)

data <- read.table(in_f, header=TRUE, row.names=1, sep="\t", quote="")
data.cl <- c(rep(1, param_G1), rep(2, param_G2))
tcc <- new("TCC", data, data.cl)

tcc <- calcNormFactors(tcc, norm.method="tmm", test.method="edger", iteration=3, FDR=0.1, floorPDEG=0.05)
tcc <- estimateDE(tcc, test.method=param_DEmethod, FDR=param_FDR)
result <- getResult(tcc, sort=FALSE)

tmp <- cbind(rownames(tcc$count), tcc$count, result)
write.table(tmp, out_f1, sep="\t", append=F, quote=F, row.names=F)

png(out_f2, pointsize=13, width=param_fig[1], height=param_fig[2])
plot(tcc, FDR=param_FDR, xlim=c(-3, 13), ylim=c(-10, 10))
legend("bottomright", c(paste("DEG(FDR<", param_FDR, ")", sep=""), "non-DEG"), col=c("magenta", "black"), pch=20)
dev.off()

2群間比較・1反復のデータの場合

1列目から遺伝子名, サンプルA_rep1, サンプルB_rep1のカウント数が並んでいるデータ(count.txt)を想定しています。

in_f <- "count.txt"
out_f1 <- "non_vs_snf_iDEGES_DE.txt"
out_f2 <- "non_vs_snf_iDEGES_DE_FDR.png"
param_G1 <- 1
param_G2 <- 1

library(TCC)

data <- read.table(in_f, header=TRUE, row.names=1, sep="\t", quote="")
data.cl <- c(rep(1, param_G1), rep(2, param_G2))
tcc <- new("TCC", data, data.cl)
tcc <- calcNormFactors(tcc, iteration=3)
normalized <- getNormalizedData(tcc)
tmp <- cbind(rownames(normalized), normalized)
write.table(tmp, out_f1, sep="\t", append=F, quote=F, row.names=F)

png(out_f2, pointsize=13, width=param_fig[1], height=param_fig[2])
plot(tcc, FDR=param_FDR, xlim=c(-3, 13), ylim=c(-10, 10))
legend("bottomright", c(paste("DEG(FDR<", param_FDR, ")", sep=""), "non-DEG"), col=c("magenta", "black"), pch=20)
dev.off()

3群間比較・3反復のデータの場合の場合

1列目から遺伝子名, サンプルA_rep1, サンプルA_rep2, サンプルA_rep3, サンプルB_rep1, サンプルB_rep2, サンプルB_rep3, サンプルC_rep1, サンプルC_rep2, サンプルC_rep3のカウント数が並んでいるデータ(count.txt)を想定しています。もし、4群間を比較したい場合は、param_G4 <- 3を追加し、data.cl の行にも追加することで、解析は可能です。

in_f <- "count.txt"
out_f <- "hoge1.txt"
param_G1 <- 3
param_G2 <- 3
param_G3 <- 3

library(TCC)

data <- read.table(in_f, header=TRUE, row.names=1, sep="\t", quote="")
data.cl <- c(rep(1, param_G1), rep(2, param_G2), rep(3, param_G3))
tcc <- new("TCC", data, data.cl)
tcc <- calcNormFactors(tcc, norm.method="tmm", test.method="edger", iteration=3, FDR=0.1, floorPDEG=0.05)
normalized <- getNormalizedData(tcc)

tmp <- cbind(rownames(data), normalized)
write.table(tmp, out_f, sep="\t", append=F, quote=F, row.names=F)

Reference

1. 東大・門田さんのHP