0

そこで、脳領域ごとの遺伝子発現値を含むデータフレームに 2 つの .csv ファイルを読み込んでいます。どちらも次のように見えます。

   reg1   reg2     reg3   reg4   reg5    reg6    reg6    reg7     reg8
1  0.2036 0.6400  -2.2647 0.5345 0.9812 -0.0390 -1.4486 -0.1340  -1.9157
2 -0.0308 0.2991  -1.8990 0.7063 0.8807 -0.0568 -1.5150 -0.8000  -1.5115
3 -0.0242 0.6358  -2.1675 0.8432 0.8149 -0.0955 -1.3789 -0.2266  -1.5513
4 -0.8990 0.5613  -2.0289 0.8544 1.1060 -0.3846 -1.4800 -0.0717  -1.4325
5 -0.6438 0.8719  -2.2708 0.6790 1.0250 -0.3035 -1.1125 -0.3746  -1.5219
6 -0.7945 0.5943  -1.8478 0.5457 0.7989 -0.1093 -2.3242  0.2995  -1.8066

以下を使用して p 値で発現率を比較し、領域ごとに遺伝子の t 検定を実行します (領域は両方の .csv ファイルで同じです)。

geneA = read.csv("geneA.csv")
geneB = read.csv("geneB.csv")
test.result = mapply(t.test, geneA, geneB)

#store p-values in a column with regions in a separate column:
p.values = stack(mapply(function(x, y) t.test(x,y)$p.value, geneA, geneB))

#order p-values
ordered = p.values[with(p.values, order(values)),]

ここで、これらの p 値でやりたいことは、次のようなマップに表示することです (p 値は構成されています)。

                  geneB

                 reg1 reg2 reg3 reg4 ...
     reg1        .02  .02  .01  .12
     reg2        .03  .03  .05  .02

遺伝子A reg3 .01 .05 .05 .05 reg4 .01 .06 .21 .02 ...などなど。助言がありますか？グラフィカルなものが好ましいです。

4

1 に答える 1

4

このような？

m1 <- as.matrix(read.table(text="   reg1   reg2     reg3   reg4   reg5    reg6    reg6    reg7     reg8
1  0.2036 0.6400  -2.2647 0.5345 0.9812 -0.0390 -1.4486 -0.1340  -1.9157
2 -0.0308 0.2991  -1.8990 0.7063 0.8807 -0.0568 -1.5150 -0.8000  -1.5115
3 -0.0242 0.6358  -2.1675 0.8432 0.8149 -0.0955 -1.3789 -0.2266  -1.5513
4 -0.8990 0.5613  -2.0289 0.8544 1.1060 -0.3846 -1.4800 -0.0717  -1.4325
5 -0.6438 0.8719  -2.2708 0.6790 1.0250 -0.3035 -1.1125 -0.3746  -1.5219
6 -0.7945 0.5943  -1.8478 0.5457 0.7989 -0.1093 -2.3242  0.2995  -1.8066",header=TRUE))

m2 <- matrix(rnorm(54),nrow=6)

fun <- Vectorize(function(i,j) t.test(m1[,i],m2[,j])$p.value)
res <- outer(1:9,1:9,FUN = "fun")

image(1:9,1:9,res,axes=FALSE,xlab="m1",ylab="m2")
axis(1, at = 1:9,labels=colnames(m1))
axis(2, at = 1:9,labels=colnames(m1))

ここに画像の説明を入力

編集

ggplot2 を使用したプロットは次のとおりです。

colnames(res) <- colnames(m1)
res <-as.data.frame(res)
res$group <- colnames(m1)

library(reshape2)
res <- melt(res,id="group")

library(ggplot2)
p <- ggplot(res, aes(x=group, y=variable)) + 
  geom_tile(aes(fill = value), colour =   "yellow") + 
  scale_fill_gradient(low = "yellow", high = "red", name="p-value") +
  geom_text(aes(label=format(value,digits=2))) +
  labs(x="m1",y="m2")

print(p)

ここに画像の説明を入力

于 2013-05-10T16:06:05.790 に答える