第85步时间序列建模实战：CNN回归建模

基于WIN10的64位系统演示

一、写在前面

这一期，我们介绍CNN回归。

同样，这里使用这个数据：

《PLoS One》2015年一篇题目为《Comparison of Two Hybrid Models for Forecasting the Incidence of Hemorrhagic Fever with Renal Syndrome in Jiangsu Province, China》文章的公开数据做演示。数据为江苏省2004年1月至2012年12月肾综合症出血热月发病率。运用2004年1月至2011年12月的数据预测2012年12个月的发病率数据。

二、CNN回归

（1）原理

卷积神经网络（CNN）最初是为图像识别和处理而设计的，但它们已经被证明对于各种类型的序列数据，包括时间序列，也是有效的。以下是一些关于CNN在时间序列预测中应用的原理：

（a）局部感受野:

-CNN的关键特点是它的局部感受野，这意味着每个卷积核只查看输入数据的一个小部分。

-对于时间序列，这意味着CNN可以捕获和学习模式中的短期依赖关系或周期性。

-这类似于在时间序列分析中使用滑动窗口来捕获短期模式。

（b）参数共享:

-在CNN中，卷积核的权重在输入的所有部分上都是共享的。

-这意味着网络可以在时间序列的任何位置都识别出相同的模式，增加了其泛化能力。

（c）多尺度特征捕获:

-通过使用多个卷积层和池化层，CNN能够在不同的时间尺度上捕获模式。

-这使得它们能够捕获长期和短期的时间序列依赖关系。

（d）堆叠结构:

多层的CNN结构使得网络可以学习时间序列中的复杂和抽象的模式。例如，第一层可能会捕获简单的趋势或周期性，而更深层的网络可能会捕获更复杂的季节性模式或其他非线性关系。

（e）自动特征学习:

-传统的时间序列分析方法通常需要手动选择和构造特征。

-使用CNN，网络可以自动从原始数据中学习和提取相关特征，这通常导致更好的性能和更少的手工工作。

（f）时间序列的结构化特征:

-和图像数据一样，时间序列数据也具有结构性。例如，过去的观察结果通常影响未来的观察结果。

-CNN利用这种结构性，通过卷积操作从数据中提取局部和全局的时间模式。

总之，虽然CNN最初是为图像设计的，但它们在处理序列数据，特别是时间序列数据时，已经显示出了很强的潜力。这是因为它们可以自动从数据中学习重要的特征，捕获多种尺度的模式，并适应时间序列中的短期和长期依赖关系。

（2）单步滚动预测

import pandas as pd
import numpy as np
from sklearn.metrics import mean_absolute_error, mean_squared_error
from tensorflow.python.keras.models import Sequential
from tensorflow.python.keras.layers import Dense, Conv1D, Flatten, MaxPooling1D
from tensorflow.python.keras.optimizers import adam_v2# 读取数据
data = pd.read_csv('data.csv')# 将时间列转换为日期格式
data['time'] = pd.to_datetime(data['time'], format='%b-%y')# 创建滞后期特征
lag_period = 6
for i in range(lag_period, 0, -1):data[f'lag_{i}'] = data['incidence'].shift(lag_period - i + 1)# 删除包含 NaN 的行
data = data.dropna().reset_index(drop=True)# 划分训练集和验证集
train_data = data[(data['time'] >= '2004-01-01') & (data['time'] <= '2011-12-31')]
validation_data = data[(data['time'] >= '2012-01-01') & (data['time'] <= '2012-12-31')]# 定义特征和目标变量
X_train = train_data[['lag_1', 'lag_2', 'lag_3', 'lag_4', 'lag_5', 'lag_6']].values
y_train = train_data['incidence'].values
X_validation = validation_data[['lag_1', 'lag_2', 'lag_3', 'lag_4', 'lag_5', 'lag_6']].values
y_validation = validation_data['incidence'].values# 对于CNN，我们需要将输入数据重塑为3D格式 [samples, timesteps, features]
X_train = X_train.reshape(X_train.shape[0], X_train.shape[1], 1)
X_validation = X_validation.reshape(X_validation.shape[0], X_validation.shape[1], 1)# 构建CNN模型
model = Sequential()
model.add(Conv1D(filters=64, kernel_size=2, activation='relu', input_shape=(X_train.shape[1], 1)))
model.add(MaxPooling1D(pool_size=2))
model.add(Flatten())
model.add(Dense(50, activation='relu'))
model.add(Dense(1))model.compile(optimizer=adam_v2.Adam(learning_rate=0.001), loss='mse')# 训练模型
history = model.fit(X_train, y_train, epochs=200, batch_size=32, validation_data=(X_validation, y_validation), verbose=0)# 单步滚动预测函数
def rolling_forecast(model, initial_features, n_forecasts):forecasts = []current_features = initial_features.copy()for i in range(n_forecasts):# 使用当前的特征进行预测forecast = model.predict(current_features.reshape(1, len(current_features), 1)).flatten()[0]forecasts.append(forecast)# 更新特征，用新的预测值替换最旧的特征current_features = np.roll(current_features, shift=-1)current_features[-1] = forecastreturn np.array(forecasts)# 使用训练集的最后6个数据点作为初始特征
initial_features = X_train[-1].flatten()# 使用单步滚动预测方法预测验证集
y_validation_pred = rolling_forecast(model, initial_features, len(X_validation))# 计算训练集上的MAE, MAPE, MSE 和 RMSE
mae_train = mean_absolute_error(y_train, model.predict(X_train).flatten())
mape_train = np.mean(np.abs((y_train - model.predict(X_train).flatten()) / y_train))
mse_train = mean_squared_error(y_train, model.predict(X_train).flatten())
rmse_train = np.sqrt(mse_train)# 计算验证集上的MAE, MAPE, MSE 和 RMSE
mae_validation = mean_absolute_error(y_validation, y_validation_pred)
mape_validation = np.mean(np.abs((y_validation - y_validation_pred) / y_validation))
mse_validation = mean_squared_error(y_validation, y_validation_pred)
rmse_validation = np.sqrt(mse_validation)print("验证集：", mae_validation, mape_validation, mse_validation, rmse_validation)
print("训练集：", mae_train, mape_train, mse_train, rmse_train)

看结果：

（3）多步滚动预测-vol. 1

import pandas as pd
import numpy as np
from sklearn.metrics import mean_absolute_error, mean_squared_error
from tensorflow.python.keras.models import Sequential
from tensorflow.python.keras.layers import Dense, Conv1D, Flatten, MaxPooling1D
from tensorflow.python.keras.optimizers import adam_v2# 读取数据
data = pd.read_csv('data.csv')
data['time'] = pd.to_datetime(data['time'], format='%b-%y')n = 6
m = 2# 创建滞后期特征
for i in range(n, 0, -1):data[f'lag_{i}'] = data['incidence'].shift(n - i + 1)data = data.dropna().reset_index(drop=True)train_data = data[(data['time'] >= '2004-01-01') & (data['time'] <= '2011-12-31')]
validation_data = data[(data['time'] >= '2012-01-01') & (data['time'] <= '2012-12-31')]# 准备训练数据
X_train = []
y_train = []for i in range(len(train_data) - n - m + 1):X_train.append(train_data.iloc[i+n-1][[f'lag_{j}' for j in range(1, n+1)]].values)y_train.append(train_data.iloc[i+n:i+n+m]['incidence'].values)X_train = np.array(X_train)
y_train = np.array(y_train)
X_train = X_train.astype(np.float32)
y_train = y_train.astype(np.float32)# 为CNN准备数据
X_train = X_train.reshape(X_train.shape[0], X_train.shape[1], 1)# 构建CNN模型
model = Sequential()
model.add(Conv1D(filters=64, kernel_size=2, activation='relu', input_shape=(X_train.shape[1], 1)))
model.add(MaxPooling1D(pool_size=2))
model.add(Flatten())
model.add(Dense(50, activation='relu'))
model.add(Dense(m))model.compile(optimizer=adam_v2.Adam(learning_rate=0.001), loss='mse')# 训练模型
model.fit(X_train, y_train, epochs=200, batch_size=32, verbose=0)def cnn_rolling_forecast(data, model, n, m):y_pred = []for i in range(len(data) - n):input_data = data.iloc[i+n-1][[f'lag_{j}' for j in range(1, n+1)]].values.astype(np.float32).reshape(1, n, 1)pred = model.predict(input_data)y_pred.extend(pred[0])# Handle overlapping predictions by averagingfor i in range(1, m):for j in range(len(y_pred) - i):y_pred[j+i] = (y_pred[j+i] + y_pred[j]) / 2return np.array(y_pred)# Predict for train_data and validation_data
y_train_pred_cnn = cnn_rolling_forecast(train_data, model, n, m)[:len(y_train)]
y_validation_pred_cnn = cnn_rolling_forecast(validation_data, model, n, m)[:len(validation_data) - n]# Calculate performance metrics for train_data
mae_train = mean_absolute_error(train_data['incidence'].values[n:len(y_train_pred_cnn)+n], y_train_pred_cnn)
mape_train = np.mean(np.abs((train_data['incidence'].values[n:len(y_train_pred_cnn)+n] - y_train_pred_cnn) / train_data['incidence'].values[n:len(y_train_pred_cnn)+n]))
mse_train = mean_squared_error(train_data['incidence'].values[n:len(y_train_pred_cnn)+n], y_train_pred_cnn)
rmse_train = np.sqrt(mse_train)# Calculate performance metrics for validation_data
mae_validation = mean_absolute_error(validation_data['incidence'].values[n:len(y_validation_pred_cnn)+n], y_validation_pred_cnn)
mape_validation = np.mean(np.abs((validation_data['incidence'].values[n:len(y_validation_pred_cnn)+n] - y_validation_pred_cnn) / validation_data['incidence'].values[n:len(y_validation_pred_cnn)+n]))
mse_validation = mean_squared_error(validation_data['incidence'].values[n:len(y_validation_pred_cnn)+n], y_validation_pred_cnn)
rmse_validation = np.sqrt(mse_validation)print("训练集：", mae_train, mape_train, mse_train, rmse_train)
print("验证集：", mae_validation, mape_validation, mse_validation, rmse_validation)

结果：

（4）多步滚动预测-vol. 2

import pandas as pd
import numpy as np
from sklearn.model_selection import train_test_split
from sklearn.metrics import mean_absolute_error, mean_squared_error
from tensorflow.python.keras.models import Sequential
from tensorflow.python.keras.layers import Dense, Conv1D, Flatten, MaxPooling1D
from tensorflow.python.keras.optimizers import adam_v2# Loading and preprocessing the data
data = pd.read_csv('data.csv')
data['time'] = pd.to_datetime(data['time'], format='%b-%y')n = 6  # 使用前6个数据点
m = 2  # 预测接下来的2个数据点# 创建滞后期特征
for i in range(n, 0, -1):data[f'lag_{i}'] = data['incidence'].shift(n - i + 1)data = data.dropna().reset_index(drop=True)train_data = data[(data['time'] >= '2004-01-01') & (data['time'] <= '2011-12-31')]
validation_data = data[(data['time'] >= '2012-01-01') & (data['time'] <= '2012-12-31')]# 只对X_train、y_train、X_validation取奇数行
X_train = train_data[[f'lag_{i}' for i in range(1, n+1)]].iloc[::2].reset_index(drop=True).values
X_train = X_train.reshape(X_train.shape[0], X_train.shape[1], 1)  # Reshape for CNN# 创建m个目标变量
y_train_list = [train_data['incidence'].shift(-i) for i in range(m)]
y_train = pd.concat(y_train_list, axis=1)
y_train.columns = [f'target_{i+1}' for i in range(m)]
y_train = y_train.iloc[::2].reset_index(drop=True).dropna().values[:, 0]  # Only take the first column for simplicityX_validation = validation_data[[f'lag_{i}' for i in range(1, n+1)]].iloc[::2].reset_index(drop=True).values
X_validation = X_validation.reshape(X_validation.shape[0], X_validation.shape[1], 1)  # Reshape for CNNy_validation = validation_data['incidence'].values# Building the CNN model
model = Sequential()
model.add(Conv1D(filters=64, kernel_size=2, activation='relu', input_shape=(X_train.shape[1], 1)))
model.add(MaxPooling1D(pool_size=2))
model.add(Flatten())
model.add(Dense(50, activation='relu'))
model.add(Dense(1))optimizer = adam_v2.Adam(learning_rate=0.001)
model.compile(optimizer=optimizer, loss='mse')# Train the model
model.fit(X_train, y_train, epochs=200, batch_size=32, verbose=0)# Predict on validation set
y_validation_pred = model.predict(X_validation).flatten()# Compute metrics for validation set
mae_validation = mean_absolute_error(y_validation[:len(y_validation_pred)], y_validation_pred)
mape_validation = np.mean(np.abs((y_validation[:len(y_validation_pred)] - y_validation_pred) / y_validation[:len(y_validation_pred)]))
mse_validation = mean_squared_error(y_validation[:len(y_validation_pred)], y_validation_pred)
rmse_validation = np.sqrt(mse_validation)# Predict on training set
y_train_pred = model.predict(X_train).flatten()# Compute metrics for training set
mae_train = mean_absolute_error(y_train, y_train_pred)
mape_train = np.mean(np.abs((y_train - y_train_pred) / y_train))
mse_train = mean_squared_error(y_train, y_train_pred)
rmse_train = np.sqrt(mse_train)print("验证集：", mae_validation, mape_validation, mse_validation, rmse_validation)
print("训练集：", mae_train, mape_train, mse_train, rmse_train)

结果：

（5）多步滚动预测-vol. 3

import pandas as pd
import numpy as np
from sklearn.metrics import mean_absolute_error, mean_squared_error
from tensorflow.python.keras.models import Sequential
from tensorflow.python.keras.layers import Dense, Conv1D, Flatten, MaxPooling1D
from tensorflow.python.keras.optimizers import adam_v2# 数据读取和预处理
data = pd.read_csv('data.csv')
data_y = pd.read_csv('data.csv')
data['time'] = pd.to_datetime(data['time'], format='%b-%y')
data_y['time'] = pd.to_datetime(data_y['time'], format='%b-%y')n = 6for i in range(n, 0, -1):data[f'lag_{i}'] = data['incidence'].shift(n - i + 1)data = data.dropna().reset_index(drop=True)
train_data = data[(data['time'] >= '2004-01-01') & (data['time'] <= '2011-12-31')]
X_train = train_data[[f'lag_{i}' for i in range(1, n+1)]]
m = 3X_train_list = []
y_train_list = []for i in range(m):X_temp = X_trainy_temp = data_y['incidence'].iloc[n + i:len(data_y) - m + 1 + i]X_train_list.append(X_temp)y_train_list.append(y_temp)for i in range(m):X_train_list[i] = X_train_list[i].iloc[:-(m-1)].valuesX_train_list[i] = X_train_list[i].reshape(X_train_list[i].shape[0], X_train_list[i].shape[1], 1)  # Reshape for CNNy_train_list[i] = y_train_list[i].iloc[:len(X_train_list[i])].values# 模型训练
models = []
for i in range(m):# Build CNN modelmodel = Sequential()model.add(Conv1D(filters=64, kernel_size=2, activation='relu', input_shape=(X_train_list[i].shape[1], 1)))model.add(MaxPooling1D(pool_size=2))model.add(Flatten())model.add(Dense(50, activation='relu'))model.add(Dense(1))optimizer = adam_v2.Adam(learning_rate=0.001)model.compile(optimizer=optimizer, loss='mse')model.fit(X_train_list[i], y_train_list[i], epochs=200, batch_size=32, verbose=0)models.append(model)validation_start_time = train_data['time'].iloc[-1] + pd.DateOffset(months=1)
validation_data = data[data['time'] >= validation_start_time]
X_validation = validation_data[[f'lag_{i}' for i in range(1, n+1)]].values
X_validation = X_validation.reshape(X_validation.shape[0], X_validation.shape[1], 1)  # Reshape for CNNy_validation_pred_list = [model.predict(X_validation) for model in models]
y_train_pred_list = [model.predict(X_train_list[i]) for i, model in enumerate(models)]def concatenate_predictions(pred_list):concatenated = []for j in range(len(pred_list[0])):for i in range(m):concatenated.append(pred_list[i][j])return concatenatedy_validation_pred = np.array(concatenate_predictions(y_validation_pred_list))[:len(validation_data['incidence'])]
y_train_pred = np.array(concatenate_predictions(y_train_pred_list))[:len(train_data['incidence']) - m + 1]
y_validation_pred = y_validation_pred.flatten()
y_train_pred = y_train_pred.flatten()mae_validation = mean_absolute_error(validation_data['incidence'], y_validation_pred)
mape_validation = np.mean(np.abs((validation_data['incidence'] - y_validation_pred) / validation_data['incidence']))
mse_validation = mean_squared_error(validation_data['incidence'], y_validation_pred)
rmse_validation = np.sqrt(mse_validation)mae_train = mean_absolute_error(train_data['incidence'][:-(m-1)], y_train_pred)
mape_train = np.mean(np.abs((train_data['incidence'][:-(m-1)] - y_train_pred) / train_data['incidence'][:-(m-1)]))
mse_train = mean_squared_error(train_data['incidence'][:-(m-1)], y_train_pred)
rmse_train = np.sqrt(mse_train)print("验证集：", mae_validation, mape_validation, mse_validation, rmse_validation)
print("训练集：", mae_train, mape_train, mse_train, rmse_train)

结果：