logo标签列表

使用ggplot2绘制DRC图。

rggplot2drc
6

我正在尝试使用ggplot2重新制作drc图。以下是我的第一次尝试（提供了MWE）。然而，我的ggplot2与基本的R绘图有点不同。我想知道我是否漏掉了什么？

library(drc)
chickweed.m1 <- drm(count~start+end, data = chickweed, fct = LL.3(), type = "event")

plot(chickweed.m1, xlab = "Time (hours)", ylab = "Proportion germinated", 
xlim=c(0, 340), ylim=c(0, 0.25), log="", lwd=2, cex=1.2)

enter image description here

library(data.table)
dt1 <- data.table(chickweed)

dt1Means1 <- dt1[, .(Germinated=mean(count)/200), by=.(start)]
dt1Means2 <- dt1Means1[, .(start=start, Germinated=cumsum(Germinated))]
dt1Means  <- data.table(dt1Means2[start!=0], Pred=predict(object=chickweed.m1))

library(ggplot2)
ggplot(data= dt1Means, mapping=aes(x=start, y=Germinated)) + 
    geom_point() +
    geom_line(aes(y = Pred)) +
    lims(y=c(0, 0.25)) +
    theme_bw()

这里输入图片描述

编辑过的内容

我按照这里给出的方法（做了一些更改）进行了操作。

1
你可以随时查看getAnywhere(plot.drc)，以了解在基本绘图中如何计算数据。 - Andy W
你想在ggplot中改变什么？两个图表中的数据似乎是相同的。 - Andrew Brēza
@MYaseen208：那就给他赏金吧。没有必要再保持它的开放状态了。 - IRTFM
@AndrewBrēza：请查看dww的答案。 - MYaseen208
当然，@42-，悬赏金是给dww的。 - MYaseen208
3个回答
9

简单答案请跳转至最后一段。本答案的其余部分是为了展示我得出答案的过程。

查看drc:::plot.drc的代码，我们可以看到最终一行通过隐式返回一个data.frame retData。

function (x, ..., add = FALSE, level = NULL, type = c("average", 
                                                      "all", "bars", "none", "obs", "confidence"), broken = FALSE, 
          bp, bcontrol = NULL, conName = NULL, axes = TRUE, gridsize = 100, 
          log = "x", xtsty, xttrim = TRUE, xt = NULL, xtlab = NULL, 
          xlab, xlim, yt = NULL, ytlab = NULL, ylab, ylim, cex, cex.axis = 1, 
          col = FALSE, lty, pch, legend, legendText, legendPos, cex.legend = 1, 
          normal = FALSE, normRef = 1, confidence.level = 0.95) 
{
  # ...lot of lines omitted...
  invisible(retData)
}

retData包含拟合模型线的坐标，因此我们可以使用它来绘制与plot.drc相同的模型

pl <- plot(chickweed.m1, xlab = "Time (hours)", ylab = "Proportion germinated", 
        xlim=c(0, 340), ylim=c(0, 0.25), log="", lwd=2, cex=1.2)
names(pl) <- c("x", "y")
ggplot(data= dt1Means, mapping=aes(x=start, y=Germinated)) + 
  geom_point() +
  geom_line(data=pl, aes(x=x, y = y)) +
  lims(y=c(0, 0.25)) +
  theme_bw()

enter image description here

这与您使用 predict(object=chickweed.m1) 在 ggplot 中创建的版本相同。 因此，不同之处不在于模型线，而在于数据点的绘制位置。 我们可以通过将函数的最后一行从 invisible(retData) 更改为 list(retData, plotPoints) 来从 drc:::plot.drc 导出数据点。 为方便起见，我将整个 drc:::plot.drc 代码复制到了一个新函数中。请注意，如果您希望重复此步骤，则需要调用 drcplot 中的一些函数，这些函数未在 drc 命名空间中导出，因此必须在对函数 parFct、addAxes、brokenAxis 和 makeLegend 的所有调用之前添加前缀 drc:::。

drcplot <- function (x, ..., add = FALSE, level = NULL, type = c("average", 
                                                      "all", "bars", "none", "obs", "confidence"), broken = FALSE, 
          bp, bcontrol = NULL, conName = NULL, axes = TRUE, gridsize = 100, 
          log = "x", xtsty, xttrim = TRUE, xt = NULL, xtlab = NULL, 
          xlab, xlim, yt = NULL, ytlab = NULL, ylab, ylim, cex, cex.axis = 1, 
          col = FALSE, lty, pch, legend, legendText, legendPos, cex.legend = 1, 
          normal = FALSE, normRef = 1, confidence.level = 0.95) 
{
  # ...lot of lines omitted...
  list(retData, plotPoints)
}

并且使用您的数据运行此程序。

pl <- drcplot(chickweed.m1, xlab = "Time (hours)", ylab = "Proportion germinated", 
          xlim=c(0, 340), ylim=c(0, 0.25), log="", lwd=2, cex=1.2)

germ.points <- as.data.frame(pl[[2]])
drc.fit <- as.data.frame(pl[[1]])
names(germ.points) <- c("x", "y")
names(drc.fit) <- c("x", "y")

现在，使用ggplot2绘制这些图形可以得到您想要的结果。
ggplot(data= dt1Means, mapping=aes(x=start, y=Germinated)) + 
  geom_point(data=germ.points, aes(x=x, y = y)) +
  geom_line(data=drc.fit, aes(x=x, y = y)) +
  lims(y=c(0, 0.25)) +
  theme_bw()

在此输入图片描述

最后，将此图的数据点值（germ.points）与您原始的ggplot图（dt1Means）进行比较，可以看出差异的原因。您在dt1Means中计算的点相对于plot.drc中的时间段提前了一个时间段。换句话说，plot.drc将事件指派给它们发生的时间段的结束时间，而您将萌发事件指派给它们发生的时间间隔的开始时间。您可以通过简单地调整来解决这个问题，例如使用

dt1 <- data.table(chickweed)
dt1[, Germinated := mean(count)/200, by=start]
dt1[, cum_Germinated := cumsum(Germinated)]
dt1[, Pred := c(predict(object=chickweed.m1), NA)]  # Note that the final time period which ends at `Inf` can not be predicted by the model, therefore added `NA` in the final row

ggplot(data= dt1, mapping=aes(x=end, y=cum_Germinated)) + 
  geom_point() +
  geom_line(aes(y = Pred)) +
  lims(y=c(0, 0.25)) +
  theme_bw()
感谢@dww提供非常有帮助的答案。根据您的回答，我理解dt1Means2 <- dt1Means1[, .(start=start, Germinated=c(0,cumsum(Germinated)))]就足够了。然而，它并没有给我所需的输出。请问您有什么想法吗？ - MYaseen208
仅使用 dt1Means2 <- dt1Means1[, .(start=start, Germinated=c(0,cumsum(Germinated)))] 并不是解决方案。它会给出以下信息 In as.data.table.list(jval) : Item 1 is of size 35 but maximum size is 36 (recycled leaving a remainder of 1 items)。 - MYaseen208
@MYaseen208，这个消息只是一个警告，而不是错误。你可以安全地忽略它。话虽如此，这只是我能够做出的最小更改，以获得你想要的正确图形 - 但并不是最整洁的方式。更好的方法是使用结束时间而不是开始时间来正确构建你的数据（这就是我在答案中建议的）。我已经进行了编辑，展示了如何做到这一点。 - dww
我正在苦恼如何获取curveid = Temp时的 ggplot2 版本。如果您能查看此问题，我将不胜感激。谢谢。 - MYaseen208
2
从@dww的回答中得到灵感后，我不得不对原始代码进行两个小改变。只需将start!=0替换为end!=Inf即可。
dt1Means1 <- dt1[, .(Germinated=mean(count)/200), by=.(start, end)]
dt1Means  <- data.table(dt1Means2[start!=0], Pred=predict(object=chickweed.m1))

给出正确的图表。
1
我很喜欢dww提出的解决方案。我可以建议对此解决方案进行概括。通过将以下行添加到自己编写的drc:::plotdrc()版本中，您可以将解决方案推广。该函数接受drc:::plotdrc()函数的输入，但输出具有与原始函数默认基本绘图输出相同的ggplot对象规格。
只需用以下内容替换invisible(retData, plotPoints) 
result <- list(retData, plotPoints) 
points <- as.data.frame(result[[2]])
drc.fit <- as.data.frame(result[[1]]) 
names(points) <- c("x", "y")
names(drc.fit) <- c("x", "y")` 

gg_plot <- ggplot2::ggplot(data=points, aes(x = x, y = y)) + 
geom_point() +
geom_line(data=drc.fit, aes(x = x, y = y)) +
scale_x_continuous(trans='log10', limits = xlim) +
ylab(ylab) +
xlab(xlab) +
lims(y = ylim) +
theme_bw()
return(gg_plot)`
网页内容由stack overflow 提供, 点击上面的
可以查看英文原文，
原文链接