4
$\begingroup$

I have a dataframe with a group variable GRP (ranging from 1-100) and an X and Y for each one. I'd like to get a list of the regression intercepts and slopes for lm(Y~X) within each group. The intercepts and slopes don't need to be in the same dataframe.

Any suggestions? R beginner here, so simplicity would be great!

$\endgroup$
8
$\begingroup$

data.table also has great tools for solving problems such as this:

library(data.table)
set.seed(1)
dat <- data.table(x=runif(100), y=runif(100), grp=rep(1:2,50))
dat[,coef(lm(y~x)),by=grp]

The first row in each group is the intercept, and the second row is the coefficient:

     grp         V1
[1,]   1  0.5991761
[2,]   1 -0.1350489
[3,]   2  0.4401174
[4,]   2  0.1400153

If you'de rather have a wide data.frame, that just takes a little more specification:

dat[,list(intercept=coef(lm(y~x))[1], coef=coef(lm(y~x))[2]),by=grp]
     grp intercept       coef
[1,]   1 0.5991761 -0.1350489
[2,]   2 0.4401174  0.1400153

I wonder if there's a more succinct way to specify this second method...

$\endgroup$
  • $\begingroup$ Thanks very much! This worked great. Just wondering: in the wide format code, what is the purpose of the first comma (right before "list")? $\endgroup$ – Lyra Feb 16 '12 at 0:51
  • $\begingroup$ @Lyra that's how you tell data.table what function you want it to use to operate on the columns. $\endgroup$ – Zach Feb 16 '12 at 2:11
11
$\begingroup$

The responses by @Henry and @Zach both work, but I think the most straight-forward way to do what you want is to use lmList in the nlme package:

dat <- data.frame(
  GRP = sample(c("A","B","C"), 100, replace=TRUE), 
  X = runif(100), 
  Y = runif(100)
)
require(nlme)
lmList(Y ~ X | GRP, data=dat)
$\endgroup$
  • 2
    $\begingroup$ That is cool - I had no idea that existed. $\endgroup$ – Matt Parker Feb 16 '12 at 4:03
  • $\begingroup$ smiling, I'm working with data with many companies and their financial data. Like Lyra, I'd like to regress select columns against each other within each company and report the coefficients for each company. I tried using your lmList method, but instead of getting unique coefficients for each company, I get the coefficients as if I regressed all the data. Do you know why this is? $\endgroup$ – user2303635 May 29 '13 at 4:47
  • $\begingroup$ @user2303635 It's difficult to say why that's happening without knowing more. I suggest you ask a new question providing the details of your issue (maybe on stackoverflow?). The answers to this question might also help you. $\endgroup$ – smillig May 29 '13 at 13:45
5
$\begingroup$

Adapting from help("by"), this example may meet your needs

mydf <- data.frame( GRP = rep(c("A","B","C"), each=100), X = rep(1:100,3), 
                    Y = rep(c(2,4,8),each=100) + 
                        rep(c(4,2,1),each=100) * rep(1:100,3) + rnorm(300))   
by(mydf, mydf$GRP, function(z) lm(Y ~ X, data = z))
$\endgroup$
0
$\begingroup$

If you use package "tidyr" you could do the following

library(data.table)
library(tidyr)
set.seed(1)
dat <- data.table(x=runif(100), y=runif(100), grp=rep(1:2,50))

ncoefs <- 1
dat <- dat[, coef( lm(y ~ x) ), by = grp]
dat[, est := rep( c("intercept", "coef"), .N/(ncoefs + 1)) ]
dat <- dat %>% spread(est, V1)

The result is

   grp       coef intercept
1:   1 -0.1350489 0.5991761
2:   2  0.1400153 0.4401174

This method is easy to scale up and must be faster than estimating the regression for each coefficient.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.