从斜率和截距绘制回归（格子或ggplot2）-Java 学习之路

我有一个微阵列数据集，我在其上进行了 limma lmFit() 测试 . 如果你没有't heard of it before, it'这是一个功能强大的线性模型包，可以测试大于20k基因的差异基因表达 . 您可以从模型中提取每个基因的斜率和截距 .

My problem is: given a table of slope and intercept values, how do I match a plot (I don't mind either ggplot2's geom_abline, lattice's panel.abline, or an alternative if necessary) with its corresponding slope and intercept?

我的表（称之为“slopeInt”）截取为第1列，斜率为第2列，并且具有与基因名称对应的行名称 . 他们的名字看起来像这样：

"202586_at"   "202769_at"   "203201_at"   "214970_s_at" "219155_at"

这些名称与另一个表（“数据”）中的基因名称相匹配，其中包含有关我的样本的一些详细信息（我有24个具有不同ID和时间/治疗组合的样本）和基因表达值 .

它是长格式的，基因名称（如上所示）每24行重复一次（相同基因的不同表达水平，对于我的每个样本）：

ID    Time  Treatment  Gene_name    Gene_exp
...   ...   ...        ...          ...

我有八个基因我有兴趣绘制，我的 Data$Gene_name 中的名字与我的 slopeInt 表的行名相匹配 . 我也可以将两个表合并在一起，这不是问题 . 但我尝试了以下两种方法，给我带有适当回归的每个基因的图形图，但无济于事：

使用 ggplot2 ：

ggplot(Data, aes(x = Time, y = Gene_exp, group = Time, color = Treatment)) +
facet_wrap(~ Gene_name, scales = "free_x") +
geom_point() +
geom_abline(intercept = Intercept, slope = Time), data = slopeInt) +
theme(panel.grid.major.y = element_blank())`

并且还使用 Lattice ：

xyplot(Gene_exp ~ Time| Gene_name, Data, 
   jitter.data = T,
   panel = function(...){
     panel.xyplot(...)
     panel.abline(a = slopeInt[,1], b = slopeInt[,2])},
   layout = c(4, 2))

我在实际的 geom_abline() 和 panel.abline() 参数中尝试了多种其他方法，包括一些for循环，但我没有R经验，我无法使它工作..我也可以使用宽格式的数据文件（单独的列对于每个基因） .

任何帮助和进一步的方向将不胜感激！

以下是可重现示例的一些代码：

Data <- data.frame(
 ID = rep(1:24, 8),
 Time = (rep(rep(c(1, 2, 4, 24), each = 3), 8)),
 Treatment = rep(rep(c("control", "smoking"), each = 12), 8),
 Gene_name = rep(c("202586_at", "202769_at", "203201_at", "214970_s_at",
   "219155_at", "220165_at", "224483_s_at", "227559_at"), each = 24),
 Gene_exp = rnorm(192))

slopeInt <- data.frame(
 Intercept = rnorm(8),
 Slope = rnorm(8))


row.names(slopeInt) <- c("202586_at", "202769_at", "203201_at",
"214970_s_at", "219155_at", "220165_at", "224483_s_at", "227559_at")

2 回答

1
有格子，这应该工作
```
xyplot(Gene_exp ~ Time| Gene_name, Data, slopeInt=slopeInt,
   jitter.data = T,
   panel = function(..., slopeInt){
     panel.xyplot(...)
     grp <- trellis.last.object()$condlevels[[1]][which.packet()]
     panel.abline(a = slopeInt[grp,1], b = slopeInt[grp,2])
   },
   layout = c(4, 2)
)
```
在生成样本数据之前使用 set.seed(15) 得到以下图表

"trick"这里是使用 trellis.last.object()$condlevels 来确定我们当前在哪个条件块 . 然后我们使用该信息从我们现在通过参数传入的附加数据中提取正确的斜率信息 . 我认为有一种更优雅的方法来确定条件变量的当前值，但是如果有的话我现在还记不住了 .
回复于 2024-04-28T02:20:14+08:00

如果您将 Gene_name 指定为 slopeInt 中的列，则它可以正常工作[据我所知您希望它] . 另请注意ggplot调用的其他一些更改 .

slopeInt$Gene_name <- rownames(slopeInt)
ggplot(Data, aes(x = Time, y = Gene_exp, color = Treatment)) +
   facet_wrap(~ Gene_name, scales = "free_x") +
   geom_point() +
   geom_abline(aes(intercept = Intercept, slope = Slope), data = slopeInt) +
   theme(panel.grid.major.y = element_blank())

回复于 2024-04-28T02:20:14+08:00

从斜率和截距绘制回归（格子或ggplot2）

2 回答

相关问题