简介这份资源是面向计算机、电子信息工程、数学等专业学生及地震数据分析初学者的MATLAB代码包围绕LSTM神经网络实现地震震级预测与地震数据可视化分析可用于课程设计、期末大作业和毕业设计等场景。压缩包共64个文件约1.96MB包含6个m脚本、3个py辅助脚本、3个csv数据文件以及31张png图表和2个fig图形文件另有txt、json、md等配置与说明文档覆盖数据获取、预处理、模型训练、仿真预测和结果可视化等环节。代码采用参数化编程参数修改方便注释详细、思路清晰并附赠案例数据可直接运行便于观察LSTM在地震时间序列上的预测效果。目前已有49人学习适合希望将神经网络理论落地到地震数据分析实践、快速搭建实验流程的读者参考。1. 从一份 MATLAB 地震震级预测代码说起它到底能跑出什么地震震级预测这件事圈外人觉得玄圈内人知道它更像一个典型的时序回归问题把前若干次地震事件的特征序列喂给网络让它输出下一次事件的震级估计。这份基于 LSTM 神经网络的地震震级预测与地震数据分析 MATLAB 代码解决的正是这条链路——从原始地震目录的清洗、特征构造到 LSTM 网络的搭建、训练、预测再到误差评估与可视化一整套流程都落在 MATLAB 里。它适合两类人一类是地球物理、防灾减灾方向的研究生和工程师手头有地震目录数据想快速验证时序模型能不能用另一类是做 MATLAB 神经网络练手的人想找一个比数字识别更有物理背景的回归任务。代码本身不依赖深度学习工具箱之外的重型框架常见做法是直接在 MATLAB 里跑通改改数据路径就能复现。2. LSTM 为什么适合震级序列原理、选型与数据准备2.1 震级序列的时序特性与 LSTM 的匹配逻辑地震事件按时间排列后本质上是一条不等间隔的时间序列。传统 BP 神经网络或者前馈神经网络处理这类数据时会把每个样本当成独立同分布的向量前一次地震对后一次的影响被彻底抹掉。而震级序列里恰恰存在这种跨事件的依赖余震序列、震群活动、应力积累释放都让相邻事件之间存在统计上的关联。LSTM 通过遗忘门、输入门、输出门三个结构控制细胞状态的更新能在较长的时间跨度上保留有用信息这正是它比普通 RNN 更适合震级预测的原因。普通 RNN 在反向传播时容易梯度消失序列一长就学不到早期信息LSTM 的门控机制缓解了这个问题。选 LSTM 而不是 Transformer是因为这份代码面向的是中小规模地震目录样本量通常几千到几万条Transformer 的自注意力机制在这种数据量下容易过拟合而 LSTM 的参数效率和训练稳定性更占优势。2.2 数据清洗与特征工程把地震目录变成网络能吃的矩阵拿到地震目录后第一步不是直接喂网络而是做清洗和特征构造。常见的地震目录包含发震时刻、经度、纬度、深度、震级这几个字段有些还带震源机制解。原始数据里经常有重复记录、震级缺失、时间格式不统一的问题。下面这段代码是我一般会先跑的预处理脚本把目录读进来、去重、按时间排序再构造滑动窗口样本。% 读取地震目录假设为 CSV 格式列顺序时间, 经度, 纬度, 深度, 震级 raw readtable(earthquake_catalog.csv); raw.Properties.VariableNames {time,lon,lat,depth,mag}; % 时间列转 datetime去除无效行 raw.time datetime(raw.time, InputFormat, yyyy-MM-dd HH:mm:ss); raw rmmissing(raw); % 按时间排序并去重同一时刻同一位置视为重复 raw sortrows(raw, time); [~, ia] unique(raw(:, {time,lon,lat}), rows, stable); raw raw(ia, :); % 构造滑动窗口用前 windowSize 次事件预测下一次震级 windowSize 10; numSamples height(raw) - windowSize; X zeros(numSamples, windowSize, 3); % 特征经度、纬度、深度 Y zeros(numSamples, 1); % 标签震级 for i 1:numSamples idx i:(i windowSize - 1); X(i, :, 1) raw.lon(idx); X(i, :, 2) raw.lat(idx); X(i, :, 3) raw.depth(idx); Y(i) raw.mag(i windowSize); end % 归一化避免量纲差异导致训练震荡 X (X - mean(X, [1 2])) ./ std(X, 0, [1 2]); Y (Y - mean(Y)) / std(Y);这段代码的逻辑分三层先去重排序保证时间顺序正确再用滑动窗口把序列切成监督学习样本最后做归一化。参数上windowSize控制用多少次历史事件预测下一次设太小捕捉不到长程依赖设太大样本数骤减且引入噪声我一般从 5 到 15 之间试。特征维度这里只用了经纬度和深度实际项目中可以把震级本身也作为历史特征加进去变成 4 维输入。注意归一化要对训练集和测试集分别做或者用训练集的均值和标准差去变换测试集否则会引入未来信息泄露。2.3 数据集划分与训练集构造的边界划分训练集和测试集时不能随机打乱。地震序列有时间顺序随机划分会让未来事件的信息泄露到训练集里评估结果虚高。正确做法是按时间切分比如前 80% 做训练后 20% 做测试。如果要做交叉验证也得用时间序列交叉验证而不是普通的 K 折。这份代码里如果用的是随机划分建议改成时间切分否则预测精度看着漂亮实际部署就翻车。3. LSTM 网络搭建与训练层参数、训练选项与回调3.1 网络结构定义层数、隐藏单元与输出层MATLAB 里搭 LSTM 回归网络核心是sequenceInputLayer、lstmLayer、fullyConnectedLayer和regressionLayer四个层。下面是一个可以直接跑的骨架。% 定义 LSTM 网络结构 inputSize 3; % 输入特征数经度、纬度、深度 numHiddenUnits 64; % LSTM 隐藏单元数 outputSize 1; % 输出震级 layers [ sequenceInputLayer(inputSize, Name, input) lstmLayer(numHiddenUnits, OutputMode, last, Name, lstm) dropoutLayer(0.2, Name, dropout) fullyConnectedLayer(outputSize, Name, fc) regressionLayer(Name, regression) ]; % 训练选项 options trainingOptions(adam, ... MaxEpochs, 100, ... MiniBatchSize, 32, ... InitialLearnRate, 0.005, ... LearnRateSchedule, piecewise, ... LearnRateDropPeriod, 30, ... LearnRateDropFactor, 0.5, ... GradientThreshold, 1, ... ValidationData, {XVal, YVal}, ... ValidationFrequency, 10, ... Shuffle, never, ... Verbose, false, ... Plots, training-progress); % 训练网络 net trainNetwork(XTrain, YTrain, layers, options);lstmLayer的OutputMode设为last表示只取最后一个时间步的输出用于回归这是序列到单值预测的标准做法。numHiddenUnits从 32 到 128 都有人用数据量小就取小一点否则容易过拟合。dropoutLayer的丢弃率 0.2 是经验值数据噪声大可以加到 0.3。GradientThreshold设 1 是为了防止梯度爆炸LSTM 训练时梯度爆炸比梯度消失更常见。Shuffle必须设成never因为序列样本之间有顺序依赖打乱会破坏时序结构。3.2 训练过程监控与早停策略训练时最怕两件事过拟合和欠拟合。过拟合的表现是训练损失持续下降但验证损失开始上升欠拟合则是两者都居高不下。MATLAB 的training-progress窗口能实时看到两条曲线但更稳妥的做法是加一个早停回调。常见做法是用ValidationData加ValidationPatience参数当验证损失连续若干轮不下降就停止训练。options trainingOptions(adam, ... MaxEpochs, 200, ... MiniBatchSize, 32, ... ValidationData, {XVal, YVal}, ... ValidationFrequency, 10, ... ValidationPatience, 15, ... Shuffle, never, ... Plots, training-progress);ValidationPatience设 15 表示验证损失 15 次评估不改善就停这个值太小会过早停止太大浪费训练时间。我一般先用 10 到 20 之间试看训练曲线再调。另外MiniBatchSize对 LSTM 影响很大太小导致训练震荡太大显存吃紧且泛化变差32 或 64 是比较稳的选择。3.3 预测与反归一化把网络输出还原成震级训练完网络后预测结果还是归一化后的值必须反归一化才能和真实震级对比。% 预测 YPred predict(net, XTest); % 反归一化 YPred YPred * std(YTrain) mean(YTrain); YTestActual YTest * std(YTrain) mean(YTrain); % 计算评估指标 rmse sqrt(mean((YPred - YTestActual).^2)); mae mean(abs(YPred - YTestActual)); r2 1 - sum((YTestActual - YPred).^2) / sum((YTestActual - mean(YTestActual)).^2); fprintf(RMSE: %.4f, MAE: %.4f, R2: %.4f\n, rmse, mae, r2);反归一化用的均值和标准差必须来自训练集不能用测试集的统计量否则评估结果不可信。RMSE 和 MAE 衡量绝对误差R2 衡量拟合优度。震级预测里 RMSE 能到 0.3 以下就算不错具体看数据质量。如果 R2 是负数说明模型还不如直接拿均值预测得回头检查数据泄露或者网络结构问题。4. 避坑与排查这份代码最容易翻车的五个地方4.1 现象训练损失正常下降验证损失剧烈震荡原因通常是MiniBatchSize太小或者学习率太高。LSTM 的梯度对批量大小敏感批量太小导致每次更新方向差异大。解决方法是把MiniBatchSize调到 64 或 128同时把InitialLearnRate从 0.01 降到 0.005 甚至 0.001观察验证曲线是否平滑。4.2 现象预测结果全部接近训练集震级均值这是典型的欠拟合网络没学到任何时序模式。原因可能是隐藏单元太少、训练轮数不够或者输入特征本身没有区分度。先把numHiddenUnits从 64 加到 128MaxEpochs加到 200如果还是这样检查特征构造是不是把震级信息漏掉了。另一个常见原因是归一化做错了比如对整个数据集归一化而不是对训练集导致测试集分布偏移。4.3 现象MATLAB 报错「Invalid training data. Responses must be a vector」这个报错一般是标签维度不对。LSTM 回归的标签应该是 N×1 的向量如果从 table 里取出来的是 N×1 的 table 列需要转成 double 数组。用Y double(raw.mag(idx))强制转换或者检查Y是不是被意外转成了 cell 数组。4.4 现象训练速度极慢一个 epoch 要跑十几分钟先看是不是用了 CPU 训练。MATLAB 的深度学习工具箱支持 GPU 加速如果有 NVIDIA 显卡在trainingOptions里加ExecutionEnvironment, gpu。另外sequenceInputLayer的输入如果是 cell 数组而不是三维数组会走慢速路径。确保X是numSamples × windowSize × numFeatures的三维数组而不是 cell。4.5 现象换一份地震目录后代码直接跑不通不同来源的地震目录列名、时间格式、震级标度都不一样。中国地震台网的数据震级可能是面波震级 Ms美国 USGS 用的是矩震级 Mw两者不能混用。换数据后先检查列名映射和时间解析格式再确认震级标度是否一致。如果目录里有负震级或者空震级rmmissing之后样本数可能骤减需要重新评估窗口大小是否还合适。5. 进阶技巧用贝叶斯优化调 LSTM 超参数并做滚动预测5.1 用 Experiment Manager 或 bayesopt 自动搜参手动调numHiddenUnits、InitialLearnRate、MiniBatchSize这三个参数很费时间。MATLAB 提供了bayesopt函数可以定义一个目标函数把验证集 RMSE 作为优化目标自动搜索参数组合。% 定义超参数搜索空间 params [ optimizableVariable(numHiddenUnits, [32, 128], Type, integer) optimizableVariable(InitialLearnRate, [1e-4, 1e-2], Transform, log) optimizableVariable(MiniBatchSize, [16, 128], Type, integer) ]; % 目标函数训练网络并返回验证集 RMSE objFcn (p) trainAndEvaluate(p, XTrain, YTrain, XVal, YVal); % 贝叶斯优化最多评估 30 次 results bayesopt(objFcn, params, ... MaxObjectiveEvaluations, 30, ... IsObjectiveDeterministic, false, ... Verbose, 1);trainAndEvaluate需要自己写内部根据传入的参数搭网络、训练、算验证 RMSE。IsObjectiveDeterministic设 false 是因为神经网络训练有随机性同一组参数跑两次结果可能不同。贝叶斯优化比网格搜索省时间30 次评估通常能找到比手动调更好的组合。注意搜索空间别设太宽InitialLearnRate用对数变换numHiddenUnits用整数类型。5.2 滚动预测模拟真实部署场景单次划分测试集只能看一个时间点的表现真实场景里模型是滚动使用的。滚动预测的做法是用前 N 次事件预测第 N1 次然后把第 N1 次的真实值加入历史窗口再预测第 N2 次如此循环。% 滚动预测 numRolling 50; YPredRolling zeros(numRolling, 1); YActualRolling zeros(numRolling, 1); currentWindow XTest(1, :, :); for i 1:numRolling pred predict(net, currentWindow); YPredRolling(i) pred * std(YTrain) mean(YTrain); YActualRolling(i) YTestActual(i); % 更新窗口去掉最早的事件加入新预测或用真实值 if i numRolling currentWindow [currentWindow(:, 2:end, :), XTest(i1, end, :)]; end end % 绘制滚动预测对比 figure; plot(YActualRolling, b-, LineWidth, 1.5); hold on; plot(YPredRolling, r--, LineWidth, 1.5); legend(实际震级, 预测震级); xlabel(滚动步数); ylabel(震级); title(LSTM 滚动预测对比); grid on;滚动预测里有个关键选择更新窗口时用预测值还是真实值。用真实值叫 teacher forcing评估偏乐观用预测值叫 free-running更接近实际部署但误差会累积。我一般两种都跑对比 RMSE 差距差距大说明模型对误差累积敏感实际用的时候要加滑动平均或者定期用真实值校正。5.3 一个我踩过的坑有次我用一份包含余震序列的目录训练测试集 RMSE 只有 0.15看着特别好。后来发现余震序列里相邻事件震级高度相关模型其实是在抄上一次的震级根本没学到物理规律。换成主震之间的事件重新训练后RMSE 回到 0.35 左右这才是真实水平。从那以后我每次做震级预测都强制把余震和主震分开评估并且检查模型输出和上一次震级的相关系数如果高得离谱基本就是抄答案。希望这份代码和这些经验帮到你少走点弯路。本文还有配套的精品资源点击获取