[R] manova question
John Fox
jfox at mcmaster.ca
Wed Mar 23 02:21:32 CET 2011
Dear Ranjan,
Looking at your data, studentgroup is a numeric variable, so your call to lm() produces a multivariate regression, not a one-way MANOVA. Here's a correction:
> morel$studentgroup <- as.factor(morel$studentgroup)
> Manova( lm( cbind(aptitude, mathematics, language, generalknowledge)
+ ~ studentgroup , data = morel), test="Wilks")
Type II MANOVA Tests: Wilks test statistic
Df test stat approx F num Df den Df Pr(>F)
studentgroup 2 0.54345 6.7736 8 152 1.384e-07 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Best,
John
--------------------------------
John Fox
Senator William McMaster
Professor of Social Statistics
Department of Sociology
McMaster University
Hamilton, Ontario, Canada
http://socserv.mcmaster.ca/jfox
> -----Original Message-----
> From: r-help-bounces at r-project.org [mailto:r-help-bounces at r-project.org]
> On Behalf Of Ranjan Maitra
> Sent: March-22-11 8:36 PM
> To: 'R-help'
> Subject: Re: [R] manova question
>
> Dear John, Peter and others,
>
> So, I now have a query at an even more elementary level and that is
> regarding my results from anova.mlm() not matching the car package's
> Manova(). Specifically, I have been trying the following out with regard
> to a simple one-way MANOVA setup. So, I try out the following using R:
>
> ******* R code *******
>
> morel <- read.table(file =
> "http://www.public.iastate.edu/~maitra/stat501/datasets/morel.dat",
> col.names = c("studentgroup", "aptitude", "mathematics", "language",
> "generalknowledge"))
>
> morel[,1] <- as.factor(morel[,1])
> fit <- anova.mlm(as.matrix(morel[,-1]) ~ morel[,1])
>
> summary(fit, test="Wilks")
>
>
>
> *** providing the output ***
>
>
> Df Wilks approx F num Df den Df Pr(>F)
> morel[, 1] 2 0.54345 6.7736 8 152 1.384e-07 ***
> Residuals 79
>
> *** end of output
>
>
>
> The above is correct, also by doing the calculations "by hand".
>
> Then, I use the car package, following the help function on Anova() and
> do the following:
>
>
> ******* R code ********
>
> morel <- read.table(file =
> "http://www.public.iastate.edu/~maitra/stat501/datasets/morel.dat",
> col.names=c("studentgroup", "aptitude", "mathematics", "language",
> "generalknowledge"))
>
> library(car)
> fit1 <- Manova( lm( cbind(aptitude, mathematics, language,
> generalknowledge) ~ studentgroup , data = morel)) summary(fit1, test =
> "Wilks")
>
>
> ****** providing the output *****
>
>
> Type II MANOVA Tests:
>
> Sum of squares and products for error:
> aptitude mathematics language generalknowledge
> aptitude 78506.68 13976.5677 11041.9434 4330.1304
> mathematics 13976.57 16040.3996 3979.9528 -416.4845
> language 11041.94 3979.9528 6035.6132 -372.8491
> generalknowledge 4330.13 -416.4845 -372.8491 7097.9562
>
> ------------------------------------------
>
> Term: studentgroup
>
> Sum of squares and products for the hypothesis:
> aptitude mathematics language generalknowledge
> aptitude 1129.7271 996.0542 237.54441 -880.4353
> mathematics 996.0542 878.1980 209.43741 -776.2594
> language 237.5444 209.4374 49.94777 -185.1266
> generalknowledge -880.4353 -776.2594 -185.12655 686.1536
>
> Multivariate Test: studentgroup
> Df test stat approx F num Df den Df Pr(>F)
> studentgroup 1 0.8620544 3.080378 4 77 0.020805 *
> ---
> Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
>
>
> ***** end of output.
>
>
>
> Which is very different from the previous results. So what am I doing
> wrong here? Same issues arise with the other tests also (Pillai, Roy,
> Hotelling-Lawley, etc).
>
> Many thanks and best wishes,
> Ranjan
>
>
>
>
>
>
>
>
>
>
>
>
> On Sun, 20 Mar 2011 19:29:41 -0500 John Fox <jfox at mcmaster.ca> wrote:
>
> > Dear Peter and Ranjan,
> >
> > In addition to Anova(), linearHypothesis() in the car package handles
> > multivariate linear models, including those for repeated measures.
> >
> > Best,
> > John
> >
> > --------------------------------
> > John Fox
> > Senator William McMaster
> > Professor of Social Statistics
> > Department of Sociology
> > McMaster University
> > Hamilton, Ontario, Canada
> > http://socserv.mcmaster.ca/jfox
> >
> >
> >
> >
> > > -----Original Message-----
> > > From: r-help-bounces at r-project.org
> > > [mailto:r-help-bounces at r-project.org]
> > > On Behalf Of peter dalgaard
> > > Sent: March-20-11 6:50 PM
> > > To: Ranjan Maitra
> > > Cc: R-help
> > > Subject: Re: [R] manova question
> > >
> > >
> > > On Mar 20, 2011, at 21:05 , Ranjan Maitra wrote:
> > >
> > > > Dear friends,
> > > >
> > > > Sorry for this somewhat generically titled posting but I had a
> > > > question with using contrasts in a manova context. So here is my
> > > question:
> > > >
> > > > Suppose I am interested in doing inference on \beta in the case of
> > > > the model given by:
> > > >
> > > > Y = X %*% \beta + e
> > > >
> > > > where Y is a n x p matrix of observations, X is a n x m design
> > > > matrix, \beta is m x p matrix of parameters, and e is a
> > > > normally-distributed random matrix with mean zero and independent
> > > > rows, each having dispersion matrix given by \Sigma. Then, I know
> > > > (I think) how to perform MANOVA. Specifically, I use:
> > > >
> > > > fit <- manova(Y ~ X)
> > > >
> > > > and
> > > >
> > > > summary(fit) will allow me to perform appropriate inference on
> beta.
> > > >
> > > > Now, suppose I am interested in doing inference on C %*% \beta %*%
> > > > M (say testing whether this is equal to zero) with C and M being q
> > > > x m and p x r matrices, respectively (with q, r both being no more
> > > > than p), then can this be done using the manova object from the
> above? How?
> > > > If not, is there an efficient way to do this?
> > >
> > > Check out anova.mlm(), it does most of this sort of testing. Not
> > > quite the "C %*% ..." bit because the linear model code is not
> > > really built to handle linear constraints, but rather compare nested
> > > models, each specified using a set of betas. (So you usually test
> > > whether a subset of betas is zero).
> > >
> > > Also check out the "car" package. Its Anova() function does some
> > > similar stuff.
> > >
> > > If noone has done so already, I wouldn't think it to be very hard to
> > > implement the general case. Most of the bits are there already.
> > >
> > > --
> > > Peter Dalgaard
> > > Center for Statistics, Copenhagen Business School Solbjerg Plads 3,
> > > 2000 Frederiksberg, Denmark
> > > Phone: (+45)38153501
> > > Email: pd.mes at cbs.dk Priv: PDalgd at gmail.com
> > >
> > > ______________________________________________
> > > R-help at r-project.org 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.
> >
>
> ______________________________________________
> R-help at r-project.org 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.
More information about the R-help
mailing list