Приведение полиномиальной модели к данным в R
Я прочитал ответы на этот question, и они очень полезны, но мне нужна помощь, особенно в 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)
Я хочу подгонять модель к этим данным, чтобы y = f(x)
. Я хочу, чтобы это была модель полинома 3-го порядка.
Как я могу это сделать в R?
Кроме того, может ли R помочь мне найти подходящую модель?
Ответы
Ответ 1
Чтобы получить полином третьего порядка по x (x ^ 3), вы можете сделать
lm(y ~ x + I(x^2) + I(x^3))
или
lm(y ~ poly(x, 3, raw=TRUE))
Вы можете поместить многочлен 10-го порядка и получить почти идеальную подгонку, но должны ли вы?
EDIT:
poly (x, 3), вероятно, лучший выбор (см. ниже @hadley).
Ответ 2
Какая модель является "лучшей подходящей моделью", зависит от того, что вы подразумеваете под "лучшим". У 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')
Какая из этих моделей лучшая? аргументы могут быть сделаны для любого из них (но я бы не хотел использовать фиолетовый для интерполяции).
Ответ 3
Что касается вопроса "может ли R помочь мне найти лучшую модель подгонки", возможно, есть функция, предполагающая, что вы можете указать набор моделей для тестирования, но это было бы хорошим первым подходом к набору n-1 степени:
polyfit <- function(i) x <- AIC(lm(y~poly(x,i)))
as.integer(optimize(polyfit,interval = c(1,length(x)-1))$minimum)
Примечания
-
Действительность этого подхода будет зависеть от ваших целей, предположений optimize()
и AIC()
, и если AIC является критерием, который вы хотите использовать,
-
polyfit()
может не иметь ни одного минимума. проверьте это следующим образом:
for (i in 2:length(x)-1) print(polyfit(i))
-
Я использовал функцию as.integer()
, потому что мне не ясно, как я буду интерпретировать нецелый многочлен.
-
для тестирования произвольного набора математических уравнений рассмотрим программу 'Eureqa, рассмотренную Эндрю Гельманом здесь
Обновление
Также см. функцию stepAIC
(в пакете MASS) для автоматизации выбора модели.
Ответ 4
Самый простой способ найти наилучшее соответствие R - это код модели:
lm.1 <- lm(y ~ x + I(x^2) + I(x^3) + I(x^4) + ...)
После использования ступенчатой регрессии AIC
lm.s <- step(lm.1)