lm 中的多项式回归（2 种方法）

Question

我正在尝试使用 R 进行多项式回归lm，我有

Time=c(1980:2016)
y= rnorm(length(Time))

我用了 ：

reg=lm(y~poly(Time,2))
round(reg$coefficients,3)

这使 ：

 (Intercept) poly(Time, 2)1 poly(Time, 2)2 
        -0.110         -1.298          0.172

和 ：

Time2=Time^2
reg2=lm(y~Time+Time2)
round(reg2$coefficients,3)

它给

(Intercept)        Time       Time2 
   1146.590      -1.128       0.000 

哪里有问题？

Answer 1

默认情况下poly()使用正交多项式表示，它在数值上更稳定。raw=TRUE如果您想匹配朴素的表示，则可以使用。

set.seed(101)
dd <- data.frame(Time=c(1980:2016),
                y=rnorm(2016-1980+1))
(c1 <- coef(lm(y~Time+I(Time^2),dd)))
##   (Intercept)          Time     I(Time^2) 
##  6.684138e+03 -6.686392e+00  1.672101e-03 
(c2 <- coef(lm(y~poly(Time,2),dd)))
##    (Intercept) poly(Time, 2)1 poly(Time, 2)2 
##    -0.04713527    -0.30359154     1.03594479 
c3 <- coef(lm(y~poly(Time,2,raw=TRUE),dd))
all.equal(unname(c1),unname(c3))  ## TRUE

您会注意到，原始多项式给出的截距和斜率指的是Time=0 时的预期值，如果您使用 CE（Common Era=AD）年的时间，则这有点荒谬。

如果您想要可解释性和数值稳定性，您可以通过将时间变量居中来获得合理的折衷（在观察期开始时将时间设置为零也是合理的）。

dd$cTime <- dd$Time-mean(dd$Time)
c4 <- coef(lm(y~poly(cTime,2,raw=TRUE),dd))
unname(c4)
## [1] -0.237754839 -0.004674513  0.001672101