[R] ordinary polynomial coefficients from orthogonal polynomials?

[Frank E Harrell Jr]

Prof Brian Ripley wrote:
> On Tue, 14 Jun 2005, James Salsman wrote:
> 
> 
>>How can ordinary polynomial coefficients be calculated
>>from an orthogonal polynomial fit?
> 
> 
> Why would you want to do that?  predict() is perfectly happy with an
> orthogonal polynomial fit and the `ordinary polynomial coefficients' are 
> rather badly determined in your example since the design matrix has a very 
> high condition number.

Brian - I don't fully see the relevance of the high condition number 
nowadays unless the predictor has a really bad origin.  Orthogonal 
polynomials are a mess for most people to deal with.

Frank

> 
> 
>>I'm trying to do something like find a,b,c,d from
>> lm(billions ~ a+b*decade+c*decade^2+d*decade^3)
>>but that gives:  "Error in eval(expr, envir, enclos) :
>>Object "a" not found"
> 
> 
> You could use
> 
> lm(billions ~ decade + I(decade^2) + I(decade^3))
> 
> except that will be numerically inaccurate, since
> 
> 
>>m <- model.matrix(~ decade + I(decade^2) + I(decade^3))
>>kappa(m)
> 
> [1] 3.506454e+16
> 
> 
> 
> 
>>>decade <- c(1950, 1960, 1970, 1980, 1990)
>>>billions <- c(3.5, 5, 7.5, 13, 40)
>>># source: http://www.ipcc.ch/present/graphics/2001syr/large/08.17.jpg
>>>
>>>pm <- lm(billions ~ poly(decade, 3))
>>>
>>>plot(decade, billions, xlim=c(1950,2050), ylim=c(0,1000),
>>
>>main="average yearly inflation-adjusted dollar cost of extreme weather
>>events worldwide")
>>
>>>curve(predict(pm, data.frame(decade=x)), add=TRUE)
>>># output: http://www.bovik.org/storms.gif
>>>
>>>summary(pm)
>>
>>Call:
>>lm(formula = billions ~ poly(decade, 3))
>>
>>Residuals:
>>      1       2       3       4       5
>> 0.2357 -0.9429  1.4143 -0.9429  0.2357
>>
>>Coefficients:
>>                 Estimate Std. Error t value Pr(>|t|)
>>(Intercept)        13.800      0.882  15.647   0.0406 *
>>poly(decade, 3)1   25.614      1.972  12.988   0.0489 *
>>poly(decade, 3)2   14.432      1.972   7.318   0.0865 .
>>poly(decade, 3)3    6.483      1.972   3.287   0.1880
>>---
>>Signif. codes:  0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1
>>
>>Residual standard error: 1.972 on 1 degrees of freedom
>>Multiple R-Squared: 0.9957,     Adjusted R-squared: 0.9829
>>F-statistic: 77.68 on 3 and 1 DF,  p-value: 0.08317
>>
>>
>>>pm
>>
>>Call:
>>lm(formula = billions ~ poly(decade, 3))
>>
>>Coefficients:
>>     (Intercept)  poly(decade, 3)1  poly(decade, 3)2  poly(decade, 3)3
>>          13.800            25.614            14.432             6.483
>>
>>______________________________________________
>>R-help at stat.math.ethz.ch mailing list
>>https://stat.ethz.ch/mailman/listinfo/r-help
>>PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
>>
> 
> 


-- 
Frank E Harrell Jr   Professor and Chair           School of Medicine
                      Department of Biostatistics   Vanderbilt University