ARTICLE DETAIL

资讯详情

深耕网站视觉设计与运营推广的一线实战洞察。

SARIMA时间序列预测实战:季节性参数定阶与Python实现教程

SARIMA时间序列预测实战:季节性参数定阶与Python实现教程 简介面向计算机、电子信息、数学等专业学生的Python时间序列预测课程设计与毕业设计资源内含基于SARIMA模型的完整实现源码与配套数据集。代码采用参数化编程关键行均有注释几乎做到一行一注释方便初学者理解建模流程并灵活调整参数。数据文件以CSV和Excel格式提供包含焦作等多组实验数据可直接用于模型验证与结果复现。压缩包共5个文件包括1个主程序、3个CSV数据文件和1个Excel辅助表整体仅61KB轻量便捷源码逻辑清晰覆盖数据读取、差分检验、季节性分解、SARIMA建模与预测等核心步骤配合详细中文注释可帮助读者快速构建可运行的时间序列预测方案。目前已有448人学习尤其适合需要完整可运行代码支撑课程报告或毕设实验的同学。1. SARIMA 时间序列预测它凭什么适合大多数周期性业务想复现一个 Python 实现 SARIMA 时间序列预测完整源码和数据的项目多半是被月度报表、库存计划和流量预警逼来的序列有趋势、有周期性波动线性回归不够用LSTM 又嫌数据量不够。SARIMA 是这种场景下最稳的起点。它在差分自回归滑动平均模型的基础上扩出一组季节参数用一组可解释的数字描述“这个月与去年同月的关系”同时给出点预测和区间预测。这份笔记按“拿到源码就能跑”的标准拆全过程数据怎么准备、六组参数怎么定、statsmodels 怎么调通、参数怎么搜、哪几个坑一定会翻车。适合想在一周内把预测模型接进报表的数据分析师和 Python 工程师。2. 把 SARIMA 的六个参数拆开(p,d,q)(P,D,Q,s) 与业务周期的对应关系2.1 从 ARIMA 到 SARIMA季节项到底加了什么ARIMA 本身处理的是“去掉趋势后的平稳序列”它的三个参数 p、d、q 分别对应自回归阶数、差分次数、移动平均阶数。自回归项描述“当前值受过去 p 个时刻自身取值的影响”差分项负责把非平稳趋势消掉移动平均项描述“当前值受过去 q 个时刻预测误差的影响”。这套东西在月度销量、日活用户这类带周期波动的数据上有一个明显短板它只建模“昨天和今天的关系”建不了“这个月和去年同月的关系”。SARIMA 的改进就是多出一组季节参数 (P, D, Q, s)。P、D、Q 的含义与 p、d、q 类似但作用对象是“同一个季节位次的历史值”s 是季节周期长度月数据通常取 12周数据取 7日数据要看你业务的自然周期。模型的实际做法是先对序列做 d 阶普通差分和 D 阶季节差分再分别对差分后的序列和季节滞后序列做自回归与移动平均拟合。写成 statsmodels 的调用形式就是两个元组model SARIMAX( data, order(p, d, q), # 非季节部分 seasonal_order(P, D, Q, s) # 季节部分 )逻辑说明order 对应非季节参数seasonal_order 对应季节参数两者叠加在一个状态空间模型里联合估计。日常分析时最容易理解的方式是先把序列画出一年以上的趋势图肉眼看清有没有年度周期再决定要不要启用季节项没有明显周期的数据硬加季节项只会让参数估计方差变大。参数说明p 和 P 的取值范围建议从 0 到 2 开始试不要一上来就给大阶数d 和 D 绝大部分业务数据 0 或 1 就够D 取 2 的情况非常少强行增大差分次数会造成过差分预测结果变成一条不动的直线。s 的含义是“一个完整周期有多少个观测点”月度数据 12季度数据 4日数据按业务定这一步定了之后整个模型的形态就定了一半。2.2 定阶依据ADF 检验、ACF/PACF 图的读法定阶是 SARIMA 里最像“手艺活”的部分因为理论上 ACF 和 PACF 图能给出参考但真实业务数据的噪声常常让截尾和拖尾的边界很模糊。我的固定动作是先用 ADF 检验解决 d 的问题再画 ACF 和 PACF 解决 p、q 以及 P、Q 的初值。ADF 检验直接看 p 值p 值小于 0.05 表明序列平稳d 可以取 0p 值大于 0.05 就先做一阶差分再把差分后的序列重新跑一遍 ADF。这个步骤执行顺序有讲究不要上来就差分因为差分本身会损失信息能用 0 阶就不用 1 阶。观察对象图形特征一般选择ACF 拖尾、PACF 截尾PACF 在第 p 阶后快速落到置信区间内自回归项取 pACF 截尾、PACF 拖尾ACF 在第 q 阶后快速落到置信区间内移动平均项取 q季节滞后处s、2s的 ACF/PACF 明显突出季节性显著P 或 Q 至少取 1ACF 和 PACF 的“截尾”到底从第几阶算起实践中常常有争议所以我把这两个图定位成“给网格搜索提供范围”而不是直接判定唯一答案。比如 PACF 在 1 阶、2 阶处都显著就把 p 的范围设成 0 到 2让网格搜索去比 AIC而不是拍脑袋认定 p 就是 1。2.3 季节周期 s 与业务是绑定关系不是数学参数s 的取值最容易忽略业务约束。月度数据固定写 12 没问题但日数据就有讲究电商日销有“周内周期性”周一到周日波动明显s7有些平台流量有“工作日与周末”两态但自然周只有 7 天强行取 14 或 30 通常没有业务含义。季度数据 s4小时级数据 s24 或 s24×7需要先问清楚排班和促销节奏。我踩过最典型的一次某饮料数据按周汇总业务周期其实来自天气和节假日没有严格的 12 个月周期硬套 s12 的结果是季节项估计不稳定预测区间被拉得极宽。后来把粒度改成月、s 保持 12换了个数据口径才正常。所以拿到数据先画两个完整周期的曲线数一数峰和谷的间隔再写 s这比任何检验都直接。3. 用 Python 跑通 SARIMA 全流程源码结构、数据准备与第一版预测3.1 工程目录与数据集格式别把时间列读成字符串“完整源码和数据”要落地先解决文件组织问题。我建议一个最小工程目录所有代码和输出不散落sarima_demo/ ├── data/ │ └── sales.csv ├── src/ │ └── sarima_pipeline.py └── output/ ├── diagnostics.png ├── forecast.png └── metrics.csvdata 目录放原始数据src 放脚本output 放图表和评估指标。这样反复调参时每次运行结果都留在 output不会把临时文件混进源码。CSV 的格式是最容易出事的点。SARIMA 需要的输入本质上是一个按时间顺序排列的单变量序列所以 CSV 只要两列date 和 value。date 必须能被 pandas 解析成日期value 是数值。我经常看到有人把日期读进来成了字符串后续画图和建模全部报错问题就出在 read_csv 少了一个 parse_dates 参数。先把数据加载这部分固定下来import pandas as pd import numpy as np # 生成一份带趋势、季节和噪声的模拟数据方便先跑通流程 np.random.seed(42) dates pd.date_range(start2020-01-01, periods72, freqMS) trend np.linspace(100, 160, 72) seasonal 12 * np.sin(2 * np.pi * np.arange(72) / 12) noise np.random.normal(0, 3, 72) value trend seasonal noise pd.DataFrame({date: dates, value: value}).to_csv(data/sales.csv, indexFalse) # 正式加载关键在 parse_dates 和 index_col df pd.read_csv(data/sales.csv, parse_dates[date], index_coldate) df df.asfreq(MS) df df.dropna()逻辑说明先造一份 6 年月度数据用于演示包含线性趋势、12 个月周期的正弦季节项和随机噪声读者如果手头没有业务数据可以直接切到这份模拟数据验证代码。正式加载后必须做 asfreq(MS)这一步给索引显式指定“月起始”频率SARIMAX 做季节差分时依赖这个频率属性缺了它后面会报错。参数说明parse_dates 让 date 列变成 DatetimeIndexindex_col 把它设为索引asfreq 的频率码要按数据粒度换月度数据用 MS日度数据用 D周度数据用 W-MON。新版 pandas 里旧的月结束频率缩写 M 有弃用警告统一用 MS 最稳妥。dropna 是防止 asfreq 补出来的间隙空值直接参与建模。3.2 数据探索与平稳性检验的固定动作模型跑通之前先做两个检查序列有没有明显季节差分阶数该取多少。这两件事不确认后面所有预测都是盲人摸象。绘图时注意中文字体设置否则图里一片方框检查效率极低import matplotlib.pyplot as plt from statsmodels.tsa.stattools import adfuller plt.rcParams[font.sans-serif] [Microsoft YaHei, SimHei, Noto Sans CJK SC] plt.rcParams[axes.unicode_minus] False df[value].plot(figsize(12, 4), titleRaw Series) plt.tight_layout() plt.savefig(output/raw.png, dpi150) adf_result adfuller(df[value].dropna()) print(ADF p-value:, adf_result[1])逻辑说明先画出原始序列图看趋势和季节是否肉眼可见再对原序列做 ADF 检验。这里 adfuller 返回的是一个元组第 2 个元素就是 p 值p 值小于 0.05 说明不需要差分。如果 p 值大于 0.05就对序列做一阶差分并重新检验。差分操作要生成一个新列不要在原始列上直接覆盖因为你后面还要做预测反变换保留原始数据能少踩不少坑df[value_diff] df[value].diff() diff_pvalue adfuller(df[value_diff].dropna())[1] print(1st diff p-value:, diff_pvalue)逻辑说明一阶差分后 p 值通常会显著变小。如果还没有先检查数据里是否混入大量重复值和异常值而不是继续加差分次数。3.3 模型训练与诊断statsmodels 最小可运行代码定阶的严谨做法放到第 4 章这里先用一组能跑通的参数把流程串起来。以月度数据为例最常用的起步组合是 order(1,1,1)seasonal_order(1,1,1,12)。这两个 (1,1,1) 分别代表“自回归 1 阶、差分 1 阶、移动平均 1 阶”和“季节自回归 1 阶、季节差分 1 阶、季节移动平均 1 阶”是经验里覆盖大多数业务序列的起始点。from statsmodels.tsa.statespace.sarimax import SARIMAX model SARIMAX( df[value], order(1, 1, 1), seasonal_order(1, 1, 1, 12), enforce_stationarityFalse, enforce_invertibilityFalse ) result model.fit(dispFalse, maxiter200) result.summary()逻辑说明这里用 SARIMAX 而不是老的 SARIMA 类是因为它支持更完整的推断和置信区间输出。enforce_stationarity 和 enforce_invertibility 设为 False目的是放宽参数空间限制避免优化器在边界上收敛失败。对很多真实业务数据默认 True 会偶尔报“似然不收敛”的错放宽之后反而稳定。参数说明maxiter200 是状态空间模型优化迭代上限数据长、参数多时可以提到 500。dispFalse 关闭迭代过程打印保持控制台干净。result.summary() 输出的表格里有每个参数的标准误和 p 值p 值过大的参数可以考虑从模型里去掉但这不是绝对标准别单看显著性裁员还要结合 AIC 一起看。训练完必须看诊断图这是判断残差是否还有可利用信息的核心步骤result.plot_diagnostics(figsize(12, 8)) plt.tight_layout() plt.savefig(output/diagnostics.png, dpi150)逻辑说明诊断图的四个子图分别展示标准化残差序列、残差直方图、QQ 图和残差自相关。重点看两处右上角 QQ 图近似一条直线说明残差接近正态右下角自相关全部落在蓝色置信区间内说明残差里没有剩下的模式可以提取。任何一张图出现明显结构都意味着当前的 p、q 阶数没有充分拟合数据。3.4 预测与评估把 MAE、RMSE 和置信区间一次算出来训练完模型就要回答“到底准不准”。我会先把数据切成训练集和测试集模拟“用过去预测未来”的真实场景再算误差指标。测试集大小通常取一个完整季节周期月度数据取 12日度数据取 7 或 30。from sklearn.metrics import mean_absolute_error, mean_squared_error n_test 12 train df.iloc[:-n_test].copy() test df.iloc[-n_test:].copy() train.index.freq MS model SARIMAX( train[value], order(1, 1, 1), seasonal_order(1, 1, 1, 12), enforce_stationarityFalse, enforce_invertibilityFalse ) result model.fit(dispFalse) forecast_obj result.get_forecast(stepsn_test) forecast_mean forecast_obj.predicted_mean forecast_ci forecast_obj.conf_int(alpha0.05) mae mean_absolute_error(test[value], forecast_mean) rmse mean_squared_error(test[value], forecast_mean, squaredFalse) print(fMAE{mae:.2f}, RMSE{rmse:.2f})逻辑说明train/test 切分后重新拟合模型get_forecast 输出的是未来 n_test 步的预测均值与 95% 置信区间。MAE 看平均绝对偏差RMSE 放大较大误差两者结合能反映真实业务损失。最后把预测结果连同置信区间画出来这张图可以直接贴进周报plt.figure(figsize(12, 5)) plt.plot(train.index, train[value], labeltrain) plt.plot(test.index, test[value], labeltest, colorgreen) plt.plot(forecast_index, forecast_mean, labelforecast, colorred) plt.fill_between( forecast_index, forecast_ci.iloc[:, 0], forecast_ci.iloc[:, 1], colorred, alpha0.2, label95% CI ) plt.legend() plt.tight_layout() plt.savefig(output/forecast.png, dpi150)逻辑说明forecast_index 就是测试集的时间索引。填充区域是置信区间区域越宽代表模型越不确定极端情况下区间宽到失去业务意义就需要回到参数优化步骤重新搜一次。4. 参数寻优这步棋AIC、网格搜索与 auto_arima 的取舍4.1 依赖 AIC 的两个前提多数人没验证过AIC 在 SARIMA 里被当成“模型选择标准”但它本质上是一个相对指标只在同一个数据集上的不同候选模型之间做比较才有意义。AIC 权衡的是拟合优度与参数数量的平衡它完全不关心业务损失函数缺货损失和库存积压成本在 AIC 眼里是一样的。所以 AIC 能帮你缩小候选范围不能替你拍板。我见过最典型的翻车现场是某个零售项目 AIC 最小的模型是 (2,1,3)(1,1,2,12)参数多、拟合漂亮但一到节假日预测就剧烈震荡因为高阶移动平均项把噪声当成了信号。后来在验证集上比 MAE才发现最稳的反而是 (1,1,1)(0,1,1,12)。所以我的习惯是先用 AIC 挑前 5 个候选再逐一在验证集上算滚动预测误差最后结合业务确定模型。4.2 网格搜索的实现自己写参数表最可靠自己写网格搜索的优势是可复现、可控制、能加业务约束。以下代码固定 s12搜索全部 18×8144 个组合对数据量不大的场景完全够用import itertools p_range range(0, 3) d_range range(0, 2) q_range range(0, 3) P_range range(0, 2) D_range range(0, 2) Q_range range(0, 2) pdq_list list(itertools.product(p_range, d_range, q_range)) seasonal_pdq_list [(P, D, Q, 12) for P, D, Q, _ in itertools.product(P_range, D_range, Q_range, range(1))] results [] for order in pdq_list: for seasonal_order in seasonal_pdq_list: try: model SARIMAX( train[value], orderorder, seasonal_orderseasonal_order, enforce_stationarityFalse, enforce_invertibilityFalse ) result model.fit(dispFalse) results.append((result.aic, order, seasonal_order)) except Exception as e: print(ffail: {order} {seasonal_order}: {e}) results.sort(keylambda x: x[0]) print(results[:5])逻辑说明itertools.product 生成笛卡尔积两层循环遍历全部非季节与季节参数组合。try/except 跳过不收敛或海森矩阵奇异的组合保留可复现的报错信息。排序后输出 AIC 最小的前 5 个组合作为候选模型列表。参数说明d_range 只取 0 和 1D_range 只取 0 和 1这是刻意的限制。把 d 或 D 放到 2 会让搜索空间膨胀而且绝大多数业务序列用不到二次差分。p、q 与 P、Q 的搜索范围也控制在 0 到 2因为更大阶数对月度数据基本是过拟合。如果你的数据明显更长、周期更复杂再逐步扩大范围。网格搜索运行时要注意时间和资源。144 个组合每个模型拟合几秒到十几秒总时长可能要十几分钟这是正常现象不要以为卡死了。数据量大时先固定 d、D只扫 p、q、P、Q能把组合数砍掉一半以上。4.3 auto_arima 的坑随机搜索导致结果不可复现pmdarima 的 auto_arima 在时间序列预测的 Python 生态里很常用它号称自动搜索参数。实际用起来有个硬伤不固定随机源时两次运行可能给出完全不同的参数组合因为它的 stepwise 搜索路径里有一部分随机行为。复现不出结果意味着你没法向同事证明“为什么选这组参数”。from pmdarima import auto_arima auto_model auto_arima( train[value], seasonalTrue, m12, stepwiseTrue, information_criterionaic, random_state42, n_jobs1 ) print(auto_model.order, auto_model.seasonal_order)逻辑说明这里显式固定 random_state42n_jobs1 避免并行计算带来额外随机性。即便如此我仍然建议把 auto_arima 的结果当作网格搜索的初始化参考而不是直接采纳。它的搜索空间设定和步长策略与你手写的网格不完全一致实际效果常常不如针对性搜索。auto_arima 挑出的组合如果没有业务解释比如季节差分 D0 而你的数据肉眼可见强周期我一般会直接否定这个结果回归到手写网格。自动工具能帮你缩小范围但没法替你理解业务周期。5. SARIMA 实战避坑五个让预测翻车的细节5.1 报错 ValueError日期索引必须带 freq现象模型.fit() 或 asfreq() 抛出“ValueError: Given an endogenous array… frequency must be set”或者干脆警告“A date index has been provided, but it has no associated frequency”。原因read_csv 只生成了 DatetimeIndex却没有给索引分配频率属性。SARIMAX 需要知道观测间隔是否恒定、季节周期 s 对应多少个观测点频率缺失时状态空间方程无法确定时间步长。解决在建模前显式指定频率。月度数据用 df df.asfreq(MS)日度数据用 df df.asfreq(D)。如果索引已经存在重复日期或时间戳乱序先 df df.sort_index() 再 asfreq。这是整个流程里出现频率最高、解决最容易的报错几乎每个新手都会撞一次。5.2 预测曲线变成一条“直线”看不到任何波动现象模型训练没有报错预测结果也正常输出但预测值几乎等于一条平移的水平线把季节波动全部抹平了。原因最常见是过差分即 d 或 D 取得过大。差分次数过多会消除序列中真实的信息成分模型只能学到均值回归。另一种可能是季节周期 s 写错模型把季节项拟合成了噪声。解决先砍差分阶数把 (p,d,q) 里的 d 从 2 降到 1seasonal_order 里的 D 从 1 降到 0重新拟合看 AIC 和预测形态有没有改善。同时回头用第 3 章的诊断图重点看残差自相关图里是否还有明显的季节峰值有就说明季节项没有充分建模。5.3 预测值整体滞后真实值一整拍现象预测曲线形状和真实值一致但整体右移了一到两个时间点像延迟回放。原因这类现象在短期预测里很常见尤其是自回归阶数 p 或季节自回归阶数 P 偏小、模型过度依赖前几个观测点时。本质上模型学到的是“跟随上一期”而不是“提前转身”。某些业务里这是可以接受的但库存和安全库存场景下滞后一拍意味着补货永远晚一周。解决增加 p 或 P 的阶数让模型看到更长的历史窗口如果序列带强趋势先确认 d1 是否真正落实最后检查是不是 s 写错。不要在评估时用“把预测曲线往前平移一个点”来掩盖问题那是给老板看图的把戏不是预测模型的解法。5.4 auto_arima 两次运行结果不一样现象同一份数据、同一段代码第一次跑出 (1,0,1)(0,1,1,12)第二次跑出 (2,0,0)(1,0,0,12)AIC 还差不多。原因前面提过pmdarima 的搜索路径里有随机成分且并行模式 n_jobs 默认行为不够确定。这不是模型问题而是搜索算法的复现性问题。解决显式固定 random_state设置 n_jobs1并且在代码注释里记录搜索范围、信息准则和所有关键参数。更好一点的做法是弃用 auto_arima 作为最终决策工具只把它拿来做候选范围参考最终参数以手写网格搜索结果为准。5.5 中文乱码与“黑匣子”式的画图失败现象plot_diagnostics 和预测图画出来后标题、图例全是方框或者 Matplotlib 直接报字体缺失警告。原因Matplotlib 默认字体不支持 CJK 字符这是环境问题不是模型问题。很多教程没写这一段导致你在复现时第一张图就看不懂。解决绘图前统一设置字体列表我的习惯是同时列三个候选字体覆盖面更广。axes.unicode_minusFalse 保证负号正常显示。如果是无图形界面的服务器记得把 matplotlib 后端设为 Agg 并把图片保存到文件再通过文件查看。6. 让预测接得住业务滚动回测、置信区间与我的调参习惯网格搜索选出的“最优”参数并不能直接证明它在未来表现稳定。一次性预测 12 个月前 3 个月的误差会拖低整段平均误差掩盖后几个月的问题。我常用的验证手段是滚动回测每预测一步就把真实值追加到历史里重新拟合再预测下一步。这更接近真实上线后的行为也更容易暴露模型在拐点处的迟钝。history train[value].copy() preds [] for i in range(len(test)): model SARIMAX( history, orderbest_order, seasonal_orderbest_seasonal_order, enforce_stationarityFalse, enforce_invertibilityFalse ) res model.fit(dispFalse) pred res.forecast(1).iloc[0] preds.append(pred) history pd.concat([history, test[value].iloc[[i]]]) if history.index.freq is None: history.index.freq MS mae_rolling mean_absolute_error(test[value], preds) print(fRolling MAE{mae_rolling:.2f})逻辑说明每次循环只用当前历史做一步预测然后立刻把测试集里的真实值并入历史。这里的 freq 检查是必须的concat 之后索引频率可能变成 None下一轮建模会爆 5.1 节那个错。滚动回测的结果如果比一次性预测差很多说明模型稳定性有隐患不要直接上线。置信区间应该被当作业务参数而不只是画图装饰。get_forecast 返回的 conf_int 在 alpha0.05 时给出 95% 区间区间宽度反映模型在预测远端时的不确定性。做库存计划时我会取区间上沿作为安全库存参考而不是用预测均值做流量预算时反而看区间下沿避免过度配置资源。数据量小或 s 较长时区间宽得吓人是正常的这恰好说明你需要更长的历史或更简单的模型。这套流程跑熟的标志是拿到一份新数据先画图、再定 d 和 D、再搜 p 和 q、最后滚动验证。我通常会把每一次搜索的参数、AIC、MAE、数据时间范围记在一个 csv 里时间久了就能沉淀出一张“哪些参数组合在哪些场景下有效”的经验表。预测模型的本质不是黑匣子出数字而是把业务周期翻译成可调参数再把不确定性摆到明面上。这个习惯救过我很多次希望帮到你。本文还有配套的精品资源点击获取
返回列表