Цикл линейной регрессии и коэффициенты экономии

Это часть набора данных (с именем «ME1»), который я использую (все переменные числовые):

   Year  AgeR   rateM
1  1751 -1.0 0.241104596
2  1751 -0.9 0.036093609
3  1751 -0.8 0.011623734
4  1751 -0.7 0.006670552
5  1751 -0.6 0.006610552
6  1751 -0.5 0.008510828
7  1751 -0.4 0.009344041
8  1751 -0.3 0.011729740
9  1751 -0.2 0.010988005
10 1751 -0.1 0.015896107
11 1751  0.0 0.018190140
12 1751  0.1 0.024588340
13 1751  0.2 0.029801362
14 1751  0.3 0.044515912
15 1751  0.4 0.055240354
16 1751  0.5 0.088476758
17 1751  0.6 0.119045309
18 1751  0.7 0.167866571
19 1751  0.8 0.239244825
20 1751  0.9 0.329683010
21 1751  1.0 0.472448318

Я хочу использовать линейную модель и сохранить коэффициенты следующим образом:

male<-lm(ME1$rateM~exp(AgeR))
summary(male)
coeff <- summary(male)$coefficients[2]

Проблема в том, что мне нужно повторять эту процедуру для каждого года (с 1751 по 2014) и я хочу сохранить все коэффициенты в один набор данных следующим образом:

Year coeff
1751 0.1556977
1752 0.0966664
...
2014 0.0420914

Я не знаю, нужно ли мне использовать цикл for, lapply или что-то еще. Кто-нибудь может мне помочь?


person DanielUp    schedule 11.01.2016    source источник


Ответы (3)


Есть несколько способов сделать это. Во-первых, мы создаем некоторые сгенерированные данные для иллюстрации:

set.seed(123)
dat <- expand.grid(year=2000:2010, AgeR=seq(-1,1,0.1))
dat$value <- rnorm(nrow(dat))

Мы можем начать с base-R. Мы разделяем наши данные по годам, подгоняем модель и извлекаем наш коэффициент. Затем связываем все вместе.

res <- do.call(rbind,lapply(split(dat, dat$year),function(x){
  fit <- lm(value~exp(AgeR), data=x)
  res <- data.frame(year=unique(x$year),coeff=coef(fit)[2])
  res
}))

Мы можем сделать то же самое, используя data.table:

library(data.table)


res2 <- setDT(dat)[,.(coeff=coef(lm(value~exp(AgeR)))[2]),year]
res2
person Heroka    schedule 11.01.2016
comment
Пожалуйста. Это легко расширяется пользовательской функцией для извлечения большего количества или разных данных из каждой модели. - person Heroka; 11.01.2016

Пакет broom здесь тоже может быть очень полезен.

library(dplyr)
library(broom)

mtcars %>%
  group_by(gear) %>%
  do(tidy(lm(mpg ~ am + cyl, data = .)))
person Benjamin    schedule 11.01.2016

Пакет nlme предлагает для этого удобную функцию:

library(nlme)
coef(lmList(value ~ exp(AgeR)| year, data=dat))[,2, drop = FALSE]
person Roland    schedule 11.01.2016