2014-12-12 2 views
2

Я хочу использовать dplyr для группировки data.frame, линейных регрессий и сохранения остатков в качестве столбца в исходном, негруппированном data.frame.сохранить остатки с `dplyr`

Вот пример

> iris %>% 
    select(Sepal.Length, Sepal.Width) %>% 
    group_by(Species) %>% 
    do(mod = lm(Sepal.Length ~ Sepal.Width, data=.)) %>% 

Возвраты:

 Species  mod 
1  setosa <S3:lm> 
2 versicolor <S3:lm> 
3 virginica <S3:lm> 

Вместо этого, я хотел бы оригинальный data.frame с новой колонки, содержащей остатков.

Например,

Sepal.Length Sepal.Width resid 
1 5.1   3.5 0.04428474 
2 4.9   3.0 0.18952960 
3 4.7   3.2 -0.14856834 
4 4.6   3.1 -0.17951937 
5 5.0   3.6 -0.12476423 
6 5.4   3.9 0.06808885 

ответ

5

Я приспособил пример из http://jimhester.github.io/plyrToDplyr/.

r <- iris %>% 
    group_by(Species) %>% 
    do(model = lm(Sepal.Length ~ Sepal.Width, data=.)) %>% 
    do((function(mod) { 
    data.frame(resid = residuals(mod$model)) 
    })(.)) 

corrected <- cbind(iris, r) 

обновление Другой способ заключается в использовании функции augment в broom пакете:

r <- iris %>% 
    group_by(Species) %>% 
    do(augment(lm(Sepal.Length ~ Sepal.Width, data=.)) 

который возвращает:

Source: local data frame [150 x 10] 
Groups: Species 

    Species Sepal.Length Sepal.Width .fitted .se.fit  .resid  .hat 
1 setosa   5.1   3.5 5.055715 0.03435031 0.04428474 0.02073628 
2 setosa   4.9   3.0 4.710470 0.05117134 0.18952960 0.04601750 
3 setosa   4.7   3.2 4.848568 0.03947370 -0.14856834 0.02738325 
4 setosa   4.6   3.1 4.779519 0.04480537 -0.17951937 0.03528008 
5 setosa   5.0   3.6 5.124764 0.03710984 -0.12476423 0.02420180 
... 
+1

(Я понимаю, что происходит, но я никогда бы не подумал об этом сам. Почему, например, мне нужна анонимная функция во втором 'do', но не в первом? –

0

Поскольку вы работать точно такой же регрессии для каждой группе, вам может быть проще просто определить вашу модель регрессии как function(), а затем выполнить его для каждой группы, используя mutate.

model<- function(y,x){ 
    a<- y + x 
    if(length(which(!is.na(a))) <= 2 ){ 
    return(rep(NA, length(a))) 
    } else { 
    m<- lm(y ~ x, na.action = na.exclude) 
    return(residuals(m)) 
    } 
} 

Обратите внимание, что первая часть этой функции застраховаться от каких-либо сообщений об ошибках выскакивают в случае, если ваш регрессия работать на группе с менее чем ноль степеней свободы (это может быть в случае, если у вас есть a dataframe с несколькими переменными группировки со многими levels или многочисленными независимыми переменными для вашей регрессии (например, lm(y~ x1 + x2)), и не может позволить себе проверять каждый из них при достаточных наблюдениях, не связанных с NA).

Так что ваш пример можно переписать следующим образом:

iris %>% group_by(Species) %>% 
    mutate(resid = model(Sepal.Length,Sepal.Width)) %>% 
    select(Sepal.Length,Sepal.Width,resid) 

Который должен уступить:.

Species Sepal.Length Sepal.Width  resid 
    <fctr>  <dbl>  <dbl>  <dbl> 
1 setosa   5.1   3.5 0.04428474 
2 setosa   4.9   3.0 0.18952960 
3 setosa   4.7   3.2 -0.14856834 
4 setosa   4.6   3.1 -0.17951937 
5 setosa   5.0   3.6 -0.12476423 
6 setosa   5.4   3.9 0.06808885 

Этот метод не должен быть вычислительно сильно отличается от того, с помощью augment() (я имел использовать оба метода на наборах данных, содержащих несколько сотен миллионов наблюдений, и полагать, что не было существенной разницы в скорости по сравнению с использованием функции do()).

Кроме того, обратите внимание, что опуская na.action = na.exclude или с использованием m$residuals вместо residuals(m), приведет к исключению строк, которые имеют Nas (упали до оценки) из выходного вектора невязок. Соответственно, соответствующий вектор не будет иметь достаточного количества length(), чтобы быть объединенным с набором данных, и может появиться некоторое сообщение об ошибке.

Смежные вопросы