R多元预测和准确性领先一步
R multivariate one step ahead forecasts and accuracy
我想使用 R 比较两个预测模型的 RMSE(均方根误差)。第一个模型使用 1966 年到 2000 年的估计值来预测 2001 年,然后使用 1966 年到 2001 年的估计值来预测 2002 年,以此类推直到 2015 年。第二个模型使用 1991 年到 2000 年的估计值来预测 2001 年,然后使用 1992 年到 2001 年的估计值预测 2002 年,依此类推,直到 2015 年。这个问题让我很困惑,非常感谢任何帮助。
DF <- data.frame(YEAR=1966:2015, TEMP=rnorm(50), PRESSURE=rnorm(50), RAINFALL=rnorm(50))
lmod <- lm(TEMP ~ PRESSURE + RAINFALL, data = DF)
rmse <- function(error) sqrt(mean(error^2))
rmse(lmod$residuals)
你可以循环它:
方法一:
pred1<-numeric(0)
rmse1<-numeric(0)
for(i in 1:15){
DF.train1<-DF[DF$YEAR < 2000+i,]
DF.test1<-DF[DF$YEAR == 2000+i,]
lmod1 <- lm(TEMP ~ PRESSURE + RAINFALL, data = DF.train1)
pred1[i]<- predict(lmod1, newdata = DF.test1)
rmse1[i]<-sqrt(mean((DF.test1$TEMP-pred1[i])^2))
}
pred1
rmse1
mean(rmse1)
方法二:
pred2<-numeric(0)
rmse2<-numeric(0)
for(i in 1:15){
DF.train2<-DF[DF$YEAR < 2000+i & DF$YEAR > 1989+i,]
DF.test2<-DF[DF$YEAR == 2000+i,]
lmod2 <- lm(TEMP ~ PRESSURE + RAINFALL, data = DF.train2)
pred2[i]<- predict(lmod2, newdata = DF.test2)
rmse2[i]<-sqrt(mean((DF.test2$TEMP-pred2[i])^2))
}
pred2
rmse2
mean(rmse2)
比较rmse1
和rmse2
的各个组成部分,以及它们各自的手段应该是有用的。向量 pred1
和 pred2
包含各自方法对每年 (2001-2015) 的单独 TEMP
预测。
编辑:现在应该可以使用了,方法 2 的训练时间跨度为 10 年。此外,我将 RMSE 视为 this 文章中为预测变量定义的 MSE 的平方根。
这是另一个解决方案,其中模拟在一个函数中。
此解决方案的好处是可以轻松修改模型规格。
例如,如果您想尝试使用范围为 15 年而不是 10 年的 model2
,只需修改函数中的输入 (range = 15
)。这也使您可以进行光敏性分析。
compare_models <- function(DF, start = 1966, end = 2000, range = 10)
{
require(hydroGOF)
for (i in (end+1):tail(DF$YEAR)[6])
{
# model1
lmod_1 = lm(TEMP ~ PRESSURE + RAINFALL, data = DF[DF$YEAR >= start & DF$YEAR < i,])
DF$model1_sim[DF$YEAR == i] <- predict(lmod_1, newdata = DF[DF$YEAR == i,])
# model2
lmod_2 = lm(TEMP ~ PRESSURE + RAINFALL, data = DF[DF$YEAR >= i-range & DF$YEAR < i,])
DF$model2_sim[DF$YEAR == i] <- predict(lmod_2, newdata = DF[DF$YEAR == i,])
}
return(DF)
}
我使用 hydroGOF
包来计算 rmse
和 NSE
,这是模型效率的常用指标(参见 Nash and Sutcliffe, 1970,目前有 11528 次引用)。
output = compare_models(DF)
require(hydroGOF) # compute RMSE and NSE
# RMSE
rmse(output$model1_sim,output$TEMP)
rmse(output$model2_sim,output$TEMP)
# Nash-Sutcliffe efficiency
NSE(output$model1_sim,output$TEMP, na.rm = T)
NSE(output$model2_sim,output$TEMP, na.rm = T)
还有一个简单的 simulated/observed 图来寻找模型预测:
# melting data for plot
output_melt = melt(output[,c("TEMP", "model1_sim", "model2_sim")], id = "TEMP")
# Plot
ggplot(output_melt, aes(x = TEMP, y = value, color = variable)) +
theme_bw() + geom_point() + geom_abline(slope = 1, intercept = 0) +
xlim(-2,2) + ylim(-2,2) + xlab("Measured") + ylab("Simulated")
这是另一个解决方案:
year <- 2000
time.frame <- 35
train.models <- function(year, time.frame) {
predictions <- sapply(year:(max(df$YEAR)-1),
function(year) {
lmod <- lm(TEMP ~ PRESSURE + RAINFALL, DF,
subset = with(DF, YEAR %in% (year - time.frame + 1):year))
pred <- predict(lmod, newdata = DF[DF$YEAR == (year + 1),])
names(pred) <- year + 1
return (pred)
})
return (predictions)
}
models1 <- train.models(2000, 35)
models2 <- train.models(2001, 10)
rmse(models1 - DF$TEMP[DF$YEAR %in% names(models1)])
rmse(models2 - DF$TEMP[DF$YEAR %in% names(models2)])
我想使用 R 比较两个预测模型的 RMSE(均方根误差)。第一个模型使用 1966 年到 2000 年的估计值来预测 2001 年,然后使用 1966 年到 2001 年的估计值来预测 2002 年,以此类推直到 2015 年。第二个模型使用 1991 年到 2000 年的估计值来预测 2001 年,然后使用 1992 年到 2001 年的估计值预测 2002 年,依此类推,直到 2015 年。这个问题让我很困惑,非常感谢任何帮助。
DF <- data.frame(YEAR=1966:2015, TEMP=rnorm(50), PRESSURE=rnorm(50), RAINFALL=rnorm(50))
lmod <- lm(TEMP ~ PRESSURE + RAINFALL, data = DF)
rmse <- function(error) sqrt(mean(error^2))
rmse(lmod$residuals)
你可以循环它:
方法一:
pred1<-numeric(0)
rmse1<-numeric(0)
for(i in 1:15){
DF.train1<-DF[DF$YEAR < 2000+i,]
DF.test1<-DF[DF$YEAR == 2000+i,]
lmod1 <- lm(TEMP ~ PRESSURE + RAINFALL, data = DF.train1)
pred1[i]<- predict(lmod1, newdata = DF.test1)
rmse1[i]<-sqrt(mean((DF.test1$TEMP-pred1[i])^2))
}
pred1
rmse1
mean(rmse1)
方法二:
pred2<-numeric(0)
rmse2<-numeric(0)
for(i in 1:15){
DF.train2<-DF[DF$YEAR < 2000+i & DF$YEAR > 1989+i,]
DF.test2<-DF[DF$YEAR == 2000+i,]
lmod2 <- lm(TEMP ~ PRESSURE + RAINFALL, data = DF.train2)
pred2[i]<- predict(lmod2, newdata = DF.test2)
rmse2[i]<-sqrt(mean((DF.test2$TEMP-pred2[i])^2))
}
pred2
rmse2
mean(rmse2)
比较rmse1
和rmse2
的各个组成部分,以及它们各自的手段应该是有用的。向量 pred1
和 pred2
包含各自方法对每年 (2001-2015) 的单独 TEMP
预测。
编辑:现在应该可以使用了,方法 2 的训练时间跨度为 10 年。此外,我将 RMSE 视为 this 文章中为预测变量定义的 MSE 的平方根。
这是另一个解决方案,其中模拟在一个函数中。
此解决方案的好处是可以轻松修改模型规格。
例如,如果您想尝试使用范围为 15 年而不是 10 年的 model2
,只需修改函数中的输入 (range = 15
)。这也使您可以进行光敏性分析。
compare_models <- function(DF, start = 1966, end = 2000, range = 10)
{
require(hydroGOF)
for (i in (end+1):tail(DF$YEAR)[6])
{
# model1
lmod_1 = lm(TEMP ~ PRESSURE + RAINFALL, data = DF[DF$YEAR >= start & DF$YEAR < i,])
DF$model1_sim[DF$YEAR == i] <- predict(lmod_1, newdata = DF[DF$YEAR == i,])
# model2
lmod_2 = lm(TEMP ~ PRESSURE + RAINFALL, data = DF[DF$YEAR >= i-range & DF$YEAR < i,])
DF$model2_sim[DF$YEAR == i] <- predict(lmod_2, newdata = DF[DF$YEAR == i,])
}
return(DF)
}
我使用 hydroGOF
包来计算 rmse
和 NSE
,这是模型效率的常用指标(参见 Nash and Sutcliffe, 1970,目前有 11528 次引用)。
output = compare_models(DF)
require(hydroGOF) # compute RMSE and NSE
# RMSE
rmse(output$model1_sim,output$TEMP)
rmse(output$model2_sim,output$TEMP)
# Nash-Sutcliffe efficiency
NSE(output$model1_sim,output$TEMP, na.rm = T)
NSE(output$model2_sim,output$TEMP, na.rm = T)
还有一个简单的 simulated/observed 图来寻找模型预测:
# melting data for plot
output_melt = melt(output[,c("TEMP", "model1_sim", "model2_sim")], id = "TEMP")
# Plot
ggplot(output_melt, aes(x = TEMP, y = value, color = variable)) +
theme_bw() + geom_point() + geom_abline(slope = 1, intercept = 0) +
xlim(-2,2) + ylim(-2,2) + xlab("Measured") + ylab("Simulated")
这是另一个解决方案:
year <- 2000
time.frame <- 35
train.models <- function(year, time.frame) {
predictions <- sapply(year:(max(df$YEAR)-1),
function(year) {
lmod <- lm(TEMP ~ PRESSURE + RAINFALL, DF,
subset = with(DF, YEAR %in% (year - time.frame + 1):year))
pred <- predict(lmod, newdata = DF[DF$YEAR == (year + 1),])
names(pred) <- year + 1
return (pred)
})
return (predictions)
}
models1 <- train.models(2000, 35)
models2 <- train.models(2001, 10)
rmse(models1 - DF$TEMP[DF$YEAR %in% names(models1)])
rmse(models2 - DF$TEMP[DF$YEAR %in% names(models2)])