Re: [R] Test for equality of coefficients in multivariate multipleregression

From: John Fox <jfox_at_mcmaster.ca>
Date: Wed 19 Jul 2006 - 21:56:13 EST


Dear Berwin,

Simply stacking the problems and treating the resulting observations as independent will give you the correct coefficients, but incorrect coefficient variances and artificially zero covariances.

The approach that I suggested originally -- testing a linear hypothesis using the coefficient estimates and covariances from the multivariate linear model -- seems simple enough. For example, to test that all three coefficients are the same across the two equations,

 b <- as.vector(coef(x.mlm))  

 V <- vcov(x.mlm)  

 L <- c(1, 0, 0,-1, 0, 0,
        0, 1, 0, 0,-1, 0,
        0, 0, 1, 0, 0,-1)

 L <- matrix(L, nrow=3, byrow=TRUE)  

 t(L %*% b) %*% (L %*% V %*% t(L)) %*% (L %*% b)

The test statistic is chi-square with 3 df under the null hypothesis. (Note: not checked carefully.)

(BTW, it's a bit unclear to me how much of this exchange was on r-help, but I'm copying to r-help since at least one of Ulrich's messages referring to alternative approaches appeared there. I hope that's OK.)

Regards,
 John



John Fox
Department of Sociology
McMaster University
Hamilton, Ontario
Canada L8S 4M4
905-525-9140x23604
http://socserv.mcmaster.ca/jfox

> -----Original Message-----
> From: Berwin A Turlach
> [mailto:berwin@bossiaea.maths.uwa.edu.au] On Behalf Of Berwin
> A Turlach
> Sent: Tuesday, July 18, 2006 9:28 PM
> To: Andrew Robinson
> Cc: Ulrich Keller; John Fox
> Subject: Re: [R] Test for equality of coefficients in
> multivariate multipleregression
>
> G'day all,
>
> >>>>> "AR" == Andrew Robinson <A.Robinson@ms.unimelb.edu.au> writes:

>
> AR> I suggest that you try to rewrite the model system into a
> AR> single mixed-effects model, [...] Why a mixed-effect
> model, wouldn't a fixed effect be o.k. too?
>
> Something like:
>
> > DF<-data.frame(x1=rep(c(0,1),each=50),x2=rep(c(0,1),50))
> > tmp<-rnorm(100)
> > DF$y1<-tmp+DF$x1*.5+DF$x2*.3+rnorm(100,0,.5)
> > DF$y2<-tmp+DF$x1*.5+DF$x2*.7+rnorm(100,0,.5)

> > x.mlm<-lm(cbind(y1,y2)~x1+x2,data=DF)
> > coef(x.mlm)
> y1 y2
> (Intercept) -0.08885266 -0.05749196
> x1 0.33749086 0.60395258
> x2 0.72017894 1.11932077
>
>
> > DF2 <- with(DF, data.frame(y=c(y1,y2)))
> > DF2$x11 <- with(DF, c(x1, rep(0,100)))
> > DF2$x21 <- with(DF, c(x2, rep(0,100)))
> > DF2$x12 <- with(DF, c(rep(0,100), x1))
> > DF2$x22 <- with(DF, c(rep(0,100), x2))
> > DF2$x1 <- with(DF, c(x1, x1))
> > DF2$wh <- rep(c(0,1), each=100)
> > fm1 <- lm(y~wh + x11 + x21 + x12 + x22, DF2)
> > fm1
>
> Call:
> lm(formula = y ~ wh + x11 + x21 + x12 + x22, data = DF2)
>
> Coefficients:
> (Intercept) wh x11 x21
> x12 x22
> -0.08885 0.03136 0.33749 0.72018
> 0.60395 1.11932
>
> > fm2 <- lm(y~wh + x1 + x21 + x22, DF2)
> > anova(fm2,fm1)

> Analysis of Variance Table
>
> Model 1: y ~ wh + x1 + x21 + x22
> Model 2: y ~ wh + x11 + x21 + x12 + x22
> Res.Df RSS Df Sum of Sq F Pr(>F)
> 1 195 246.919
> 2 194 246.031 1 0.888 0.6998 0.4039
>
>
> Cheers,
>
> Berwin
>
> ========================== Full address ============================
> Berwin A Turlach Tel.: +61 (8) 6488 3338
> (secr)
> School of Mathematics and Statistics +61 (8) 6488 3383
> (self)
> The University of Western Australia FAX : +61 (8) 6488 1028
> 35 Stirling Highway
> Crawley WA 6009 e-mail: berwin@maths.uwa.edu.au
> Australia http://www.maths.uwa.edu.au/~berwin
>



R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Received on Wed Jul 19 22:03:33 2006

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.1.8, at Thu 20 Jul 2006 - 00:17:54 EST.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help. Please read the posting guide before posting to the list.