将名称保存为运行众多单变量回归模型的函数的列表

Question

我正在使用包含大约 40 个变量的数据集构建逻辑回归模型。在构建这些类型的模型时，我使用的第一步是使用 DV 单变量运行每个变量（Hosmer、Lemeshow 和 Sturdivant，2013 年）。我已经构建了一个函数来为我执行此操作并返回每个函数的 p 值。

Fit Univariate logistic regression model for each covariate
uni.log2 <- function(x) {
  log.mod2 <- glm(Renewf ~ x, data = dt.train2, family = binomial())
  return(coef(summary(log.mod2))[,4]) #get p-values only
}

然后我将此函数应用于我的 dt 中的每个选定列

#apply function to selected IV's 
apply(X = dt.train2[c(3:16)], MARGIN = 2, FUN = uni.log2)

我想做的下一步是筛选这些变量的阈值为 p

有谁知道如何做到这一点？

我可以使用以下代码设置阈值并从多变量模型中复制名称列表：

threshold <- 0.001
signif_form <- as.formula(paste("Renewf ~ 
",paste(names(which((summary(log.mod2)$coefficients[2: 
(nrow(summary(log.mod2)$coefficients)), 4] < threshold) == TRUE)), collapse 
= "+")))

但是，我不知道如何从一系列单变量回归模型中粘贴名称。如果有人知道如何做到这一点，我将不胜感激。

提前谢谢你！

Answer 1

如果您在查看@BenBolker 提供的链接后仍想使用此方法（可能还有其他关于stepwise regression 和statistical significance 危险的资源）...

以下代码将为每个回归中的自变量返回一个 p 值向量。我使用了内置的mtcars 数据框进行说明。

library(tidyverse)
library(broom)

pvals = sapply(names(mtcars)[names(mtcars) != "vs"], function(x) {
  glm(paste("vs ~ ", x), data=mtcars, family=binomial) %>% 
    tidy %>% 
    filter(term==x) %>% pull(p.value)
})

pvals

        mpg         cyl        disp          hp        drat          wt        qsec          am 
0.006590045 0.001917098 0.002453817 0.012340143 0.021777872 0.008672977 0.008813419 0.343628917 
       gear        carb 
0.250981095 0.004157666

上面的代码使用管道运算符 (%>%) 将函数链接在一起。使用glm 创建模型后，tidy 将系数和 p 值作为数据框返回：

glm(vs ~ mpg, data=mtcars, family=binomial) %>% 
    tidy

         term   estimate std.error statistic     p.value
1 (Intercept) -8.8330726  3.162274 -2.793266 0.005217877
2         mpg  0.4304135  0.158422  2.716880 0.006590045

然后filter 和pull 函数为正在考虑的特定变量选择 p 值：

glm(vs ~ mpg, data=mtcars, family=binomial) %>% 
  tidy %>% filter(term=="mpg") %>% pull(p.value)

[1] 0.006590045

将整个内容包装在 sapply 中会返回一个 p 值的命名向量，其中名称是每个单变量回归中的自变量。

仅返回低于 p 值阈值的元素：

pvals[pvals < 0.25]

        mpg         cyl        disp          hp        drat          wt        qsec        carb 
0.006590045 0.001917098 0.002453817 0.012340143 0.021777872 0.008672977 0.008813419 0.004157666

如果你只想要满足阈值标准的变量名称：

names(pvals[pvals < 0.25])

直接返回低于 p 值阈值的元素：

pvals = sapply(names(mtcars)[names(mtcars) != "vs"], function(x) {
  glm(paste("vs ~ ", x), data=mtcars, family=binomial) %>% 
    tidy %>% 
    filter(term==x) %>% pull(p.value)
}) %>% .[. < 0.25]

最后，将其打包成一个函数，返回所需变量的名称：

select_vars = function(DV, data, threshold) {
  sapply(names(data)[names(data) != DV], function(x) {
    glm(paste(DV, " ~ ", x), data=data, family=binomial) %>% 
      tidy %>% 
      filter(term==x) %>% pull(p.value)
  }) %>% .[. < threshold] %>% names
}

select_vars("vs", mtcars, 0.25)

[1] "mpg"  "cyl"  "disp" "hp"   "drat" "wt"   "qsec" "carb"

select_vars("Species", iris %>% filter(Species %in% c("versicolor","virginica")), 0.001)

[1] "Sepal.Length" "Petal.Length" "Petal.Width"