I haven't used poly() to be able to check out the error you reported. However
in the days of accurate statistical computing, I don't see any advantage that
justifies the
complexity of orthogonal polynomials. Fitting ordinary polynomials works fine
although
splines are better. Note that ordinary polynomials have highly collinear terms
but this is
of no consequence as this does not affect predicted values, standard errors of
predicted values, or multiple degree of freedom tests of association.
In the Design library I have a little function pol for regular polynomials.
You could easily
strip this function to come up with something like the following.
Speed of the code could be improved
pol <- function(x, poly.degree) {
if(poly.degree < 2) stop('poly.degree must be >= 2')
xd <- matrix(0,nrow=length(x),ncol=poly.degree)
nam <- as.character(substitute(x))
name <- character(poly.degree)
name[1] <- nam
xd[,1] <- x
for(j in 2:poly.degree) {
name[j] <- paste(nam,"^",j,sep="")
xd[,j] <- x^j }
dimnames(xd) <- list(names(x), name)
xd
}
. -Frank Harrell
Daniel Erdin wrote:
> Hello,
> I'm working with Splus 4.5 on Windows NT. I have some data with polynomial
> dependencies. As I've understood, the use of poly() has statistical
> advantages in the modelling of polynomial dependencies.
> In the Splus-books I have (for version 3.3, 1995), it is shown, that for a
> model
> lmfit1_lm(y~poly(x1,3))
> the coefficients corresponding to the untransformed values of the
> polynomials of x can be obtained through
> poly.transform(poly(x1,3),lmfit1)
> This works fine.
>
> But I have more than the covariable x1 in my model and if I try
> lmfit2_lm(y~poly(x1,3)+x2+x3)
> I get no result with the following function
> poly.transform(poly(x1,3),lmfit2)
> just an error message
> Error in poly.transform(........): wrong length for coefs
>
> I understand the reason for the error message, but how do I have to use
> poly.transform() to get what I want? Or is there an other solution to get
> the estimates of the coefficients for the not transformed polynomials? At
> last, I have the same problem in lme-models and would be interested if the
> solution in lme is the same.
>
> I've searched the archives (but maybe not enough...). I'll report the
> result to the list.
> Thanks
> Daniel
>
> ***************************************************
> Daniel Erdin
> Swiss Federal Institute of Technology
> Institute of Animal Science
> Animal Breeding
> Research Station Chamau
> CH-6331 Huenenberg
> Phone: +41 41 780 10 47
> Fax: +41 41 780 43 83
> e-mail: erdin@inw.agrl.ethz.ch
> homepage: http://www.tz.inw.agrl.ethz.ch/~erdin/
> -----------------------------------------------------------------------
> 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
--
Frank E Harrell Jr
Professor of Biostatistics and Statistics
Division of Biostatistics and Epidemiology
Department of Health Evaluation Sciences
University of Virginia School of Medicine
http://hesweb1.med.virginia.edu/biostat
-----------------------------------------------------------------------
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
|