在 ggplot2 中绘制 drc 模型; seq( ) 问题

Plotting drc model in ggplot2; issue with seq( )

我的模型在 ggplot2 中绘制时不会继续朝向渐近线,尽管它在 R 基础图形中是这样。在 ggplot2 中,它停在 X 轴上的某些点(下面包含图像),我 90% 确定这与 seq().

有关

我正在使用 predict() 从 drm(dose-response 包)logit 模型转换数据。在基础图形中,sigmoidal 曲线看起来很棒,在 ggplot2 中,没那么多:

library(drc)
library(ggplot2)

将数据拟合到 logit 模型:

mod1 <- drm(probability ~ (dose), weights = total, data = mydata1,     type ="binomial", fct = LL.2())

plot(mod1,broken=FALSE,type="all",add=FALSE, col= "purple", xlim=c(0, 10000))

基本图形 2 参数 logit 图像:

使用作者演示(链接如下)中的代码提取模型数据,我有:

newdata1 <-expand.grid(dose=exp(seq(log(0.5),log(100),length=200)))

pm1<- predict(mod1, newdata=newdata1,interval="confidence")

newdata1$p1 <-pm1[,1]

newdata1$pmin1 <-pm1[,2]

newdata1$pmax1 <- pm1[,3]

最后是 ggplot2 图形:

p1 <- ggplot(mydata1, aes(x=dose01,y=probability))+
  geom_point()+
   geom_ribbon(data=newdata1, aes(x=dose,y=p1,
 ymin=pmin1,ymax=pmax1),alpha=0.2,color="blue",fill="pink") +
   geom_step(data=newdata1, aes(x=dose,y=p1))+
  coord_trans(x="log") +  #creates logline for x axis
  xlab("dose")+ylab("response")

图片 1&2 和 3&4 显示了我的绘图中的差异,具体取决于序列:

当以下用于seq()时,图形被拉过数据! (seq(log(0.5),**log(10000)**,length=200)))(图片 1 和 2)

尽管研究了,我还是不明白seq()我的剧情怎么了?

似乎seq中的第一项定义了下限,那么第三项定义的是什么?您可以在图像 3 和 4 中看到这一点 - 图形很不错;我有点混淆了这个问题,但它仍然没有继续走向 infin。这是一个小问题,因为我将 co-plotting 8 个 logit 模型。

对于使用 ggplot2 绘制 drc/drm 模型有困难的任何人,以下帖子非常有帮助:搜索 绘制-dose-response-曲线-with-ggplot2-and-drc
这个标题:plotting-dose-response-curves-with-ggplot2-and-drc

我已经按照DRC的作者的说明做了,可以在他的文章的支持信息中找到——上面使用了部分代码。文章标题:Dose-Response 使用 R 进行分析,Christopher Ritz。 PlosOne.

数据:

> dput(mydata1)

structure(list(dose = c(25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 
25L, 25L, 25L, 25L, 25L, 25L, 25L, 75L, 75L, 75L, 75L, 75L, 75L, 
75L, 75L, 75L, 75L, 75L, 75L, 75L, 75L, 75L, 100L, 100L, 100L, 
100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 100L, 
100L, 150L, 150L, 150L, 150L, 150L, 150L, 150L, 150L, 150L, 150L, 
150L, 150L, 150L, 150L, 150L, 200L, 200L, 200L, 200L, 200L, 200L, 
200L, 200L, 200L, 200L, 200L, 200L, 200L, 200L, 200L), total = c(25L, 
25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 
25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 
25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 
25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 
25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 
25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L, 25L), affected = c(2L, 
0L, 0L, 0L, 0L, 1L, 1L, 2L, 2L, 4L, 1L, 1L, 4L, 0L, 10L, 0L, 
1L, 0L, 1L, 0L, 3L, 0L, 4L, 2L, 0L, 2L, 0L, 1L, 2L, 3L, 2L, 0L, 
2L, 0L, 0L, 4L, 0L, 1L, 2L, 3L, 0L, 21L, 1L, 3L, 1L, 2L, 7L, 
0L, 0L, 0L, 0L, 8L, 7L, 3L, 7L, 2L, 2L, 10L, 3L, 4L, 0L, 7L, 
0L, 3L, 3L, 20L, 25L, 22L, 23L, 22L, 18L, 14L, 20L, 20L, 21L), 
    probability = c(0.08, 0, 0, 0, 0, 0.04, 0.04, 0.08, 0.08, 
    0.16, 0.04, 0.04, 0.16, 0, 0.4, 0, 0.04, 0, 0.04, 0, 0.12, 
    0, 0.16, 0.08, 0, 0.08, 0, 0.04, 0.08, 0.12, 0.08, 0, 0.08, 
    0, 0, 0.16, 0, 0.04, 0.08, 0.12, 0, 0.84, 0.04, 0.12, 0.04, 
    0.08, 0.28, 0, 0, 0, 0, 0.32, 0.28, 0.12, 0.28, 0.08, 0.08, 
    0.4, 0.12, 0.16, 0, 0.28, 0, 0.12, 0.12, 0.8, 1, 0.88, 0.92, 
    0.88, 0.72, 0.56, 0.8, 0.8, 0.84)), .Names = c("dose", "total", 
"affected", "probability"), row.names = c(NA, -75L), class = "data.frame")

请参阅 ?seq 以了解 seq 的条款的作用。 seq(log(.5), log(1000), length=200) 使 200 个数字均匀分布,从 log(0.5) 到 log(1000)。如果你不命名第三个参数(或指定by=XYZ),它是数字之间的间距。

因此,当您执行 seq(log(0.5), log(1000), length=200) 时,它将计算从 log(0.5) 到 log(1000) 的 200 个点的拟合度。

我认为 "go out to infinity" 你的意思是你希望线条离开情节的边缘,而不是像你的链接图片那样停在边缘之前。默认情况下,ggplot 将尝试确保您绘制的所有内容都适合您的绘图,因此它会将轴延伸到数据范围之外。

如果要限制的话就用+ xlim(c(lower, upper)).

(我在安装 drc 时遇到问题,因此无法重现您的示例;这是一个玩具)

x = seq(0.5, 100, length=200)
df <- data.frame(x=x, y=x^2)
ggplot(df, aes(x=x, y=y)) + geom_line() + coord_trans(x="log")

在上面,线延伸到 100(如预期的那样)并且限制超出了一点。如果我想让线触及图的边缘,那么我可以(例如)将 x 限制精确地剪裁在 100 - 使用 limx 参数 coord_trans:

ggplot(df, aes(x=x, y=y)) + geom_line() + coord_trans(x="log", limx=c(0.5, 100))

因此,当您绘制模型时,请确定 x 轴的边界,并确保您根据这些值预测所有模型。然后将 x 限制限制在这些范围内,这些行将显示为 "go to infinity".

您将 dose 值误认为是 predict()aes(x) 值:

log10000 <- exp(seq(log(0.5), log(10000), length=200))
log1000 <- exp(seq(log(0.5), log(1000), length=200))

log10000df <- as.data.frame(cbind(dose = log10000, predict(mod1, data.frame(dose = log10000), interval="confidence")))
log1000df <- as.data.frame(cbind(dose = log1000, predict(mod1, data.frame(dose = log1000), interval="confidence")))

 ## a common part
p0 <- ggplot(mydata1, aes(x = dose, y = probability)) +
  geom_point() + coord_trans(x="log") + 
  xlab("dose") + ylab("response") + xlim(0.5, 10001)

p10000 <- p0 + geom_line(data = log10000df, aes(x = dose, y = Prediction)) +
  geom_ribbon(data = log10000df, aes(x = dose, y = Prediction, ymin = Lower, ymax = Upper),
              alpha = 0.2, color = "blue", fill = "pink")
  
p1000 <- p0 + geom_line(data = log1000df, aes(x = dose, y = Prediction)) +
  geom_ribbon(data = log1000df, aes(x = dose, y = Prediction, ymin = Lower, ymax = Upper),
              alpha = 0.2, color = "blue", fill = "pink")