基于信号分解和多种深度学习结合的上证指数预测模型

大家好，我是带我去滑雪！

为了给投资者提供更准确的投资建议、帮助政府和监管部门更好地制定相关政策，维护市场稳定，本文对股民情绪和上证指数之间的关系进行更深入的研究，并结合信号分解、优化算法和深度学习对上证指数进行预测，以期更好地理解股市的运行规律，为股民提供一定的参考。

选取了2014年12月17日至2024年4月17日的共计2269个交易日(覆盖牛市、熊市及调整期多个市场周期)的上证指数数据，建立了基于信号分解和多种深度学习结合的预测模型。在建模过程中，针对变分模态分解(VMD)参数难以选择的问题，我们利用鲸鱼捕食者算法(MPA)寻优能力强的特点自适应选取VMD的关键参数，使VMD分解效果最大化；针对预测精度问题，我们使用LSTM有效传递和表达经MPA-VMD分解后的上证指数信息，使用双向门控循环单元(BiGRU)模型通过正向和反向传播进行信息的交互和整合，添加注意力机制(Attention)为特征提供更具针对性的权重分配，最大程度上提高预测模型的效率和性能。

最终建立了基于上证指数的VMD-LSTM-BiGRU-Attention预测模型，该模型的拟合优度达到了0.98，而均方根误差(RMSE)、平均绝对误差(MAE)和平均绝对百分比误差(MAPE)也明显优于本文用于比较的其他模型，说明本文模型可为上证指数提供更加精准的预测。不仅如此，在相同VMD-LSTM-BiGRU-Attention模型的条件下，加入情绪指数与上证指数结合进行预测，比单独使用上证指数进行预测的准确率更高。这很好的体现了本文研究股民情绪指数和上证指数之间的关联性是有实际意义的，实现了基于股民情绪与上证指数混合模型的股票走向预测且性能优秀，可用于上证指数的预测分析，对股民及相关研究领域有一定的参考价值和指导作用。

下面开始代码实战：

（1）数据展示

（2）MPA算法优化的VMD分解

（3）预测模型

（4）模型评价指标

（1）数据展示

选取东方财富网2014年12月17日至2024年4月17日的上证指数数据，共计2269个样本，将其划分为训练集和测试集，其中训练集为数据前80%，共计1815个，测试集为数据的后20%，共计454个。利用python将上证指数数据可视化，如图16所示。其中黄色实线为训练集数据，粉色虚线为测试集数据。

（2）MPA算法优化的VMD分解

import numpy as np
import pandas as pd
from matplotlib import pyplot as plt
import MPA
from scipy.signal import hilbert
from vmdpy import VMD
f= pd.read_csv('E:\工作\硕士\上证指数.csv',encoding="ANSI")
f=f.收盘
# VMD参数
tau = 0.            # noise-tolerance (no strict fidelity enforcement)               # 3 modes  
DC = 0             # no DC part imposed  
init = 1           # initialize omegas uniformly  
tol = 1e-7 
#MPA设置参数
num_particles = 30 #种群数量
MaxIter = 1000 #最大迭代次数
dim = 10 #维度
lb = -10*np.ones([dim, 1]) #下边界
ub = 10*np.ones([dim, 1])#上边界

# 计算每个IMF分量的包络熵
import numpy as np
from scipy.signal import hilbert
from scipy.stats import kurtosisdef calculate_entropy(imf):# 计算希尔伯特变换的包络线env = np.abs(hilbert(imf))# 归一化包络线env_norm = env / np.max(env)# 计算归一化包络线的概率分布p = env_norm / np.sum(env_norm)# 计算包络熵entropy = -np.sum(p * np.log2(p))# 计算峭度kurt = kurtosis(imf)# 返回包络熵和峭度的和return entropy + kurt# 定义适应度函数，即最大包络熵
def fitness_func(x):if x[1] < 0:x[1] = np.random.choice([3, 12])u, u_hat, omega = VMD(f, int(round(x[0])), tau, int(round(x[1])), DC, init, tol)num_modes = u.shape[0]entropy = np.zeros(num_modes)for i in range(num_modes):entropy[i] = calculate_entropy(u[i,:])# 找到最小的包络熵对应的模态min_entropy_index = np.argmin(entropy)min_entropy_mode = u[min_entropy_index]#     print("最小包络熵对应的模态：", min_entropy_index)# x为VMD参数向量# signal为要分解的信号# 分解信号并计算最大包络熵# 返回最大包络熵值return entropy[min_entropy_index]def mpa_optimization(num_particles, dim,fitness_func, MaxIter):# 初始化海洋捕食者位置和速度predators_pos = np.zeros((num_particles, num_dimensions))for i in range(num_particles):predators_pos[i, 0] = np.random.uniform(500, 3000)predators_pos[i, 1] = np.random.uniform(3, 12)predators_vel = np.zeros((num_particles, num_dimensions))# 记录每个海洋捕食者的最佳位置和适应度值predators_best_pos = np.copy(predators_pos)predators_best_fit = np.zeros(num_particles)# 记录整个群体的最佳位置和适应度值global_best_pos = np.zeros(num_dimensions)global_best_fit = float('inf')# 迭代更新for i in range(max_iter):# 计算每个海洋捕食者的适应度值predators_fitness = np.array([fitness_func(p) for p in predators_pos])# 更新每个海洋捕食者的最佳位置和适应度值for j in range(num_particles):if predators_fitness[j] < predators_best_fit[j]:predators_best_fit[j] = predators_fitness[j]predators_best_pos[j] = np.copy(predators_pos[j])# 更新整个群体的最佳位置和适应度值global_best_idx = np.argmin(predators_fitness)if predators_fitness[global_best_idx] < global_best_fit:global_best_fit = predators_fitness[global_best_idx]global_best_pos = np.copy(predators_pos[global_best_idx])# 更新每个海洋捕食者的速度和位置for j in range(num_particles):# 计算新速度r1 = np.random.rand(num_dimensions)r2 = np.random.rand(num_dimensions)cognitive_vel = 2.0 * r1 * (predators_best_pos[j] - predators_pos[j])social_vel = 2.0 * r2 * (global_best_pos - predators_pos[j])predators_vel[j] = predators_vel[j] + cognitive_vel + social_vel# 更新位置predators_pos[j] = predators_pos[j] + predators_vel[j]# 记录每次迭代的global_best_pos和global_best_fitglobal_best_pos_list.append(global_best_pos)global_best_fit_list.append(global_best_fit)print("第：" + str(i) + '次迭代')# 返回全局最优# 返回全局最优位置和适应度值return global_best_pos, global_best_fit

# 初始化空列表用于存储每次迭代的global_best_pos和global_best_fit
global_best_pos_list = []
global_best_fit_list = []
# 使用PSO算法优化VMD参数
num_particles = 2
num_dimensions = 2 # 假设有2个VMD参数
max_iter =30
best_pos, best_fit = mpa_optimization(num_particles, dim,fitness_func, MaxIter=max_iter)# 输出结果
print("Best VMD parameters:", best_pos)
print("Best fitness value:", best_fit)

得到的最优参数组合和适应度变化曲线如下：

得到的分解序列和中心模态分别为如下所示：

在分解结果中，IMF1显示了一种相对平滑的趋势，可能代表了原始序列中的最低频成分或趋势项。而IMF2显示了略微增加的频率和振幅，但依旧比较平滑，可能捕捉到了次低频的周期性变化。IMF3-IMF7则展示了随着序号增加频率逐渐增高的成分。随着IMF序号的增加，我们可以看到更高频的振荡，这表明它们捕获了原始上证指数序列中的更快变化部分。可认为经MPA优化后的VMD在处理非线性和非平稳信号时具有较好的性能，它能够提取出不同频率范围的振动模态，并具有一定的抗噪能力和平滑性，适用于本文上证指数的分解。

（3）预测模型

根据以上分解过程，VMD将原始上证指数序列分解成多个代表不同频率成分的模态。通过对这些模态分别进行预测，最终可以实现上证指数的预测。本文之后就将基于VMD分解和深度学习结合的方法对上证指数进行预测。为了选取最适合本文上证指数预测的VMD-深度学习模型，在这里我们做了大量实验。首先将上证指数数据集以8：2的比例划分训练集和测试集，用训练集分别训练LSTM单一模型、VMD-LSTM、VMD-LSTM-BiGRU、VMD-LSTM-BiGRU-Attention。碍于篇幅原因，这里以VMD-LSTM-BiGRU-Attention模型为例，展示各个子模态的预测情况如下所示：

通过上面可以发现，VMD-LSTM-BiGRU-Attention模型的各个子模态的预测效果均较好，预测模型能够较好地跟踪和预测多个频率层面上的序列数据。经训练后分别使用了LSTM单一模型、VMD-LSTM、VMD-LSTM-BiGRU、VMD-LSTM-BiGRU-Attention对测试集进行拟合预测，并与实际值进行对比，绘制了不同模型预测效果对比图如下图所示：

部分代码：

from keras import backend as K
from keras.layers import Layerclass Embedding(Layer):def __init__(self, vocab_size, model_dim, **kwargs):self._vocab_size = vocab_sizeself._model_dim = model_dimsuper(Embedding, self).__init__(**kwargs)def build(self, input_shape):self.embeddings = self.add_weight(shape=(self._vocab_size, self._model_dim),initializer='glorot_uniform',name="embeddings")super(Embedding, self).build(input_shape)def call(self, inputs):if K.dtype(inputs) != 'int32':inputs = K.cast(inputs, 'int32')embeddings = K.gather(self.embeddings, inputs)embeddings *= self._model_dim ** 0.5 # Scalereturn embeddingsdef compute_output_shape(self, input_shape):return input_shape + (self._model_dim,)class PositionEncoding(Layer):def __init__(self, model_dim, **kwargs):self._model_dim = model_dimsuper(PositionEncoding, self).__init__(**kwargs)def call(self, inputs):seq_length = inputs.shape[1]position_encodings = np.zeros((seq_length, self._model_dim))for pos in range(seq_length):for i in range(self._model_dim):position_encodings[pos, i] = pos / np.power(10000, (i-i%2) / self._model_dim)position_encodings[:, 0::2] = np.sin(position_encodings[:, 0::2]) # 2iposition_encodings[:, 1::2] = np.cos(position_encodings[:, 1::2]) # 2i+1position_encodings = K.cast(position_encodings, 'float32')return position_encodingsdef compute_output_shape(self, input_shape):return input_shape
class Add(Layer):def __init__(self, **kwargs):super(Add, self).__init__(**kwargs)def call(self, inputs):input_a, input_b = inputsreturn input_a + input_bdef compute_output_shape(self, input_shape):return input_shape[0]class ScaledDotProductAttention(Layer):def __init__(self, masking=True, future=False, dropout_rate=0., **kwargs):self._masking = maskingself._future = futureself._dropout_rate = dropout_rateself._masking_num = -2**32+1super(ScaledDotProductAttention, self).__init__(**kwargs)def mask(self, inputs, masks):masks = K.cast(masks, 'float32')masks = K.tile(masks, [K.shape(inputs)[0] // K.shape(masks)[0], 1])masks = K.expand_dims(masks, 1)outputs = inputs + masks * self._masking_numreturn outputsdef future_mask(self, inputs):diag_vals = tf.ones_like(inputs[0, :, :])tril = tf.linalg.LinearOperatorLowerTriangular(diag_vals).to_dense()  future_masks = tf.tile(tf.expand_dims(tril, 0), [tf.shape(inputs)[0], 1, 1])paddings = tf.ones_like(future_masks) * self._masking_numoutputs = tf.where(tf.equal(future_masks, 0), paddings, inputs)return outputsdef call(self, inputs):if self._masking:assert len(inputs) == 4, "inputs should be set [queries, keys, values, masks]."queries, keys, values, masks = inputselse:assert len(inputs) == 3, "inputs should be set [queries, keys, values]."queries, keys, values = inputsif K.dtype(queries) != 'float32':  queries = K.cast(queries, 'float32')if K.dtype(keys) != 'float32':  keys = K.cast(keys, 'float32')if K.dtype(values) != 'float32':  values = K.cast(values, 'float32')matmul = K.batch_dot(queries, tf.transpose(keys, [0, 2, 1])) # MatMulscaled_matmul = matmul / int(queries.shape[-1]) ** 0.5  # Scaleif self._masking:scaled_matmul = self.mask(scaled_matmul, masks) # Mask(opt.)if self._future:scaled_matmul = self.future_mask(scaled_matmul)softmax_out = K.softmax(scaled_matmul) # SoftMax# Dropoutout = K.dropout(softmax_out, self._dropout_rate)outputs = K.batch_dot(out, values)return outputsdef compute_output_shape(self, input_shape):return input_shapeclass MultiHeadAttention(Layer):def __init__(self, n_heads, head_dim, dropout_rate=.1, masking=True, future=False, trainable=True, **kwargs):self._n_heads = n_headsself._head_dim = head_dimself._dropout_rate = dropout_rateself._masking = maskingself._future = futureself._trainable = trainablesuper(MultiHeadAttention, self).__init__(**kwargs)def build(self, input_shape):self._weights_queries = self.add_weight(shape=(input_shape[0][-1], self._n_heads * self._head_dim),initializer='glorot_uniform',trainable=self._trainable,name='weights_queries')self._weights_keys = self.add_weight(shape=(input_shape[1][-1], self._n_heads * self._head_dim),initializer='glorot_uniform',trainable=self._trainable,name='weights_keys')self._weights_values = self.add_weight(shape=(input_shape[2][-1], self._n_heads * self._head_dim),initializer='glorot_uniform',trainable=self._trainable,name='weights_values')super(MultiHeadAttention, self).build(input_shape)def call(self, inputs):if self._masking:assert len(inputs) == 4, "inputs should be set [queries, keys, values, masks]."queries, keys, values, masks = inputselse:assert len(inputs) == 3, "inputs should be set [queries, keys, values]."queries, keys, values = inputsqueries_linear = K.dot(queries, self._weights_queries) keys_linear = K.dot(keys, self._weights_keys)values_linear = K.dot(values, self._weights_values)queries_multi_heads = tf.concat(tf.split(queries_linear, self._n_heads, axis=2), axis=0)keys_multi_heads = tf.concat(tf.split(keys_linear, self._n_heads, axis=2), axis=0)values_multi_heads = tf.concat(tf.split(values_linear, self._n_heads, axis=2), axis=0)if self._masking:att_inputs = [queries_multi_heads, keys_multi_heads, values_multi_heads, masks]else:att_inputs = [queries_multi_heads, keys_multi_heads, values_multi_heads]attention = ScaledDotProductAttention(masking=self._masking, future=self._future, dropout_rate=self._dropout_rate)att_out = attention(att_inputs)outputs = tf.concat(tf.split(att_out, self._n_heads, axis=0), axis=2)return outputsdef compute_output_shape(self, input_shape):return input_shape

def build_model(X_train,mode='LSTM',hidden_dim=[32,16]):set_my_seed()if mode=='RNN':#RNNmodel = Sequential()model.add(SimpleRNN(hidden_dim[0],return_sequences=True, input_shape=(X_train.shape[-2],X_train.shape[-1])))model.add(SimpleRNN(hidden_dim[1]))  model.add(Dense(1))elif mode=='MLP':model = Sequential()model.add(Dense(hidden_dim[0],activation='relu',input_shape=(X_train.shape[-1],)))model.add(Dense(hidden_dim[1],activation='relu'))model.add(Dense(1))elif mode=='LSTM':# LSTMmodel = Sequential()model.add(LSTM(hidden_dim[0],return_sequences=True, input_shape=(X_train.shape[-2],X_train.shape[-1])))model.add(LSTM(hidden_dim[1]))model.add(Dense(1))elif mode=='GRU':#GRUmodel = Sequential()model.add(GRU(hidden_dim[0],return_sequences=True, input_shape=(X_train.shape[-2],X_train.shape[-1])))model.add(GRU(hidden_dim[1]))model.add(Dense(1))elif mode=='CNN':#一维卷积model = Sequential()model.add(Conv1D(hidden_dim[0],17,activation='relu',input_shape=(X_train.shape[-2],X_train.shape[-1])))model.add(GlobalAveragePooling1D())model.add(Flatten())model.add(Dense(hidden_dim[1],activation='relu'))model.add(Dense(1))elif mode=='CNN+LSTM': model = Sequential()model.add(Conv1D(filters=hidden_dim[0], kernel_size=3, padding="same",activation="relu"))model.add(MaxPooling1D(pool_size=2))model.add(LSTM(hidden_dim[1]))model.add(Dense(1))elif mode=='BiLSTM':model = Sequential()model.add(Bidirectional(LSTM(hidden_dim[0],return_sequences=True, input_shape=(X_train.shape[-2],X_train.shape[-1]))))model.add(Bidirectional(LSTM(hidden_dim[1])))model.add(Dense(1))elif mode=='BiGRU':model = Sequential()model.add(Bidirectional(GRU(hidden_dim[0],return_sequences=True, input_shape=(X_train.shape[-2],X_train.shape[-1]))))model.add(Bidirectional(GRU(hidden_dim[1])))model.add(Dense(1))elif mode=='BiLSTM+Attention':inputs = Input(name='inputs',shape=[X_train.shape[-2],X_train.shape[-1]], dtype='float64')attention_probs = Dense(32, activation='softmax', name='attention_vec')(inputs)attention_mul =  Multiply()([inputs, attention_probs])mlp = Dense(64)(attention_mul) #原始的全连接gru=Bidirectional(LSTM(32))(mlp)mlp = Dense(16,activation='relu')(gru)output = Dense(1)(mlp)model = Model(inputs=[inputs], outputs=output)elif mode=='Attention':inputs = Input(name='inputs',shape=[X_train.shape[-2],X_train.shape[-1]], dtype='float32')attention_probs = Dense(hidden_dim[0], activation='softmax', name='attention_vec')(inputs)attention_mul =  Multiply()([inputs, attention_probs])mlp = Dense(hidden_dim[1])(attention_mul) #原始的全连接fla=Flatten()(mlp)output = Dense(1)(fla)model = Model(inputs=[inputs], outputs=output)  elif mode=='BiGRU+Attention':inputs = Input(name='inputs',shape=[X_train.shape[-2],X_train.shape[-1]], dtype='float64')attention_probs = Dense(32, activation='softmax', name='attention_vec')(inputs)attention_mul =  Multiply()([inputs, attention_probs])mlp = Dense(64)(attention_mul) #原始的全连接gru=Bidirectional(GRU(32))(mlp)mlp = Dense(16,activation='relu')(gru)output = Dense(1)(mlp)model = Model(inputs=[inputs], outputs=output)elif mode=='MultiHeadAttention':         inputs = Input(shape=[X_train.shape[-2],X_train.shape[-1]], name="inputs")#masks = Input(shape=(X_train.shape[-2],), name='masks')encodings = PositionEncoding(X_train.shape[-2])(inputs)encodings = Add()([inputs, encodings])x = MultiHeadAttention(8, hidden_dim[0],masking=False)([encodings, encodings, encodings])x = GlobalAveragePooling1D()(x)x = Dropout(0.2)(x)x = Dense(hidden_dim[1], activation='relu')(x)outputs = Dense(1)(x)model = Model(inputs=[inputs], outputs=outputs)elif mode=='BiGRU+MAttention':    inputs = Input(shape=[X_train.shape[-2],X_train.shape[-1]], name="inputs")encodings = PositionEncoding(X_train.shape[-2])(inputs)encodings = Add()([inputs, encodings])x = MultiHeadAttention(8, hidden_dim[0],masking=False)([encodings, encodings, encodings])
#         x = GlobalAveragePooling1D()(x)x = Bidirectional(GRU(32))(x)x = Dropout(0.2)(x)output = Dense(1)(x)model = Model(inputs=[inputs], outputs=output)#     elif mode=='BiGRU+Attention':
#         inputs = Input(name='inputs',shape=[max_words,], dtype='float64')
#         x = Embedding(top_words, input_length=max_words, output_dim=embed_dim)(inputs)
#         x = Bidirectional(GRU(32,return_sequences=True))(x)
#         x = MultiHeadAttention(2, key_dim=embed_dim)(x,x,x)
#         x = Bidirectional(GRU(32))(x)
#         x = Dropout(0.2)(x)
#         output = Dense(num_labels, activation='softmax')(x)
#         model = Model(inputs=[inputs], outputs=output)model.compile(optimizer='Adam', loss='mse',metrics=[tf.keras.metrics.RootMeanSquaredError(),"mape","mae"])return modeldef plot_loss(hist,imfname):plt.subplots(1,4,figsize=(16,2), dpi=600)for i,key in enumerate(hist.history.keys()):n=int(str('14')+str(i+1))plt.subplot(n)plt.plot(hist.history[key], 'k', label=f'Training {key}')plt.title(f'{imfname} Training {key}')plt.xlabel('Epochs')plt.ylabel(key)plt.legend()plt.tight_layout()plt.show()