如何在R中嵌套和循环步骤函数，然后输出到单独的Excel表中

Question

我正在尝试通过多个国家/地区（澳元、欧元和新加坡元）的 Excel 工作表循环

step

中的

R

函数。我已将这些工作表合并到一个数据框中，并尝试使用

nest

然后使用

step

对每个国家/地区进行逐步回归。但是，我无法将循环集成到

step

函数中。 循环

step

函数后，我想将

coefficients

数据（包括 t-stats、p-value 等）输出到每个国家/地区标记的单独 Excel 工作表中。 有人可以帮忙吗？我已经包含了我的循环尝试、无循环的

step

代码以及

dput

中的数据。

### Loop with error
dfraw %>% group_by(dfraw, 'Country') %>%
          nest() %>%
          mutate(both_model = map(data, ~lm(VarY ~ ., data = dfraw))) %>%
          step(both_model, direction = "both") %>%
          group_by(Country) %>%
          summary(both_model)

### Functional Step code
both_model <- lm(VarY ~ ., data = dfSGD)
both_model <- step(both_model, direction = "both")
summary(both_model)
Out <- createWorkbook()
addWorksheet(Out, "SGD")
writeData(Out, sheet = "SGD", x = both_model$coefficients)
saveWorkbook(Out, "output.xlsx")

### dput(dfraw)
structure(list(Country = c("AUD", "AUD", "AUD", "AUD", "AUD", 
"AUD", "AUD", "AUD", "AUD", "AUD", "EUR", "EUR", "EUR", "EUR", 
"EUR", "EUR", "EUR", "EUR", "EUR", "EUR", "SGD", "SGD", "SGD", 
"SGD", "SGD", "SGD", "SGD", "SGD", "SGD", "SGD"), VarY = c(-0.0244845360824741, 
-0.0624174372523117, 0.0452624163437831, 0.0530749789385003, 
0.00624000000000002, -0.0295754491970107, 0.0327707684745209, 
-0.0182452800253847, 0.025048480930834, 0.0425666088601608, 0.0141000829416644, 
0.00154489276626668, 0.0621540695036749, 0.00905518537502115, 
-0.0209109380291228, 0.0147859922178988, -0.0318677573278799, 
-0.0293082203837354, -0.024208903799075, -0.0178405500836277, 
-0.0188187608569775, -0.0269326121253097, 0.0527939257325898, 
0.0383738835848475, -0.0163586791881248, 0.000606244316459614, 
-0.0250029554320841, -0.0177659080352996, -0.00352907144923342, 
0.0195835545331209), VarX1 = c(-0.000207994636117426, -0.14021214097687, 
0.019987000188707, 0.0905973411305785, 0.0606271485403591, 0.0460107119360711, 
0.0197995343444735, -0.0253560472301735, 0.045979341688096, 0.0432113473174549, 
-0.000207994636117426, -0.14021214097687, 0.019987000188707, 
0.0905973411305785, 0.0606271485403591, 0.0460107119360711, 0.0197995343444735, 
-0.0253560472301735, 0.045979341688096, 0.0432113473174549, -0.000207994636117426, 
-0.14021214097687, 0.019987000188707, 0.0905973411305785, 0.0606271485403591, 
0.0460107119360711, 0.0197995343444735, -0.0253560472301735, 
0.045979341688096, 0.0432113473174549), VarX2 = c(-0.04, 0.390000000000001, 
-1.015, -0.149999999999999, 0.234999999999999, 0.0650000000000004, 
0.0750000000000002, 0.22, -0.0449999999999999, 0.115, -0.115, 
-0.241, -0.252, -0.102, -0.0249999999999999, -0.0840000000000001, 
-0.47202, 0.28388, -0.21801, -0.2159, 0.43, 1.02, -0.919999999999999, 
-0.65, 0.0700000000000003, -0.63, -0.02, 0, -0.44, -0.29), VarX3 = c(0.0199999999999996, 
-0.734999999999999, 0.229, 0.220000000000001, 0.0739999999999998, 
-0.0760000000000005, 0.00800000000000001, 0.151, -0.04, -0.0739999999999998, 
0.146, -0.117999999999999, -0.234, -0.024, 0.266, 0.101999999999999, 
0.31286, 0.25693, -0.0321999999999996, 0.236699999999999, -0.0819999999999999, 
-0.778, -0.116, 1.075, -0.151, 0.073999999999999, 0.133000000000001, 
0.545999999999999, 0.175000000000001, 0.346), VarX4 = c(-0.0602317110633257, 
-0.0633604956143422, 0.0775846309806723, -0.0297417548946277, 
-0.0720518811692337, -0.0339499828708462, -0.00436185680343282, 
-0.0348429263427839, 0.0958382153499828, 0.0404280217883632, 
-0.0602317110633257, -0.0633604956143422, 0.0775846309806723, 
-0.0297417548946277, -0.0720518811692337, -0.0339499828708462, 
-0.00436185680343282, -0.0348429263427839, 0.0958382153499828, 
0.0404280217883632, -0.0602317110633257, -0.0633604956143422, 
0.0775846309806723, -0.0297417548946277, -0.0720518811692337, 
-0.0339499828708462, -0.00436185680343282, -0.0348429263427839, 
0.0958382153499828, 0.0404280217883632)), class = c("tbl_df", 
"tbl", "data.frame"), row.names = c(NA, -30L))

Answer 1

正如您在

lm()

调用中包装

map()

一样，您需要包装

step()

中的

map()

部分也是如此。这是它的样子：

library(tidyverse)

res <- 
  dfraw |> 
  # If you want to run one model for each country you onle need to
  # `group_by()` `Country`.
  group_by(Country) |>
  nest() |>
  mutate(both_model = map(data, ~lm(VarY ~ ., data = dfraw)),
         step_model = map(both_model, step, direction = "both", trace = 0))

summary()

对于以编程方式评估模型结果没有那么有用。 您可以使用

broom::glance()

和

broom::tidy()

代替：

library(broom)

eval_res <- 
  res |> 
  mutate(quality = map(step_model, glance),
         estimates = map(step_model, tidy)) 

eval_res |> 
  select(Country, quality) |> 
  unnest(quality)
#> # A tibble: 3 × 13
#> # Groups:   Country [3]
#>   Country r.squared adj.r.squared  sigma statistic p.value    df logLik   AIC
#>   <chr>       <dbl>         <dbl>  <dbl>     <dbl>   <dbl> <dbl>  <dbl> <dbl>
#> 1 AUD         0.313         0.262 0.0270      6.14 0.00632     2   67.4 -127.
#> 2 EUR         0.313         0.262 0.0270      6.14 0.00632     2   67.4 -127.
#> 3 SGD         0.313         0.262 0.0270      6.14 0.00632     2   67.4 -127.
#> # ℹ 4 more variables: BIC <dbl>, deviance <dbl>, df.residual <int>, nobs <int>

eval_res |> 
  select(Country, estimates) |> 
  unnest(estimates)
#> # A tibble: 9 × 6
#> # Groups:   Country [3]
#>   Country term        estimate std.error statistic p.value
#>   <chr>   <chr>          <dbl>     <dbl>     <dbl>   <dbl>
#> 1 AUD     (Intercept) -0.00349   0.00515    -0.679  0.503 
#> 2 AUD     VarX1        0.131     0.0912      1.43   0.163 
#> 3 AUD     VarX2       -0.0306    0.0138     -2.21   0.0354
#> 4 EUR     (Intercept) -0.00349   0.00515    -0.679  0.503 
#> 5 EUR     VarX1        0.131     0.0912      1.43   0.163 
#> 6 EUR     VarX2       -0.0306    0.0138     -2.21   0.0354
#> 7 SGD     (Intercept) -0.00349   0.00515    -0.679  0.503 
#> 8 SGD     VarX1        0.131     0.0912      1.43   0.163 
#> 9 SGD     VarX2       -0.0306    0.0138     -2.21   0.0354