☰
R语言实现LSTM时间序列预测:从数据预处理到多步预测实战
2026/10/3 8:05:58 网站建设 项目流程

上周有个做水文预报的朋友问我:R里到底能不能跑长短期记忆模型?他说自己在网上搜了一圈,中文教程清一色是Python,偶尔看到几个R的版本也是零散代码,根本跑不通。我说你把数据发我,我现场搭一个给你看。结果半小时后,他问我为什么以前没人告诉他R也能这么干净利落地做LSTM。

这个概念确实被低估了。R的深度学习生态不像Python那样高调,但keras包和torch包都维护得很好,尤其是配合tidyverse做时间序列数据预处理,写起来比Python的pandas还要顺手。这篇内容面向想用R做LSTM但不知道从哪下手的读者,无论你手中是水文径流数据、气象观测数据、电网负荷数据还是销售序列,思路完全一致:把时间序列改造成LSTM能吃的三维数组,然后建模、训练、评估、滚动预测。我会把每一步的原理和坑都讲透,而不是丢一段能跑的代码就完事。

1. 为什么“用R跑LSTM”是个被低估的组合

1.1 R不是不能做深度学习,只是入口和Python不一样

很多人对R做深度学习的印象停留在“装个包试试,跑不了就换Python”,这个印象其实来自两个误区。第一,R的早期深度学习包确实不成熟,keras之前有个叫deepnet的包,功能非常原始,支持的网络层类型少、训练困难,用过的人都会留下心理阴影。第二,当时R和Python的互操作不如现在顺畅,跑一次模型要来回导数据,体验很差。

但现在的格局完全变了。RStudio团队维护的keras包是Python Keras的完整封装,底层由reticulate调用TensorFlow执行。换句话说,你写的model <- keras_model_sequential()本质上就是在操控一个Python环境里的Keras后端,两者的数学实现、GPU调用方式、权重格式完全一致。真正的区别只在API的语法习惯上:Python里写LSTM(...),R里写layer_lstm(...),函数名带了个layer_前缀,仅此而已。

另一个被忽略的点是,R在数据预处理阶段往往比Python更顺手。时间序列建模第一步是清洗、对齐、构造滞后特征,这个过程R有dplyr、tidyr、slider这些专用工具,管道操作写起来非常紧凑。比如你要给一个日径流序列同时构造12天的滞后降水、滞后气温和滞后径流,Python可能要写循环或者反复merge,R里几行mutate就能完成。这也是我在实际项目中倾向把整条链路都留在R里的原因——不需要在R和Python之间来回切换。

1.2 keras和torch,R语言里两条不同的路

如果你愿意深入R的深度学习生态,会发现其实有两条路线可选。第一条是keras包,高层API,封装程度高,适合快速建模调参,也是这篇博文的主线。第二条是torch包,它是PyTorch在R中的移植版,使用方式和Python的PyTorch几乎一一对应,灵活性更强,代价是你要自己写训练循环、自己管梯度,对新手来说门槛高不少。

我的建议很直接:如果你要做的是表格型数据、时间序列单步或多步预测,keras足够了,它的fit()函数把训练循环、批次打乱、梯度更新全部封装好,出结果快,出问题也好定位。torch则更适合你在跑一些论文里的特殊网络结构、需要自定义损失函数里做复杂控制流的场景。两者的底层性能差别不大,因为torch和TensorFlow在CPU/GPU上都不会成为瓶颈,真正的瓶颈往往在数据读取和预处理上。

选择keras还有一个现实原因:文档和案例数量更多。RStudio官网的keras文档更新很及时,Stack Overflow上针对R keras的问答量也更大。对中文用户来说,Python Keras的资料同样可以作为参考,因为函数逻辑基本一一对应,只是换个名字。这是torch包暂时比不了的。

1.3 环境安装的完整路径(含GPU确认)

R里安装keras的流程比大多数人想象中简单,但有几个细节需要注意。基础部分两行搞定:

install.packages("keras") library(keras)

装完并不是马上就能用,还要安装TensorFlow后端。官方推荐直接用install_keras(),这个函数会为R自动创建一个独立的conda环境,在环境里装入Python和TensorFlow,好处是它不会污染系统里已有的Python环境。

install_keras()

这一步下载的东西比较多,耗时会根据网络情况从几分钟到十几分钟不等。装完之后我建议先确认环境是否完整:

library(keras) tensorflow::tf_config()

如果输出里显示了TensorFlow版本和设备列表,说明安装成功。看到“Device: CPU”就是CPU模式,看到“Device: GPU: 0”说明TensorFlow已经识别到了你的显卡。需要注意一点,GPU版的TensorFlow需要你提前装好CUDA和cuDNN,版本匹配问题比较麻烦。如果你只是想跑通流程或者数据量不大,CPU完全够用,不必一上来就折腾GPU。

安装阶段最常见的报错是reticulate找不到Python环境,通常发生在你用了多个conda环境、又手动指定过Python路径的情况下。处理方法就是在R里重新指向正确的Python:

reticulate::use_condaenv("r-reticulate", required = TRUE)

这个命令会告诉reticulate使用名为r-reticulate的环境,也就是install_keras()创建的那个。跑完这条再执行tensorflow::tf_config(),一般就能恢复正常。

2. 从原始径流序列到“能塞进LSTM”的三维数组

2.1 时间序列变成监督学习样本:滞后窗口

LSTM不是那种把一整条序列扔进去就能直接预测下一个值的黑盒子。它学习的最小单元是一个“样本”,这个样本由两个部分组成:过去lookback个时间步的特征序列,以及对应的目标值。这一步在深度学习里叫做把无监督的时间序列改造成监督学习的表格形式,是整个过程最容易出错也最关键的环节。

以水文径流预报为例,假设你的原始数据是每天一行,包含三列:日降水量、日均气温、日径流量。你想用过去12天的数据预测第13天的径流量,那么一个样本就是过去12天每天的三维特征组成的矩阵,形状是(12, 3),标签是第13天的径流量这一个数值。整个样本集就可以看成一个大的三维数组:样本数量 × 时间步数 × 特征数量。

在R里构造滞后窗口,我以前用循环写,效率很低,后来改用purrr和dplyr的组合,干净又直观:

library(tidyverse) # 模拟一份日径流数据,真实项目替换为自己的数据即可 set.seed(2024) n <- 1500 daily_data <- tibble( date = seq.Date(as.Date("2018-01-01"), by = "day", length.out = n), rainfall = rgamma(n, shape = 1.2, rate = 0.8), temp = sin(seq(0, 4 * pi, length.out = n)) * 8 + 15 + rnorm(n, 0, 1), discharge = 40 + 12 * rainfall + rnorm(n, 0, 3) ) lookback <- 12 lagged_data <- map_dfc(1:lookback, function(k) { daily_data %>% arrange(date) %>% transmute( !!paste0("rainfall_lag", k) := lag(rainfall, k), !!paste0("temp_lag", k) := lag(temp, k), !!paste0("discharge_lag", k) := lag(discharge, k) ) }) lagged_data <- bind_cols(daily_data, lagged_data) %>% filter(!is.na(discharge_lag12))

这里有一个非常容易踩的坑:构造标签的时候,目标必须是“未来第13天”的径流,而不是行内已有的discharge。具体来说,当某一行代表的是第t天的滞后特征组合时,你对应的y应该是第t+1天的实际观测值。用dplyr实现就一句话:

model_df <- lagged_data %>% mutate(target = lead(discharge, 1)) %>% filter(!is.na(target))

如果你的purpose是直接多步预测,要预测未来7天的径流,那就构造7个标签列,等于把模型输出维度从1改成7。这个在后面第5节展开讲。

2.2 归一化必须在切分数据之后再做

时间序列深度学习模型训练前几乎必须做归一化。LSTM内部使用sigmoid和tanh作为激活函数,它们的输出范围分别是(0,1)和(-1,1),如果输入特征数值跨度太大,比如气温在零下20度到零上40度之间波动,径流可能在10到500之间跳跃,不归一化会让梯度更新非常不稳定,模型很难收敛。

但归一化有一个特别容易被忽略的顺序问题:必须在切分训练集和测试集之后,用“训练集”的均值和标准差去变换所有数据。如果你把整条序列放在一起计算均值,再分别切分训练集和测试集,等于让测试集的信息提前泄露到了训练过程中。这会让模型评估结果偏乐观,但在真实部署时,未来数据是不可见的,模型表现会大打折扣。

R里用caret包或者base R的scale()都能做。我习惯手动提取统计量,方便后面反归一化:

feature_cols <- grep("_lag", names(model_df), value = TRUE) train_size <- floor(nrow(model_df) * 0.8) train_df <- model_df[1:train_size, ] test_df <- model_df[(train_size + 1):nrow(model_df), ] train_mean <- apply(train_df[, feature_cols], 2, mean) train_sd <- apply(train_df[, feature_cols], 2, sd) train_scaled <- as.data.frame(scale(train_df[, feature_cols], center = train_mean, scale = train_sd)) test_scaled <- as.data.frame(scale(test_df[, feature_cols], center = train_mean, scale = train_sd)) y_train <- train_df$target y_test <- test_df$target

这段代码里的train_mean和train_sd保存下来很重要,预测阶段做反归一化要原样使用。有很多人训练时很开心,等到预测完成把结果画出来发现量级不对,才意识到忘了保存归一化参数,只能重新跑一遍流程,纯属浪费时间。

2.3 三维数组的两个容易翻车的维度坑

数据表准备好了,接下来要把它从二维表格改造成LSTM需要的三维数组,形状应该是(样本数, 滞后步数, 特征数)。这个转换看起来简单,实际动手时有两个坑非常常见。

第一个坑是特征顺序和滞后顺序的混淆。假设你有3个特征、12个滞后步,在宽表格里列的排列方式可能是先放完所有特征的lag1,再放所有特征的lag2,也可能先放rainfall的全部12个滞后,再放temp的全部12个滞后。无论哪种排列,转换成三维数组时都必须明确地按“时间步”重组,否则模型读到的“一个时间步”的数据混杂了不同时间点的信息,预测结果会莫名其妙。

我用的方法是用嵌套列表逐个构造矩阵,看起来代码多一点,但逻辑绝对不出错:

make_3d_array <- function(scaled_df, n_features, lookback) { feature_names <- colnames(scaled_df) var_names <- unique(sub("_lag\\d+$", "", feature_names)) samples <- vector("list", nrow(scaled_df)) for (i in seq_len(nrow(scaled_df))) { mat <- matrix(0, nrow = lookback, ncol = n_features) for (j in seq_len(n_features)) { col_pattern <- paste0(var_names[j], "_lag") mat[, j] <- as.numeric(scaled_df[i, feature_names[grep(col_pattern, feature_names)]]) } samples[[i]] <- mat } array(unlist(samples), dim = c(length(samples), lookback, n_features)) } x_train <- make_3d_array(train_scaled, n_features = 3, lookback = lookback) x_test <- make_3d_array(test_scaled, n_features = 3, lookback = lookback)

第二个坑是维度顺序本身。R的数组填充方式和Python的NumPy略有不同,如果你用array(..., dim = c(n, lookback, n_features))来装数据,默认是按列填充的,这会导致时间步的位置被错误地展开了。上面这个逐样本构造再统一unlist的方法,最终维度一定正确,因为每个样本都已经是(时间步, 特征)的矩阵。建议构造完打印一下维度确认:

dim(x_train) # 应该是 c(样本数, 12, 3)

只要这一步确认无误,后面的模型结构就有了稳的基础。

3. 搭一个可用的LSTM模型:R中的建模与训练参数

3.1 模型结构怎么选:几层LSTM、多少units

在R里用keras搭LSTM模型,语法上非常简洁。一个最常用的回归结构是两层LSTM加一个全连接输出层:

library(keras) model <- keras_model_sequential() %>% layer_lstm(units = 64, return_sequences = TRUE, input_shape = c(lookback, 3)) %>% layer_lstm(units = 32) %>% layer_dense(units = 1) model

为什么选两层LSTM?一层LSTM只能学到单一时间尺度上的依赖关系,两层结构让上层网络有机会在底层提取的特征之上再抽象一层更高阶的模式。比如底层可能学到“连续三天强降水后径流开始爬升”这样的子模式,上层再去学习这些子模式之间的组合关系。如果你的数据更复杂,可以加到三层,但层数增加意味着参数量暴增,训练时间变长,过拟合风险也更大。

units的选择没有绝对标准,经验上看64起步是性价比比较高的选择。units太少,LSTM的记忆容量不够,复杂时间依赖学不出来;units太多,模型容量过剩,训练集上拟合得很好,测试集上一塌糊涂。我遇到过的绝大多数水文和气象序列预测项目,两层LSTM加64和32个units的组合都够用,不需要一上来就堆128甚至256。

return_sequences = TRUE这个参数是新手最容易搞混的地方。它在第一个LSTM层里设置为TRUE,表示把每个时间步的隐藏状态都输出给下一层;如果设成FALSE,则只输出最后一个时间步的隐藏状态。如果你打算堆叠两层LSTM,前面层的return_sequences必须为TRUE,否则第二层LSTM拿到的输入形状不对。最后一层LSTM通常设为FALSE,因为它后面只连接一个全连接层,不需要保留完整的时间步输出。

3.2 compile和fit里的关键参数解析

模型结构定义好之后,下一步是编译和训练。compile阶段要指定三件事:优化器、损失函数、评估指标。

model %>% compile( optimizer = optimizer_adam(learning_rate = 0.001), loss = "mse", metrics = c("mae") )

优化器我几乎固定用Adam,它对学习率不那么敏感,自带动量机制,在多数时间序列任务上收敛速度比传统SGD快很多。learning_rate设置为0.001是最常用的起点,如果发现loss下降过慢,可以考虑调到0.003;如果发现loss震荡不收敛,就调小到0.0003。

损失函数用mse(均方误差),适用于绝大多数回归问题。如果你的数据里有明显的异常峰值,改用mae(平均绝对误差)会更稳健一些。metrics参数里加个mae只是为了方便看训练过程的直观反馈,不影响模型本身的优化目标。

真正需要仔细处理的是fit()函数里的参数:

history <- model %>% fit( x = x_train, y = y_train, validation_data = list(x_test, y_test), epochs = 100, batch_size = 32, callbacks = list( callback_early_stopping(patience = 10, restore_best_weights = TRUE), callback_reduce_lr_on_plateau(factor = 0.5, patience = 5, min_lr = 0.00001) ) )

batch_size控制每次梯度更新使用的样本数。32是最常见的默认值,如果你的数据量很小,比如只有几千个样本,用16会更稳定。数据量很大、内存吃紧的时候,可以增加到64或128,但要注意batch_size太大会让模型趋向于收敛到尖锐的局部最小值,泛化能力反而下降。

epochs我通常会设一个偏大的数字,比如100或200,但不会真的让它跑完。因为配合了early stopping,模型在验证集上连续多轮不改善就会自动停止,并把权重恢复到验证集表现最好的那个时刻。这也是训练深度学习模型最重要的一个习惯:不要手动盯loss曲线盯着盯到过拟合才反应过来要停止。

3.3 用回调函数避免训练过程失控

回调函数是keras里仅次于fit()本身的实用功能。上面代码里已经用到了两个最关键的:early stopping和reduce learning rate on plateau。

callback_early_stopping的意义在于自动判断“什么时候停”。它的patience参数表示连续多少个epoch验证loss没有改善就停止训练。10是一个比较安全的值,太小容易在验证集短暂波动时误停,太大又增加了无谓的训练时间。restore_best_weights = TRUE保证训练结束后模型参数是验证集上表现最好的那一份,而不是最后一次迭代的参数。

callback_reduce_lr_on_plateau解决的是“训练卡住”的问题。深度模型在训练过程中经常出现loss降到某个平台后怎么都降不下去,这时候把学习率减半,往往能突破瓶颈继续下降。patience设为5,意思是验证集连续5轮没有提升就自动把学习率缩为原来的0.5倍,最低降到0.00001,防止学习率被减到几乎为0导致完全无法学习。

这两个回调配合起来,大部分训练问题都能被自动化解。我见过很多手写训练循环的人,每跑一个模型就要守着看loss曲线,用keras根本不需要这样,把回调设好,让它自己跑完就行。等训练结束,用plot(history)看一眼loss曲线,确认训练集和验证集的下降趋势基本一致,就可以进入下一步评估了。

4. 训练过程中的翻车现场与排查顺序

4.1 验证集划分的隐藏坑:keras默认shuffle会打乱时序

这是我踩过最深的一个坑,也是我在帮别人排查代码时遇到频率最高的问题。keras的fit()函数默认在每个epoch开始时对训练集执行shuffle,这个设计本来是为了打破样本间的顺序相关性,让梯度更新更稳定。但对时间序列来说,这恰恰可能是有害的。

问题出在验证集和训练集的划分逻辑上。如果你用validation_split = 0.2,keras内部会把样本集中最后20%的数据划出来当验证集,然后对前面80%的训练数据做shuffle。这看起来没问题,但LSTM的训练样本本身是重叠窗口构造出来的:第t个样本包含第t天到第t+12天的数据,第t+1个样本包含第t+1天到第t+13天的数据,两个样本有11天是完全重叠的。对这样的重叠样本做随机打乱,等于让模型在一个epoch里用高度重复的信息做跳跃式更新,训练曲线会变得很毛糙,验证集表现也不稳定。

更隐蔽的风险是,如果你的数据构造过程中不小心混入了未来信息,shuffle会把这个问题完全掩盖掉。我自己在早期项目里就干过这种事:构造滞后特征时没删干净空值,导致某个样本的输入里包含了目标值本身的信息。由于shuffle破坏了时间顺序,模型居然在验证集上表现得“非常完美”,直到实际部署才发现效果一塌糊涂。

稳妥的做法是不依赖validation_split,手动按时间顺序切分好验证集,然后传进validation_data参数。这样验证集永远是时间上靠后的那段数据,也和真实预测场景一致。我在前面的代码里已经这样做了,这是时间序列LSTM训练必须遵守的纪律。

4.2 Loss变NaN、loss不下降、验证loss飙升

训练过程中最常见的三个异常,按出现频率排分别是loss变成NaN、loss完全不下降、验证loss先降后飙升。

Loss变成NaN基本是数值稳定性的问题,排查顺序很固定。第一步检查输入数据里有没有NA或Inf,这个在时间序列里特别容易发生——滞后特征构造完,前12行全是NA,如果你忘了filter掉,NaN会顺着计算图一直蔓延到loss。第二步检查归一化参数,如果某个特征的方差为0,比如某列全是同一个值,scale()之后会出现除以0,产生Inf,后患无穷。第三步才是检查学习率,学习率过大会让梯度爆炸,数值溢出变成NaN。前两步在数据层,最后一步在模型层,按这个顺序查,大概率十几分钟内定位。

Loss完全不下降的情况,先看基线是不是就没做对。如果你的y值范围很大,比如径流在几百的量级,但模型最后输出层的activation是默认的linear,初始预测值可能是随机的,mse损失巨大。我遇到过一次怎么也降不下去,后来打印了pred和y的分布,发现数据没归一化,网络输入输出尺度差了三个数量级,梯度根本传不回去。把y也做一次归一化,或者至少在compile里用mae看趋势,问题立刻就明朗了。

验证loss先降后飙升是过拟合的典型信号,说明模型开始“背”训练数据了。处理手段优先级从低到高排列:先加early stopping防止它继续恶化,然后减小units或者加dropout层,最后才是减少epochs或者增大数据量。Dropout是LSTM过拟合最有效的正则化手段之一,R里用法很简单:

model <- keras_model_sequential() %>% layer_lstm(units = 64, return_sequences = TRUE, input_shape = c(lookback, 3), dropout = 0.2) %>% layer_lstm(units = 32, dropout = 0.2) %>% layer_dense(units = 1)

注意dropout参数可以直接写在layer_lstm()里,写法比Python版本还简洁。一般来说dropout从0.2开始试,不要一上来就设0.5,太高会让模型欠拟合。

4.3 预测结果全是均值/全一样怎么办

模型训练完了,predict()也跑出来了,拿过预测结果一看,所有预测值都挤在一个极小的范围内,几乎是同一个常数。这个问题在新手项目里出现频率特别高,而且让人非常慌张。

本质原因通常是模型退化成了一种“保守策略”:当输入信号对输出的影响没有学透彻时,模型发现预测值接近训练集的均值,损失也不会增加太多,于是它选择了一个在数学上“安全”的常数输出。这背后对应的实际原因有三个。

第一个是特征信息不够,模型从输入里看不到和target有明显相关性的模式。水文预报就是典型例子:如果只看气温预测明天的径流,模型大概只能学到“输出均值”,因为径流和气温的关系本来就弱。第二是lookback窗口太短,LSTM看不到足够长的历史来捕捉事件的持续影响。径流对降水的响应可能滞后好几天,你只给它看过去3天,信息当然不够。第三是模型容量太小,units太少导致它没能力拟合复杂函数。

排查手段也很直接:把输入的滞后特征和target之间的相关性打出来,看看哪些特征真正有预测力;增加lookback、增加units;再不行就把batch_size调小一点。另外一个容易被忽略的检查点是,预测前一定要确认你给predict()的输入做了和训练阶段完全一样的归一化,如果特征尺度和训练时不一致,模型输出会理所当然地变成一个奇怪的常数。

5. 评估和多步预测:从“预测明天”到“预测未来七天”

5.1 除了RMSE,水文预测还要看NSE

评估一个LSTM模型,不能只看训练集loss,要看它在中长期预测中的实际价值。很多从传统统计模型转过来的朋友习惯用RMSE和MAE,这两个指标当然没问题,但在水文等领域有一个指标更受认可:纳什效率系数NSE。

NSE的计算公式是:

NSE = 1 - sum((obs - sim)^2) / sum((obs - mean(obs))^2)

它的含义可以理解成“你的模型预测误差相比直接用历史均值当预测,好了多少”。如果NSE等于1,说明预测完全准确;等于0,说明和直接用平均值没区别;小于0,说明你的模型还不如摆烂预测均值。

在R里算NSE非常简单:

nse <- 1 - sum((y_test - pred)^2) / sum((y_test - mean(y_test))^2)

实际项目里,径流预报模型NSE超过0.8通常被认为表现优秀,0.6到0.8之间是可以接受的水平,低于0.6说明模型还有不小改进空间。我建议同时报告RMSE和NSE,因为它们回答的问题不同:RMSE告诉你平均误差是多少个立方米每秒,NSE告诉你这个误差相对基线好多少。两个都看,你对模型的判断才不会偏。以下是我常用的一套评估汇总:

指标计算方式适用场景好坏参考
RMSEsqrt(mean((obs - sim)^2))常规回归误差与观测均值对比才有意义
MAEmean(abs(obs - sim))对异常值不敏感场景越小越好
NSE1 - SS_res / SS_tot水文模型>0.8优秀,0.6-0.8可用
KGE1 - sqrt((r-1)^2 + (a-1)^2 + (b-1)^2)综合评估越接近1越好

KGE是近年在水文领域越来越受推荐的指标,它同时考虑相关性、偏差和变异性三个方面,比NSE更全面。R里计算KGE可以自己写函数,也可以用hydroGOF包,一行代码搞定。

5.2 递归预测与直接多步预测的取舍

如果你只有一日预报的需求,上面训练好的单输出模型就已经能用了。但流域防汛调度、水库发电计划这些真实场景通常要求未来7天甚至15天的径流预报。这时候你必须面对多步预测的问题。

多步预测最朴素的做法是递归预测:模型一次只预测明天,然后把“明天的预测值”当作已知的输入,去预测后天,再预测大后天,这样一步步滚下去。实现起来很简单:

recursive_forecast <- function(model, initial_input, steps = 7) { preds <- numeric(steps) current_input <- initial_input for (i in seq_len(steps)) { pred <- model %>% predict(array(current_input, dim = c(1, lookback, 3))) preds[i] <- pred[1] # 把当前预测值滚动进输入序列,同时推进所有特征 new_row <- current_input new_row[1:(lookback - 1), ] <- current_input[2:lookback, ] new_row[lookback, ] <- pred[1] # 实际项目里需要把特征也同步推进 current_input <- new_row } preds }

递归预测的致命弱点是误差累积。第1天的预测误差会作为输入影响第2天的预测,第2天的误差又影响第3天,整个误差随预测步长滚雪球。尤其是径流这种对初始条件敏感的过程,预测到第5天之后精度会明显恶化。

另一个思路是直接多步预测:把模型输出维度从1改成H,让模型一次性输出未来H天的值。R里的实现只需要改两处:

horizon <- 7 # 构造多步标签 model_df <- model_df %>% mutate(across(1:horizon, ~ lead(discharge, .x), .names = "target_{.x}")) %>% filter(!is.na(target_7)) # 模型输出维度改为 horizon model <- keras_model_sequential() %>% layer_lstm(units = 64, return_sequences = TRUE, input_shape = c(lookback, 3)) %>% layer_lstm(units = 32) %>% layer_dense(units = horizon)

直接多步避免了误差累积,但代价是模型要学习“从同一个输入序列同时输出7个不同时间尺度的预测”,这个学习任务难度更高。实际项目中两种方案可以都跑一遍对比验证集上的NSE,没有哪一种是永远更好的。我在防汛预报项目里的经验是:预测未来3天用递归即可,预测未来7天以上时直接多步往往更稳。

5.3 反归一化与结果的实用化输出

预测完成后必须把结果从归一化的尺度还原回真实的物理尺度,否则画出来的图没有任何意义。反归一化其实就是归一化的逆运算,利用第2节保存的train_mean和train_sd:

# 假设你对target也做了归一化 target_mean <- mean(train_df$target) target_sd <- sd(train_df$target) pred_ori <- pred_rnn * target_sd + target_mean obs_ori <- y_test * target_sd + target_mean

如果你没有对target归一化,只是在训练前减均值,那反归一化也只需要加回均值。关键在于保存的参数要和你实际做的变换一致,别搞混。

反归一化之后,我习惯把预测值和观测值放进一个tibble里,方便后续画图和计算指标:

result_df <- tibble( date = test_df$date[1:length(pred_ori)], observed = obs_ori, predicted = pred_ori ) library(ggplot2) ggplot(result_df, aes(x = date)) + geom_line(aes(y = observed, color = "观测值")) + geom_line(aes(y = predicted, color = "预测值")) + labs(y = "径流量 (m3/s)", color = NULL) + theme_minimal()

画图不是走形式。时间序列预测的误差很容易集中出现在某个特定的时间模式里,比如汛期的洪峰、融雪时段的上涨段。只看平均指标看不出这些问题,把预测和观测曲线叠在一起,哪里高估、哪里低估、洪峰相位有没有滞后,一眼就能看出来。这也是我在每个项目里都会保留的最后一个步骤。

另外,在实际生产环境里,模型预测的输入数据需要实时更新。你不能每次预测都重新构造整个lagged_data,应该维护一个实时更新的滚动窗口:新的一天数据到达后,把旧的滞后特征平移一格,塞入新数据,再调用predict()。这个逻辑和递归预测里的滚动输入是一样的,建议提前封装成函数,不要每次都在生产脚本里复制粘贴。

我在实际项目中更习惯把整个流程封装成一个“训练函数加预测函数”的组合,因为水文预报往往有几十个站点,每个站点都要独立建模,手动逐个跑一遍完全不现实。封装好之后只需要循环每个站点传入数据,自动完成滞后构造、归一化、建模、评估,最后输出一个汇总表。这个经验同样适合气象、电网负荷等批量时间序列预测场景。继续往下扩展,LSTM模型的下一步可以考虑加入时间特征编码,比如把月份和日序数用正弦/余弦变换后作为额外特征输入,模型在季节性明显的数据上往往能再涨一截精度。

需要专业的网站建设服务?

联系我们获取免费的网站建设咨询和方案报价,让我们帮助您实现业务目标

立即咨询