Прогнозы на 1 шаг вперед с пакетом dynlm R

11

Я установил модель с несколькими независимыми переменными, одна из которых - задержка зависимой переменной, используя пакет dynlm.

Предполагая, что у меня есть прогнозы на 1 шаг вперед для моих независимых переменных, как мне получить прогнозы на 1 шаг вперед для моих зависимых переменных?

Вот пример:

library(dynlm)

y<-arima.sim(model=list(ar=c(.9)),n=10) #Create AR(1) dependant variable
A<-rnorm(10) #Create independant variables
B<-rnorm(10)
C<-rnorm(10)
y<-y+.5*A+.2*B-.3*C #Add relationship to independant variables 
data=cbind(y,A,B,C)

#Fit linear model
model<-dynlm(y~A+B+C+L(y,1),data=data)

#Forecast
A<-c(A,rnorm(1)) #Assume we already have 1-step forecasts for A,B,C
B<-c(B,rnorm(1))
C<-c(C,rnorm(1))
y=window(y,end=end(y)+c(1,0),extend=TRUE)
newdata<-cbind(y,A,B,C)
predict(model,newdata)

А вот пример использования пакета dyn, который работает.

library(dyn)

#Fit linear model
model<-dyn$lm(y~A+B+C+lag(y,-1),data=data)

#Forecast
predict(model,newdata)the dyn packages, which works:
Zach
источник
Использование только dynlmпакета не обеспечит прогнозы для ваших зависимых переменных. Для предоставления прогнозов для ваших зависимых переменных потребуется модель для их объяснения и, возможно, дополнительные данные. Я предлагаю вам прочитать кое-что о многомерной регрессии, например, «Прикладной многомерный статистический анализ» Джонсона и Вихерна. или курс по прогнозированию: duke.edu/~rnau/411home.htm
deps_stats
1
@deps_stats Зависимая переменная - это то, что я хочу прогнозировать. Я предполагаю, что у меня уже есть прогнозы для моих независимых переменных. В моем примере кода y - это зависимая переменная, которую я пытаюсь прогнозировать, а A, B, C - независимые переменные, для которых у меня уже есть прогнозы. Если вы запустите пример кода, который я разместил, вы поймете природу моей проблемы.
Зак
@ Зак: Хороший рейтинг Kaggle! (Я щелкнул по ссылке в вашем профиле при наведении курсора мыши)
Хью Перкинс

Ответы:

13

Поздравляем, вы нашли ошибку. Прогнозирование dynlmс новыми данными нарушается, если используются лаговые переменные. Чтобы понять, зачем смотреть на вывод

predict(model)
predict(model,newdata=data)

Результаты должны быть одинаковыми, но это не так. Без newdataаргументов predictфункция в основном захватывает modelэлемент из dynlmвывода. С newdataаргументом predictпытается сформировать новую модель матрицы из newdata. Так как это включает в себя формулу синтаксического анализа, которая передается, dynlmи формула имеет функцию L, которая определена только внутри функции dynlm, формируется неправильная матрица модели. Если вы попытаетесь отладить, вы увидите, что в случае newdataаргумента lagged зависимая переменная не отстает .

Что вы можете сделать, так это отстать от зависимой переменной и включить ее в newdata. Вот код, иллюстрирующий этот подход. Я использую set.seedтак, чтобы это было легко воспроизводимо.

library(dynlm)

set.seed(1)
y<-arima.sim(model=list(ar=c(.9)),n=10) #Create AR(1) dependant variable
A<-rnorm(10) #Create independant variables
B<-rnorm(10)
C<-rnorm(10)
y<-y+.5*A+.2*B-.3*C #Add relationship to independant variables 
data=cbind(y,A,B,C)

#Fit linear model
model<-dynlm(y~A+B+C+L(y,1),data=data)

Вот ошибочное поведение:

> predict(model)
       2        3        4        5        6        7        8        9       10 
3.500667 2.411196 2.627915 2.813815 2.468595 1.733852 2.114553 1.423225 1.470738 
> predict(model,newdata=data)
        1         2         3         4         5         6         7         8         9        10 
2.1628335 3.7063579 2.9781417 2.1374301 3.2582376 1.9534558 1.3670995 2.4547626 0.8448223 1.8762437 

Форма newdata

#Forecast fix.
A<-c(A,rnorm(1)) #Assume we already have 1-step forecasts for A,B,C
B<-c(B,rnorm(1))
C<-c(C,rnorm(1))

newdata<-ts(cbind(A,B,C),start=start(y),freq=frequency(y))

newdata<-cbind(lag(y,-1),newdata)
colnames(newdata) <- c("y","A","B","C")

Сравните прогноз с подгонкой модели:

> predict(model)
       2        3        4        5        6        7        8        9       10 
3.500667 2.411196 2.627915 2.813815 2.468595 1.733852 2.114553 1.423225 1.470738 
> predict(model,newdata=newdata)
       1        2        3        4        5        6        7        8        9       10       11 
      NA 3.500667 2.411196 2.627915 2.813815 2.468595 1.733852 2.114553 1.423225 1.470738 1.102367 

Как видно из исторических данных, прогноз совпадает, и последний элемент содержит прогноз на 1 шаг вперед.

mpiktas
источник
Как вы можете справиться со случаем, когда у вас есть два лага в одной формуле? lag(y,-1)+lag(y,-2)?
Хью Перкинс,
1
Ну, тогда это решение не работает. Вам нужно написать свою собственную функцию прогнозирования.
mpiktas
Ах, это то, что я сделал на самом деле :-P
Хью Перкинс
1
Вы рассматривали возможность представить его авторам dynlm? Это странная ситуация, которую нельзя предсказать, используя dynlm.
mpiktas
Хммм, вы говорите, что они не собираются волшебным образом отслеживать переполнение стека и исправлять ошибки? Я думаю, это, вероятно, правда!
Хью Перкинс,
2

Следуя запросу @ md-azimul-haque, я просмотрел свой 4-летний исходный код и нашел следующую функцию с соответствующим именем. Не уверены, что это то, что ищет @ md-azimul-haque?

# pass in training data, test data,
# it will step through one by one
# need to give dependent var name, so that it can make this into a timeseries
predictDyn <- function( model, train, test, dependentvarname ) {
    Ntrain <- nrow(train)
    Ntest <- nrow(test)
    # can't rbind ts's apparently, so convert to numeric first
    train[,dependentvarname] <- as.numeric(train[,dependentvarname])
    test[,dependentvarname] <- NA
    testtraindata <- rbind( train, test )
    testtraindata[,dependentvarname] <- ts( as.numeric( testtraindata[,dependentvarname] ) )
    for( i in 1:Ntest ) {
       cat("predicting i",i,"of",Ntest,"\n")
       result <- predict(model,newdata=testtraindata,subset=1:(Ntrain+i-1))
       testtraindata[Ntrain+i,dependentvarname] <- result[Ntrain + i + 1 - start(result)][1]
    }
    testtraindata <- testtraindata[(Ntrain+1):(Ntrain + Ntest),dependentvarname]
    names(testtraindata) <- 1:Ntest
    return( testtraindata )
}
Хью Перкинс
источник