ARTICLE DETAIL

资讯详情

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

SARIMA时间序列预测Python源码拆解与实战避坑指南

SARIMA时间序列预测Python源码拆解与实战避坑指南 简介这份Python实现SARIMA时间序列预测的资源面向计算机、电子信息工程、数学等专业的大学生课程设计、期末大作业与毕业设计也适合刚入门时间序列分析的小白。资源包含完整源码与配套数据代码采用参数化编程关键参数可方便更改编程思路清晰几乎每一行都有注释便于理解SARIMA建模、参数寻优与预测流程。压缩包共5个文件包括1个Python脚本、3个CSV数据文件和1个XLSX表格分别用于模型实现、历史序列数据存储与辅助查看整体大小仅61KB轻量易用。目前已有448人学习下载。作者为资深算法工程师长期从事Python与Matlab算法仿真读者可借助这份源码快速上手SARIMA实战节省调参和排错时间并在此基础上扩展自己的预测任务。1. SARIMA时间序列预测这份源码能帮你解决什么问题时间序列预测在气象、交通、能源和电商场景里几乎是标配需求而SARIMA季节性差分自回归移动平均模型是其中最“能打”的统计学模型之一。今天拆的这份Python实现SARIMA时间序列预测源码自带完整数据集和保姆级注释适用于课程设计、期末大作业和毕业设计。数据包里包含焦作市的多组时间序列焦作.csv、焦作全.csv以及A.xlsx、A.csv程序本身是参数化编程换一组数据只要改文件路径和几个模型参数就能跑。它解决的是“拿到一份带季节性趋势的时间序列如何用Python完成建模、预测、可视化”这件事适合刚入门Python但需要快速产出结果的计算机、数学、电子信息专业学生也适合想用SARIMA做基线模型的量化分析新手。2. 数据准备与平稳性检验从原始表格到可建模的时间序列拿到SARIMA.zip解压之后里面除了SARIMA.py主程序还有五种数据文件。我一般先把数据读进来做初步探查再做平稳性处理。很多人在这一步直接跳过去拟合模型结果预测曲线一条直线回头还以为是模型问题。2.1 数据加载pandas读取Excel和CSV的差异import pandas as pd import numpy as np # 读取CSV文件 - 焦作市数据 df_jz pd.read_csv(焦作.csv, parse_dates[date], index_coldate) print(df_jz.head()) print(df_jz.info()) # 读取Excel文件 - 注意sheet_name参数 df_a pd.read_excel(A.xlsx, sheet_name0, parse_dates[date], index_coldate) print(df_a.head())这段代码先把日期列解析为datetime类型并设为索引这是时间序列建模的第一步。parse_dates参数让pandas自动识别日期字符串index_col指定哪一列作为行索引。读取Excel时要特别留意sheet_name有时候数据不在第一个Sheet里sheet_name0表示读取第一个工作表也可以传工作表名字。我发现焦作.csv和焦作全.csv的区别在于时间跨度不同前者数据量小适合快速测试后者适合做完整训练。2.2 平稳性检验与差分ADF检验与季节差分SARIMA模型要求输入序列是平稳的也就是说均值、方差不随时间变化。做ADF检验是最常见的做法from statsmodels.tsa.stattools import adfuller # 对原始序列做ADF检验 result adfuller(df_jz[value].dropna()) print(fADF统计量: {result[0]:.4f}) print(fp值: {result[1]:.4f}) for key, val in result[4].items(): print(f{key}: {val:.4f}) # p值 0.05 说明序列非平稳需要差分 if result[1] 0.05: df_jz_diff df_jz[value].diff().dropna() result_diff adfuller(df_jz_diff) print(f一阶差分后p值: {result_diff[1]:.4f})ADF检验的p值小于0.05表示拒绝原假设序列是平稳的。如果原始序列不平稳先做一阶差分如果数据有明显的季节波动——比如每个月呈现相同走势——往往还需要季节差分。SARIMA的数学形式是SARIMA(p,d,q)(P,D,Q)[S]这里的d是普通差分阶数D是季节差分阶数S是周期长度。焦作这份数据如果做月度观测S一般取12如果是日度数据且存在周规律S取7。2.3 数据划分训练集与测试集的切分策略# 按时间顺序切分不能随机打乱 train_size int(len(df_jz) * 0.8) train_data df_jz.iloc[:train_size] test_data df_jz.iloc[train_size:] print(f训练集: {train_data.index[0]} 至 {train_data.index[-1]}, 共{len(train_data)}条) print(f测试集: {test_data.index[0]} 至 {test_data.index[-1]}, 共{len(test_data)}条)时间序列的切分跟普通机器学习不一样必须保留时间顺序用前80%的数据做训练后20%做测试。这份源码里的参数化设计可以让切分比例直接通过变量调整我一般习惯把train_size写成int(len(data) * 0.8)而不是写死一个数字这样换数据时不用改两处。3. SARIMA模型构建从ACF/PACF到statsmodels拟合这一章是源码的核心价值所在。SARIMA的难点不在调用库而在确定参数组合(p,d,q)(P,D,Q)[S]。源码的注释几乎一行一解释对理解模型内部逻辑帮助很大。3.1 自相关图与偏自相关图初步确定p和qimport matplotlib.pyplot as plt from statsmodels.graphics.tsaplots import plot_acf, plot_pacf # 绘制差分后序列的ACF和PACF图 fig, axes plt.subplots(1, 2, figsize(14, 4)) plot_acf(df_jz_diff, axaxes[0], lags24, titleACF - 差分后序列) plot_pacf(df_jz_diff, axaxes[1], lags24, titlePACF - 差分后序列) plt.tight_layout() plt.savefig(acf_pacf.png, dpi150) plt.show()ACF图显示自相关系数随滞后阶数的衰减模式PACF图显示偏自相关系数。经验法则是如果ACF拖尾、PACF截尾p由PACF的显著滞后阶数决定如果ACF截尾、PACF拖尾q由ACF的显著滞后阶数决定。但实际数据往往不这么干净我见过很多人对着图纠结半天。源码的做法是把这两个图先画出来然后用AIC/BIC在候选参数里筛选这样从图形得到的只是初值范围p∈[0,2]、q∈[0,2]最终以信息准则为准。3.2 季节性参数S的确定月度数据与周周期# 数据重采样观察季节性 monthly_mean df_jz[value].resample(M).mean() monthly_mean.plot(markero, figsize(10, 3)) plt.title(月度均值走势 - 观察季节性) plt.savefig(monthly_pattern.png, dpi150) plt.show() # 周期S月度数据通常S12日数据观察是否有7天周期 # SARIMA(p,d,q)(P,D,Q)[S] 中的S就在这里确定 S 12 # 月度周期确定S的方式有两种一是业务知识比如气温数据有明确的年度周期S12二是画重采样图看波动规律。源码里将S作为参数单独定义这是参数化编程的典型做法——所有可能变化的量全部提到文件头部只需要改一处就能适配不同的数据周期。焦作全.csv如果跨度超过两年你会看到明显的年度重复形态。3.3 SARIMAX模型拟合与预测核心代码拆解from statsmodels.tsa.statespace.sarimax import SARIMAX import warnings warnings.filterwarnings(ignore) # 模型参数 - 修改这里即可适配不同数据 p, d, q 1, 1, 1 # 非季节性参数 P, D, Q 1, 1, 1, 12 # 季节性参数最后一位是S # 注意这里是P, D, Q, S四个值 model SARIMAX( train_data[value], order(p, d, q), seasonal_order(P, D, Q, S), enforce_stationarityFalse, enforce_invertibilityFalse ) # 拟合模型 - disp0关闭迭代信息输出 model_fit model.fit(disp0) print(model_fit.summary())order(p,d,q)对应非季节性的AR、差分、MA阶数seasonal_order(P,D,Q,S)对应季节性部分。enforce_stationarityFalse和enforce_invertibilityFalse是给优化器松绑避免参数估计过程中因为边界约束导致收敛失败。源码这里注释得非常细几乎每一行都有解释对入门者来说这些参数的含义比调用本身更难理解。拟合完成后打印summary重点看P|z|这一列如果某个参数的p值很大说明该项可能不显著可以在下一轮调参中把它降为0。预测部分# 训练集内预测拟合值 in_sample_pred model_fit.predict(starttrain_data.index[0], endtrain_data.index[-1]) # 测试集外预测 - 动态预测 n_test len(test_data) forecast model_fit.forecast(stepsn_test) # 对比预测与真实值 plt.figure(figsize(12, 5)) plt.plot(train_data.index, train_data[value], label训练集真实值, color#2C3E50) plt.plot(test_data.index, test_data[value], label测试集真实值, color#27AE60) plt.plot(test_data.index, forecast, labelSARIMA预测值, color#E74C3C, linestyle--) plt.legend() plt.title(SARIMA模型预测效果对比) plt.savefig(forecast_compare.png, dpi150) plt.show()这里有个关键区别predict用于样本内拟合forecast用于样本外预测。新手最容易混的就是这两个方法。源码里先调用predict看训练集拟合效果——如果拟合值偏离严重说明模型结构有问题——再调用forecast产生真正的未来预测。测试集上的对比图直观展示了模型效果预测值如果紧贴真实值曲线这份大作业基本就能过关了。4. 模型评估与残差诊断不只是看RMSE完成了拟合和预测只是第一步模型好不好还要靠评估指标和残差分析说话。源码里同时实现了定量指标和图形化诊断这是很多课程设计里容易忽略但老师爱挑刺的部分。4.1 定量评估MAE、RMSE与MAPEfrom sklearn.metrics import mean_absolute_error, mean_squared_error # 计算评估指标 mae mean_absolute_error(test_data[value], forecast) rmse np.sqrt(mean_squared_error(test_data[value], forecast)) # MAPE - 平均绝对百分比误差 mape np.mean(np.abs((test_data[value] - forecast) / test_data[value])) * 100 print(fMAE: {mae:.4f}) print(fRMSE: {rmse:.4f}) print(fMAPE: {mape:.2f}%)RMSE对大的预测误差更敏感MAE反映平均绝对偏差MAPE则以百分比形式呈现误差幅度。实际写报告时这三个指标放一个表里就够了。我发现MAPE有一个坑当真实值接近0时MAPE会变得异常大甚至无穷如果数据里有接近零的观测值建议改用对称MAPE或直接报告MAE和RMSE。4.2 残差白噪声检验Ljung-Box与正态性from statsmodels.stats.diagnostic import acorr_ljungbox # 计算残差 residuals model_fit.resid # Ljung-Box检验p值0.05说明残差是白噪声 lb_test acorr_ljungbox(residuals, lags[12, 24], return_dfTrue) print(lb_test) # 残差QQ图 - 观察是否接近正态分布 from scipy import stats stats.probplot(residuals, distnorm, plotplt) plt.title(残差Q-Q图) plt.savefig(residual_qq.png, dpi150) plt.show()残差诊断的逻辑是好的模型已经把时间序列中的规律全部提取干净剩下的残差应该是白噪声。Ljung-Box检验的p值大于0.05表示没有自相关性符合白噪声假设Q-Q图上的点近似落在直线上说明残差接近正态分布。如果残差还有明显的自相关结构说明模型阶数不够或者季节性部分没提干净需要回头调参。源码里这段的注释详细到每行都有说明照着敲一遍基本就懂了这个检验流程。4.3 AIC/BIC对比在候选模型里选最优# 候选参数组合 candidates [ (1, 1, 1, 1, 1, 1, 12), (0, 1, 1, 1, 1, 1, 12), (1, 1, 0, 1, 1, 1, 12), (1, 1, 1, 0, 1, 1, 12), ] results [] for params in candidates: p, d, q, P, D, Q, S params model SARIMAX( train_data[value], order(p, d, q), seasonal_order(P, D, Q, S), enforce_stationarityFalse, enforce_invertibilityFalse ) fit model.fit(disp0) results.append({ params: f({p},{d},{q})({P},{D},{Q})[{S}], AIC: fit.aic, BIC: fit.bic }) # 打印对比结果 for r in results: print(f参数: {r[params]}, AIC: {r[AIC]:.2f}, BIC: {r[BIC]:.2f})AIC和BIC都是越小越好但BIC对参数数量的惩罚更重。如果两个模型的AIC接近我倾向于选参数更少的那个如果差异超过10才认为复杂模型有显著优势。源码把这组候选循环写在函数里换数据集时只需要扩展candidates列表。这个自动化搜索思路可以极大提高调参效率不用手动一个个试。5. SARIMA实战避坑指南五个高频踩坑点这一章写的是我拆这份源码时遇到的实际问题以及帮别人调试 SARIMA 模型时反复出现的翻车现场每一条都是真金白银的教训。5.1 日期索引不连续导致预测错位现象训练集效果不错但测试集预测曲线整体向右平移了一截看着像“延迟了一个周期”。 原因原始数据有缺失日期pandas的索引不连续forecast(stepsn)按索引位置生成未来点但画图时索引对不上看起来就像预测滞后。 解决建模前先重采样补齐缺失日期# 将日期索引补齐为连续日序列 df_full df_jz[value].asfreq(D) df_full df_full.interpolate(methodlinear)这里asfreq(D)把索引重排为连续的日序列缺失值用线性插值填充。对月度数据用asfreq(M)对小时数据用asfreq(H)。血泪经验时间序列必须先看索引是否连续再谈建模否则后面所有诊断图都可能误导。5.2 ADF检验与差分次数过度的矛盾现象ADF检验显示一阶差分后p值仍大于0.05于是继续做二阶差分模型预测结果却不如一阶差分模型。 原因ADF检验对样本长度敏感短序列检验功效不足过度差分会损失信息量导致预测方差增大。 解决一阶差分后如果p值在0.05到0.15之间先用肉眼观察差分序列是否有明显的趋势残留。如果没有可见趋势就用一阶差分继续建模不要机械地追求p值小于0.05。我一般同时观察差分序列的均值和方差是否稳定结合业务判断而不是只信统计检验。5.3 网格搜索参数时忽略P和D现象用itertools.product跑参数搜索AIC最低的组合跑出来预测却是直线。 原因候选组合里P和D全部为0模型退化为SARIMA的普通ARIMA自然无法捕捉季节性。 解决至少保留季节性部分的一个非零分量。如果数据有明显季节性(P,D,Q)三个值不能全为0。更稳妥的做法是把季节差分D固定为1只搜索P和Q的取值组合。5.4 预测结果方差偏小但均值偏移现象预测曲线很平滑甚至接近直线整体均值比真实值低了一截RMSE看着还行但图很难看。 原因forecast(stepsn)是点预测只会给出条件期望当模型被过度差分化或季节性参数设置错误时预测趋于平均水平无法跟随波动幅度。 解决先检查model_fit.resid中是否有明显的周期分量残留如果残差ACF图在滞后12处显著说明季节项没提干净。另一个技巧是用get_forecast获取置信区间把区间画出来观察不确定性范围。# 带置信区间的预测 forecast_result model_fit.get_forecast(stepsn_test) forecast_mean forecast_result.predicted_mean conf_int forecast_result.conf_int(alpha0.05)conf_int(alpha0.05)返回95%置信区间的上下界。画图时把上下界填充在预测曲线周围比单独画一条预测线信息量大得多报告中加分效果好。5.5 中文路径和中文列名的编码问题现象pd.read_csv(焦作.csv)报错UnicodeDecodeError或者读取成功但列名显示乱码。 原因Windows下CSV文件常以GBK编码保存pandas默认按UTF-8解析。 解决显式指定编码参数df_jz pd.read_csv(焦作.csv, encodinggbk, parse_dates[date], index_coldate)如果encodinggbk还是失败就再加一个encodinggb18030这个编码覆盖字符范围更广。Excel文件一般不需要处理编码问题但列名如果带中文建议在加载后统一重命名成英文避免后续画图时字体警告。6. 滚动预测与模型更新把SARIMA用出生产感的技巧模型调试到测试集上指标不错不代表这份源码就只能交作业。最后一个实用技巧是滚动预测——每预测一个点就把真实值并入训练集重新拟合模拟在线更新的场景。# 滚动预测每次预测1步然后更新训练集 history list(train_data[value].values) test_values list(test_data[value].values) rolling_pred [] for t in range(len(test_values)): # 用当前history重新拟合SARIMA model SARIMAX( history, order(1, 1, 1), seasonal_order(1, 1, 1, 12), enforce_stationarityFalse, enforce_invertibilityFalse ) fit model.fit(dispFalse) # 预测1步 yhat fit.forecast(steps1)[0] rolling_pred.append(yhat) # 把真实值并入history作为下一步的输入 history.append(test_values[t]) # 计算滚动预测的RMSE rolling_rmse np.sqrt(mean_squared_error(test_values, rolling_pred)) print(f滚动预测RMSE: {rolling_rmse:.4f}) # 对比一次性预测与滚动预测 plt.figure(figsize(12, 4)) plt.plot(test_data.index, test_values, label真实值, color#27AE60) plt.plot(test_data.index, forecast, label一次性预测, color#E74C3C, linestyle--) plt.plot(test_data.index, rolling_pred, label滚动预测, color#2C3E50, linestyle:) plt.legend() plt.savefig(rolling_forecast.png, dpi150) plt.show()滚动预测的核心思想是“预测完一步就要吸取真实反馈”。这种做法在业务上更贴近真实使用场景——预测明天的量到了明天看真实值再预测后天。一次性的多步预测误差会随步长累积滚动预测每一步都基于最新的真实数据所以RMSE通常更低。需要注意滚动预测的代价是计算量加大每步都要重新拟合对大样本数据来说速度会明显变慢。如果追求效率可以改为每5步或每10步重拟合一次用最近的真实值更新历史窗口牺牲一点精度换取速度。从那以后我每次做时间序列预测都会强制走一遍完整流程先查索引连续性再做ADF检验画ACF/PACF初定参数AIC搜索细化最后一定补一组滚动预测作为对照——这个习惯帮我避免了很多“训练集完美、测试集翻车”的尴尬。希望这份SARIMA源码的拆解能帮到你少走一些我走过的弯路。本文还有配套的精品资源点击获取
返回列表