基于SARIMA模型的新能源汽车产销量分析

创作时间:

作者:

@小白创作中心

基于SARIMA模型的新能源汽车产销量分析

引用

CSDN

https://m.blog.csdn.net/m0_73962294/article/details/138876620

本文通过一个实际案例，介绍了如何使用季节性自回归积分滑动平均模型（SARIMA）对新能源汽车的月产销量进行预测分析。通过2015-2023年的历史数据，预测了2024年的月产销量，展示了SARIMA模型在处理具有季节性变化的时间序列数据中的应用。

SARIMA模型简介

SARIMA模型是ARIMA模型的扩展，专门用于处理具有季节性变化的时间序列数据。它在ARIMA模型的基础上添加了季节性成分，可以更好地捕捉季节性变化的影响。季节性不仅局限于自然季节的周期性变化，任何类型的周期性都可以用SARIMA模型来处理，例如一天内的小时变化（24）、一周内的星期变化（7）、一个月内的日变化（30）、一年内的月变化（12）等。

数据准备

本案例使用2015年1月1日至2024年3月1日的新能源汽车月产销量数据。数据读入后，需要进行一些预处理，包括将时间列转换为周期索引格式。

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt 
import seaborn as sns
from tqdm import tqdm  
from pyecharts.globals import CurrentConfig, NotebookType
CurrentConfig.NOTEBOOK_TYPE = NotebookType.JUPYTER_LAB
import pyecharts.options as opts
from pyecharts.charts import Line
from statsmodels.tsa.statespace.sarimax import SARIMAX
from sklearn.metrics import mean_squared_error

# 读取数据
df1 = pd.read_excel("2014-2024年新能源汽车产销量（月）.xlsx", sheet_name = 0)

# 将时间列转换为日期格式
def date(para):
    if type(para) == int:
        delta = pd.Timedelta(str(int(para))+'days')
        time = pd.to_datetime('1899-12-30') + delta
        return time
    else:
        return para
df1['时间'] = df1['时间'].apply(date)
df1.set_index('时间', inplace=True)
df1.index = pd.DatetimeIndex(df1.index).to_period('M') # 设置周期频率

数据探索

检测空缺值

通过绘制产销量曲线，可以直观地查看数据的完整性和趋势。如果发现数据中存在空缺值，可以考虑使用SARIMA模型进行填补。

# 绘制产销量曲线
x_axis_data = [str(period) for period in df1.index.tolist()]
yields = (
    Line(init_opts=opts.InitOpts(theme=ThemeType.LIGHT))
    .add_xaxis(x_axis_data)
    .add_yaxis(series_name="月产量", 
               is_smooth=True,
               y_axis=df1['产量(万辆)'].tolist(), 
               symbol_size=10)
    .add_yaxis(series_name="月销量", 
               is_smooth=True,
               y_axis=df1['销量(万辆)'].tolist(), 
               symbol_size=10)
    .set_global_opts(title_opts=opts.TitleOpts(title="2015年-2024年3月新能源汽车的月产销量"),
                     xaxis_opts=opts.AxisOpts(type_='category', name="时间"),
                     yaxis_opts=opts.AxisOpts(type_='value', name="产销量/万辆"))
)
yields.load_javascript()
yields.render_notebook()

检测异常值

通过绘制小提琴图，可以查看产销量的分布情况。如果数据分布异常，可能需要进行异常值处理。

custom_palette = ["#ff7f0e", "#2ca02c"]
sns.violinplot(data=df1, inner="quartile", palette=custom_palette)
plt.title('产销量小提琴图')
plt.xlabel('Groups')
plt.ylabel('Values')
plt.show()

检测平稳性

通过绘制每年月销量的趋势图，可以分析数据的平稳性。如果数据存在明显的趋势或季节性变化，需要进行相应的处理。

sales_data = [df1['产量(万辆)'].tolist()[i*12: (i+1)*12] for i in range(10)]
years = [str(i) for i in range(2015, 2025)]
line = Line()
for year, sales in zip(years, sales_data):
    line.add_xaxis(range(1, 13))
    line.add_yaxis(f"{year}", sales)
line.set_global_opts(
    title_opts=opts.TitleOpts(title="月产量"),
    xaxis_opts=opts.AxisOpts(type_="category", name="月份"),
    yaxis_opts=opts.AxisOpts(name="产量(万辆)"),
    tooltip_opts=opts.TooltipOpts(trigger="axis", axis_pointer_type="cross"),
    legend_opts=opts.LegendOpts(pos_left="center", orient="horizontal")
)
line.load_javascript()
line.render_notebook()

SARIMA建模

参数选择

使用网格搜索法选择最优的SARIMA模型参数。这里我们限定p、d、q的参数范围，计算每组参数的均方根误差（RMSE），并选择误差最小的参数组。

train = df1[:108] # 选取完整的时序数据
p_range = range(0, 5)  
d_range = range(0, 3)  
q_range = range(0, 5)  
parameters = itertools.product(p_range, d_range, q_range)
best_score = float('inf')
best_params = None

for param in tqdm(parameters, total=len(p_range)*len(d_range)*len(q_range)):
    try:
        model = SARIMAX(train['产量(万辆)'], order=param, seasonal_order=(1, 1, 1, 12))
        result = model.fit()
        predicted_sales = result.predict(start=train.index[0], end=train.index[-1])
        rmse = np.sqrt(mean_squared_error(train['产量(万辆)'], predicted_sales))
        print("RMSE:", rmse)
        if rmse < best_score:
            best_score = rmse
            best_params = param
    except:
        continue
print("Best Parameters:", best_params)

模型训练与预测

使用选择的最优参数训练SARIMA模型，并对2024年的月产销量进行预测。

model = SARIMAX(train['产量(万辆)'], order=(2, 1, 4), seasonal_order=(1, 1, 1, 12))
result = model.fit()
predicted_sales = result.predict(start=predicted_index[0], end=predicted_index[-1])

# 绘制预测结果
x_axis_data = [str(period) for period in df1.index.tolist()]
predict = (
    Line(init_opts=opts.InitOpts(theme=ThemeType.LIGHT))
    .add_xaxis(x_axis_data)
    .add_yaxis(series_name="月产量", 
               is_smooth=True,
               y_axis=df3['产量(万辆)'].tolist(), 
               symbol_size=10)
    .add_yaxis(series_name="预测产量", 
               is_smooth=True,
               y_axis=predicted_sales.tolist(), 
               symbol_size=10, 
               linestyle_opts=opts.LineStyleOpts(color="red"))
    .set_global_opts(title_opts=opts.TitleOpts(title="2015年-2024年3月新能源汽车的月产销量"),
                     xaxis_opts=opts.AxisOpts(type_='category', name="时间"),
                     yaxis_opts=opts.AxisOpts(type_='value', name="产量/万辆"))
)
predict.load_javascript()
predict.render_notebook()

2024年预测结果

设定预测区间为2024年1月到2024年12月，得到2024年每月新能源汽车产量的预测结果。

predicted_2024 = result.predict(start='2024-01', end='2024-12', freq='M')
predicted_2024 = [round(i, 2) for i in predicted_2024]

x_axis_data = [f'2024-{i}' for i in range(1, 13)]
predict_2024 = (
    Line(init_opts=opts.InitOpts(theme=ThemeType.LIGHT))
    .add_xaxis(x_axis_data)
    .add_yaxis(series_name="预测产量", 
               is_smooth=True,
               y_axis=predicted_2024, 
               symbol_size=10, 
               label_opts=opts.LabelOpts(is_show=True))
    .set_global_opts(title_opts=opts.TitleOpts(title="2024年新能源汽车的月产量"),
                     xaxis_opts=opts.AxisOpts(type_='category', name="时间"),
                     yaxis_opts=opts.AxisOpts(type_='value', name="产量/万辆"))
)
predict_2024.load_javascript()