简介:面向需要构建空气质量预测模型或研究时序预测算法的学生与开发者,这份Python工程基于麻雀搜索算法(SSA)与长短时记忆网络(LSTM)实现了完整预测流程。代码按模块拆分为MLP基线、LSTM训练、SSA自动寻优LSTM超参数、优化参数回填复现及多模型对比,可直接对照运行,便于理解SSA与LSTM结合的建模思路与调参逻辑。麻雀搜索优化部分覆盖学习率、第一隐含层节点、第二隐含层节点等关键超参数,并输出适应度曲线与优化过程图,结果数据中同时保留MLP、LSTM和SSA-LSTM的预测效果,方便横向比较。压缩包共16个文件,以py源码为核心,配有mat结果数据、png优化曲线与可视化图、xls数据集及md说明文档,整体约778KB,目录清晰,适合毕业设计、课程实验或算法对比研究使用。目前已有807人学习下载,借助其中完整的参数寻优过程与结果对比,可快速复现基于SSA-LSTM的空气质量预测实验。
1. 先拆解SSA-LSTM这件事:它到底在优化什么,又解决了什么
手里有一份某城市监测站的小时级PM2.5浓度数据,要预测未来24小时,传统多项式回归跑不出非线性,LSTM能跑但超参很难定,这正是基于SSA-LSTM(麻雀搜索算法优化LSTM)的空气质量预测Python代码要解决的问题。这条路线的核心不是换个网络结构,而是用麻雀搜索算法替LSTM找一组更合适的超参数:学习率、隐含层单元数、批大小、Dropout。找到之后再训练LSTM做时序预测。适合三类人:做大气环境数据分析的研究生、做毕业设计需要“智能优化算法+深度学习”组合的学生,以及被固定超参折磨的中小型环境预报工程师。下面直接讲清楚原理、代码骨架和坑,你拿到后能换成自己的数据跑通。
2. 麻雀搜索算法为什么适合时序预测:SSA机制、选型理由与最小实现
2.1 发现者、跟随者、警戒者:SSA的三个角色怎么协作
SSA是群体智能优化算法,灵感来自麻雀觅食与反捕食。种群中每个个体代表一个候选解,对应LSTM的一组超参数。发现者负责大范围探索,位置更新偏向大步长;跟随者向当前最优位置靠拢,负责精细开发;警戒者随机出现,如果感知到危险就跳开,避免整个种群快速掉进同一个局部最优。这个协作方式很关键,因为LSTM超参空间是非凸的,学习率8e-4和9e-4可能带来完全不同的收敛路径,SSA相当于在超参空间里撒一批麻雀,让发现者探索、跟随者收割。
实际LSTM超参搜索里,我们通常把每个麻雀的位置定义成一个4维向量:学习率、LSTM单元数、批大小、Dropout。适应度函数就是“用这组超参训练一个LSTM,再在验证集上算RMSE”。和遗传算法相比,SSA不需要选择、交叉、变异三个独立算子,代码量更小;和粒子群相比,SSA多了发现者和警戒者的角色,后期跳出局部最优的概率稍高。但注意,SSA的数学性质并不保证一定优于粒子群或贝叶斯优化,它胜在容易实现、可解释,尤其适合论文里需要展示“智能算法收敛曲线”的场景。
2.2 用SSA还是用网格搜索:给定一组超参的代价对比
如果只用网格搜索,四五个超参全组合的搜索次数会迅速爆炸。比如学习率10个候选、单元数5个候选、批大小3个候选、Dropout 3个候选,全组合是450次LSTM训练,大多数项目根本跑不完。随机搜索砍到50组,但它是盲目的,不会根据上一条结果调整下一组候选。SSA有记忆,保存了全局最优位置,还通过发现者和跟随者协作,让搜索方向有反馈。当评估次数只有几十次时,SSA通常会给出比随机搜索更低的验证RMSE。
但SSA不是免费的。每一次迭代要评估所有麻雀个体,种群10、迭代10,等于额外训练100次LSTM。所以我建议把SSA当成“粗筛工具”:先在小数据上快速找到超参范围,再在这个范围附近做一次细调,而不是让它和全量数据硬碰硬。很多时候,SSA跑出来的最优超参和人工经验值很接近,但它能帮你判断“这个模型到底还能压榨多少性能”。另外,代码里一定要给所有个体加边界约束,否则学习率可能变成负数,批大小可能小于1,LSTM直接报错。
2.3 SSA优化LSTM的两种常见接法
第一种常见做法是超参数寻优,也就是把SSA放在LSTM外面,每个麻雀个体代表learning_rate、units、batch_size、dropout,适应度是该超参下LSTM在验证集上的RMSE。第二种是把LSTM的全部权重和偏置展开成一个一维向量,让SSA直接搜索权重空间。我基本只用第一种,因为一个中型LSTM的权重可能超过十万维,SSA在连续高维空间的搜索效率和精度都很难保证,就像让麻雀在高原上找一粒沙子。
论文里很多“SSA-LSTM”指的其实是超参寻优,代码也更好改。如果你要做权值初始化方向,建议把SSA的每个个体只优化最后全连接层或attention参数,而不是整网权重。网上能搜到不少S SA-LSTM示例代码,但很多是Matlab版本,转成Python时要注意维度顺序:Matlab默认列优先,Python的numpy和TensorFlow默认行优先,直接复制容易踩坑。
2.4 手写一个最小SSA更新循环(Python)
要接上后面的LSTM训练,先得跑通SSA的位置更新。下面是最小实现,省略了警戒者,只保留发现者和跟随者,目的是看清超参位置怎么移动。
import numpy as np def ssa_update(positions, sorted_idx, it, max_it): n, dim = positions.shape new_positions = positions.copy() r2 = np.random.rand() alpha = np.random.rand() # 发现者:前20%的麻雀,负责向大步长方向搜索 n_discover = int(n * 0.2) for i in range(n_discover): if r2 < 0.8: # 环境安全,在当前解附近小步探索 new_positions[i] = positions[i] * np.exp( -i / (alpha * max_it) ) else: # 发现危险,随机跳开 new_positions[i] = positions[i] + np.random.randn(dim) * 0.1 # 跟随者:后80%的麻雀向当前最优位置靠拢 best_pos = positions[sorted_idx[0]] worst_pos = positions[sorted_idx[-1]] for i in range(n_discover, n): if i < n / 2: new_positions[i] = best_pos + np.random.randn(dim) * np.abs( positions[i] - best_pos ) else: new_positions[i] = np.random.randn(dim) * np.exp( (worst_pos - positions[i]) / i**2 ) # 边界裁剪,防止学习率/单元数跑到无效区间 lb = np.array([0.0001, 8, 16, 0.1]) ub = np.array([0.01, 128, 64, 0.5]) return np.clip(new_positions, lb, ub)这个函数接收当前种群位置和按照适应度从小到大排列的序号,返回下一代位置。维度设计成4维:第一维学习率,第二维LSTM单元数,第三维批大小,第四维Dropout。批大小和单元数本质是整数,但SSA是连续优化,所以边界里用浮点,解出来之后再转int。
代码里random.randn(dim)*0.1是随机扰动,避免所有发现者在几步之内挤到一个点,这是实际使用中最常见的翻车点。如果你想要完整SSA,可以再加10%的警戒者,每次对全局最优位置施加一个小随机扰动,并让扰动幅度不超过当前位置的10%。但先跑通这个简化版,再想象加警戒者。
3. 用Python把SSA-LSTM跑起来:数据构造、适应度函数与训练主循环
3.1 滑动窗口造样本:空气监测数据怎么变成监督学习
空气质量预测最常见的输入是CSV表格,一列时间戳,后面跟着PM2.5、PM10、NO2、SO2、CO、O3以及温度、湿度、风速等。LSTM不认时间索引,只能认固定长度的二维数组,所以要把原始序列切成时间步长×特征数的窗口。我习惯用过去24小时预测下一小时PM2.5,如果做未来24小时预测,可以改成多步输出,但下面先给出单步版本。
处理顺序很有讲究:先按时间切分训练集和测试集,再只用训练集做MinMaxScaler的fit,绝不能用全量数据fit,否则测试集的信息会通过缩放器的min/max泄漏进训练过程。PM2.5浓度可能从个位数跳到300以上,归一化这一步不能省略。
import pandas as pd import numpy as np from sklearn.preprocessing import MinMaxScaler # 示例列名,实际按你的监测站点字段调整 df = pd.read_csv("air_quality.csv", parse_dates=["time"]) df = df.set_index("time").sort_index() df = df.replace([np.inf, -np.inf], np.nan) # 缺失值先按时间插值,再前后填充 df = df.interpolate(method="time").fillna(method="ffill").fillna(method="bfill") features = ["PM25", "PM10", "NO2", "SO2", "CO", "O3", "TEMP", "HUMI", "WSPD"] data = df[features].values.astype(np.float32) train_len = int(len(data) * 0.7) scaler = MinMaxScaler(feature_range=(0, 1)) scaler.fit(data[:train_len]) # 只用训练段拟合 scaled = scaler.transform(data) # 再转换全量,用于后续切片注意interpolate(method="time")要求索引是时间类型,且时间步长均匀。如果数据是小时级,这个方法能很好处理短时间缺失;如果出现连续一整天的缺失,插值会插出一条平滑的疑问曲线,最好还是用同一个监测站前后几日的均值或移除这一段。这个细节在第4章继续讲。
然后是切窗口的函数:
def make_samples(scaled, time_step=24, target_col=0): X, y = [], [] for i in range(time_step, len(scaled)): X.append(scaled[i - time_step:i, :]) # 过去24小时所有特征 y.append(scaled[i, target_col]) # 当前时刻PM2.5 return np.array(X), np.array(y) time_step = 24 X_all, y_all = make_samples(scaled, time_step) split = train_len - time_step X_train, y_train = X_all[:split], y_all[:split] X_val, y_val = X_all[split:split + int(len(X_all) * 0.1)], y_all[split:split + int(len(X_all) * 0.1)] X_test, y_test = X_all[split + int(len(X_all) * 0.1):], y_all[split + int(len(X_all) * 0.1):]这里split是训练段最后能构造窗口的位置,因为归一化只用了训练段,切分时训练样本的行号最多到train_len - time_step - 1,它的y索引恰好是train_len - 1,不会跨到测试段。验证集和测试集从train_len开始往后取。target_col=0表示PM2.5列,如果要预测PM10,把索引改成对应列。窗口大小24是经验和效率的折中,LSTM对更长时间依赖不是线性增长,不要盲目设置72步。
3.2 适应度函数:每次进化都要“真正训练一次LSTM”
SSA的适应度是整个代码里最耗时的部分。每个麻雀个体代表一组超参,适应度就是用这组超参构建LSTM、在训练集上快速训练几个epoch、在验证集上算RMSE。这里有个血泪经验:不要在适应度里用完整训练epoch,比如50甚至100,否则一次SSA迭代就好几小时。我把epoch压到3,验证集RMSE略有波动但能区分超参好坏;如果数据集比较长,可以加到5。
模型内部随机性会让相同超参训练结果不同,导致SSA排序不稳定,所以在模型构建后固定随机种子。数据量大时小epoch的噪声可以被随机种子抵消一部分,但不要完全依赖它。
import tensorflow as tf def build_lstm(units, lr, dropout, time_step, n_features): tf.random.set_seed(20240630) model = tf.keras.Sequential([ tf.keras.layers.LSTM(units, return_sequences=True, dropout=dropout), tf.keras.layers.LSTM(units, dropout=dropout), tf.keras.layers.Dense(1) ]) model.compile( optimizer=tf.keras.optimizers.Adam(learning_rate=lr), loss="mse" ) return model def fitness(params, X_train, y_train, X_val, y_val): lr = params[0] units = int(params[1]) batch_size = int(params[2]) dropout = np.clip(params[3], 0.1, 0.5) model = build_lstm(units, lr, dropout, time_step=X_train.shape[1], n_features=X_train.shape[2]) model.fit( X_train, y_train, epochs=3, batch_size=batch_size, validation_data=(X_val, y_val), verbose=0 ) pred = model.predict(X_val, verbose=0).ravel() return float(np.sqrt(np.mean((pred - y_val) ** 2)))return_sequences=True让第一层输出完整序列,第二层输出最后一个时间步,再接Dense(1)。如果数据量小,可以只保留一层LSTM,但空气质量序列通常有早晚峰值,两层LSTM更容易抓长期依赖。units边界8到128,太小欠拟合,太大过拟合且训练慢。batch_size 16到64,取值必须转int,浮点batch_size会直接报错。dropout 0.1到0.5,LSTM对dropout敏感,太高会丢记忆信息。
3.3 SSA-LSTM训练主循环:种群、迭代与最优超参
现在把第2章的SSA更新和3.2的fitness接起来。种群数建议6到10,迭代5到8,不要一上来就设50。种群8、迭代6、适应度3 epoch,相当于额外训练48次小LSTM,已经不小了。正式跑全量之前,先用500条数据试通,再上全量。
n_pop = 8 max_iter = 6 lb = np.array([0.0001, 8, 16, 0.1]) ub = np.array([0.01, 128, 64, 0.5]) # 均匀初始化种群 positions = np.random.uniform(lb, ub, (n_pop, 4)) fitness_vals = np.zeros(n_pop) for it in range(max_iter): for i in range(n_pop): fitness_vals[i] = fitness( positions[i], X_train, y_train, X_val, y_val ) # 按RMSE从小到大排序 sorted_idx = np.argsort(fitness_vals) print(f"Iter {it + 1}/{max_iter}, best RMSE: {fitness_vals[sorted_idx[0]]:.4f}") # 除最后一轮外更新种群 if it < max_iter - 1: positions = ssa_update(positions, sorted_idx, it, max_iter) best_idx = sorted_idx[0] best_params = positions[best_idx] print("Best params:", best_params)这里容易混淆的点是:sorted_idx是每轮根据当前适应度算出来的,但ssa_update内部通过positions[sorted_idx[0]]取最优位置。因为positions数组没有被排序,所以必须用argsort索引去访问原始数组,不能直接把fitness_vals排序后覆盖positions。
得到最优超参后,不要再只用原始训练集训练,因为SSA已经在验证集上挑过超参,验证集信息被用掉了。常见做法是把训练集和验证集合并,重新训练一个epoch更多、batch_size固定的最终模型,再用从未碰过的测试集去评估。
best_lr = best_params[0] best_units = int(best_params[1]) best_batch = int(best_params[2]) best_dropout = np.clip(best_params[3], 0.1, 0.5) final_model = build_lstm( best_units, best_lr, best_dropout, time_step=X_train.shape[1], n_features=X_train.shape[2] ) final_model.fit( np.concatenate([X_train, X_val], axis=0), np.concatenate([y_train, y_val], axis=0), epochs=20, batch_size=best_batch, validation_split=0.1, verbose=1 ) final_model.save("ssa_lstm_air_quality.h5")为什么训练最终模型时还要validation_split=0.1?这不是给SSA用,只是观察训练过程是否发散。如果验证损失升高,说明最优超参仍然过拟合,可以回退epoch或增大dropout。
3.4 预测、反归一化与评价指标
预测阶段用加载好的模型对测试集做前向推断。但SSA-LSTM训练时用的所有特征都被缩放到0到1之间,输出的PM2.5也是缩放后的值,必须反归一化回真实浓度微克每立方米,否则RMSE和MAE没有实际意义。
反归一化时不能直接用inverse_transform只处理单列,因为训练时的MinMaxScaler是在9列特征上fit的,反变换需要9列的二维数组。常见做法是构造一个和特征维度相同的零数组,把预测值放到目标列,反变换后取回目标列。
def inverse_scale_target(scaled_pred, scaler, target_col=0, n_features=9): dummy = np.zeros((len(scaled_pred), n_features), dtype=np.float32) dummy[:, target_col] = scaled_pred.ravel() return scaler.inverse_transform(dummy)[:, target_col] pred_scaled = final_model.predict(X_test, verbose=0).ravel() y_test_scaled = y_test.ravel() pred_true = inverse_scale_target(pred_scaled, scaler, 0, X_test.shape[2]) true_true = inverse_scale_target(y_test_scaled, scaler, 0, X_test.shape[2]) rmse = np.sqrt(np.mean((pred_true - true_true) ** 2)) r2 = 1 - np.sum((true_true - pred_true) ** 2) / np.sum( (true_true - np.mean(true_true)) ** 2 ) print(f"Test RMSE: {rmse:.2f} ug/m3, R2: {r2:.3f}")注意n_features必须和训练时的特征数一致,这里用X_test.shape[2]。如果单独只取PM2.5一列,scaler是在一列上fit的,直接inverse_transform(pred_scaled.reshape(-1,1))即可。我保留9个特征是因为空气预测通常会加气象变量,只用浓度序列会变成“自己预测自己”,缺乏外因驱动。
评价指标上,RMSE对异常高浓度更敏感,MAE反映普通误差,R²能看出模型是否比“直接取均值”更好。污染事件中PM2.5会突然从30跳到150,RMSE会被放大,所以还要记录峰值时段的误差,这个用第5章的滚动预测来缓解。
4. 实际跑空气质量预测时容易翻车的6个坑:现象、原因与排查方法
4.1 预测结果是一条平线:损失函数正常但输出不变
现象:训练loss下降,验证loss也不高,但画出测试集预测曲线是一条水平直线,几乎等于训练集目标均值。
原因:回归任务输出层误用了sigmoid,或者所有特征被归一化到0到1后目标方差过低,模型选择输出一个近似平均值来最小化MSE。另一个常见原因是有代码样例在输出层加了activation='relu'或tanh,把预测值截断在一个区间内,导致高浓度时段全部预测成同一个上限值。
解决:输出层保持无激活函数;检查目标列归一化后的方差,如果小于0.01,考虑增大时间步或加入更多解释变量;把batch_size降到32以下,防止模型训练不充分时过早收敛到均值。
4.2 反归一化后预测值整体偏移一个常数
现象:预测曲线和真实曲线走势一致,但整体低或高10到20微克每立方米。
原因:MinMaxScaler是在全量数据上fit的,或者单独用一个目标列的scaler,没有和训练时保持一致。网上很多示例代码为了省事,直接scaler = MinMaxScaler()然后scaler.fit(data),当测试集有未来更高的极值时,训练时模型的输出范围被压缩,反归一化就把未来新高拉到极值之外,产生偏移。
解决:严格按时间顺序fit训练段;保存缩放器参数;所有inverse_transform都使用同一个scaler。检查方法:打印scaler.data_min_里的PM2.5列,如果最小值是从未来数据里取的,赶紧重新预处理。
4.3 SSA迭代后期种群挤在一起,LSTM陷入重复超参
现象:适应度曲线前三轮明显下降,第四轮后完全不动,打印每个个体的超参发现所有行几乎一模一样。
原因:SSA发现者更新公式随迭代次数按指数缩小,后期步长几乎为0,跟随者快速贴到最优位置,种群多样性丢失。这和LSTM本身无关,是SSA的早熟收敛。
解决:给发现者的搜索步长加一个下限,比如max(0.1 * (1 - it / max_iter), 0.02);每轮更新后随机替换10%的个体为边界内的新随机位置;更简单的是把迭代压到6次以内,不让种群有时间完全聚合。第2章代码里省略警戒者,就是为了避免完整的预警机制在小数据上引发震荡,先求稳定再求精密。
4.4 空气质量数据有缺失和传感器漂移,训练时loss变成NaN
现象:训练到第二个epoch,loss突然变成nan,或者预测结果全是nan。
原因:原始数据里有NaN或无穷值,插值只做了fillna(method='ffill'),没有处理开头和结尾的空洞。传感器漂移会产生连续几天同一个数值,LSTM把这种假周期当成了真实规律。
解决:训练前统一做replace([np.inf, -np.inf], np.nan),然后按时间插值,再用前后填充;对连续重复超过6小时的CO、O3列做一次rolling median平滑,消除触点式漂移;如果PM2.5负值,用clip(lower=0)截断。另外可以关掉混合精度,半精度训练在数据不干净时会放大NaN。
4.5 只看训练集损失,模型在验证集上漂移严重
现象:训练loss很低,验证RMSE反而比随机初始化还高,且每隔一段时间预测值大幅偏离。
原因:LSTM容量太大,加上空气质量序列本身有强烈的日周期和气象扰动,模型把训练集的噪声也背了下来。SSA的适应度如果用训练集而不是验证集,选出的一定是过拟合超参。
解决:在fitness函数里必须传验证集,RMSE在验证集上计算;最终模型训练时多留出10%验证,并且用EarlyStopping(monitor='val_loss', patience=3, restore_best_weights=True)。这样才能避免报告一个“看似完美”的测试成绩。
4.6 种群和迭代设得太大,一次实验跑一个通宵
现象:SSA每轮都要训练几十个LSTM,种群20、迭代20,总训练次数等于400次小LSTM训练,即使每次只有3个epoch,全量数据也要数小时,最后没跑完就被强制中断。
原因:SSA的时间开销和训练次数线性相关,很多人误以为它能像超参搜索一样并行跑多GPU,但代码里所有个体是串行fit。
解决:先在小样本上比如最近2000个小时调好SSA参数范围,再上全量;把verbose调成0,batch_size尽量调大减少step轮次;如果机器有多GPU,可以用ThreadPoolExecutor并行评估fitness,但注意每个线程要设置不同随机种子,避免所有并行任务共享全局随机状态。
5. 往工程化走的两个小技巧:多步滚动预测与残差校正
5.1 多步滚动预测:用预测值反哺输入,避免单步强制
单步模型只能做“给定过去24小时,预测下一小时”。实际业务要求预报未来24小时甚至72小时,最直接的做法是把预测结果拼回输入窗口,再用模型预测下一步,这个循环在论文里叫滚动预测。好处是输入输出维度不变,代码改动很小;坏处是误差会累计,预测第24小时基本只能依赖趋势。我一般会同时报告第1、3、6、12、24小时的RMSE,而不是只报最后一步。
def rolling_forecast(model, scaler, history, steps=24, target_col=0): current = history[-time_step:].copy() # shape (time_step, n_features) preds = [] for _ in range(steps): p = model.predict(current[np.newaxis, :, :], verbose=0)[0, 0] preds.append(p) new_row = current[-1].copy() new_row[target_col] = p # 用预测值替换PM2.5列 current = np.vstack([current[1:], new_row]) return inverse_scale_target(np.array(preds), scaler, target_col, current.shape[1])history要传测试集中的最后一个真实窗口;如果要预报明天,那最后一个窗口来自当前时间。每次循环把最旧一行丢出去,再加一行预测结果,窗口始终是24步。如果模型对初始值很敏感,第二步之后预测曲线会越来越平滑,这是正常现象,因为模型失去了新的真实输入。
5.2 残差校正:用一个小线性模型吃掉LSTM的尾部误差
滚动预测误差会累积,单纯增加LSTM参数并不能消除。我习惯再加一个残差校正环节:把LSTM在训练集上的残差作为一个小序列,训练一个线性回归或Ridge模型,用上一时刻的残差和预测值预测下一时刻残差,最终预测等于LSTM输出加残差校正。这个技巧来自时序预测里的误差修正,不需要改LSTM结构,实现成本低,能显著提高第24小时步长的精度。
from sklearn.linear_model import Ridge # 假设 train_true, train_pred 是训练集的真实值和LSTM预测值 residuals = train_true - train_pred X_res = np.column_stack([train_pred, residuals])[:-1] y_res = residuals[1:] res_model = Ridge(alpha=1.0).fit(X_res, y_res) # 测试时没有上一时刻真实残差,用LSTM预测和上一轮校正后的残差递推 corr = np.zeros_like(pred_true) for i in range(1, len(pred_true)): corr[i] = res_model.predict([[pred_true[i - 1], corr[i - 1]]])[0] final_pred = pred_true + corr注意train_true和train_pred要来自训练集,而不是验证集,否则残差校正会偷看测试信息。实际数据里残差通常有自相关性,Ridge的alpha可以稍微调大,避免校正模型把异常点学进残差。
我自己跑这段代码时,一开始把SSA迭代次数设到50,结果一晚上只跑了两轮。后来把所有SSA参数压到种群6、迭代5,再逐步放开,才找到合理的超参范围。这个思路比追求“全局最优”更实在。希望这些步骤和坑能帮到你跑通自己的空气质量预测代码。
本文还有配套的精品资源,点击获取