1

在 ggplot2 中绘制时，我的模型不会继续朝向渐近线，尽管它在 R 基础图形中会如此。在 ggplot2 中，它停在 X 轴上的某些点（如下图所示），我 90% 确定这与seq().

我正在使用predict()从 drm（剂量响应包）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 的图像：

基本图形 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：

2

当以下用于 时seq()，图形被拉出数据！ (seq(log(0.5),**log(10000)**,length=200))) （图 1 和 2）

seq()尽管研究了它，但我不明白。我的阴谋怎么了？

seq 中的第一项似乎定义了下限，但是第三项定义的是什么？您可以在图像 3 和 4 中看到这一点 - 图形很不错；我有点混淆了这个问题，但它仍然没有继续朝着infin方向发展。这是一个小问题，因为我将共同绘制 8 个 logit 模型。

对于使用 ggplot2 绘制 drc/drm 模型时遇到问题的任何人，以下帖子非常有帮助：搜索 plotting-dose-response-curves-with-ggplot2-and-drc
和此标题： plotting -dose-response-curves- with- ggplot2-and-drc

我已经按照 DRC 的作者的说明进行操作，可以在他的文章的支持信息中找到 - 上面使用了部分代码。文章标题：使用 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")
4

2 回答 2

1

您误认为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")

在此处输入图像描述 在此处输入图像描述

于 2016-10-12T05:37:52.313 回答
0

请参阅有关?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 点的适合度。

我认为“走向无穷大”是指您希望线离开情节的边缘，而不是像链接图片中那样在边缘之前停止。默认情况下，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 限制限制在这些范围内，这些线条将看起来“趋于无穷大”。

于 2016-10-11T00:08:04.027 回答