r - lm()中的poly():原始与正交之间的差异

Question

我有

library(ISLR)
attach(Wage)

# Polynomial Regression and Step Functions

fit=lm(wage~poly(age,4),data=Wage)
coef(summary(fit))

fit2=lm(wage~poly(age,4,raw=T),data=Wage)
coef(summary(fit2))

plot(age, wage)
lines(20:350, predict(fit, newdata = data.frame(age=20:350)), lwd=3, col="darkred")
lines(20:350, predict(fit2, newdata = data.frame(age=20:350)), lwd=3, col="darkred")

预测线似乎相同，但是为什么系数如此不同？您如何用 raw=T和 raw=F解释它们。

我看到 poly(...,raw=T)产生的系数与 ~age+I(age^2)+I(age^3)+I(age^4)产生的系数匹配。

如果我想使用系数“手动”获得预测(不使用 predict()函数)，我应该注意些什么？我应该如何解释 poly()中正交多项式的系数。

Answer 1

默认情况下，使用raw = FALSE，poly()计算正交多项式。它在内部用原始编码x，x ^ 2，x ^ 3，...来建立模型矩阵，然后缩放列，以使每列与先前的列正交。这不会更改拟合值，但具有的优势是，您可以看到多项式中的某个阶数是否比低阶阶数显着改善了回归。

考虑具有响应的简单cars数据，该响应会停止dist并驱动speed。从物理上讲，这应该具有二次关系，但是在这个(旧的)数据集中，二次项并不重要:

m1 <- lm(dist ~ poly(speed, 2), data = cars)
m2 <- lm(dist ~ poly(speed, 2, raw = TRUE), data = cars)

在正交编码中，您可以在 summary(m1)中获得以下系数:

                Estimate Std. Error t value Pr(>|t|)    
(Intercept)       42.980      2.146  20.026  < 2e-16 ***
poly(speed, 2)1  145.552     15.176   9.591 1.21e-12 ***
poly(speed, 2)2   22.996     15.176   1.515    0.136    

这表明存在高度显着的线性效应，而二阶不显着。后一个p值(即多项式中的最高阶之一)与原始编码中的相同:

                            Estimate Std. Error t value Pr(>|t|)
(Intercept)                  2.47014   14.81716   0.167    0.868
poly(speed, 2, raw = TRUE)1  0.91329    2.03422   0.449    0.656
poly(speed, 2, raw = TRUE)2  0.09996    0.06597   1.515    0.136

但是较低阶的p值会发生巨大变化。原因是在 m1模型中，回归变量是正交的，而它们在 m2中高度相关:

cor(model.matrix(m1)[, 2], model.matrix(m1)[, 3])
## [1] 4.686464e-17
cor(model.matrix(m2)[, 2], model.matrix(m2)[, 3])
## [1] 0.9794765

因此，在原始编码中，如果 speed保留在模型中，则只能解释 speed^2的p值。而且由于两个回归指标都高度相关，因此可以删除其中之一。但是，在正交编码中， speed^2仅捕获线性项尚未捕获的二次部分。然后清楚的是，线性部分是有效的，而二次部分则没有附加的意义。