我正在尝试使用drc
重现ggplot2
图。这是我的第一次尝试(下面给出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)
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()
编辑
给定了here,我遵循了方法(进行了一些更改)。
最佳答案
注意,您可以跳至最后一段,以获取简单的答案。该答案的其余部分记录了我如何得出该解决方案
查看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包含拟合模型线的坐标,因此我们可以使用它来ggplot 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()
该版本与您在ggplot中使用predict(object = chickweed.m1)创建的版本相同。因此,差异不在于模型线,而在于绘制数据点的位置。通过将函数的最后一行从
invisible(retData)
更改为list(retData, plotPoints)
,我们可以从drc::: plot.drc导出数据点。为了方便起见,我将drc::: plot.drc的整个代码复制到了一个新函数中。请注意,如果您希望复制此步骤,则drcplot调用的一些函数不会在drc命名空间中导出,因此drc:::
必须在对函数parFct
,addAxes
,brokenAxis
和makeLegend
的所有调用之前。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()
关于r - drc::带有ggplot2的drc图,我们在Stack Overflow上找到一个类似的问题:https://stackoverflow.com/questions/38169765/