r - 为我的逻辑回归模型生成逻辑曲线

Question

我想编写用于绘制逻辑回归模型的代码，即“S”形逻辑曲线。怎么可能做到，因为我有两个独立的协变量 ?我附上了我的数据集和我的模型的代码。先感谢您。

239 0.72    1
324.6   0.83    1
331.8   0.95    1
334.3   0.83    1
259.7   0.89    1
212.3   0.88    1
204.7   0.65    1
253.86  0.75    1
258.94  0.85    1
329.66  0.95    0
469.68  1.46    0
459.74  1.11    0
293.2   0.64    0
297.88  0.98    0
267.9   0.82    0
374.1   1.29    0
333.62  0.74    0


dat <- read.table("data.txt")
colnames(dat)<-c("press","v","gender")

# logostic regression
dat$gender <- factor(dat$gender)
mylogit<- glm(gender~press+v,data=dat,family="binomial")
summary(mylogit)

######## the code below are irrelevant to making plot, ignore if you want

mylogit$fitted.values

newdat <- data.frame(t(c(300,0.1)))
colnames(newdat)<-c("press","v")
   # this is your new dataset, we name it as "newdat"
pred <- predict(mylogit,newdata = newdat,type="response")
pred # the probability of being in class 1 will stored in this object

pred <- predict(mylogit,newdata = dat,type="response")
pred # the probability of being in class 1 will stored in this object
# accuracy
dat$pred <- 0
factor(dat$pred)
dat$pred[which(pred>0.5)] <- 1

table(dat$gender,dat$pred)

Answer 1

您有 2 个连续的非分类变量，因此逻辑曲线将是 3D 曲线。我将为您提供两种演示方式。

使用 persp产生真实3D平滑曲线的功能；

修复 v在多个值处，然后生成多个 2D 逻辑曲线(您称之为“S”形曲线)。

3D曲线

press_grid <- seq(200, 480, by = 5)
v_grid <- seq(0.6, 1.5, by = 0.1)
newdat <- data.frame(press = rep(press_grid, times = length(v_grid)), v = rep(v_grid, each = length(press_grid)))
pred <- predict.glm(mylogit, newdata = newdat, type="response")
z <- matrix(pred, length(press_grid))
persp(press_grid, v_grid, z, xlab = "pressure", ylab = "velocity", zlab = "predicted probability", main = "logistic curve (3D)", theta = 30, phi = 20)

您需要先生成一个二维网格。 newdat持有这个网格，你可以做plot(newdat)看到这个网格。然后通过调用 predict.glm(..., type = "response") 在这个网格上进行预测.结果pred是一个向量。要绘制它，请将其转换为矩阵 z ，然后调用 persp制作 3D 绘图。 xlab , ylab和 zlab是三轴的标签。参数theta和 phi用于调整您的视角。

在上面，press 的边际网格和 v是基于您原始数据的范围:range(dat$press)和 range(dat$v) .我们不会做出超出这个范围的预测太远。但即使在这个范围内，您也只有 17 个观测值。所以你仍然需要对情节持怀疑态度。

这是曲线:



二维曲线

这个玩具函数对于制作二维曲线很有用，v固定为某个级别:

curve_2D_fix_v <- function(model, v = 1, press_grid = seq(200, 480, by = 5), add = FALSE, col = "black") {
  newdat <- data.frame(press = press_grid, v = v)
  pred <- predict.glm(model, newdat, type = "response")
  if (add) lines(press_grid, pred, col = col) else {
    plot(press_grid, pred, xlab = "pressure", ylab = "predicted probability", type = "l", col = col, main = "logistic curve (2D)")
    abline(h = c(0, 0.5, 1), lty = 2, col = col)
    }
  }

如 add = FALSE , 它会打开一个新的绘图窗口；而它是 TRUE ，它在前一个窗口上绘图(但您有责任确保有这样一个窗口!)2D 绘图提供了更多信息，因为您可以在 0、0.5 和 1 处添加一条水平线。

让我们试一试:

curve_2D_fix_v(mylogit, v = 0.4, add = FALSE, col = "black")
curve_2D_fix_v(mylogit, v = 0.6, add = TRUE, col = "red")
curve_2D_fix_v(mylogit, v = 0.8, add = TRUE, col = "green")
curve_2D_fix_v(mylogit, v = 1, add = TRUE, col = "blue")
curve_2D_fix_v(mylogit, v = 1.2, add = TRUE, col = "cyan")
curve_2D_fix_v(mylogit, v = 0.4, add = TRUE, col = "yellow")

这是曲线:



讨论

在这两个图中，我们看到 gender 之间的关系(预测概率)和 v (速度)不是很强。在二维图中，几乎所有 v 的值产生相同的曲线。另一方面，press (压力)是一种强烈的影响。

回到你的模型:

> summary(mylogit)
Coefficients:
            Estimate Std. Error z value Pr(>|z|)  
(Intercept)  8.08326    4.45463   1.815   0.0696 .
press       -0.02575    0.01618  -1.591   0.1115  
v           -0.15385    4.83824  -0.032   0.9746  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

你可以看到v一点都不重要!虽然严格来说，press在 0.1 水平上也不显着。 所以这是一个非常弱的模型 .我建议你删除变量 v并使用 press 再次创建模型作为唯一的变量。