Как создать линию прогноза для квадратичной модели
Я пытаюсь создать линию квадратичного прогноза для квадратичной модели. Я использую набор данных Auto, который поставляется с R. У меня не было проблем с созданием линии прогноза для линейной модели. Однако квадратичная модель дает сумасшедшие линии. Вот мой код.
# Linear Model
plot(Auto$horsepower, Auto$mpg,
main = "MPG versus Horsepower",
pch = 20)
lin_mod = lm(mpg ~ horsepower,
data = Auto)
lin_pred = predict(lin_mod)
lines(
Auto$horsepower, lin_pred,
col = "blue", lwd = 2
)
# The Quadratic model
Auto$horsepower2 = Auto$horsepower^2
quad_model = lm(mpg ~ horsepower2,
data = Auto)
quad_pred = predict(quad_model)
lines(
Auto$horsepower,
quad_pred,
col = "red", lwd = 2
)
Я на 99% уверен, что проблема заключается в функции прогнозирования. Почему я не могу построить аккуратно выглядящую кривую квадратичного прогноза? Следующий код, который я попробовал, не работает - может быть это связано ?:
quad_pred = predict(quad_model, data.frame(horsepower = Auto$horsepower))
Спасибо!
Ответы 3
Проблема в том, что значения x-axis не отсортированы. Неважно, если бы это была линейная модель, но было бы заметно, если бы она была полиномиальной. Я создал новый отсортированный набор данных, и он отлично работает:
library(ISLR) # To load data Auto
# Linear Model
plot(Auto$horsepower, Auto$mpg,
main = "MPG versus Horsepower",
pch = 20)
lin_mod = lm(mpg ~ horsepower,
data = Auto)
lin_pred = predict(lin_mod)
lines(
Auto$horsepower, lin_pred,
col = "blue", lwd = 2
)
# The Quadratic model
Auto$horsepower2 = Auto$horsepower^2
# Sorting Auto by horsepower2
Auto2 <- Auto[order(Auto$horsepower2), ]
quad_model = lm(mpg ~ horsepower2,
data = Auto2)
quad_pred = predict(quad_model)
lines(
Auto2$horsepower,
quad_pred,
col = "red", lwd = 2
)
Скорее всего, вы получите ошибку Error in eval(predvars, data, env) : object 'horsepower2' not found. Вы должны изменить его на это: quad_pred = predict(quad_model, data.frame(horsepower2 = sort(Auto$horsepower^2)))
Один из вариантов - создать последовательность значений x, для которой вы хотите построить подходящие линии. Это может быть полезно, если в ваших данных есть «пробел» или если вы хотите построить соответствующие линии за пределами диапазона переменных x.
# load dataset; if necessary run install.packages("ISLR")
data(Auto, package = "ISLR")
# since only 2 variables at issue, use short names
mpg <- Auto$mpg
hp <- Auto$horsepower
# fit linear and quadratic models
lmod <- lm(mpg ~ hp)
qmod <- lm(mpg ~ hp + I(hp^2))
# plot the data
plot(x=hp, y=mpg, pch=20)
# use predict() to find coordinates of points to plot
x_coords <- seq(from=floor(min(hp)), to=ceiling(max(hp)), by=1)
y_coords_lmod <- predict(lmod, newdata=data.frame(hp=x_coords))
y_coords_qmod <- predict(qmod, newdata=data.frame(hp=x_coords))
# alternatively, calculate this manually using the fitted coefficients
y_coords_lmod <- coef(lmod)[1] + coef(lmod)[2]*x_coords
y_coords_qmod <- coef(qmod)[1] + coef(qmod)[2]*x_coords + coef(qmod)[3]*x_coords^2
# add the fitted lines to the plot
points(x=x_coords, y=y_coords_lmod, type = "l", col = "blue")
points(x=x_coords, y=y_coords_qmod, type = "l", col = "red")
Спасибо! Это выглядит довольно интуитивно. Мне нравится, как вы вручную создали последовательность. Я думаю, это поможет сгладить линию / кривую прогноза.
В качестве альтернативы, используя ggplot2:
ggplot(Auto, aes(x = horsepower, y = mpg)) + geom_point() +
stat_smooth(aes(x = horsepower, y = mpg), method = "lm", formula = y ~ x, colour = "red") +
stat_smooth(aes(x = horsepower, y = mpg), method = "lm", formula = y ~ poly(x, 2), colour = "blue")
Спасибо! Я ценю показ альтернативного метода на ggplot. Выглядит очень аккуратно. Мне очень нравится диапазон интервалов предсказания.
Спасибо! Похоже, вы создали полностью новый фрейм данных только для сортировки квадратов значений мощности в лошадиных силах. Почему я не могу создать отсортированные значения в функции прогнозирования? Например, как я упоминал выше, я получаю сообщение об ошибке при запуске строки:
predict(quad_model, data.frame(horsepower = sort(Auto$horsepower))). Я мог бы поклясться, что использование функции data.frame () в качестве аргумента в функции прогнозирования позволит управлять желаемыми значениями x для прогнозирования - следовательно, сортировать их там, а не создавать фрейм данных Auto2 в целом.