[R] Anova in 'car': "SSPE apparently deficient rank"

[John Fox]

Dear Peter and Colleen,

I think that Peter realizes this, but what Anova() does in this case is
equivalent to the MANOVA 

> anova(multmodel, M = ~ Afac*Bfac, X = ~Afac + Bfac, idata=poke.idata)
Error in anova.mlm(multmodel, M = ~Afac * Bfac, X = ~Afac + Bfac, idata =
poke.idata) : 
  residuals have rank 4 < 6

which in turn is equivalent to

> anova(multmodel, idata=poke.idata, X=~Afac+Bfac)
Error in anova.mlm(multmodel, idata = poke.idata, X = ~Afac + Bfac) : 
  residuals have rank 4 < 6

both of which fail for the same reason that Anova() does: Because the
within-subject interaction has 6 df and there are just 5 subjects, the
residual SSP matrix, say SSPE, is of rank 4. The hypothesis of no
interaction has (3 - 1)*(4 - 1) = 6 df, and thus the response-transformation
matrix for this hypothesis, say P, has 6 columns. The error SSP matrix for
the interaction, t(P) %*% SSPE %*% P, is also therefore of rank 4 < 6.

I believe that under these circumstances, it's possible to do the univariate
F-tests but not the multivariate repeated-measures ANOVA. Since Anova()
always computes the multivariate tests, however, I don't see a way around
the problem without entirely changing how Anova() gets the univariate tests.

What's unclear to me is whether the full data set really has just 5
subjects.

Regards,
 John

--------------------------------
John Fox
Senator William McMaster 
  Professor of Social Statistics
Department of Sociology
McMaster University
Hamilton, Ontario, Canada
web: 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: January-03-10 5:32 AM
> To: Colleen F. Moore
> Cc: r-help at r-project.org
> Subject: Re: [R] Anova in 'car': "SSPE apparently deficient rank"
> 
> Colleen F. Moore wrote:
> > I have design with two repeated-measures factor, and no grouping
> > factor. I can analyze the dataset successfully in other software,
> > including my legacy DOS version BMDP, and R's 'aov' function. I would
> > like to use 'Anova' in 'car' in order to obtain the sphericity tests
> > and the H-F corrected p-values.   I do not believe the data are truly
> > deficient in rank. I followed the methods for this kind of analysis
> > outlined in Bennett's excellent handouts for his Psychology 710 course
>
http://www.psychology.mcmaster.ca/bennett/psy710/lectures/maxwell_chp12.pdf
> >     I am trying to convert my own similar course to R for my students
> > for next fall.  I have been successful at analyzing a segment of the
> > data as a 2-way repeated measures design.
> >
> > Here is my code:
> >  > your.data=read.table(pipe("pbpaste"),header=T)
> >  > your.data
> >    partic A1B1 A1B2 A1B3 A1B4 A2B1 A2B2 A2B3 A2B4 A3B1 A3B2 A3B3 A3B4
> > 1     p1    1    1    2    3    1    2    4    7    1    3    7   10
> > 2     p2    2    2    3    3    2    2    5    6    2    4    6    9
> > 3     p3    1    2    2    3    2    3    2    6    1    4    7    9
> > 4     p4    1    1    2    2    1    2    3    6    2    3    8   10
> > 5     p5    2    2    3    3    2    3    5    7    2    3    7    9
> >  > attach(your.data)
> >  > multmodel=lm(cbind(A1B1, A1B2, A1B3, A1B4, A2B1, A2B2, A2B3, A2B4,
> > A3B1, A3B2, A3B3, A3B4)~1)
> >  > poke.idata=read.table(pipe("pbpaste"),header=T)
> >  > poke.idata
> >     Afac Bfac
> > 1    A1   B1
> > 2    A1   B2
> > 3    A1   B3
> > 4    A1   B4
> > 5    A2   B1
> > 6    A2   B2
> > 7    A2   B3
> > 8    A2   B4
> > 9    A3   B1
> > 10   A3   B2
> > 11   A3   B3
> > 12   A3   B4
> >  > attach(poke.idata)
> >  >
> > pokeAnova
> > =Anova(multmodel,idata=poke.idata,idesign=~Afac*Bfac,type="III")
> > Error in linear.hypothesis.mlm(mod, hyp.matrix, SSPE = SSPE, idata =
> > idata,  :
> >    The error SSP matrix is apparently of deficient rank = 4 < 6
> >
> > Thanks for any help or advice. And thanks for the 'car' package, which
> > is a great asset to my course. I'm just stuck on this one example.
> >
> 
> Hmm, this does seem to work with regular anova.mlm:
> 
>  > anova(multmodel, idata=poke.idata, X=~Afac+Bfac,test="Sph")
> Analysis of Variance Table
> 
> 
> Contrasts orthogonal to
> ~Afac + Bfac
> 
> Greenhouse-Geisser epsilon: 0.2880
> Huynh-Feldt epsilon:        0.4871
> 
>              Df     F num Df den Df    Pr(>F)     G-G Pr     H-F Pr
> (Intercept)  1 36.67      6     24 6.164e-11 2.5249e-04 3.3530e-06
> Residuals    4
> 
> 
> As far as I recall, the epsilon corrections do not have a formal
> requirement of a nonsingular SSD of the relevant contrast. Not sure
> about the accuracy of the F probabilities in such cases, though.
> 
> --
>     O__  ---- Peter Dalgaard             Øster Farimagsgade 5, Entr.B
>    c/ /'_ --- Dept. of Biostatistics     PO Box 2099, 1014 Cph. K
>   (*) \(*) -- University of Copenhagen   Denmark      Ph:  (+45) 35327918
> ~~~~~~~~~~ - (p.dalgaard at biostat.ku.dk)              FAX: (+45) 35327907
> 
> ______________________________________________
> 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.