tidyverse の list-columns データ構造を使用して、データ フレームの行によって異なるさまざまなモデル式を適合させる最良の方法は何ですか?
In R for Data Science で、Hadley は list-columns データ構造を使用して多くのモデルを簡単に適合させる方法の素晴らしい例を示しています ( http://r4ds.had.co.nz/many-models.html#gapminder )。式がわずかに異なる多くのモデルを適合させる方法を見つけようとしています。彼の元の例から改作された以下の例では、各大陸に異なるモデルを適合させる最良の方法は何ですか?
library(gapminder)
library(dplyr)
library(tidyr)
library(purrr)
library(broom)
by_continent <- gapminder %>%
group_by(continent) %>%
nest()
by_continent <- by_continent %>%
mutate(model = map(data, ~lm(lifeExp ~ year, data = .)))
by_continent %>%
mutate(glance=map(model, glance)) %>%
unnest(glance, .drop=T)
## A tibble: 5 × 12
# continent r.squared adj.r.squared sigma statistic p.value df
# <fctr> <dbl> <dbl> <dbl> <dbl> <dbl> <int>
#1 Asia 0.4356350 0.4342026 8.9244419 304.1298 6.922751e-51 2
#2 Europe 0.4984659 0.4970649 3.8530964 355.8099 1.344184e-55 2
#3 Africa 0.2987543 0.2976269 7.6685811 264.9929 6.780085e-50 2
#4 Americas 0.4626467 0.4608435 6.8618439 256.5699 4.354220e-42 2
#5 Oceania 0.9540678 0.9519800 0.8317499 456.9671 3.299327e-16 2
## ... with 5 more variables: logLik <dbl>, AIC <dbl>, BIC <dbl>,
## deviance <dbl>, df.residual <int>
by_continent を反復処理することでそれを実行できることはわかっています (すべての大陸の各モデルを推定するため、効率的ではありません:
formulae <- list(
Asia=~lm(lifeExp ~ year, data = .),
Europe=~lm(lifeExp ~ year + pop, data = .),
Africa=~lm(lifeExp ~ year + gdpPercap, data = .),
Americas=~lm(lifeExp ~ year - 1, data = .),
Oceania=~lm(lifeExp ~ year + pop + gdpPercap, data = .)
)
for (i in 1:nrow(by_continent)) {
by_continent$model[[i]] <- map(by_continent$data, formulae[[i]])[[i]]
}
by_continent %>%
mutate(glance=map(model, glance)) %>%
unnest(glance, .drop=T)
## A tibble: 5 × 12
# continent r.squared adj.r.squared sigma statistic p.value df
# <fctr> <dbl> <dbl> <dbl> <dbl> <dbl> <int>
#1 Asia 0.4356350 0.4342026 8.9244419 304.1298 6.922751e-51 2
#2 Europe 0.4984677 0.4956580 3.8584819 177.4093 3.186760e-54 3
#3 Africa 0.4160797 0.4141991 7.0033542 221.2506 2.836552e-73 3
#4 Americas 0.9812082 0.9811453 8.9703814 15612.1901 4.227928e-260 1
#5 Oceania 0.9733268 0.9693258 0.6647653 243.2719 6.662577e-16 4
## ... with 5 more variables: logLik <dbl>, AIC <dbl>, BIC <dbl>,
## deviance <dbl>, df.residual <int>
しかし、ベースRのループに戻らずにこれを行うことは可能ですか(そして、必要のないモデルのフィッティングを回避します)?
私が試したのは次のようなものです：
by_continent <- by_continent %>%
left_join(tibble::enframe(formulae, name="continent", value="formula"))
by_continent %>%
mutate(model=map2(data, formula, est_model))
しかし、機能する est_model 関数を思い付くことができないようです。機能しないこの関数 (h/t: https://gist.github.com/multidis/8138757 ) を試しました:
est_model <- function(data, formula, ...) {
mc <- match.call()
m <- match(c("formula","data"), names(mc), 0L)
mf <- mc[c(1L, m)]
mf[[1L]] <- as.name("model.frame")
mf <- eval(mf, parent.frame())
data.st <- data.frame(mf)
return(data.st)
}
(確かに、これは不自然な例です。私の実際のケースは、データの主要な独立変数が欠落している観測がかなりあるということです。そのため、1 つのモデルを完全な観測のすべての変数に適合させ、別のモデルを完全な観測の変数のサブセットのみに適合させたいと考えています。残りの観察。）
アップデート
動作する est_model 関数を思いつきました (ただし、おそらく効率的ではありません)。
est_model <- function(data, formula, ...) {
map(list(data), formula, ...)[[1]]
}
by_continent <- by_continent %>%
mutate(model=map2(data, formula, est_model))
by_continent %>%
mutate(glance=map(model, glance)) %>%
unnest(glance, .drop=T)
## A tibble: 5 × 12
# continent r.squared adj.r.squared sigma statistic p.value df
# <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <int>
#1 Asia 0.4356350 0.4342026 8.9244419 304.1298 6.922751e-51 2
#2 Europe 0.4984677 0.4956580 3.8584819 177.4093 3.186760e-54 3
#3 Africa 0.4160797 0.4141991 7.0033542 221.2506 2.836552e-73 3
#4 Americas 0.9812082 0.9811453 8.9703814 15612.1901 4.227928e-260 1
#5 Oceania 0.9733268 0.9693258 0.6647653 243.2719 6.662577e-16 4
## ... with 5 more variables: logLik <dbl>, AIC <dbl>, BIC <dbl>, deviance <dbl>,
## df.residual <int>