[R] How to compare two regression line slopes
Ben Bolker
bolker at ufl.edu
Tue Jan 27 22:52:43 CET 2009
Etienne Toffin <etoffin <at> ulb.ac.be> writes:
>
> I've made a research about how to compare two regression line slopes
> (of y versus x for 2 groups, "group" being a factor ) using R.
>
> I knew the method based on the following statement :
> t = (b1 - b2) / sb1,b2
> where b1 and b2 are the two slope coefficients and sb1,b2 the pooled
> standard error of the slope (b)
>
> However, I also found a procedure in Wonnacott & Wonnacott, that is
> based on the use of a mute variable D that will have a binary value
> according to the group to which a given point belongs (group : D=0;
> group 2: D=1). Then the equation that is computed is as follow:
> y = b0 + b1.x + D.b2.x
>
> which can be computed in R with:
> > fit <- lm(y ~ group + x + x:group)
> where y is the response of the 2 groups.
> The p-value of x:group gives the probability for the two slopes to be
> different, and the estimated values of parameters are these of both
> populations.
>
> These two methods have already been described in the mailing list but
> not confronted and discussed.
> So, my questions are:
> - are these methods different ?
> - which one should be preferentially used ?
>
I think you're perfectly clear.
These procedures are identical: the first has the virtue
of being very mechanical and transparent, but the
second is much easier (and easier to extend, e.g. to
multiple groups):
dat <- data.frame(x=rep(1:3,2),y=rep(1:3,2)+rnorm(6),
f=factor(rep(1:2,each=3)))
test1 <- function(dat) {
fits <- lapply(split(dat,dat$f),lm,formula=y~x)
sums <- lapply(fits,summary)
coefs <- lapply(sums,coef)
db <- coefs[[2]]["x","Estimate"]-coefs[[1]]["x","Estimate"]
sd <- sqrt(sum(sapply(coefs,function(x) x["x","Std. Error"])^2))
df <- sum(sapply(fits,"[[","df.residual"))
td <- db/sd
c(est=db,sd=sd,tstat=td,prt=2*pt(-abs(td),df))
}
test2 <- function(dat) {
fit <- lm(y~x*f,data=dat)
coef(summary(fit))["x:f2",]
}
rbind(test1(dat),test2(dat))
est sd tstat prt
[1,] 0.5333555 1.382019 0.3859249 0.7367364
[2,] 0.5333555 1.382019 0.3859249 0.7367364
>
More information about the R-help
mailing list