（sigmoid）曲线拟合r中的glm

Question

我希望可视化两个类别（发射机）的响应变量、检测概率 (P.det) 和预测变量（距离）之间的关系，显示误差线并通过平均数据点绘制（S 形）曲线。

数据集是这样的：

df <- structure(list(distance = c(50L, 100L, 200L, 300L, 400L, 50L, 
100L, 200L, 300L, 400L), Transmitter = structure(c(1L, 1L, 1L, 
1L, 1L, 2L, 2L, 2L, 2L, 2L), .Label = c("CT", "PT"), class = "factor"), 
    P.det = c(0.918209097, 0.88375438, 0.709288774, 0.534977488, 
    0.341724516, 0.828123952, 0.822201191, 0.543289433, 0.352886247, 
    0.10082457), st.error = c(0.01261614, 0.014990469, 0.024136478, 
    0.027311169, 0.026941438, 0.018663591, 0.019420587, 0.02754911, 
    0.026809247, 0.017041264), ly = c(0.905592958, 0.868763911, 
    0.685152295, 0.50766632, 0.314783078, 0.809460361, 0.802780604, 
    0.515740323, 0.326077, 0.083783306), uy = c(0.930825237, 
    0.898744849, 0.733425252, 0.562288657, 0.368665955, 0.846787544, 
    0.841621778, 0.570838544, 0.379695494, 0.117865833), Valid.detections = c(18, 
    12.5472973, 8.608108108, 4.287162162, 2.158783784, 12.46959459, 
    7.956081081, 4.550675676, 1.682432432, 0.39527027), False.detections = c(0.388513514, 
    0.550675676, 0.368243243, 0.263513514, 0.131756757, 0.533783784, 
    0.385135135, 0.277027027, 0.182432432, 0.14527027)), .Names = c("distance", 
"Transmitter", "P.det", "st.error", "ly", "uy", "Valid.detections", 
"False.detections"), class = "data.frame", row.names = c(NA, 
-10L))

我设法完成了前两部分，但卡在了最后一部分。用误差线绘制图形的代码：

library(lattice)
    prepanel.ci <- function(x, y, ly, uy, subscripts, ...)
    {
      x <- as.numeric(x)
      ly <- as.numeric(ly[subscripts])
      uy <- as.numeric(uy[subscripts])
      list(ylim = range(y, uy, ly, finite = TRUE))
    }

    panel.ci <- function(x, y, ly, uy, subscripts, pch = 16, ...)
    {
      x <- as.numeric(x)
      y <- as.numeric(y)
      ly <- as.numeric(ly[subscripts])
      uy <- as.numeric(uy[subscripts])
      panel.arrows(x, ly, x, uy, col = "black",
                   length = 0.25, unit = "native",
                   angle = 90, code = 3)
      panel.xyplot(x, y, pch = pch, ...)
    }

xyplot(P.det~distance, type=c("p","g"),
       ylim=c(0,1),
       ylab="Detection probability", xlab="Distance (m)", 
       group=Transmitter,
       data=df,
       ly = df$ly,
       uy = df$uy,
       prepanel = prepanel.ci,
       panel = panel.superpose,
       panel.groups = panel.ci,
       col=c(1,1),
       layout=c(1,1),
       between=list(x=2),
       scales=list(x=list(alternating=c(1,1), tck=c(1,0)),y=list(alternating=c(1,1), tck=c(1,0))), # ticks inside = tck=c(-1,0)
       aspect=1,
       main="Detection probability vs distance per transmitter type",
)

我之所以在标题中注明“glm”，是因为数据分析是使用 lme4 包的二项式 glm() 进行的。

我注意到另一个与我的线程相似的线程：find the intersection of abline with fitted curve，但不同之处在于，虽然我的图表也是基于每 1 x 1 y，但我的 glm 是基于每 x y 的多个。因此，在此线程中遵循相同的代码会返回一个错误，指出长度不相等。它似乎也不适用于“xyplot”。

谢谢

Answer 1

在 Lattice 中，您可以将平滑添加到自定义面板功能中。您可以将其更改为包含

panel.ci <- function(x, y, ly, uy, subscripts, type="p", pch = 16, ...)
{
  x <- as.numeric(x)
  y <- as.numeric(y)
  ly <- as.numeric(ly[subscripts])
  uy <- as.numeric(uy[subscripts])
  panel.arrows(x, ly, x, uy, col = "black",
               length = 0.25, unit = "native",
               angle = 90, code = 3)
  panel.xyplot(x, y, pch = pch, type=type, ...)

  #calculate smooth curve
  gg <- glm(y~x, family="binomial")
  gx <- seq(min(x), max(x), length.out=100)
  panel.xyplot(gx, predict(gg, data.frame(x=gx), type="response"),
      pch = pch, type="l", ...)
}

这里我们自己做glm并在情节上绘制响应。

Answer 2

使用ggplot 相当简单：

library(ggplot2)
ggplot(data = df, aes(x = distance, y = P.det, colour = Transmitter)) +
  geom_pointrange(aes(ymin = P.det - st.error, ymax = P.det + st.error)) +
  geom_smooth(method = "glm", family = binomial, se = FALSE)

关于glmwarning 消息，参见例如here。