David,
2.3.1 is a bit old to be reporting bugs -- we do ask people to check that their
problem is still present in a contemporary version of R. However, your data do
still give the same output in R 2.7.2 (which is not current, but was current
less than a year ago).
I've tidied up the code to remove all the weird characters:
booth<-read.table(tmp<-textConnection("COUNT REWGRP COMMIT SATIS STAY
1 1 16 19 18
2 1 18 15 17
3 1 18 14 14
4 1 16 20 10
5 1 15 13 17
6 1 12 15 11
7 2 16 20 13
8 2 18 14 16
9 2 13 10 14
10 2 17 13 19
11 2 14 18 15
12 2 19 16 18
13 3 20 18 16
14 3 18 15 19
15 3 13 14 17
16 3 12 16 15
17 3 16 17 18
18 3 14 19 15
"),header=TRUE)
fit<-manova(cbind(COMMIT,SATIS,STAY)~REWGRP,data=booth)
Now, as to the question of whether this is a bug. You don't give the SAS
answers that you are happy with, just the R answers. This makes it a lot more
difficult.
It's possible that there is a bug in manova(), but another possibility, since
you are concerned about degrees of freedom, and based on the last three letters
of the name of your predictor variable, is that you wanted
fit2<-manova(cbind(COMMIT,SATIS,STAY)~factor(REWGRP),data=booth)
summary(fit2, test="Pillai")
Df Pillai approx F num Df den Df Pr(>F)
factor(REWGRP) 2 0.28342 0.77049 6 28 0.5995
Residuals 15
summary(fit2, test="Roy")
Df Roy approx F num Df den Df Pr(>F)
factor(REWGRP) 2 0.31963 1.49159 3 14 0.2599
Residuals 15
summary(fit2, test="Hotelling")
Df Hotelling-Lawley approx F num Df den Df Pr(>F)
factor(REWGRP) 2 0.36261 0.72521 6 24 0.6336
Residuals 15
Googling for the variable names and SAS, MANOVA found some programs in which
REWGRP was specified as a CLASS variable, ie, a factor.
Also
http://my.safaribooksonline.com/9781590474174/ch11lev1sec3
has what might be the output of this code. The test statistics all match the
ones in R, but the p-values are slightly different except for Wilks' lambda.
So, it looks as though at least you need to specify that the variable is a
factor. I will have to leave the question of whether the p-values are correct
to someone with more knowledge of MANOVA. It does seem from the documentation
that agreement with SAS is intended at least for the Pillai trace and Roy's
largest root.
We do appreciate bug reports, but it shouldn't be necessary to do all this work
to find out what you think the correct answer is.
-thomas
On Mon, 16 Mar 2009 dvdbo...@cs.com wrote:
Hi.? There appears to be a bug in R function manova.? My friend and I both ran
it the same way as shown below (his run) with the shown data set. His results
are shown below. we both got the same results.? I was running with R 2.3.1. I'm
not sure what version he used.
Thanks very much,
David Booth
Kent State University
-----Original Message-----
From: dvdbo...@cs.com
To: kb...@ilstu.edu
Sent: Sun, 15 Mar 2009 7:01 pm
Subject: Re: MANOVA Data
Ken,
Did you notice that Wilks, Roy, etc p-values are all the same?? Pillai is
almost the SAS result.? Can't figure it out.? I'll submit a bug report. What's
Velleman going to talk about?? Thanks for looking at the R.
Best,
Dave
-----Original Message-----
From: Ken Berk <kb...@ilstu.edu>
To: dvdbo...@cs.com
Sent: Sun, 15 Mar 2009 3:45 pm
Subject: Re: Fwd: MANOVA Data
At 08:07 PM 3/5/2009, you wrote:
Hi Ken,
I've run the attached data set ( a one way MANOVA ex. from the SAS manual
chapter on MANOVA) in both SAS and R and I don't get the same
results.? Do you have any suggestions about how I can find out
what's going on?
Thanks,
Dave
-----Original Message-----
From: dvdbo...@cs.com
To: dvdbo...@aol.com
Sent: Thu, 5 Mar 2009 5:06 pm
Subject: MANOVA Data
Email message sent from CompuServe - visit us today at
http://www.cs.com
Email message sent from CompuServe - visit us today at
http://www.cs.com
Hello, David
My R results are clearly crap, as shown below.
The degrees of freedom are clearly wrong, as is apparent when looking at
the univariate anovas.
SAS gives the correct answers.
I don't know what to do about R.
Ken
COUNT??? REWGRP??? COMMIT???
SATIS??? STAY
1????
1????????
16???????
19?????? 18
2????
1????????
18???????
15?????? 17
3????
1????????
18???????
14?????? 14
4????
1????????
16???????
20?????? 10
5????
1????????
15???????
13?????? 17
6????
1????????
12???????
15?????? 11
7????
2????????
16???????
20?????? 13
8????
2????????
18???????
14?????? 16
9????
2????????
13???????
10?????? 14
10??? 2????????
17???????
13?????? 19
11??? 2????????
14???????
18?????? 15
12??? 2????????
19???????
16?????? 18
13??? 3????????
20???????
18?????? 16
14??? 3????????
18???????
15?????? 19
15??? 3????????
13???????
14?????? 17
16??? 3????????
12???????
16?????? 15
17??? 3????????
16???????
17?????? 18
18??? 3????????
14???????
19?????? 15
attach(booth)
Y <- cbind(COMMIT, SATIS, STAY)
fit <- manova(Y ~ REWGRP)
summary(fit, test="Pillai")
????????? Df? Pillai
approx F num Df den Df Pr(>F)
REWGRP???? 1 0.22731?
1.37283????? 3???? 14
0.2918
Residuals
16?????????????????????????????????????
summary(fit, test="Wilks")
????????? Df??
Wilks approx F num Df den Df Pr(>F)
REWGRP???? 1 0.77269?
1.37283????? 3???? 14
0.2918
Residuals
16?????????????????????????????????????
summary(fit, test="Hotelling-Lawley")
????????? Df
Hotelling-Lawley approx F num Df den Df Pr(>F)
REWGRP????
1????????? 0.29418?
1.37283????? 3???? 14
0.2918
Residuals
16??????????????????????????????????????????????
summary(fit, test="Roy")
?????????
Df???? Roy approx F num Df den Df Pr(>F)
REWGRP???? 1 0.29418?
1.37283????? 3???? 14
0.2918
Residuals
16?????????????????????????????????????
summary(fit)
????????? Df? Pillai
approx F num Df den Df Pr(>F)
REWGRP???? 1 0.22731?
1.37283????? 3???? 14
0.2918
Residuals
16?????????????????????????????????????
summary.aov(fit)
?Response COMMIT :
???????????
Df? Sum Sq Mean Sq F value Pr(>F)
REWGRP?????? 1??
0.333?? 0.333? 0.0532 0.8204
Residuals?? 16 100.167??
6.260??????????????
?Response SATIS :
???????????
Df? Sum Sq Mean Sq F value Pr(>F)
REWGRP?????? 1??
0.750?? 0.750? 0.0945 0.7625
Residuals?? 16 127.028??
7.939??????????????
?Response STAY :
??????????? Df Sum
Sq Mean Sq F value Pr(>F)
REWGRP?????? 1 14.083? 14.083?
2.3013 0.1488
Residuals?? 16 97.917??
6.120??????????????
Email message sent from CompuServe - visit us today at http://www.cs.com
________________________________________________________________________
Email message sent from CompuServe - visit us today at http://www.cs.com
[[alternative HTML version deleted]]
______________________________________________
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
Thomas Lumley Assoc. Professor, Biostatistics
tlum...@u.washington.edu University of Washington, Seattle
______________________________________________
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
