将多项式模型拟合为R中的数据

Meh*_*lar 80 r curve-fitting data-analysis polynomial-math

我已经阅读了这个问题的答案并且它们非常有用,但我需要特别是在R中提供帮助.

我在R中有一个示例数据集,如下所示:

x <- c(32,64,96,118,126,144,152.5,158)  
y <- c(99.5,104.8,108.5,100,86,64,35.3,15)
Run Code Online (Sandbox Code Playgroud)

我想为这些数据拟合一个模型y = f(x).我希望它是一个三阶多项式模型.

我怎么能在R?

另外,R可以帮我找到最合适的模型吗?

Gre*_*reg 91

要获得x(x ^ 3)中的三阶多项式,您可以这样做

lm(y ~ x + I(x^2) + I(x^3))
Run Code Online (Sandbox Code Playgroud)

要么

lm(y ~ poly(x, 3, raw=TRUE))
Run Code Online (Sandbox Code Playgroud)

你可以拟合10阶多项式并得到近乎完美的拟合,但是你应该这样做吗?

编辑:poly(x,3)可能是更好的选择(参见下面的@hadley).

  • 你为什么用`raw = T`?最好使用不相关的变量. (9认同)
  • 是在问"你应该".样本数据只有8个点.这里的自由度非常低.当然,现实生活中的数据可能还要多得多. (6认同)
  • 这取决于您对"最佳模型"的定义.给出最大R ^ 2(10阶多项式将)的模型不一定是"最佳"模型.您需要合理选择模型中的术语.你可以通过很多参数得到近乎完美的拟合,但是模型没有任何预测能力,除了通过这些点绘制最佳拟合线之外的其他任何东西都是无用的. (4认同)
  • 我这样做是为了得到与`lm(y~x + I(x ^ 2)+ I(x ^ 3))`相同的结果.也许不是最优的,只是给同一个目的的两个手段. (2认同)

Gre*_*now 41

哪种型号是"最合适的型号"取决于您所说的"最佳".R有工具可以提供帮助,但您需要提供"最佳"的定义,以便在它们之间进行选择.请考虑以下示例数据和代码:

x <- 1:10
y <- x + c(-0.5,0.5)

plot(x,y, xlim=c(0,11), ylim=c(-1,12))

fit1 <- lm( y~offset(x) -1 )
fit2 <- lm( y~x )
fit3 <- lm( y~poly(x,3) )
fit4 <- lm( y~poly(x,9) )
library(splines)
fit5 <- lm( y~ns(x, 3) )
fit6 <- lm( y~ns(x, 9) )

fit7 <- lm( y ~ x + cos(x*pi) )

xx <- seq(0,11, length.out=250)
lines(xx, predict(fit1, data.frame(x=xx)), col='blue')
lines(xx, predict(fit2, data.frame(x=xx)), col='green')
lines(xx, predict(fit3, data.frame(x=xx)), col='red')
lines(xx, predict(fit4, data.frame(x=xx)), col='purple')
lines(xx, predict(fit5, data.frame(x=xx)), col='orange')
lines(xx, predict(fit6, data.frame(x=xx)), col='grey')
lines(xx, predict(fit7, data.frame(x=xx)), col='black')
Run Code Online (Sandbox Code Playgroud)

哪种型号最好?可以为它们中的任何一个做出参数(但是我不想使用紫色参数进行插值).


Dav*_*uer 15

关于'可以帮助我找到最合适的模型'的问题,可能有一个函数可以做到这一点,假设您可以说明要测试的模型集,但这对于n-1集合来说这将是一个很好的第一种方法.度多项式:

polyfit <- function(i) x <- AIC(lm(y~poly(x,i)))
as.integer(optimize(polyfit,interval = c(1,length(x)-1))$minimum)
Run Code Online (Sandbox Code Playgroud)

笔记

  • 这种方法的有效性将取决于你的目标,假设optimize()AIC()如果AIC的是,你要使用的标准,

  • polyfit()可能没有一个最低限度.检查这个类似于:

    for (i in 2:length(x)-1) print(polyfit(i))
    
    Run Code Online (Sandbox Code Playgroud)
  • 我使用了这个as.integer()函数,因为我不清楚如何解释非整数多项式.

  • 为了测试任意一组数学方程式,请考虑安德鲁·格尔曼(Andrew Gelman)在此审查的"Eureqa"计划

更新

另请参阅stepAIC功能(在MASS包中)以自动选择模型.


小智 5

在R中找到最佳拟合的最简单方法是将模型编码为:

lm.1 <- lm(y ~ x + I(x^2) + I(x^3) + I(x^4) + ...)
Run Code Online (Sandbox Code Playgroud)

使用降压AIC回归后

lm.s <- step(lm.1)
Run Code Online (Sandbox Code Playgroud)

  • 使用"I(x ^ 2)"等不能给出适当的正交多项式拟合. (4认同)