Подбор полиномиальной модели к данным в R

86

Я прочитал ответы на этот вопрос, и они очень полезны, но мне нужна помощь, особенно в 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 помочь мне найти наиболее подходящую модель?

Мехпер К. Палавузлар
источник

Ответы:

102

Чтобы получить полином третьего порядка от 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 ниже).

Грег
источник
6
правильно спрашивает "вы должны". В данных выборки всего 8 точек. Степени свободы здесь довольно низкие. Конечно, реальных данных может быть намного больше.
JD Long
1
Спасибо за Ваш ответ. Как насчет того, чтобы R подобрал наиболее подходящую модель? Есть ли для этого какие-то функции?
Mehper C. Palavuzlar
5
Это зависит от вашего определения «лучшей модели». Модель, которая дает вам наибольшее R ^ 2 (что дает многочлен 10-го порядка), не обязательно является «лучшей» моделью. Условия в вашей модели должны быть выбраны разумно. Вы можете получить почти идеальное совпадение с множеством параметров, но модель не будет иметь предсказательной силы и будет бесполезна для чего-либо, кроме построения наилучшей линии через точки.
Грег,
11
Почему вы используете raw = T? Лучше использовать некоррелированные переменные.
Хэдли
2
Я сделал это, чтобы получить те же результаты, что и lm(y ~ x + I(x^2) + I(x^3)). Возможно, не оптимально, просто предоставление двух средств для достижения одной цели.
Грег,
46

Какая модель является «наиболее подходящей», зависит от того, что вы подразумеваете под «наилучшей». В 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')

Какая из этих моделей лучшая? аргументы могут быть приведены для любого из них (но я, например, не хотел бы использовать фиолетовый для интерполяции).

Грег Сноу
источник
15

Что касается вопроса `` может ли 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) для автоматизации выбора модели.

Дэвид Лебауэр
источник
Как я могу связать Eurequa с R?
adam.888
@ adam.888 отличный вопрос - я не знаю ответа, но вы можете опубликовать его отдельно. Последний пункт был небольшим отступлением.
Дэвид Лебауэр
Примечание. AIC - это информационный критерий Акаике , который награждает близкое соответствие и штрафует большее количество параметров модели, что было показано как оптимальное в различных смыслах. en.wikipedia.org/wiki/Akaike_information_criterion
Евгений Сергеев
4

Самый простой способ найти наилучшее соответствие в R - это закодировать модель как:

lm.1 <- lm(y ~ x + I(x^2) + I(x^3) + I(x^4) + ...)

После использования понижающей регрессии AIC

lm.s <- step(lm.1)
Мэтью Фидлер
источник
6
Использование I(x^2)и т.д. не дает подходящих ортогональных многочленов для подбора.
Брайан Диггс,