为什么在R中建立线性多项式模型时我() "AsIs“是必需的？

Question

当使用线性多项式模型或函数I()时，我试图理解R中的poly基函数的作用。当我计算模型时

1. q + q^2
2. q + I(q^2)
3. poly(q, 2)

我有不同的答案。

下面是一个示例：

set.seed(20)
q <- seq(from=0, to=20, by=0.1)
y <- 500 + .1 * (q-5)^2
noise <- rnorm(length(q), mean=10, sd=80)
noisy.y <- y + noise
model3 <- lm(noisy.y ~ poly(q,2))
model1 <- lm(noisy.y ~ q + I(q^2))
model2 <- lm(noisy.y ~ q + q^2)
I(q^2)==I(q)^2
I(q^2)==q^2

summary(model1)
summary(model2)
summary(model3)

这是输出：

> summary(model1)

Call:
lm(formula = noisy.y ~ q + I(q^2))

Residuals:
     Min       1Q   Median       3Q      Max 
-211.592  -50.609    4.742   61.983  165.792 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) 489.3723    16.5982  29.483   <2e-16 ***
q             5.0560     3.8344   1.319    0.189    
I(q^2)       -0.1530     0.1856  -0.824    0.411    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 79.22 on 198 degrees of freedom
Multiple R-squared:  0.02451,   Adjusted R-squared:  0.01466 
F-statistic: 2.488 on 2 and 198 DF,  p-value: 0.08568

> summary(model2)

Call:
lm(formula = noisy.y ~ q + q^2)

Residuals:
    Min      1Q  Median      3Q     Max 
-219.96  -54.42    3.30   61.06  170.79 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) 499.5209    11.1252  44.900   <2e-16 ***
q             1.9961     0.9623   2.074   0.0393 *  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 79.16 on 199 degrees of freedom
Multiple R-squared:  0.02117,   Adjusted R-squared:  0.01625 
F-statistic: 4.303 on 1 and 199 DF,  p-value: 0.03933

> summary(model3)

Call:
lm(formula = noisy.y ~ poly(q, 2))

Residuals:
     Min       1Q   Median       3Q      Max 
-211.592  -50.609    4.742   61.983  165.792 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  519.482      5.588  92.966   <2e-16 ***
poly(q, 2)1  164.202     79.222   2.073   0.0395 *  
poly(q, 2)2  -65.314     79.222  -0.824   0.4107    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 79.22 on 198 degrees of freedom
Multiple R-squared:  0.02451,   Adjusted R-squared:  0.01466 
F-statistic: 2.488 on 2 and 198 DF,  p-value: 0.08568

为什么在R中做多项式模型时I()是必需的。

另外，poly函数没有给出与q + I(q^2)相同的结果，这是否正常？

Answer 1

R中的公式语法在?formula帮助页中进行了描述。^符号还没有被赋予通常的乘法幂的含义。相反，它用于指数基础上所有项之间的相互作用。例如

y ~ (a+b)^2

是相同的

y ~ a + b + a:b

但如果你这么做

y ~ a + b^2
y ~ a + b    # same as above, no way to "interact" b with itself.

该插入符号将只包含b术语，因为它不能包括与自身的交互。因此，公式中的^和*与乘法无关，就像+实际上并不意味着通常意义上的变量的加法。

如果您想要^2的“通常”定义，您需要将其放在as is函数中。否则，它根本不适合平方项。

默认情况下，poly()函数返回帮助页面中描述的正交多项式。这有助于减少协变量中的共线性。但是，如果您不需要正交版本，只需要“原始”多项式项，那么只需将raw=TRUE传递给您的poly调用即可。例如

lm(noisy.y ~ poly(q,2, raw=TRUE))

将返回与model1相同的估计值。