在 R 中，如何通过最小二乘法估计多项式函数？答案

【问题标题】：In R, how to estimate a polynomial function by Least Square method?在 R 中，如何通过最小二乘法估计多项式函数？
【发布时间】：2019-11-30 07:23:57
【问题描述】：

n = 50
set.seed(100)
x = matrix(runif(n, -2, 2), nrow=n)
y = 2 + 0.75*sin(x) - 0.75*cos(x) + rnorm(n, 0, 0.2)

在 R 中，我想用最小二乘法估计上面的多项式函数。

这意味着我想知道 γ0、γ1、γ2 和 γ3 的估计值。

而且我也想知道这个估计函数下的MSE。

我用过这个

summary(lm(y ~ x+ x^2+ x^3))

但只要得到this output:

Call:
lm(formula = y ~ x + (x^2) + (x^3))

Residuals:
    Min      1Q  Median      3Q     Max 
-0.66448 -0.22251 -0.07694  0.20647  0.79429

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  1.54972    0.04761   32.55    <2e-16 ***
x            0.65279    0.04633   14.09    <2e-16 ***   
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.3319 on 48 degrees of freedom
Multiple R-squared:  0.8053,    Adjusted R-squared:  0.8013
F-statistic:  198.6 on 1 and 48 DF,  p-value: < 2.2e-16

请告诉我在 R 中，我可以使用什么函数或包来做到这一点。

谢谢。

【问题讨论】：

试试summary(lm(y ~ x + I(x^2)+ I(x^3))) 或summary(lm(y ~ poly(x, 3, raw=TRUE))。
谢谢，为什么只在 (x^2) 和 (x^3) 上加上 ''I'' 就可以工作？

标签： r statistics least-squares

【解决方案1】：

您需要将 I(.) 包裹在多项式周围，以指示 R 应将 ^ 运算符用作算术运算符而不是公式运算符。

让我们创建一些示例数据来说明。

set.seed(42)
x <- runif(100)
y <- x + x^2 + x^3 + rnorm(100)

你想要的是这个：

lm(y ~ x + I(x^2) + I(x^3))$coe
# (Intercept)           x      I(x^2)      I(x^3) 
#  -0.4069207   3.4580770  -4.0516060   4.3227674

或：

lm(y ~ poly(x, 3, raw=TRUE))$coe
# (Intercept) poly(x, 3, raw = TRUE)1 poly(x, 3, raw = TRUE)2 poly(x, 3, raw = TRUE)3 
#  -0.4069207               3.4580770              -4.0516060               4.3227674

没有I(.) 结果是不同的，

lm(y ~ x + x^2 + x^3)$coe
# (Intercept)           x 
#  -0.6149136   3.3590395

因为出于某种原因，R 将论坛解释为

lm(y ~ x)$coe
# (Intercept)           x 
#  -0.6149136   3.3590395 

all(lm(y ~ x)$coe - lm(y ~ x + x^2 + x^3)$coe == 0)
# [1] TRUE

但是，请考虑：

lm(y ~ I(x + x^2 + x^3))$coe
# (Intercept) I(x + x^2 + x^3) 
#   -0.191611         1.141912

创建变量z

z <- x + x^2 + x^3

给出相同的结果：

lm(y ~ z)$coe
# (Intercept)           z 
#   -0.191611    1.141912

或者，更明确地说，假设我们想计算x 与a 的交互，我们想使用* 作为公式运算符：

a <- runif(100)
lm(y ~ x*a)$coe
# (Intercept)           x           a         x:a 
# -0.71920356  3.42008631  0.19180499 -0.07049342

当我们使用*作为算术运算符时，

lm(y ~ I(x*a))$coe
# (Intercept)    I(x * a) 
#   0.5317547   2.7361742

结果是一样的：

xa <- x*a

lm(y ~ xa)$coe
# (Intercept)          xa 
#   0.5317547   2.7361742

【讨论】：