On May 6, 2011, at 12:31 PM, David Winsemius wrote: > On May 6, 2011, at 11:35 AM, Ben Haller wrote: > >> Hi all! I'm getting a model fit from glm() (a binary logistic regression >> fit, but I don't think that's important) for a formula that contains powers >> of the explanatory variable up to fourth. So the fit looks something like >> this (typing into mail; the actual fit code is complicated because it >> involves step-down and so forth): >> >> x_sq <- x * x >> x_cb <- x * x * x >> x_qt <- x * x * x * x >> model <- glm(outcome ~ x + x_sq + x_cb + x_qt, binomial) > > It might have been easier with: > > model <- glm(outcome ~ poly(x, 4) , binomial)
Interesting. The actual model I'm fitting has lots more terms, and needs to be able to be stepped down; sometimes some of the terms of the polynomials will get dropped while others get retained, for example. But more to the point, poly() seems to be doing something tricky involving "orthogonal polynomials" that I don't understand. I don't think I want whatever that is; I want my original variable x, plus its powers. For example, if I do: x < runif(10) poly(x, 3) the columns I get are not x, x^2, x^3; they are something else. So the poly() fit is not equivalent to my fit. > Since the authors of confint might have been expecting a poly() formulation > the results might be more reliable. I'm just using lm() but I think the > conclusions are more general: > > ... > > Coefficients: > Estimate Std. Error t value Pr(>|t|) > (Intercept) 0.45737 0.52499 0.871 0.3893 > x -0.75989 1.15080 -0.660 0.5131 > x2 1.30987 0.67330 1.945 0.0594 . > x3 -0.03559 0.11058 -0.322 0.7494 > > ... > > Coefficients: > Estimate Std. Error t value Pr(>|t|) > (Intercept) 5.4271 0.1434 37.839 < 2e-16 *** > poly(x, 3)1 30.0235 0.9184 32.692 < 2e-16 *** > poly(x, 3)2 8.7823 0.9184 9.563 1.53e-11 *** > poly(x, 3)3 -0.2956 0.9184 -0.322 0.749 Here, in your illustration, is underscored what I mean. Whatever these orthogonal polynomial terms are that you're using, they are clearly not the original x, x^2 and x^3, and they're giving you a different fit than those do. I probably ought to learn about this technique, since it looks interesting; but for my purposes I need the fit to actually be in terms of x, since x is my explanatory variable. And the fit I'm getting is highly significant (all terms < 2e-16), so the lack of fit problem you're illustrating does not seem to apply to my case. Or maybe I'm totally misunderstanding your point...? :-> Thanks! Ben Haller McGill University http://biology.mcgill.ca/grad/ben/ ______________________________________________ R-help@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.
