迭代预测动态模型
我已经编写了一个函数来迭代预测使用 dyn 包构建的模型,并且我想要一些关于它的反馈。有更好的方法吗?有人为 dyn 类(或 dynlm 类)编写了规范的“预测”方法,还是我正在冒险进入未知领域?
ipredict <-function(model, newdata, interval = "none",
level = 0.95, na.action = na.pass, weights = 1) {
P<-predict(model,newdata=newdata,interval=interval,
level=level,na.action=na.action,weights=weights)
for (i in seq(1,dim(newdata)[1])) {
if (is.na(newdata[i])) {
if (interval=="none") {
P[i]<-predict(model,newdata=newdata,interval=interval,
level=level,na.action=na.action,weights=weights)[i]
newdata[i]<-P[i]
}
else{
P[i,]<-predict(model,newdata=newdata,interval=interval,
level=level,na.action=na.action,weights=weights)[i,]
newdata[i]<-P[i,1]
}
}
}
P_end<-end(P)[1]*frequency(P)+(end(P)[2]-1) #Convert (time,period) to decimal time
P<-window(P,end=P_end-1*frequency(P)) #Drop last observation, which is NA
return(P)
}
用法示例:
library(dyn)
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)
Error<-rnorm(10)
y<-y+.5*A+.2*B-.3*C+.1*Error #Add relationship to independant variables
data=cbind(y,A,B,C)
#Fit linear model
model.dyn<-dyn$lm(y~A+B+C+lag(y,-1),data=data)
summary(model.dyn)
#Forecast linear model
A<-c(A,rnorm(5))
B<-c(B,rnorm(5))
C<-c(C,rnorm(5))
y=window(y,end=end(y)+c(5,0),extend=TRUE)
newdata<-cbind(y,A,B,C)
P1<-ipredict(model.dyn,newdata)
P2<-ipredict(model.dyn,newdata,interval="prediction")
#Plot
plot(y)
lines(P1,col=2)
I've written a function to iteratively forecast models built using the package dyn, and I'd like some feedback on it. Is there a better way to do this? Has someone written canonical "forecast" methods for the dyn class (or dynlm class), or am I venturing into uncharted territory here?
ipredict <-function(model, newdata, interval = "none",
level = 0.95, na.action = na.pass, weights = 1) {
P<-predict(model,newdata=newdata,interval=interval,
level=level,na.action=na.action,weights=weights)
for (i in seq(1,dim(newdata)[1])) {
if (is.na(newdata[i])) {
if (interval=="none") {
P[i]<-predict(model,newdata=newdata,interval=interval,
level=level,na.action=na.action,weights=weights)[i]
newdata[i]<-P[i]
}
else{
P[i,]<-predict(model,newdata=newdata,interval=interval,
level=level,na.action=na.action,weights=weights)[i,]
newdata[i]<-P[i,1]
}
}
}
P_end<-end(P)[1]*frequency(P)+(end(P)[2]-1) #Convert (time,period) to decimal time
P<-window(P,end=P_end-1*frequency(P)) #Drop last observation, which is NA
return(P)
}
Example usage:
library(dyn)
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)
Error<-rnorm(10)
y<-y+.5*A+.2*B-.3*C+.1*Error #Add relationship to independant variables
data=cbind(y,A,B,C)
#Fit linear model
model.dyn<-dyn$lm(y~A+B+C+lag(y,-1),data=data)
summary(model.dyn)
#Forecast linear model
A<-c(A,rnorm(5))
B<-c(B,rnorm(5))
C<-c(C,rnorm(5))
y=window(y,end=end(y)+c(5,0),extend=TRUE)
newdata<-cbind(y,A,B,C)
P1<-ipredict(model.dyn,newdata)
P2<-ipredict(model.dyn,newdata,interval="prediction")
#Plot
plot(y)
lines(P1,col=2)
如果你对这篇内容有疑问,欢迎到本站社区发帖提问 参与讨论,获取更多帮助,或者扫码二维码加入 Web 技术交流群。
绑定邮箱获取回复消息
由于您还没有绑定你的真实邮箱,如果其他用户或者作者回复了您的评论,将不能在第一时间通知您!
发布评论
评论(1)
R 核心中的
predict.Arima
有n.ahead
参数来预测前进n
步,看来这就是您正在寻找的for 与 dyn 结合使用,但predict.dyn
目前不支持该功能。要获得这种效果,必须像您所做的那样迭代调用dyn$whatever
。predict.Arima
in the core of R has then.ahead
argument to forecastn
steps ahead and it seems that that is what you are looking for in conjunction with dyn butpredict.dyn
does not currently support that functionality. To get that effect one must iteratively calldyn$whatever
as you are doing.