I fitted a several models on a single variable (a) as a mod1. Output mod1 contains regression fitting list of the 5 models. I did several operation in the list of models. Now I want to unlist mod1 into a single regression model output like fit1, fit2,...fit5 etc.
func <-function(z){
fit1 <- lm( y~ x + z )
fit2 <- lm( y~x + I(z^2))
fit3 <- lm( y~poly(x,3) + z)
fit4 <- lm( y~ns(x, 3) + z)
fit5 <- lm( y~ns(x, 9) + z)
return(list(fit1, fit2, fit3, fit4, fit5))
}
mod1 <- func(data$a)
test <- unlist(mod1, recursive = TRUE, use.names = TRUE)
When I used unlist(), output turned very long un-understandable strings. Even when I tried this following, it is still a list of 1 model.
fit1 <- mod1[1]
Does anyone have any idea about how to separate models one by one from a list of regression model output?
Thanks in advance!
I'd recommend a combination of using broom to tidy the output of lm, and standard list methods lapply and [[ indexing to work with lists.
Read more about working with lists in R here.
library(splines)
# create some example data
d <- data.frame(x = rnorm(100, 0, 1),
y = rnorm(100, 0, 1),
z = rnorm(100, 0, 1))
# function to fit 5 models
func <-function(d){
fit1 <- lm( y~ x + z, data = d)
fit2 <- lm( y~x + I(z^2), data = d)
fit3 <- lm( y~poly(x,3) + z, data = d)
fit4 <- lm( y~ns(x, 3) + z, data = d)
fit5 <- lm( y~ns(x, 9) + z, data = d)
# store models in a list
l <- list(fit1, fit2, fit3, fit4, fit5)
# name the models
names(l) <- paste0("fit", 1:5)
return(l)
}
# run the function
mods <- func(d)
Access each element of the list with double brackets [[
mods[[1]]
Call:
lm(formula = y ~ x + z, data = d)
Coefficients:
(Intercept) x z
0.03339 -0.05128 -0.15288
mods[[2]]
Call:
lm(formula = y ~ x + I(z^2), data = d)
Coefficients:
(Intercept) x I(z^2)
0.01739 -0.04490 0.01258
Use broom to "tidy" model output
library(broom)
tidy(mods[[1]])
# A tibble: 3 x 5
term estimate std.error statistic p.value
<chr> <dbl> <dbl> <dbl> <dbl>
1 (Intercept) 0.0334 0.0985 0.339 0.735
2 x -0.0513 0.103 -0.499 0.619
3 z -0.153 0.102 -1.50 0.138
Use lapply (or purrr::map) to tidy the list of model output.
tidy_mods <- lapply(mods, tidy)
# add names to each data frame and combine into one big data frame
for(i in 1:length(tidy_mods)) tidy_mods[[i]]$mod <- names(tidy_mods[i])
do.call(rbind.data.frame, tidy_mods)
# A tibble: 27 x 6
term estimate std.error statistic p.value mod
* <chr> <dbl> <dbl> <dbl> <dbl> <chr>
1 (Intercept) 0.0334 0.0985 0.339 0.735 fit1
2 x -0.0513 0.103 -0.499 0.619 fit1
3 z -0.153 0.102 -1.50 0.138 fit1
4 (Intercept) 0.0174 0.130 0.134 0.894 fit2
5 x -0.0449 0.105 -0.429 0.669 fit2
6 I(z^2) 0.0126 0.0894 0.141 0.888 fit2
7 (Intercept) 0.0309 0.0975 0.317 0.752 fit3
8 poly(x, 3)1 -0.493 0.975 -0.505 0.614 fit3
9 poly(x, 3)2 -0.569 0.975 -0.584 0.561 fit3
10 poly(x, 3)3 1.78 0.976 1.83 0.0709 fit3