r-循环数据帧中的列,将模型beta和SE提取到新的数据帧中



我之前提出了一个类似的问题,这个问题用简单线性回归的嵌套循环得到了充分的回答(谢谢!(

这一次,我正在寻找一种方法,在数据帧中的列上循环,并将输出保存在新的数据帧中,其中新输出数据帧的每一行代表循环的一次迭代(即我们正在循环的原始数据帧的列(。

这是我迄今为止的尝试:

library(tidyverse)
set.seed(42)
n <- 1000
dat <- data.frame(id=1:n,
q=runif(n, min=45, max=85),
r=runif(n, min=2.4, max=6.0),
s=runif(n, min=24, max=60),
t=runif(n, min=0.28, max=1.73))
vd <- runif(n, min=15, max=125)
my_models <- list()
x <- 1
for (i in 1:colnames(dat)){ 
model <- lm(paste(dat[[i]], "~", vd, "-1"))
my_models[[x]]<-data.frame(ModelNo=i,coefficients(summary(model)))
x <- x+1
} 
output <- do.call(rbind,my_models)
reg.tbl <- output %>% 
mutate(beta=(round(exp(Estimate), digits=3)),
lower=(round(exp(Estimate-1.96*Std..Error), digits=3)),
upper=(round(exp(Estimate+1.96*Std..Error), digits=3))) %>% 
select(-Estimate, -Std..Error, -t.value, -Pr...t..) %>% 
mutate(coef=paste(beta, " ", "(", lower, " ", "to", " ", upper, ")")) %>% 
mutate(outcome=if_else(ModelNo==1, "MI", "Stroke")) %>% 
select(-ModelNo, -beta, -lower, -upper)

我没有将以前的解决方案应用于我的当前数据的原因是,我有40多个变量要运行回归,并且将每个变量保存为一个对象,然后按照以前的解决方法保存到一个"预测列表"中是耗时的,并且无法实现循环某些内容以自动分析的目的。这正是我18个月前转到R的原因。

我们自然非常感谢任何援助。

我不明白为什么要将最终结果系数和CI作为字符串,但下面的代码会产生这样的结果。

  • 我已经更改了for循环,问题不起作用
  • 回归公式是用reformulate而不是用paste创建的
  • 在管道中,数字不是四舍五入的,sprintf负责处理这一点
suppressPackageStartupMessages(library(tidyverse))
set.seed(42)
n <- 1000
dat <- data.frame(id=1:n,
q=runif(n, min=45, max=85),
r=runif(n, min=2.4, max=6.0),
s=runif(n, min=24, max=60),
t=runif(n, min=0.28, max=1.73))
vd <- runif(n, min=15, max=125)
my_models <- vector("list", length = ncol(dat) - 1L)
for(i in seq_along(names(dat))[-1]){
resp <- names(dat)[i]
fmla <- reformulate(termlabels = c(-1, "vd"), response = resp)
model <- lm(fmla, data = dat)
smry <- summary(model)
my_models[[i - 1L]] <- data.frame(ModelNo = i - 1L, coefficients(smry), check.names = FALSE)
} 
output <- do.call(rbind, my_models)
row.names(output) <- NULL
output %>% 
mutate(beta = exp(Estimate),
lower = exp(Estimate - 1.96*`Std. Error`),
upper = exp(Estimate + 1.96*`Std. Error`)) %>% 
select(-Estimate, -`Std. Error`, -`t value`, -`Pr(>|t|)`) %>% 
mutate(coef = sprintf("%.03g (%.03g to %.03g)", beta, lower, upper)) %>% 
mutate(outcome = if_else(ModelNo == 1, "MI", "Stroke")) %>% 
select(-ModelNo, -beta, -lower, -upper)
#>                  coef outcome
#> 1 2.13 (2.08 to 2.18)      MI
#> 2 1.05 (1.05 to 1.05)  Stroke
#> 3 1.64 (1.61 to 1.66)  Stroke
#> 4 1.01 (1.01 to 1.01)  Stroke

创建于2022-03-07由reprex包(v2.0.1(


编辑

这是另一种方式。起初,它看起来更复杂,但它简化了主循环,使其成为一个简单的purrr::map调用。

编写一个回归函数reg_fun,处理所有问题,回归和置信区间下限和上限的计算
然后让它返回一个S3类对象子类"data.frame"。我已经将这个自定义类称为"Sandro",请随意更改
自定义类具有自己的print方法。像这样,回归的返回值仍然是数字,可以使用标准提取运算符提取它们,但它们将以所需的格式打印。

然后用一个简单得多的管道调用purrr::map_dfr

suppressPackageStartupMessages(library(tidyverse))
set.seed(42)
n <- 1000
dat <- data.frame(id=1:n,
q=runif(n, min=45, max=85),
r=runif(n, min=2.4, max=6.0),
s=runif(n, min=24, max=60),
t=runif(n, min=0.28, max=1.73))
vd <- runif(n, min=15, max=125)
reg_fun <- function(y, x, data, conf = 0.95){
alpha <- qnorm(1 - (1 - conf)/2)
fit <- lm(y ~ 0 + x, data = data)
smry <- summary(fit)
estimate <- coef(smry)[, 1]
se <- coef(smry)[, 2]
out <- data.frame(
beta = exp(estimate),
lower = exp(estimate - alpha * se),
upper = exp(estimate + alpha * se)
)
class(out) <- c("Sandro", class(out))
out
}
# custom print method
print.Sandro <- function(x, digits = 3){
outcome <- rep("Stroke", nrow(x))
outcome[1] <- "MI"
fmt <- paste0("%.0", digits, "g")
fmt <- paste0(fmt, " (", fmt, " to ", fmt, ")")
out <- data.frame(
coef = sprintf(fmt, x[["beta"]], x[["lower"]], x[["upper"]]),
outcome = outcome
)
print.data.frame(out)
}
dat %>%
select(-id) %>%
map_dfr(reg_fun, x = vd, data = dat)
#>                  coef outcome
#> 1 2.13 (2.08 to 2.18)      MI
#> 2 1.05 (1.05 to 1.05)  Stroke
#> 3 1.64 (1.61 to 1.66)  Stroke
#> 4 1.01 (1.01 to 1.01)  Stroke

创建于2022-03-07由reprex包(v2.0.1(


编辑2

要同时打印响应,请将print方法和管道更改为以下内容。

print.Sandro <- function(x, digits = 3){
fmt <- paste0("%.0", digits, "g")
fmt <- paste0(fmt, " (", fmt, " to ", fmt, ")")
out <- data.frame(
coef = sprintf(fmt, x[["beta"]], x[["lower"]], x[["upper"]])
)
print.data.frame(out)
}
dat %>%
select(-id) %>%
map_dfr(reg_fun, x = vd, data = dat) %>%
mutate(outcome = names(dat)[-1]) %>% 
relocate(outcome, .before = "beta")

要查看带有系数、下限和上限的返回data.frame,请以对as.data.frame的调用结束管道。

dat %>%
select(-id) %>%
map_dfr(reg_fun, x = vd, data = dat) %>%
mutate(outcome = names(dat)[-1]) %>% 
relocate(outcome, .before = "beta") %>%
as.data.frame()
#>   outcome     beta    lower    upper
#> 1       q 2.128118 2.079911 2.177443
#> 2       r 1.050421 1.048759 1.052085
#> 3       s 1.638222 1.612149 1.664716
#> 4       t 1.012015 1.011534 1.012497

相关内容

最新更新