[R] Decompose numerical factor into orthog. poly parts

Bill.Venables@cmis.csiro.au Bill.Venables at cmis.csiro.au
Mon Sep 30 09:26:47 CEST 2002


Multiple degrees of freedom terms can be split up in the anova table using

lsmod <- aov(time ~ spinner + order + d)  # using aov( ) makes things a bit
easier.
summary(lsmod, split = list(.....))

ie the splitting is done using summary.aov and the split= argument allows
the split to be controlled in an orderly way.  The details are just too long
to fit into this margin, though...

Bill Venables.

>  -----Original Message-----
> From: 	John Maindonald [mailto:john.maindonald at anu.edu.au] 
> Sent:	Monday, September 30, 2002 3:58 PM
> To:	Murray Jorgensen
> Cc:	r-help at stat.math.ethz.ch
> Subject:	Re: [R] Decompose numerical factor into orthog. poly parts
> 
> The following perhaps ought to work:
> 
> d <- ordered(diam[treat])
> lsmod <- lm(time ~ spinner + order + d)
> anova(lsmod)
> 
> Unfortunately, it gives exactly the same result as before.
> The following, which is inelegant, does split up the orthogonal
> terms:
> 
> d <- ordered(diam[treat])
> p3 <- model.matrix(~ -1 + d)
> lsmod <- lm(time ~ spinner + order + p3[,1]+p3[,2]+p3[,3])
> anova(lsmod)
> 
> Note that orthogonal polynomial contrasts are the default for
> ordered factors.
> 
> John Maindonald
> 
> >Consider the following analysis of a class experiment done as a Latin
> Square:
> >
> >  > spinner <- gl(4,4,16,label=c("Murray","Angela","Shasha","Stephen"))
> >>  order <- gl(4,1,16)
> >>  treat <- scan()
> >1:  1 2 4 3
> >5:  4 3 1 2
> >9:  3 4 2 1
> >13:  2 1 3 4
> >17:
> >Read 16 items
> >>  coin <- factor(treat,label=c("5c","10c","20c","50c"))
> >>  time <- scan()
> >1: 12.0 13.3 13.0 8.9
> >5: 9.7 7.2 7.3 8.6
> >9: 9.2 9.0 7.6 7.4
> >13: 10.1 5.7 7.5 8.9
> >17:
> >Read 16 items
> >>  diam <- c(19,24,28,33)
> >>  treat
> >  [1] 1 2 4 3 4 3 1 2 3 4 2 1 2 1 3 4
> >>  diam[treat]
> >  [1] 19 24 33 28 33 28 19 24 28 33 24 19 24 19 28 33
> >>  d <- diam[treat]
> >>
> >>  lsmod <- lm(time ~ spinner + order + poly(d,3))
> >  > anova(lsmod)
> >Analysis of Variance Table
> >
> >Response: time
> >            Df Sum Sq Mean Sq F value   Pr(>F)  
> >spinner     3 39.367  13.122  9.9476 0.009594 **
> >order       3  7.587   2.529  1.9172 0.228001  
> >poly(d, 3)  3 14.208   4.736  3.5900 0.085602 .
> >Residuals   6  7.915   1.319                   
> >---
> >Signif. codes:  0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1
> >
> >I'm a bit disappointed: I thought that the poly(d,3) sum of squares would
> >be decomposed into its three parts. How would I do that?
> >
> >(Usual apologies for mental laziness.)
> >
> >
> >
> >Dr Murray Jorgensen      http://www.stats.waikato.ac.nz/Staff/maj.html
> >Department of Statistics, University of Waikato, Hamilton, New Zealand
> >Email: maj at waikato.ac.nz                            Fax +64-7 838 4155
> >Phone  +64-7 838 4773 wk    +64 7 849 6486 home     Mobile 021 395 862
> >
> >-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
> .-.-.-
> >r-help mailing list -- Read
> http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
> >Send "info", "help", or "[un]subscribe"
> >(in the "body", not the subject !)  To: r-help-request at stat.math.ethz.ch
> >_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._
> ._._._
> 
> -- 
> John Maindonald             email: john.maindonald at anu.edu.au
> phone : +61 2 (6125)3473    fax  : +61 2(6125)5549
> Centre for Bioinformation Science, Room 1194,
> John Dedman Mathematical Sciences Building (Building 27)
> Australian National University, Canberra ACT 0200.
> -.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.
> -.-.-
> r-help mailing list -- Read
> http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
> Send "info", "help", or "[un]subscribe"
> (in the "body", not the subject !)  To: r-help-request at stat.math.ethz.ch
> _._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._.
> _._._
-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !)  To: r-help-request at stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._



More information about the R-help mailing list