Прогнозирование значений в линейной модели OLS с индексированной независимой переменной (годы)

Это один из тех вопросов, для ответа на который, вероятно, есть миллион способов, которые сделают фактический ответ неуместным, но упрямство мешает...

При попытке понять применение временных рядов становится ясно, что удаление тренда данных делает прогнозирование будущих значений неправдоподобным. Например, используя набор данных gtemp из пакета astsa, необходимо учитывать тенденцию к росту в последние десятилетия:

введите здесь описание изображения

Таким образом, я получаю модель ARIMA (правильную или неправильную) для данных без тренда, которая позволяет мне «прогнозировать» на 10 лет вперед:

fit = arima(gtemp, order = c(4, 1, 1))
pred = predict(fit, n.ahead = 10)

и оценщик тренда OLS, основанный на значениях с 1950 года:

gtemp1 = window(gtemp, start = 1950, end = 2009)
fit2 = lm(gtemp1 ~ I(1950:2009))

Проблема заключается в том, как использовать predict() для получения оценочного значения для линейной части модели в следующие 10 лет.

Если я запускаю predict(fit2, data.frame(I(2010:2019))), я получаю 60 значений, которые я получил бы при запуске predict(fit2), плюс сообщение об ошибке: 'newdata' had 10 rows but variables found have 60 rows.


r time-series predict lm
person Antoni Parellada    schedule 15.08.2016    source источник

Ответы (1)


arrow_upward
1
arrow_downward

Тебе нужно:

dat <- data.frame(year = 1950:2009, gtemp1 = as.numeric(gtemp1))
fit2 <- lm(gtemp1 ~ year, data = dat)
unname( predict(fit2, newdata = data.frame(year = 2010:2019)) )
# [1] 0.4928475 0.5037277 0.5146079 0.5254882 0.5363684 0.5472487 0.5581289
# [8] 0.5690092 0.5798894 0.5907697

В качестве альтернативы, если вы не хотите использовать аргумент data в lm, вам нужно:

year <- 1950:2009
fit2 <- lm(gtemp1 ~ year)
unname( predict(fit2, newdata = data.frame(year = 2010:2019)) )
# [1] 0.4928475 0.5037277 0.5146079 0.5254882 0.5363684 0.5472487 0.5581289
# [8] 0.5690092 0.5798894 0.5907697

Почему исходный код не работает

Когда вы выполняете fit2 <- lm(gtemp1 ~ I(1950:2009)), lm предполагает наличие ковариаты с именем I(1950:2009):

attr(fit2$terms, "term.labels")  ## names of covariates
# [1] "I(1950:2009)"

Когда вы сделаете прогноз позже, predict попытается найти переменную в вашем новом фрейме данных, называемую I(1950:2009). Однако взгляните на имена столбцов вашего newdata:

colnames( data.frame(I(2010:2019)) )
# [1] "X2010.2019"

В результате predict.lm не может найти переменную I(1950:2009) в newdata, тогда он будет использовать фрейм внутренней модели fit2$model как newdata и по умолчанию будет возвращать подогнанные значения (что объясняет, почему вы получаете 60 значений).

person Zheyuan Li    schedule 15.08.2016
comment
Вы знаете, как я мог бы добавить эти предсказанные значения, основанные на OLS, к предсказанным значениям ARIMA, чтобы получить что-то разумное? - person Antoni Parellada; 15.08.2016
comment
Я только начинаю пытаться понять временные ряды. Я видел, что вы можете отказаться от тренда, выполнив что-то вроде model=lm(y ~ I(1950:2009)); st.y = gtemp - predict(model) (не совсем то же самое, что и OP). Итак, я пытался теперь вместо вычитания добавить результаты прогнозирования () ARIMA к прогнозу OLS. Но я получаю огромный шаг вверх. - person Antoni Parellada; 15.08.2016
comment
Что касается вашего последнего комментария, pred будет fit = arima(gtemp, order = c(4, 1, 1)); pred = predict(fit, n.ahead = 10)? Отказаться от ОЛС полностью? - person Antoni Parellada; 15.08.2016
comment
Объявить вклад OLS ошибочным? - person Antoni Parellada; 15.08.2016
comment
Да, (fit = arima(gtemp, order = c(4, 1, 1)));pred = predict(fit, n.ahead = 10);pred$pred = as.ts(cumsum(pred$pred)); ts.plot(gtemp, pred$pred, lty = c(1,3), col=c(5,2))... Некрасивое зрелище... Я уверен, что по-царски все испортил. - person Antoni Parellada; 15.08.2016