Python实现灰狼优化算法（GWO）+SVR的组合预测模型

1.导包

import csv
import numpy as np
from sklearn.preprocessing import StandardScaler
from sklearn.model_selection import train_test_split
import matplotlib.pyplot as plt
from sklearn.metrics import explained_variance_score
from sklearn.svm import SVR 
from sklearn.metrics import mean_squared_error 
from sklearn.metrics import accuracy_score
from sklearn.metrics import confusion_matrix
from sklearn.metrics import classification_report
from sklearn.metrics import explained_variance_score
from sklearn import metrics
from sklearn.metrics import mean_absolute_error # 平方绝对误差
import random
from keras.models import Sequential
from keras.layers import Dense
from keras.layers import LSTM
from datetime import datetime 
from pandas import Series
from pandas import concat
from pandas import read_csv
from pandas import DataFrame

2.导数据

#特征
feature=[]
#目标
target=[]

csv_file = csv.reader(open('东营海水数据.csv'))
for content in csv_file:
    content=list(map(float,content))
    if len(content)!=0:
        feature.append(content[1:4])
        target.append(content[0:1])
        
targets=[]
for i in target:
    targets.append(i[0])

# 标准化转换
scaler = StandardScaler()
# 训练标准化对象
scaler.fit(feature)
# 转换数据集
feature= scaler.transform(feature)

fig = plt.gcf()
fig.set_size_inches(20, 5)
plt.plot(targets)
plt.xlim(0,1350)
plt.xlabel('Serial Number')
plt.ylabel('Dissolved Oxygen')

 3.看效果

看一眼裸SVR效果如何，参数我平时都是，1、0.1、10来设置的。

from sklearn.svm import SVR 
feature_train, feature_test, target_train, target_test = train_test_split(feature, targets, test_size=0.1,random_state=10)

feature_test=feature[int(len(feature)*0.9):int(len(feature))]
target_test=targets[int(len(targets)*0.9):int(len(targets))]

model_svr = SVR(C=1,epsilon=0.1,gamma=10)
model_svr.fit(feature_train,target_train)
predict_results=model_svr.predict(feature_test)

plt.plot(target_test)#测试数组
plt.plot(predict_results)#测试数组
fig = plt.gcf()
fig.set_size_inches(20, 5)
plt.xlim(0,135)
plt.xlabel('Serial Number')
plt.ylabel('Dissolved Oxygen')
SVR1=predict_results

print("MSE:",mean_squared_error(target_test,predict_results))
print("R2 = ",metrics.r2_score(target_test,predict_results)) # R2
print("MAE = ",mean_absolute_error(target_test,predict_results))

4.跑GWO

忘了是哪位大佬的代码了。

import numpy  as np
import math
import random
import os
 
fitness=[]
# Function
def target_function():
    return
 
# Function: Initialize Variables
def initial_position(pack_size = 5, min_values = [-5,-5], max_values = [5,5], target_function = target_function):
    position = np.zeros((pack_size, len(min_values)+1))
    for i in range(0, pack_size):
        for j in range(0, len(min_values)):
             position[i,j] = random.uniform(min_values[j], max_values[j])
        position[i,-1] = target_function(position[i,0:position.shape[1]-1])
    return position
 
# Function: Initialize Alpha
def alpha_position(dimension = 2, target_function = target_function):
    alpha = np.zeros((1, dimension + 1))
    for j in range(0, dimension):
        alpha[0,j] = 0.0
    alpha[0,-1] = target_function(alpha[0,0:alpha.shape[1]-1])
    return alpha
 
# Function: Initialize Beta
def beta_position(dimension = 2, target_function = target_function):
    beta = np.zeros((1, dimension + 1))
    for j in range(0, dimension):
        beta[0,j] = 0.0
    beta[0,-1] = target_function(beta[0,0:beta.shape[1]-1])
    return beta
 
# Function: Initialize Delta
def delta_position(dimension = 2, target_function = target_function):
    delta =  np.zeros((1, dimension + 1))
    for j in range(0, dimension):
        delta[0,j] = 0.0
    delta[0,-1] = target_function(delta[0,0:delta.shape[1]-1])
    return delta
 
# Function: Updtade Pack by Fitness
def update_pack(position, alpha, beta, delta):
    updated_position = np.copy(position)
    for i in range(0, position.shape[0]):
        if (updated_position[i,-1] > alpha[0,-1]):
            alpha[0,:] = np.copy(updated_position[i,:])
        if (updated_position[i,-1] < alpha[0,-1] and updated_position[i,-1] > beta[0,-1]):
            beta[0,:] = np.copy(updated_position[i,:])
        if (updated_position[i,-1] < alpha[0,-1] and updated_position[i,-1] < beta[0,-1]  and updated_position[i,-1] > delta[0,-1]):
            delta[0,:] = np.copy(updated_position[i,:])
    return alpha, beta, delta
 
# Function: Updtade Position
def update_position(position, alpha, beta, delta, a_linear_component = 2, min_values = [-5,-5], max_values = [5,5], target_function = target_function):
    updated_position = np.copy(position)
    for i in range(0, updated_position.shape[0]):
        for j in range (0, len(min_values)):   
            r1_alpha              = int.from_bytes(os.urandom(8), byteorder = "big") / ((1 << 64) - 1)
            r2_alpha              = int.from_bytes(os.urandom(8), byteorder = "big") / ((1 << 64) - 1)
            a_alpha               = 2*a_linear_component*r1_alpha - a_linear_component
            c_alpha               = 2*r2_alpha      
            distance_alpha        = abs(c_alpha*alpha[0,j] - position[i,j]) 
            x1                    = alpha[0,j] - a_alpha*distance_alpha   
            r1_beta               = int.from_bytes(os.urandom(8), byteorder = "big") / ((1 << 64) - 1)
            r2_beta               = int.from_bytes(os.urandom(8), byteorder = "big") / ((1 << 64) - 1)
            a_beta                = 2*a_linear_component*r1_beta - a_linear_component
            c_beta                = 2*r2_beta            
            distance_beta         = abs(c_beta*beta[0,j] - position[i,j]) 
            x2                    = beta[0,j] - a_beta*distance_beta                          
            r1_delta              = int.from_bytes(os.urandom(8), byteorder = "big") / ((1 << 64) - 1)
            r2_delta              = int.from_bytes(os.urandom(8), byteorder = "big") / ((1 << 64) - 1)
            a_delta               = 2*a_linear_component*r1_delta - a_linear_component
            c_delta               = 2*r2_delta            
            distance_delta        = abs(c_delta*delta[0,j] - position[i,j]) 
            x3                    = delta[0,j] - a_delta*distance_delta                                 
            updated_position[i,j] = np.clip(((x1 + x2 + x3)/3),min_values[j],max_values[j])     
        updated_position[i,-1] = target_function(updated_position[i,0:updated_position.shape[1]-1])
    return updated_position
 
# GWO Function
def grey_wolf_optimizer(pack_size = 5, min_values = [-5,-5], max_values = [5,5], iterations = 50, target_function = target_function):    
    count    = 0
    alpha    = alpha_position(dimension = len(min_values), target_function = target_function)
    beta     = beta_position(dimension  = len(min_values), target_function = target_function)
    delta    = delta_position(dimension = len(min_values), target_function = target_function)
    position = initial_position(pack_size = pack_size, min_values = min_values, max_values = max_values, target_function = target_function)
    while (count <= iterations):      
        print("Iteration = ", count, " f(x) = ", alpha[-1])      
        a_linear_component = 2 - count*(2/iterations)
        alpha, beta, delta = update_pack(position, alpha, beta, delta)
        position           = update_position(position, alpha, beta, delta, a_linear_component = a_linear_component, min_values = min_values, max_values = max_values, target_function = target_function)       
        count              = count + 1       
        fitness.append(alpha[-1][3])  
    return alpha
 
#===============================================================================应用测试
 
def six_hump_camel_back(variables_values = [0,0,0]):
    if variables_values[0]==0:variables_values[0]=0.1
    if variables_values[1]==0:variables_values[1]=0.1
    if variables_values[2]==0:variables_values[2]=0.1
    model_svr = SVR(C=variables_values[0],epsilon=variables_values[1],gamma=variables_values[2])
    model_svr.fit(feature_train,target_train)
    predict_results=model_svr.predict(feature_test)
    return metrics.r2_score(target_test,predict_results)
 
gwo = grey_wolf_optimizer(pack_size = 15, min_values = [0,0,0], max_values = [10,10,100], iterations = 50, target_function = six_hump_camel_back)

plt.figure()
fig = plt.gcf()
fig.set_size_inches(18.5, 10.5)
plt.xlabel("iterators", size=14)
plt.ylabel("fitness", size=14)
plt.plot(fitness, color='b', linewidth=2)
plt.show()

 5.跑PSO

代码略，见我专栏

6.跑GA

不说了，更拉

7.结果对比

GWO：

PSO：

GA：