EDIT：我转而使用mtcars数据集，因为它更好地表示了我的真实数据的结构。忽略下面的虹膜示例：
我希望运行几个具有多个因变量和自变量的回归模型。
假设我正在研究虹膜数据集。
我的因变量是Petal.Length和Petal.Width。我的自变量是Sepal.Length和Sepal.Width。
我想分别为数据集中的每个物种运行以下模型(setosa、versicolor、virginica(：
lm(Petal.Length ~ Sepal.Length)
lm(Petal.Length ~ Sepal.Width)
lm(Petal.Width ~ Sepal.Length)
lm(Petal.Width ~ Sepal.Width)
我的真实数据集有更多的因变量和自变量，这就是为什么我想编写不那么重复的代码。到目前为止，这是我的尝试，我已经用Tidyr嵌套了Species列，并突变了一个包含每个物种回归数据的新列：
library(dplyr)
library(tidyr)
library(purrr)
iris %>%
nest(data = !Species) %>%
mutate(
fit = map(data, ~ lm (cbind(Petal.Length, Petal.Width) ~ ., data = .x))
)
对于每一个物种，这给出了相当的结果：
lm(Petal.Length ~ Sepal.Length + Sepal.Width)
lm(Petal.Width ~ Sepal.Length + Sepal.Width)
让我们使用mtcars数据集：
mpg cyl disp hp drat wt qsec vs am gear carb
Mazda RX4 21.0 6 160.0 110 3.90 2.620 16.46 0 1 4 4
Mazda RX4 Wag 21.0 6 160.0 110 3.90 2.875 17.02 0 1 4 4
Datsun 710 22.8 4 108.0 93 3.85 2.320 18.61 1 1 4 1
Hornet 4 Drive 21.4 6 258.0 110 3.08 3.215 19.44 1 0 3 1
Hornet Sportabout 18.7 8 360.0 175 3.15 3.440 17.02 0 0 3 2
Valiant 18.1 6 225.0 105 2.76 3.460 20.22 1 0 3 1
...
我的因变量是disp和hp。我的自变量是drat、wt和qsec。
我想分别为数据集中的每个vs和am类型运行以下模型：
lm(disp ~ drat)
lm(disp ~ wt)
lm(disp ~ qsec)
lm(hp ~ drat)
lm(hp ~ wt)
lm(hp ~ qsec)
我的真实数据集有更多的因变量和自变量，这就是为什么我想编写不那么重复的代码。到目前为止，这是我的尝试，我已经用Tidyr嵌套了vs和am列，并突变了一个新列，其中包含每个vs和am的回归数据：
library(dplyr)
library(tidyr)
library(purrr)
mtcars %>%
group_by(vs,am) %>%
nest() %>%
mutate(
fit = map(data, ~ lm (cbind(disp, hp) ~ drat + wt + qsec, data = .x))
)
# A tibble: 4 x 4
# Groups: vs, am [4]
vs am data fit
<dbl> <dbl> <list> <list>
1 0 1 <tibble [6 x 9]> <mlm>
2 1 1 <tibble [7 x 9]> <mlm>
3 1 0 <tibble [7 x 9]> <mlm>
4 0 0 <tibble [12 x 9]> <mlm>
它给出了vs和am的每个组合的以下公式的等价物：
lm(disp ~ drat + wt + qsec)
lm(hp ~ drat + wt + qsec)
但这不是我想要的。当我想要六个简单的回归时，上面的代码为每个vs和am组合提供了两个多重回归。我不知道如何在R.中生产这个
这是建模，在重塑数据后使用正确的公式，您可以在没有循环的情况下获得所需的结果。
new_data <- pivot_longer(iris, c(Sepal.Length, Sepal.Width),
names_to = 'Sepal', names_prefix = 'Sepal')
lm(cbind(Petal.Width, Petal.Length) ~ value/(Species:Sepal) - value +
Sepal/Species - Sepal + 0, new_data)
Call:
lm(formula = cbind(Petal.Width, Petal.Length) ~ value/(Species:Sepal) -
value + Sepal/Species - Sepal + 0, data = new_data)
Coefficients:
Petal.Width Petal.Length
Speciessetosa:Sepal.Length -0.17022 0.80305
Speciesversicolor:Sepal.Length 0.08326 0.18512
Speciesvirginica:Sepal.Length 1.22611 0.61047
Speciessetosa:Sepal.Width 0.02418 1.18292
Speciesversicolor:Sepal.Width 0.16691 1.93492
Speciesvirginica:Sepal.Width 0.66406 3.51090
value:Speciessetosa:Sepal.Length 0.08314 0.13163
value:Speciesversicolor:Sepal.Length 0.20936 0.68647
value:Speciesvirginica:Sepal.Length 0.12142 0.75008
value:Speciessetosa:Sepal.Width 0.06471 0.08141
value:Speciesversicolor:Sepal.Width 0.41845 0.83938
value:Speciesvirginica:Sepal.Width 0.45795 0.68632
如果你想运行for循环：
iris %>%
nest(data = !Species) %>%
summarise(model = map(c('Sepal.Length', 'Sepal.Width'),
~map(data,~lm(reformulate(.y, 'cbind(Petal.Width, Petal.Length)'),.x) %>%
coef()%>%
as_tibble(rownames = 'rn'), .y = .x)%>%
set_names(Species)%>%
bind_rows(.id='Species')))%>%
unnest(model)
# A tibble: 12 x 4
Species rn Petal.Width Petal.Length
<chr> <chr> <dbl> <dbl>
1 setosa (Intercept) -0.170 0.803
2 setosa Sepal.Length 0.0831 0.132
3 versicolor (Intercept) 0.0833 0.185
4 versicolor Sepal.Length 0.209 0.686
5 virginica (Intercept) 1.23 0.610
6 virginica Sepal.Length 0.121 0.750
7 setosa (Intercept) 0.0242 1.18
8 setosa Sepal.Width 0.0647 0.0814
9 versicolor (Intercept) 0.167 1.93
10 versicolor Sepal.Width 0.418 0.839
11 virginica (Intercept) 0.664 3.51
12 virginica Sepal.Width 0.458 0.686
编辑：在mtcars数据集的情况下，以下内容就足够了：
pivot_longer(mtcars, c(drat, wt, qsec)) %>%
unite('am_vs', c(am, vs)) %>%
lm(cbind(disp, hp)~value/(am_vs:name)-value + name/am_vs - name + 0, .)
使用map:运行
mtcars%>%
nest_by(vs, am)%>%
rowwise()%>%
summarise(vs, am, model = map(c('drat', 'wt', 'qsec'),
~lm(sprintf('cbind(disp, hp)~%s',.x), data)%>%
coef()%>%
as_tibble(rownames = 'rn')))%>%
unnest(model)