Skip to content

Bug Report Fwd: MANOVA Data (PR#13595)

2 messages · dvdbooth at cs.com, Thomas Lumley

Original
dvdbooth at cs.com
# 
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: dvdbooth at cs.com
To: kberk at 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 <kberk at ilstu.edu>

To: dvdbooth at 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: dvdbooth at cs.com


To: DvdBooth at 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
????????? Df? Pillai
approx F num Df den Df Pr(>F)


REWGRP???? 1 0.22731?
1.37283????? 3???? 14
0.2918


Residuals
16?????????????????????????????????????
????????? Df??
Wilks approx F num Df den Df Pr(>F)


REWGRP???? 1 0.77269?
1.37283????? 3???? 14
0.2918


Residuals
16?????????????????????????????????????
????????? Df
Hotelling-Lawley approx F num Df den Df Pr(>F)


REWGRP????
1????????? 0.29418?
1.37283????? 3???? 14
0.2918


Residuals
16??????????????????????????????????????????????
?????????
Df???? Roy approx F num Df den Df Pr(>F)


REWGRP???? 1 0.29418?
1.37283????? 3???? 14
0.2918


Residuals
16?????????????????????????????????????
????????? Df? Pillai
approx F num Df den Df Pr(>F)


REWGRP???? 1 0.22731?
1.37283????? 3???? 14
0.2918


Residuals
16?????????????????????????????????????
?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
Thomas Lumley
# 
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
Df  Pillai approx F num Df den Df Pr(>F)
factor(REWGRP)  2 0.28342  0.77049      6     28 0.5995
Residuals      15
Df     Roy approx F num Df den Df Pr(>F)
factor(REWGRP)  2 0.31963  1.49159      3     14 0.2599
Residuals      15
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 dvdbooth at cs.com wrote:

            

      
      
      
    

        

      
      
    

        
Thomas Lumley			Assoc. Professor, Biostatistics
tlumley at u.washington.edu	University of Washington, Seattle