Re: [S] Split-plot analysis

Bill Venables (wvenable@attunga.stats.adelaide.edu.au)
Fri, 24 Apr 1998 10:12:23 +0930


Jim Robison-Cox writes:

> My first problem: when I use the above code, I get a different std error:
> > oats.aov1<-aov(Y~N*V+Error(B/V),data=oats,qr=T)
> > se.contrast(oats.aov1,list(V=="Victory",V=="Marvellous"),data=oats)
> [,1]
> [1,] 8.09018

This is odd. Here is what happens in 3.4 on SunOS:

> oats.aov1 <- aov(Y ~ N*V + Error(B/V), data = oats, qr = T)
> se.contrast(oats.aov1, list(V=="Victory", V=="Marvellous"), data=oats)
[,1]
[1,] 7.0789

I'm puzzled about the "paste" form, too. If you omit the
unnecessary paste(B,V,N) it seems to work:

> oats.aov2 <- aov(Y ~ B + N*V + Error(paste(B,V)),
+ data = oats, qr = T)
> se.contrast(oats.aov2, list(V=="Victory", V=="Marvellous"), data=oats)
[,1]
[1,] 7.0789

but in my view this situation is unacceptable. It's a bug.

> I don't find much help on se.contrast in the manuals or in
> V&R2. ...

Yes, well, confidentially I didn't understand it fully either,
and in view of model.tables it seemed unnecessary. Moreover my
efforts to come to grips with it led me to believe it was flaky.

It may sound recalcitrant of me, but in any future edition of V&R
(a recurrent nightmare of mine) I intend to let the same rather
unhelpful lack of information on se.contrast remain. So there!

I have made some progress in another direction, though. Thanks
to a suggestion from Trevor Hastie (I think) I now have a
prototype version of dummy.coef that produces coefficients and
their variance matrices as well. This will let you find standard
errors of any linear function of the coefficients you like (even
the stupid ones, such as between main effect levels in the
presence of pertinent interactions) in an entirely natural (if
still not very easy) way. It is not yet public software since it
needs polishing, but I'm thinking about it. Again, describing
how to use it is going to be tricky and awkward.

> So I tried:
> > oats.aov2<-aov(Y~ B+ N*V + Error(B+B:V+B:V:N),data=oats,qr=T)
> > se.contrast(oats.aov2,list(V=="Victory",V=="Marvellous"),data=oats)
> [,1]
> [1,] NA
>
> Must an error term like B appear in the error statement and
> only in the error statement?

Yes. When each model formula is expanded every distinct full
term may only appear in one place, fixed or Error. (Essentially
because each term in the model must be either fixed or random.)

Also it should be more widely understood that aov(... + Error())
is very restricted in what it can do. Unless your error
structure is completely balanced (orthogonal and has equal cell
sizes at all levels) the analysis will not be fully correct. For
example just one missing value renders the analysis technically
incorrect, with *no warning or comment* from aov. (With missing
values you can repair most aspects of the analysis by the usual
"dummy covariate" trick, but that's another long story.)

The present version of lme goes some way to handling this sort of
analysis better, of course, and the beta version of the next
release promises an almost complete solution with automatic
recovery of interblock information (as it used to be called), but
for designed and fully balanced experiments a stratum by stratum
analysis is still very enlightening.

-- 
Bill Venables, Head, Dept of Statistics,    Tel.: +61 8 8303 5418
University of Adelaide,                     Fax.: +61 8 8303 3696
South AUSTRALIA.     5005.   Email: Bill.Venables@adelaide.edu.au

----------------------------------------------------------------------- This message was distributed by s-news@wubios.wustl.edu. To unsubscribe send e-mail to s-news-request@wubios.wustl.edu with the BODY of the message: unsubscribe s-news