[R] ordinary polynomial coefficients from orthogonal polynomials?

Roger D. Peng rpeng at jhsph.edu
Tue Jun 14 15:48:41 CEST 2005


I think this is covered in the FAQ:

http://cran.r-project.org/doc/FAQ/R-FAQ.html#Models

-roger

James Salsman wrote:
> How can ordinary polynomial coefficients be calculated
> from an orthogonal polynomial fit?
> 
> 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"
> 
>  > 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
> 

-- 
Roger D. Peng
http://www.biostat.jhsph.edu/~rpeng/




More information about the R-help mailing list