[R] Orthogonal Polynomials

Douglas Bates bates at stat.wisc.edu
Tue Oct 8 23:41:50 CEST 2002


"Bliese, Paul D MAJ WRAIR-Wash DC" <Paul.Bliese at NA.AMEDD.ARMY.MIL> writes:

> Looking to the wonderful statistical advice that this group can offer.
> 
> In behavioral science applications of stats, we are often introduced to
> coefficients for orthogonal polynomials that are nice integers.  For
> instance, Kirk's experimental design book presents the following
> coefficients for p=4:
> 
> Linear     -3 -1  1  3
> Quadratic   1 -1 -1  1
> Cubic      -1  3 -3  1
> 
> In R orthogonal polynomials are not integers. For instance, in R where p =4:
> 
> > poly(c(1:4),3)
>               1    2          3
> [1,] -0.6708204  0.5 -0.2236068
> [2,] -0.2236068 -0.5  0.6708204
> [3,]  0.2236068 -0.5 -0.6708204
> [4,]  0.6708204  0.5  0.2236068
> 
> Where, of course, column 1 is linear, column 2 Quadratic and 3 cubic.
> 
> My experience is that the coding scheme used in R works "better" than the
> integer scheme discussed in Kirk for many regression type analyses.
> 
> Can anyone enlighten me as to why?

I think the only difference is that the columns in the orthogonal
polynomial representation in R are scaled to have unit length.  The
rows in the table you give from Kirk's book have lengths sqrt(20), 2,
and sqrt(20) respectively so

> poly(1:4,3)*sqrt(20)
      1         2  3
[1,] -3  2.236068 -1
[2,] -1 -2.236068  3
[3,]  1 -2.236068 -3
[4,]  3  2.236068  1

gives you the first and third rows from Kirk in the first and third
columns.

Although there is some slight numerical advantage in having the
columns of a model matrix of comparable length I don't think it would
be noticeable here.





-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
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