如何使用 R 中的 purrr::map 函数从多个 tidymodels 对象中提取模型摘要？
How can I extract model summary from multiple tidymodels objects using purrr::map functions in R?
我想使用 purrr::map_* 函数从涉及线性回归方法的多个模型中提取信息。我首先创建一些随机数据集。该数据集具有三个因变量和一个自变量。
set.seed(123)
library(olsrr)
#;-)
#;-) Attaching package: 'olsrr'
#;-) The following object is masked from 'package:datasets':
#;-)
#;-) rivers
y1 <- rnorm(100) + runif(100)
y2 <- rnorm(100) + runif(100)
y3 <- rnorm(100) + runif(100)
x <- rnorm(100)
df1 <- data.frame(y1, y2, y3, x)
创建数据后，现在我正在创建三个模型。这一次，我使用基本函数 lm() 来创建模型。我正在使用 lm() 函数来证明它可以在不使用 tidymodels 框架的情况下工作。
m1 <- lm(y1 ~ x, data = df1)
m2 <- lm(y2 ~ x, data = df1)
m3 <- lm(y3 ~ x, data = df1)
现在我尝试通过 purrr::map 函数获取这些模型的摘要。我先创建一个模型列表，然后使用映射函数。
list_of_models <- list(m1, m2, m3)
model_summaries <- purrr::map(list_of_models, summary)
model_normality <- purrr::map(list_of_models, olsrr::ols_test_normality)
model_normality
#;-) [[1]]
#;-) -----------------------------------------------
#;-) Test Statistic pvalue
#;-) -----------------------------------------------
#;-) Shapiro-Wilk 0.9872 0.4517
#;-) Kolmogorov-Smirnov 0.0865 0.4423
#;-) Cramer-von Mises 8.9238 0.0000
#;-) Anderson-Darling 0.5807 0.1276
#;-) -----------------------------------------------
#;-)
#;-) [[2]]
#;-) -----------------------------------------------
#;-) Test Statistic pvalue
#;-) -----------------------------------------------
#;-) Shapiro-Wilk 0.9499 8e-04
#;-) Kolmogorov-Smirnov 0.1252 0.0869
#;-) Cramer-von Mises 9.115 0.0000
#;-) Anderson-Darling 1.6102 4e-04
#;-) -----------------------------------------------
#;-)
#;-) [[3]]
#;-) -----------------------------------------------
#;-) Test Statistic pvalue
#;-) -----------------------------------------------
#;-) Shapiro-Wilk 0.99 0.6685
#;-) Kolmogorov-Smirnov 0.0551 0.9213
#;-) Cramer-von Mises 5.8198 0.0000
#;-) Anderson-Darling 0.3632 0.4347
#;-) -----------------------------------------------
因此，我可以看到可以使用 map 函数来获取模型摘要，以及使用基函数 lm() 时的残差正态性测试。但我想使用 tidymodels 框架进行分析。
library(tidymodels)
#;-) Registered S3 method overwritten by 'tune':
#;-) method from
#;-) required_pkgs.model_spec parsnip
#;-) -- Attaching packages --------------------------------- tidymodels 0.1.4.9000 --
#;-) v broom 0.7.11 v recipes 0.1.17
#;-) v dials 0.0.10 v rsample 0.1.1
#;-) v dplyr 1.0.7 v tibble 3.1.6
#;-) v ggplot2 3.3.5 v tidyr 1.1.4
#;-) v infer 1.0.0 v tune 0.1.6
#;-) v modeldata 0.1.1 v workflows 0.2.4
#;-) v parsnip 0.1.7 v workflowsets 0.1.0
#;-) v purrr 0.3.4 v yardstick 0.0.9
#;-) -- Conflicts ----------------------------------------- tidymodels_conflicts() --
#;-) x purrr::discard() masks scales::discard()
#;-) x dplyr::filter() masks stats::filter()
#;-) x dplyr::lag() masks stats::lag()
#;-) x recipes::step() masks stats::step()
#;-) * Dig deeper into tidy modeling with R at https://www.tmwr.org
lm_spec <- linear_reg()
lm_fit1 <- fit(lm_spec, y1 ~ x, data = df1)
lm_fit2 <- fit(lm_spec, y2 ~ x, data = df1)
lm_fit3 <- fit(lm_spec, y3 ~ x, data = df1)
list_tidymodels <- c(lm_fit1, lm_fit2, lm_fit3)
get_lm_coefs <- function(x) {
x %>%
# get the lm model object
extract_fit_engine() %>%
# transform its format
tidy()
}
map(list_tidymodels, get_lm_coefs)
#;-) Error in UseMethod("extract_fit_engine"): no applicable method for 'extract_fit_engine' applied to an object of class "NULL"
以上命令报错。现在终于，我试图在不使用 map 函数的情况下单独查找每个模型的摘要。
lm_fit1 %>%
extract_fit_engine() %>%
summary()
#;-)
#;-) Call:
#;-) stats::lm(formula = y1 ~ x, data = data)
#;-)
#;-) Residuals:
#;-) Min 1Q Median 3Q Max
#;-) -2.90637 -0.57659 0.02682 0.46190 2.45272
#;-)
#;-) Coefficients:
#;-) Estimate Std. Error t value Pr(>|t|)
#;-) (Intercept) 0.59226 0.09883 5.993 3.44e-08 ***
#;-) x -0.10387 0.10697 -0.971 0.334
#;-) ---
#;-) Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#;-)
#;-) Residual standard error: 0.9748 on 98 degrees of freedom
#;-) Multiple R-squared: 0.00953, Adjusted R-squared: -0.0005771
#;-) F-statistic: 0.9429 on 1 and 98 DF, p-value: 0.3339
lm_fit1 %>%
extract_fit_engine() %>%
olsrr::ols_test_normality()
#;-) -----------------------------------------------
#;-) Test Statistic pvalue
#;-) -----------------------------------------------
#;-) Shapiro-Wilk 0.9872 0.4517
#;-) Kolmogorov-Smirnov 0.0865 0.4423
#;-) Cramer-von Mises 8.9238 0.0000
#;-) Anderson-Darling 0.5807 0.1276
#;-) -----------------------------------------------
所以，很明显只有当我将 tidymodels 与 map 函数结合使用时才会出现错误。我在这里做错了什么？
需要使用 list() 而不是 c() 创建 list_tidymodels。
set.seed(123)
library(olsrr)
library(tidymodels)
y1 <- rnorm(100) + runif(100)
y2 <- rnorm(100) + runif(100)
y3 <- rnorm(100) + runif(100)
x <- rnorm(100)
df1 <- data.frame(y1, y2, y3, x)
lm_spec <- linear_reg()
lm_fit1 <- fit(lm_spec, y1 ~ x, data = df1)
lm_fit2 <- fit(lm_spec, y2 ~ x, data = df1)
lm_fit3 <- fit(lm_spec, y3 ~ x, data = df1)
list_tidymodels <- list(lm_fit1, lm_fit2, lm_fit3) # this line used `c(lm_fit1, ...)`.
get_lm_coefs <- function(x) {
x %>%
# get the lm model object
extract_fit_engine() %>%
# transform its format
tidy()
}
map(list_tidymodels, get_lm_coefs)
#> [[1]]
#> # A tibble: 2 x 5
#> term estimate std.error statistic p.value
#> <chr> <dbl> <dbl> <dbl> <dbl>
#> 1 (Intercept) 0.592 0.0988 5.99 0.0000000344
#> 2 x -0.104 0.107 -0.971 0.334
#>
#> [[2]]
#> # A tibble: 2 x 5
#> term estimate std.error statistic p.value
#> <chr> <dbl> <dbl> <dbl> <dbl>
#> 1 (Intercept) 0.523 0.101 5.16 0.00000129
#> 2 x 0.0705 0.110 0.642 0.522
#>
#> [[3]]
#> # A tibble: 2 x 5
#> term estimate std.error statistic p.value
#> <chr> <dbl> <dbl> <dbl> <dbl>
#> 1 (Intercept) 0.443 0.114 3.87 0.000197
#> 2 x 0.0437 0.124 0.353 0.725
由 reprex package (v2.0.1)
于 2022 年 1 月 20 日创建
我想使用 purrr::map_* 函数从涉及线性回归方法的多个模型中提取信息。我首先创建一些随机数据集。该数据集具有三个因变量和一个自变量。
set.seed(123)
library(olsrr)
#;-)
#;-) Attaching package: 'olsrr'
#;-) The following object is masked from 'package:datasets':
#;-)
#;-) rivers
y1 <- rnorm(100) + runif(100)
y2 <- rnorm(100) + runif(100)
y3 <- rnorm(100) + runif(100)
x <- rnorm(100)
df1 <- data.frame(y1, y2, y3, x)
创建数据后，现在我正在创建三个模型。这一次，我使用基本函数 lm() 来创建模型。我正在使用 lm() 函数来证明它可以在不使用 tidymodels 框架的情况下工作。
m1 <- lm(y1 ~ x, data = df1)
m2 <- lm(y2 ~ x, data = df1)
m3 <- lm(y3 ~ x, data = df1)
现在我尝试通过 purrr::map 函数获取这些模型的摘要。我先创建一个模型列表，然后使用映射函数。
list_of_models <- list(m1, m2, m3)
model_summaries <- purrr::map(list_of_models, summary)
model_normality <- purrr::map(list_of_models, olsrr::ols_test_normality)
model_normality
#;-) [[1]]
#;-) -----------------------------------------------
#;-) Test Statistic pvalue
#;-) -----------------------------------------------
#;-) Shapiro-Wilk 0.9872 0.4517
#;-) Kolmogorov-Smirnov 0.0865 0.4423
#;-) Cramer-von Mises 8.9238 0.0000
#;-) Anderson-Darling 0.5807 0.1276
#;-) -----------------------------------------------
#;-)
#;-) [[2]]
#;-) -----------------------------------------------
#;-) Test Statistic pvalue
#;-) -----------------------------------------------
#;-) Shapiro-Wilk 0.9499 8e-04
#;-) Kolmogorov-Smirnov 0.1252 0.0869
#;-) Cramer-von Mises 9.115 0.0000
#;-) Anderson-Darling 1.6102 4e-04
#;-) -----------------------------------------------
#;-)
#;-) [[3]]
#;-) -----------------------------------------------
#;-) Test Statistic pvalue
#;-) -----------------------------------------------
#;-) Shapiro-Wilk 0.99 0.6685
#;-) Kolmogorov-Smirnov 0.0551 0.9213
#;-) Cramer-von Mises 5.8198 0.0000
#;-) Anderson-Darling 0.3632 0.4347
#;-) -----------------------------------------------
因此，我可以看到可以使用 map 函数来获取模型摘要，以及使用基函数 lm() 时的残差正态性测试。但我想使用 tidymodels 框架进行分析。
library(tidymodels)
#;-) Registered S3 method overwritten by 'tune':
#;-) method from
#;-) required_pkgs.model_spec parsnip
#;-) -- Attaching packages --------------------------------- tidymodels 0.1.4.9000 --
#;-) v broom 0.7.11 v recipes 0.1.17
#;-) v dials 0.0.10 v rsample 0.1.1
#;-) v dplyr 1.0.7 v tibble 3.1.6
#;-) v ggplot2 3.3.5 v tidyr 1.1.4
#;-) v infer 1.0.0 v tune 0.1.6
#;-) v modeldata 0.1.1 v workflows 0.2.4
#;-) v parsnip 0.1.7 v workflowsets 0.1.0
#;-) v purrr 0.3.4 v yardstick 0.0.9
#;-) -- Conflicts ----------------------------------------- tidymodels_conflicts() --
#;-) x purrr::discard() masks scales::discard()
#;-) x dplyr::filter() masks stats::filter()
#;-) x dplyr::lag() masks stats::lag()
#;-) x recipes::step() masks stats::step()
#;-) * Dig deeper into tidy modeling with R at https://www.tmwr.org
lm_spec <- linear_reg()
lm_fit1 <- fit(lm_spec, y1 ~ x, data = df1)
lm_fit2 <- fit(lm_spec, y2 ~ x, data = df1)
lm_fit3 <- fit(lm_spec, y3 ~ x, data = df1)
list_tidymodels <- c(lm_fit1, lm_fit2, lm_fit3)
get_lm_coefs <- function(x) {
x %>%
# get the lm model object
extract_fit_engine() %>%
# transform its format
tidy()
}
map(list_tidymodels, get_lm_coefs)
#;-) Error in UseMethod("extract_fit_engine"): no applicable method for 'extract_fit_engine' applied to an object of class "NULL"
以上命令报错。现在终于，我试图在不使用 map 函数的情况下单独查找每个模型的摘要。
lm_fit1 %>%
extract_fit_engine() %>%
summary()
#;-)
#;-) Call:
#;-) stats::lm(formula = y1 ~ x, data = data)
#;-)
#;-) Residuals:
#;-) Min 1Q Median 3Q Max
#;-) -2.90637 -0.57659 0.02682 0.46190 2.45272
#;-)
#;-) Coefficients:
#;-) Estimate Std. Error t value Pr(>|t|)
#;-) (Intercept) 0.59226 0.09883 5.993 3.44e-08 ***
#;-) x -0.10387 0.10697 -0.971 0.334
#;-) ---
#;-) Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#;-)
#;-) Residual standard error: 0.9748 on 98 degrees of freedom
#;-) Multiple R-squared: 0.00953, Adjusted R-squared: -0.0005771
#;-) F-statistic: 0.9429 on 1 and 98 DF, p-value: 0.3339
lm_fit1 %>%
extract_fit_engine() %>%
olsrr::ols_test_normality()
#;-) -----------------------------------------------
#;-) Test Statistic pvalue
#;-) -----------------------------------------------
#;-) Shapiro-Wilk 0.9872 0.4517
#;-) Kolmogorov-Smirnov 0.0865 0.4423
#;-) Cramer-von Mises 8.9238 0.0000
#;-) Anderson-Darling 0.5807 0.1276
#;-) -----------------------------------------------
所以，很明显只有当我将 tidymodels 与 map 函数结合使用时才会出现错误。我在这里做错了什么？
需要使用 list() 而不是 c() 创建 list_tidymodels。
set.seed(123)
library(olsrr)
library(tidymodels)
y1 <- rnorm(100) + runif(100)
y2 <- rnorm(100) + runif(100)
y3 <- rnorm(100) + runif(100)
x <- rnorm(100)
df1 <- data.frame(y1, y2, y3, x)
lm_spec <- linear_reg()
lm_fit1 <- fit(lm_spec, y1 ~ x, data = df1)
lm_fit2 <- fit(lm_spec, y2 ~ x, data = df1)
lm_fit3 <- fit(lm_spec, y3 ~ x, data = df1)
list_tidymodels <- list(lm_fit1, lm_fit2, lm_fit3) # this line used `c(lm_fit1, ...)`.
get_lm_coefs <- function(x) {
x %>%
# get the lm model object
extract_fit_engine() %>%
# transform its format
tidy()
}
map(list_tidymodels, get_lm_coefs)
#> [[1]]
#> # A tibble: 2 x 5
#> term estimate std.error statistic p.value
#> <chr> <dbl> <dbl> <dbl> <dbl>
#> 1 (Intercept) 0.592 0.0988 5.99 0.0000000344
#> 2 x -0.104 0.107 -0.971 0.334
#>
#> [[2]]
#> # A tibble: 2 x 5
#> term estimate std.error statistic p.value
#> <chr> <dbl> <dbl> <dbl> <dbl>
#> 1 (Intercept) 0.523 0.101 5.16 0.00000129
#> 2 x 0.0705 0.110 0.642 0.522
#>
#> [[3]]
#> # A tibble: 2 x 5
#> term estimate std.error statistic p.value
#> <chr> <dbl> <dbl> <dbl> <dbl>
#> 1 (Intercept) 0.443 0.114 3.87 0.000197
#> 2 x 0.0437 0.124 0.353 0.725
由 reprex package (v2.0.1)
于 2022 年 1 月 20 日创建