如何在 R 中使用 glm 使循环在二进制逻辑回归中包含混杂因素(协变量)?

Question

有50多个暴露变量和总共16,000个数据。
我需要分析所有暴露变量与结果之间的关联。

我想重复以下公式。

example.data <- data.frame(outcome = c(0,0,0,1,0,1,0,0,1,0),
exposure_1 = c(2.03, 2.13, 0.15, -0.14, 0.32,2.03, 2.13, 0.15, -0.14, 0.32),
exposure_2 = c(-0.11, 0.93, -1.26, -0.95, 0.24,-0.11, 0.93, -1.26, -0.95, 0.24),
age = c(20, 25, 30, 35, 40, 50, 55, 60, 65, 70),
bmi = c(20, 23, 21, 20, 25, 18, 20, 25, 26, 27))
logit_1 <- glm(outcome ~exposure_1, family = binomial, data = example.data)
logit_2 <- glm(outcome~ exposure_2 + age+ bmi, family = binomial, data = example.data)

我做了一个这样的公式。

Model1 <- function(x) {
temp <- glm(reformulate(x,response="outcome"), data=example.data, family=binomial)
c(exp(summary(temp)$coefficients[2,1]), # OR
exp(confint(temp)[2,1]), # CI 2.5 
exp(confint(temp)[2,2]), # CI 97.5
summary(temp)$coefficients[2,4], # P-value
colAUC(predict(temp, type = "response"),example.data$outcome)) #AUC
Model1.data <- as.data.frame(t(sapply(setdiff(names(example.data),"outcome"), Model1)))
}
Model2 <- function(x) {
temp <- glm(reformulate(x + age + bmi, response="outcome"), data=example.data, family=binomial)
c(exp(summary(temp)$coefficients[2,1]), # OR
exp(confint(temp)[2,1]), # CI 2.5 
exp(confint(temp)[2,2]), # CI 97.5
summary(temp)$coefficients[2,4], # P-value
colAUC(predict(temp, type = "response"),example.data$outcome)) #AUC
}
Model2.data <- as.data.frame(t(sapply(setdiff(names(example.data),"outcome"), Model2)))

但是，函数"模型 2"不起作用。
我编写的代码只运行单个二进制逻辑回归，但不能通过添加多个混杂因素来分析。

Answer 1

> 在reformulate中使用c()不+。在这两个函数中，x都采用列名的值。我将使用mtcars列名称来说明：

## Model1 works
reformulate("hp", response = "mpg")
# mpg ~ hp
## Model2 doesn't work
reformulate("hp" + wt + cyl, response = "mpg")
# Error in reformulate("hp" + wt + cyl, response = "mpg") : 
#   object 'wt' not found
## Fix it with `c()` and quoted column names
reformulate(c("hp", "wt", "cyl"), response = "mpg")
# mpg ~ hp + wt + cyl
## Showing it works with a mix of variables and quoted column names
x = "hp"
reformulate(c(x, "wt", "cyl"), response = "mpg")
# mpg ~ hp + wt + cyl

所以在你的Model2，改为reformulate(c(x, "age", "bmi"), response="outcome")

您还有其他问题 - 您在除outcome(setdiff(names(example.data),"outcome")) 之外的所有列上运行Model2，但您还应该排除bmi和age，因为它们包含在函数中。

Answer 2

所以，如果我理解正确的话，你想要一个模型outcome ~ exposure,50 个暴露变量中的每一个都有一个模型outcome ~ exposure + age + bmi......

以下是使用函数式编程的一种方法，假设所有暴露变量都命名为exposure_1、exposure_2等：

library(tidyverse)
# formulas for exposure only
lexpo1 <-
str_subset(string = names(example.data), pattern = "^exposure") %>%
as.list() %>%
modify(., ~ reformulate(termlabels = ., response = "outcome"))
# formulas for exposure + age + bmi
lexpo2 <-
lexpo1 %>%
modify(~ modelr::add_predictors(.x, ~age, ~bmi))
# put together
list_mods <-
c(lexpo1, lexpo2)
# fit the models
list_res <-
list_mods %>%
map(., ~ glm(., family = binomial, data = example.data))