CODEWITHSUNDEEP
rdplyrlm
user72716
user72716

Reputation: 273

Apply linear model formula to grouped data

I want to group my dataframe by Participant and iteratively apply a simple linear model formula, lm(Outcome ~ A, data = mydata), so that I end up with a new, separate dataframe with one coefficient for each Participant.

Here's a sample of mydata:

structure(list(Participant = c(1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
2, 2, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 
4, 4, 4, 4, 4, 4, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 6, 6, 
6, 6, 6, 6, 7, 7, 7, 7, 7, 7, 7, 7, 8, 8, 8, 9, 9, 9, 9, 9, 9, 
9, 9, 9, 10, 10, 10, 10, 11, 11, 11, 11, 12, 12, 12, 12, 12, 
12, 12, 12, 13, 13, 13, 13, 13, 13, 13, 14), Outcome = c(15, 
-4, 5, 25, 0, 3, 16, 0, 5, 0, 10, 0, 5, 0, 0, 0, 0, 9, 5, 1, 
20, 11, 8, 15, 0, 0, 13, 22, 20, 0, 0, 0, 0, 0, 0, 10, 0, 12, 
0, 0, 0, 0, 0, -12, 0, 0, 0, 0, 0, 0, 5, 9, 5, 0, 0, 10, 20, 
0, 10, 0, 0, 20, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 14, 0, 0, 
11, 12, 19, 0, 0, 10, 0, 10, -10, 0, 0, 0, 6, -13, 0, 0, 0, -4, 
0, 0, 0, 0, 0), A = c(16, 50, 9, 25, 33, 3, 23, 13, 20, 11, 21, 
20, 19, 36, 6, 22, 18, 20, 5, 6, 23, 43, 14, 46, 7, 18, 20, 78, 
35, 5, 8, 5, 18, 9, 17, 71, 18, 26, 8, 56, 45, 29, 21, 10, 14, 
15, 21, 11, 38, 26, 15, 9, 22, 20, 21, 51, 20, 29, 14, 48, 10, 
21, 9, 11, 29, 6, 21, 25, 20, 27, 29, 36, 31, 7, 27, 38, 30, 
32, 3, 43, 19, 28, 31, 33, 10, 9, 36, 45, 46, 27, 7, 21, 25, 
15, 20, 35, 23, 22, 16, 24), B = c(11, 42, 17, 26, -1, -8, 18, 
7, -25, 6, 11, 10, 14, 41, 11, 18, 23, 16, 10, 4, 47, 26, 14, 
16, 12, 23, 0, 66, 20, -3, 5, 0, 53, 17, 10, 66, 20, 14, 8, 11, 
25, 14, -6, 22, 2, -2, -29, 3, 31, 26, 10, 9, 17, -20, -19, 31, 
0, -1, -6, -2, -10, 31, -11, -29, -21, -19, 21, 25, 18, 6, 13, 
24, -31, 2, 2, 8, 3, 10, -19, 33, 5, 4, 16, 18, 10, 19, -14, 
-25, 21, 16, 20, 13, 4, 5, -8, -15, 16, 12, -1, 14)), row.names = c(2041L, 
2281L, 2521L, 2641L, 3901L, 4141L, 4201L, 4681L, 4801L, 4921L, 
161L, 241L, 321L, 361L, 401L, 481L, 1241L, 2L, 42L, 82L, 122L, 
162L, 202L, 362L, 482L, 1242L, 1562L, 1682L, 1802L, 1842L, 1922L, 
43L, 123L, 323L, 483L, 1683L, 1963L, 2042L, 2102L, 2282L, 2402L, 
2522L, 2642L, 2762L, 3482L, 3962L, 4382L, 4922L, 4982L, 5042L, 
44L, 204L, 484L, 1444L, 1564L, 1684L, 45L, 325L, 965L, 1165L, 
1445L, 1685L, 1765L, 1925L, 86L, 366L, 406L, 2043L, 2103L, 2343L, 
2523L, 2583L, 2643L, 4083L, 4323L, 4983L, 407L, 1247L, 1407L, 
1807L, 48L, 208L, 408L, 1248L, 2104L, 2164L, 2284L, 2404L, 2584L, 
2644L, 2764L, 4384L, 2045L, 2105L, 2345L, 2405L, 2645L, 2765L, 
4385L, 2046L), class = "data.frame")

And here's what my desired output would look like (with hypothetical coefficients):

   Participant   Coef
1            1   0.09
2            2   0.07
3            3   0.11
...

In the past, I've used the group_by function to group by Participant and calculate a descriptive stat (e.g., mean, median)for each. For instance, I can use the code below to create a dataframe, myMeans, with the mean Outcome for each participant: 

myMeans<- as.data.frame(mydata %>%
                           group_by(Participant) %>%
                           select(Outcome) %>%
                           summarise_each(list(mean)))
head(myCoefficients)

  Participant    Outcome
1           1  7.0454545
2           2  9.8510638
3           3 10.0652174
4           4  5.2156863
5           5  0.5319149
6           6  6.1041667

I was hoping something like this would work to create a dataframe, myCoefficients:

myCoefficients<- as.data.frame(mydata %>%
                            group_by(Participant) %>%
                              coef(lm(Outcome ~ A)))

...but it evidently does not.

Any suggestions?

Upvotes: 0

Views: 256

Answers (3)

MSR
MSR

Reputation: 2891

For a tidyverse solution, there is a similar use case in ?do. Reframing this for the current example:

library(tidyverse)

data %>% 
  group_by(Participant) %>% 
  do(mod = lm(Participant ~ A, data = .)) %>% 
  summarise(Participant = Participant, 
            coef = list(mod$coefficients)) %>% 
  unnest_wider(coef)

Note that this requires the relatively recent tidyr 1.0.0 for unnest_wider().

Upvotes: 2

Dave2e
Dave2e

Reputation: 24139

Here is a solution using sapply. 

#find the slope and intercept
intercept<-sapply(unique(mydata$Participant), function(x){
  lm(Outcome ~ A, data=mydata[mydata$Participant==x,])$coefficients[1]})
A_coefficient<-sapply(unique(mydata$Participant), function(x){
  lm(Outcome ~ A, data=mydata[mydata$Participant==x,])$coefficients[2]})

#combine results into a dataframe
answer<-data.frame(Participant=unique(mydata$Participant), intercept, A_coefficient)

 #slightly more compact coding:
fit<-sapply(unique(mydata$Participant), function(x){
  lm(Outcome ~ A, data=mydata[mydata$Participant==x,])$coefficients})

answer<-cbind(Participant=unique(mydata$Participant), as.data.frame(t(fit)))

Another reasonable option is mentioned in the comments is to use split and lapply

Upvotes: 2

G. Grothendieck
G. Grothendieck

Reputation: 270195

Try lmList. Note that the nlme package already comes with R.

library(nlme)

coef(lmList(Outcome ~ A | Participant, mydata))

giving:

   (Intercept)            A
1     8.122188 -0.079910741
2     2.111455  0.001547988
3     1.722062  0.304546146
4    -2.127148  0.164948454
5    -1.883623  0.076522166
6     2.463768  0.103024575
7     7.133361 -0.043622767
8     0.000000  0.000000000
9     1.370920  0.006923838
10    8.286374  0.081986143
11   -5.359477  0.283224401
12   -4.486884  0.143756558
13   -1.333333  0.034188034
14    0.000000           NA

Upvotes: 6

Related Questions