1.5 基于改进蛇优化VGG13SE故障诊断方法的有效性分析
本博客来源于CSDN机器鱼，未同意任何人转载。
更多内容，欢迎点击本专栏，查看更多内容。
目录
0 引言
1 改进优化算法的有效性分析
2 改进算法用于优化VGG13SE的有效性分析
3 其他对比
4 结语
0 引言
此博客是轴承故障【从数据处理到模型优化全流程】系列的最后一篇，主要是从多个角度对比所提方法的有效性。
1 改进优化算法的有效性分析
在【博客】中，我们利用tent映射与模拟退火对蛇算法进行了改进，所以首先对比改进前后并与其他算法进行对比，分别是改进蛇、蛇、粒子群、模拟退火四种算法，代码如下：
import numpy as np
import matplotlib.pyplot as plt
import math
import random
# 修正后的适应度函数
def chung_reynolds(x):
return np.sum(np.power(x, 4)) # 直接计算x的四次方和
def boundary(pop,lb,ub,dim):
for j in range(dim):
if pop[j]>ub[j] or pop[j]<lb[j]:
pop[j] =(ub[j]-lb[j])*np.random.rand()+lb[j]
return pop
'''粒子群优化算法'''
def PSO(pN,max_iter,dim,bound):
lb=bound[0]*np.ones(dim)
ub=bound[1]*np.ones(dim)
wmax = 0.8;wmin = 0.2;
c1 = 1.5;c2 = 1.5;r1 = 0.8;r2 = 0.3
#初始化
X = np.zeros((pN,dim))
V = np.zeros((pN,dim))
pbest = np.zeros((pN,dim))
gbest = np.zeros((1,dim))
p_fit = np.zeros(pN)
fit = np.inf
for i in range(pN):
for j in range(dim):
X[i][j] = (ub[j]-lb[j])*np.random.rand()+lb[j]
V[i][j] = np.random.rand()
pbest[i] = X[i].copy()
tmp = chung_reynolds(X[i,:])
p_fit[i] = tmp.copy()
if(tmp < fit):
fit = tmp.copy()
gbest = X[i].copy()
# 开始循环迭代
trace=[]
for t in range(max_iter):
w=wmax-(wmax-wmin)*t/max_iter
for i in range(pN):
V[i,:] = w*V[i,:] + c1*r1*(pbest[i] - X[i,:])+c2*r2*(gbest - X[i,:])
X[i,:] = X[i,:] + V[i,:]
# 加入变异操作，避免陷入局部最优
prob=0.5*t/max_iter +0.5#变异，随着进化代数的增加，变异几率越小
if np.random.rand()>prob:
for j in range(dim):
X[i][j] =(ub[j]-lb[j])*np.random.rand()+lb[j]
for i in range(pN): #更新gbest\pbest
X[i,:]=boundary(X[i,:],lb,ub,dim)#位置边界判断
V[i,:]=boundary(V[i,:],np.zeros(dim,),np.ones(dim,),dim)#速度边界判断
temp = chung_reynolds(X[i,:])
if(temp < p_fit[i]): #更新个体最优
p_fit[i] = temp
pbest[i,:] = X[i,:].copy()
if(p_fit[i] < fit): #更新全局最优
gbest = X[i,:].copy()
fit = p_fit[i].copy()
trace.append(fit)
return gbest ,fit,trace
'''模拟退火优化算法'''
def simulated_annealing(dim, bounds, max_iter, initial_temp, cooling_rate):
# Initialize variables
lower_bound, upper_bound = bounds
current_solution = np.random.uniform(lower_bound, upper_bound, dim)
current_fitness = chung_reynolds(current_solution)
best_solution = current_solution.copy()
best_fitness = current_fitness
fitness_history = []
temperature = initial_temp
for iteration in range(max_iter):
# Generate a new candidate solution by perturbing the current solution
new_solution = current_solution + np.random.uniform(-1, 1, dim)
new_solution = np.clip(new_solution, lower_bound, upper_bound) # Ensure bounds
new_fitness = chung_reynolds(new_solution)
# Calculate acceptance probability
delta_fitness = new_fitness - current_fitness
acceptance_probability = np.exp(-delta_fitness / temperature) if delta_fitness > 0 else 1.0
# Accept or reject the new solution
if np.random.rand() < acceptance_probability:
current_solution = new_solution
current_fitness = new_fitness
# Update the best solution if necessary
if current_fitness < best_fitness:
best_solution = current_solution.copy()
best_fitness = current_fitness
# Record fitness history
fitness_history.append(best_fitness)
# Decrease the temperature
temperature *= cooling_rate
return best_solution, best_fitness, fitness_history
def snake_optimization(search_agents_no, max_iter, dim, solution_bound):
# 参数初始化
T = 1000
T_min = 0.01
alpha = 0.9
vec_flag = [1, -1]
food_threshold = 0.25
temp_threshold = 0.6
c1 = 0.5
model_threshold = 0.6
c2 = 0.5
c3 = 2
# 种群初始化（使用numpy数组代替矩阵）
X = solution_bound[0] + np.random.rand(search_agents_no, dim) * (solution_bound[1] - solution_bound[0])
# 适应度计算
fitness = np.array([chung_reynolds(x) for x in X])
g_best_idx = np.argmin(fitness)
food = X[g_best_idx].copy()
gy_best = fitness[g_best_idx]
# 性别分组
male_number = search_agents_no // 2
female_number = search_agents_no - male_number
male = X[:male_number].copy()
female = X[male_number:].copy()
# 性别相关参数
male_fitness = fitness[:male_number]
female_fitness = fitness[male_number:]
male_best_idx = np.argmin(male_fitness)
male_best = male_fitness[male_best_idx]
male_best_solution = male[male_best_idx].copy()
female_best_idx = np.argmin(female_fitness)
female_best = female_fitness[female_best_idx]
female_best_solution = female[female_best_idx].copy()
gene_best = np.zeros(max_iter)
for t in range(max_iter):
temp = math.exp(-t / max_iter)
quantity = c1 * math.exp((t - max_iter) / max_iter)
quantity = min(quantity, 1)
new_male = male.copy()
new_female = female.copy()
if quantity < food_threshold:
# 无食物模式
for i in range(male_number):
leader = random.choice(male)
am = math.exp(-(chung_reynolds(leader) / (male_fitness[i] + 1e-30)))
direction = random.choice(vec_flag)
new_male[i] = leader + direction * c2 * am * np.random.rand(dim) * (
solution_bound[1] - solution_bound[0])
for i in range(female_number):
leader = random.choice(female)
am = math.exp(-(chung_reynolds(leader) / (female_fitness[i] + 1e-30)))
direction = random.choice(vec_flag)
new_female[i] = leader + direction * c2 * am * np.random.rand(dim) * (
solution_bound[1] - solution_bound[0])
else:
if temp > temp_threshold:
# 热态探索
for i in range(male_number):
direction = random.choice(vec_flag)
new_male[i] = food + direction * c3 * temp * np.random.rand(dim) * (food - male[i])
for i in range(female_number):
direction = random.choice(vec_flag)
new_female[i] = food + direction * c3 * temp * np.random.rand(dim) * (food - female[i])
else:
if random.random() < model_threshold:
# 战斗模式
for i in range(male_number):
fm = math.exp(-female_best / (male_fitness[i] + 1e-30))
direction = random.choice(vec_flag)
new_male[i] = male[i] + direction * c3 * fm * np.random.rand(dim) * (
quantity * male_best_solution - male[i])
for i in range(female_number):
ff = math.exp(-male_best / (female_fitness[i] + 1e-30))
direction = random.choice(vec_flag)
new_female[i] = female[i] + direction * c3 * ff * np.random.rand(dim) * (
quantity * female_best_solution - female[i])
else:
# 交配模式
for i in range(male_number):
mm = math.exp(-female_fitness[i % female_number] / (male_fitness[i] + 1e-30))
direction = random.choice(vec_flag)
new_male[i] = male[i] + direction * c3 * np.random.rand() * mm * (
quantity * female[i % female_number] - male[i])
for i in range(female_number):
mf = math.exp(-male_fitness[i % male_number] / (female_fitness[i] + 1e-30))
direction = random.choice(vec_flag)
new_female[i] = female[i] + direction * c3 * np.random.rand() * mf * (
quantity * male[i % male_number] - female[i])
if random.random() < 0.5:
worst_male = np.argmax(male_fitness)
new_male[worst_male] = solution_bound[0] + np.random.rand(dim) * (
solution_bound[1] - solution_bound[0])
worst_female = np.argmax(female_fitness)
new_female[worst_female] = solution_bound[0] + np.random.rand(dim) * (
solution_bound[1] - solution_bound[0])
# 边界处理
new_male = np.clip(new_male, solution_bound[0], solution_bound[1])
new_female = np.clip(new_female, solution_bound[0], solution_bound[1])
# 更新雄性种群
for i in range(male_number):
new_fit = chung_reynolds(new_male[i])
if new_fit < male_fitness[i]:
male[i] = new_male[i]
male_fitness[i] = new_fit
# 更新雌性种群
for i in range(female_number):
new_fit = chung_reynolds(new_female[i])
if new_fit < female_fitness[i]:
female[i] = new_female[i]
female_fitness[i] = new_fit
# 更新全局最优
current_male_best = np.min(male_fitness)
current_female_best = np.min(female_fitness)
gene_best[t] = min(current_male_best, current_female_best)
if current_male_best < current_female_best:
if current_male_best < gy_best:
gy_best = current_male_best
food[:] = male[np.argmin(male_fitness)]
else:
if current_female_best < gy_best:
gy_best = current_female_best
food[:] = female[np.argmin(female_fitness)]
return food, gy_best, gene_best
'''tent映射'''
def tent(dim):
x=np.zeros([dim])
x[0]=np.random.rand()#初始点
a=0.7 #参数a的值
for i in range(dim-1):
if x[i]<a:
x[i+1]=x[i]/a
else:
x[i+1]=(1-x[i]) / (1 - a)
return x
''' Tent种群初始化函数 '''
def initialnew(pop, dim, ub, lb):
X = np.zeros([pop, dim])
for i in range(pop):
tentvalue=tent(dim)
for j in range(dim):
X[i, j] = tentvalue[j]*(ub - lb) + lb
return X
# 修正后的邻域生成函数
def neighbor(x, solution_bound):
perturbation = np.random.uniform(-1, 1, size=x.shape)
x_new = x + perturbation
# 确保扰动后仍在边界内
x_new = np.clip(x_new, solution_bound[0], solution_bound[1])
return x_new
def sa_snake_optimization(search_agents_no, max_iter, dim, solution_bound):
# 参数初始化
T = 1000
T_min = 0.01
alpha = 0.9
vec_flag = [1, -1]
food_threshold = 0.25
temp_threshold = 0.6
c1 = 0.5
model_threshold = 0.6
c2 = 0.5
c3 = 2
# 种群初始化（使用numpy数组代替矩阵）
#X = solution_bound[0] + np.random.rand(search_agents_no, dim) * (solution_bound[1] - solution_bound[0])
X = initialnew(search_agents_no, dim, solution_bound[1] , solution_bound[0]) #初始化种群
# 适应度计算
fitness = np.array([chung_reynolds(x) for x in X])
g_best_idx = np.argmin(fitness)
food = X[g_best_idx].copy()
gy_best = fitness[g_best_idx]
# 性别分组
male_number = search_agents_no // 2
female_number = search_agents_no - male_number
male = X[:male_number].copy()
female = X[male_number:].copy()
# 性别相关参数
male_fitness = fitness[:male_number]
female_fitness = fitness[male_number:]
male_best_idx = np.argmin(male_fitness)
male_best = male_fitness[male_best_idx]
male_best_solution = male[male_best_idx].copy()
female_best_idx = np.argmin(female_fitness)
female_best = female_fitness[female_best_idx]
female_best_solution = female[female_best_idx].copy()
gene_best = np.zeros(max_iter)
for t in range(max_iter):
temp = math.exp(-t / max_iter)
quantity = c1 * math.exp((t - max_iter) / max_iter)
quantity = min(quantity, 1)
new_male = male.copy()
new_female = female.copy()
if quantity < food_threshold:
# 无食物模式
for i in range(male_number):
leader = random.choice(male)
am = math.exp(-(chung_reynolds(leader) / (male_fitness[i] + 1e-30)))
direction = random.choice(vec_flag)
new_male[i] = leader + direction * c2 * am * np.random.rand(dim) * (
solution_bound[1] - solution_bound[0])
for i in range(female_number):
leader = random.choice(female)
am = math.exp(-(chung_reynolds(leader) / (female_fitness[i] + 1e-30)))
direction = random.choice(vec_flag)
new_female[i] = leader + direction * c2 * am * np.random.rand(dim) * (
solution_bound[1] - solution_bound[0])
else:
if temp > temp_threshold:
# 热态探索
for i in range(male_number):
direction = random.choice(vec_flag)
new_male[i] = food + direction * c3 * temp * np.random.rand(dim) * (food - male[i])
for i in range(female_number):
direction = random.choice(vec_flag)
new_female[i] = food + direction * c3 * temp * np.random.rand(dim) * (food - female[i])
else:
if random.random() < model_threshold:
# 战斗模式
for i in range(male_number):
fm = math.exp(-female_best / (male_fitness[i] + 1e-30))
direction = random.choice(vec_flag)
new_male[i] = male[i] + direction * c3 * fm * np.random.rand(dim) * (
quantity * male_best_solution - male[i])
for i in range(female_number):
ff = math.exp(-male_best / (female_fitness[i] + 1e-30))
direction = random.choice(vec_flag)
new_female[i] = female[i] + direction * c3 * ff * np.random.rand(dim) * (
quantity * female_best_solution - female[i])
else:
# 交配模式
for i in range(male_number):
mm = math.exp(-female_fitness[i % female_number] / (male_fitness[i] + 1e-30))
direction = random.choice(vec_flag)
new_male[i] = male[i] + direction * c3 * np.random.rand() * mm * (
quantity * female[i % female_number] - male[i])
for i in range(female_number):
mf = math.exp(-male_fitness[i % male_number] / (female_fitness[i] + 1e-30))
direction = random.choice(vec_flag)
new_female[i] = female[i] + direction * c3 * np.random.rand() * mf * (
quantity * male[i % male_number] - female[i])
if random.random() < 0.5:
worst_male = np.argmax(male_fitness)
new_male[worst_male] = solution_bound[0] + np.random.rand(dim) * (
solution_bound[1] - solution_bound[0])
worst_female = np.argmax(female_fitness)
new_female[worst_female] = solution_bound[0] + np.random.rand(dim) * (
solution_bound[1] - solution_bound[0])
# 边界处理
new_male = np.clip(new_male, solution_bound[0], solution_bound[1])
new_female = np.clip(new_female, solution_bound[0], solution_bound[1])
# 更新雄性种群
for i in range(male_number):
new_fit = chung_reynolds(new_male[i])
if new_fit < male_fitness[i]:
male[i] = new_male[i]
male_fitness[i] = new_fit
# 更新雌性种群
for i in range(female_number):
new_fit = chung_reynolds(new_female[i])
if new_fit < female_fitness[i]:
female[i] = new_female[i]
female_fitness[i] = new_fit
# 更新全局最优
current_male_best = np.min(male_fitness)
current_female_best = np.min(female_fitness)
gene_best[t] = min(current_male_best, current_female_best)
if current_male_best < current_female_best:
if current_male_best < gy_best:
gy_best = current_male_best
food[:] = male[np.argmin(male_fitness)]
else:
if current_female_best < gy_best:
gy_best = current_female_best
food[:] = female[np.argmin(female_fitness)]
# 模拟退火
if T > T_min:
food_new = neighbor(food, solution_bound)
f_new = chung_reynolds(food_new)
delta = f_new - gy_best
if delta < 0 or np.random.rand() < math.exp(-delta / T):
food[:] = food_new
gy_best = f_new
T *= alpha
return food, gy_best, gene_best
if __name__ == "__main__":
# Problem settings
dim = 10
search_agents_no = 26
bounds = (-100, 100)
max_iter = 100
initial_temp = 1000
cooling_rate = 0.98
# Run
best_solution0, best_fitness0, fitness_history0 = simulated_annealing(
dim, bounds, max_iter, initial_temp, cooling_rate)
best_solution1, best_fitness1, fitness_history1 = PSO(
search_agents_no, max_iter, dim, bounds)
best_solution2, best_fitness2, fitness_history2 = snake_optimization(
search_agents_no, max_iter, dim, bounds
)
best_solution3, best_fitness3, fitness_history3 = sa_snake_optimization(
search_agents_no, max_iter, dim, bounds
)
print("SA最佳解决方案：", best_solution0)
print("SA最佳适应度：", best_fitness0)
print("PSO最佳解决方案：", best_solution1)
print("PSO最佳适应度：", best_fitness1)
print("SO最佳解决方案：", best_solution2)
print("SO最佳适应度：", best_fitness2)
print("SASO最佳解决方案：", best_solution3)
print("SASO最佳适应度：", best_fitness3)
# Plot optimization process
plt.figure(figsize=(10, 6))
plt.plot(fitness_history0, linewidth=2,label='SA')
plt.plot(fitness_history1, linewidth=2,label='PSO')
plt.plot(fitness_history2, linewidth=2,label='SO')
plt.plot(fitness_history3, linewidth=2,label='SASO')
plt.legend()
plt.title("Optimization Process")
plt.xlabel("Iteration")
plt.ylabel("Best Fitness")
plt.grid(True)
# plt.savefig('result/函数极值寻优适应度函数对比.jpg')
plt.show()
2 改进算法用于优化VGG13SE的有效性分析
与1一样，为了验证所提方法的有效性，将上述四种方法用于优化VGG13SE网络。代码如下
# -*- coding: utf-8 -*-
import matplotlib.pyplot as plt
import numpy as np
import math
from model import ProposedNet2D
from tensorflow.keras.losses import SparseCategoricalCrossentropy
from tensorflow.keras.optimizers import Adam,SGD
import tensorflow as tf
import math
import random
plt.rcParams['font.family'] = ['sans-serif']
plt.rcParams['font.sans-serif'] = ['SimHei']
plt.rcParams['axes.unicode_minus']=False
def HUATU(mode,trace):
plt.figure()
plt.plot(trace)
plt.title('fitness curve')
plt.xlabel('iteration')
plt.ylabel('fitness value')
plt.savefig('result/'+mode+'_cnn_fitness.png')
# plt.show()
# 修正后的适应度函数
def chung_reynolds(pop,X_train,y_train,X_valid,y_valid):
tf.random.set_seed(0)
para = [pop[i] if i==0 else int(pop[i]) for i in range(len(pop))]
print(para)
learningrate=para[0]#学习率
epochs = para[1]#迭代次数
batch_size = para[2] #BATCH_SIZE
k1_num,k1_size,k2_num,k2_size,k3_num,k3_size,k4_num,k4_size,k5_num,k5_size,fc1,fc2=para[3:]
model= ProposedNet2D(k1_num,k1_size,k2_num,k2_size,k3_num,k3_size,k4_num,k4_size,k5_num,k5_size,fc1,fc2)
#model.summary()
model.compile(optimizer=Adam(learning_rate=learningrate),
loss=SparseCategoricalCrossentropy(),
metrics=['accuracy'])
history = model.fit(X_train,y_train, epochs=epochs,
validation_data=(X_valid, y_valid),
batch_size=batch_size, verbose=0)
return 1-history.history['val_accuracy'][-1]#以最小化验证集错误率为适应度函数
def boundary(pop,lb,ub,dim):
for j in range(dim):
if pop[j]>ub[j] or pop[j]<lb[j]:
pop[j] =(ub[j]-lb[j])*np.random.rand()+lb[j]
return pop
def Optim(mode,p_train,t_trian,p_valid,t_valid):
search_agents_no=5
max_iter=10
initial_temp = 10
cooling_rate = 0.98
#一共有10个参数需要优化
#分别是 learningrate epoches batchsize k1num k1size k2num k2size k3num k3size k4num k4size k5num k5size fc1 fc2
# lb ub是各参数的寻优上下界，如学习率是0.0001-0.01的范围内找最优解，epoches是50-100的范围内找最优解，fc1是10-200的范围内找最优解
lb=[0.0001,30, 32,5 ,1,10 ,1 ,20 ,1,30,1,40,1, 10 ,10]
ub=[0.01, 100,128,10,5,20, 5, 30, 5,40,5,50,5,200,100]#
dim = len(lb)
if mode in ['SASO','PSO','SO']:
return eval(mode)(search_agents_no,max_iter,dim,lb,ub,p_train,t_trian,p_valid,t_valid)
else:
return SA(max_iter,dim,lb,ub,initial_temp,cooling_rate,p_train,t_trian,p_valid,t_valid)
'''粒子群优化算法'''
def PSO(pN,max_iter,dim,lb,ub,p_train,t_trian,p_valid,t_valid):
wmax = 0.8;wmin = 0.2;
c1 = 1.5;c2 = 1.5;r1 = 0.8;r2 = 0.3
#初始化
X = np.zeros((pN,dim))
V = np.zeros((pN,dim))
pbest = np.zeros((pN,dim))
gbest = np.zeros((1,dim))
p_fit = np.zeros(pN)
fit = np.inf
for i in range(pN):
for j in range(dim):
X[i][j] = (ub[j]-lb[j])*np.random.rand()+lb[j]
V[i][j] = np.random.rand()
pbest[i] = X[i].copy()
p_fit[i] = chung_reynolds(X[i,:],p_train,t_trian,p_valid,t_valid)
if(p_fit[i] < fit):
fit = p_fit[i].copy()
gbest = X[i].copy()
# 开始循环迭代
trace=[]
for t in range(max_iter):
w=wmax-(wmax-wmin)*t/max_iter
for i in range(pN):
V[i,:] = w*V[i,:] + c1*r1*(pbest[i] - X[i,:])+c2*r2*(gbest - X[i,:])
X[i,:] = X[i,:] + V[i,:]
# 加入变异操作，避免陷入局部最优
prob=0.5*t/max_iter +0.5#变异，随着进化代数的增加，变异几率越小
if np.random.rand()>prob:
for j in range(dim):
X[i][j] =(ub[j]-lb[j])*np.random.rand()+lb[j]
for i in range(pN): #更新gbest\pbest
X[i,:]=boundary(X[i,:],lb,ub,dim)#位置边界判断
V[i,:]=boundary(V[i,:],np.zeros(dim,),np.ones(dim,),dim)#速度边界判断
temp = chung_reynolds(X[i,:],p_train,t_trian,p_valid,t_valid)
if(temp < p_fit[i]): #更新个体最优
p_fit[i] = temp
pbest[i,:] = X[i,:].copy()
if(p_fit[i] < fit): #更新全局最优
gbest = X[i,:].copy()
fit = p_fit[i].copy()
print('第',t,'轮优化完成，最优适应度值为:',fit)
trace.append(fit)
return gbest ,fit,np.array(trace)
'''模拟退火优化算法'''
def SA(max_iter,dim,lower_bound, upper_bound, initial_temp, cooling_rate,p_train,t_trian,p_valid,t_valid):
# Initialize variables
lower_bound, upper_bound =np.array(lower_bound),np.array(upper_bound)
current_solution = lower_bound + np.random.rand(dim,) * (upper_bound - lower_bound)
current_fitness = chung_reynolds(current_solution,p_train,t_trian,p_valid,t_valid)
best_solution = current_solution.copy()
best_fitness = current_fitness
fitness_history = []
temperature = initial_temp
for iteration in range(max_iter):
# Generate a new candidate solution by perturbing the current solution
new_solution = current_solution + np.random.uniform(-1, 1, dim)
new_solution = np.clip(new_solution, lower_bound, upper_bound) # Ensure bounds
new_fitness = chung_reynolds(new_solution,p_train,t_trian,p_valid,t_valid)
# Calculate acceptance probability
delta_fitness = new_fitness - current_fitness
acceptance_probability = np.exp(-delta_fitness / temperature) if delta_fitness > 0 else 1.0
# Accept or reject the new solution
if np.random.rand() < acceptance_probability:
current_solution = new_solution
current_fitness = new_fitness
# Update the best solution if necessary
if current_fitness < best_fitness:
best_solution = current_solution.copy()
best_fitness = current_fitness
# Record fitness history
fitness_history.append(best_fitness)
print('第',iteration,'轮优化完成，最优适应度值为:',best_fitness)
# Decrease the temperature
temperature *= cooling_rate
return best_solution, best_fitness, np.array(fitness_history)
'''蛇算法优化'''
def SO(search_agents_no, max_iter, dim, lower_bound,upper_bound,p_train,t_trian,p_valid,t_valid):
# 参数初始化
vec_flag = [1, -1]
food_threshold = 0.25
temp_threshold = 0.6
c1 = 0.5
model_threshold = 0.6
c2 = 0.5
c3 = 2
lower_bound, upper_bound =np.array(lower_bound),np.array(upper_bound)
# 种群初始化（使用numpy数组代替矩阵）
X = lower_bound + np.random.rand(search_agents_no, dim) * (upper_bound - lower_bound)
# 适应度计算
fitness = np.array([chung_reynolds(x,p_train,t_trian,p_valid,t_valid) for x in X])
g_best_idx = np.argmin(fitness)
food = X[g_best_idx].copy()
gy_best = fitness[g_best_idx]
# 性别分组
male_number = search_agents_no // 2
female_number = search_agents_no - male_number
male = X[:male_number].copy()
female = X[male_number:].copy()
# 性别相关参数
male_fitness = fitness[:male_number]
female_fitness = fitness[male_number:]
male_best_idx = np.argmin(male_fitness)
male_best = male_fitness[male_best_idx]
male_best_solution = male[male_best_idx].copy()
female_best_idx = np.argmin(female_fitness)
female_best = female_fitness[female_best_idx]
female_best_solution = female[female_best_idx].copy()
gene_best = np.zeros(max_iter)
for t in range(max_iter):
temp = math.exp(-t / max_iter)
quantity = c1 * math.exp((t - max_iter) / max_iter)
quantity = min(quantity, 1)
new_male = male.copy()
new_female = female.copy()
if quantity < food_threshold:
# 无食物模式
for i in range(male_number):
leader = random.choice(male)
am = math.exp(-(chung_reynolds(leader,p_train,t_trian,p_valid,t_valid) / (male_fitness[i] + 1e-30)))
direction = random.choice(vec_flag)
new_male[i] = leader + direction * c2 * am * np.random.rand(dim) * (
upper_bound - lower_bound)
for i in range(female_number):
leader = random.choice(female)
am = math.exp(-(chung_reynolds(leader,p_train,t_trian,p_valid,t_valid) / (female_fitness[i] + 1e-30)))
direction = random.choice(vec_flag)
new_female[i] = leader + direction * c2 * am * np.random.rand(dim) * (
upper_bound - lower_bound)
else:
if temp > temp_threshold:
# 热态探索
for i in range(male_number):
direction = random.choice(vec_flag)
new_male[i] = food + direction * c3 * temp * np.random.rand(dim) * (food - male[i])
for i in range(female_number):
direction = random.choice(vec_flag)
new_female[i] = food + direction * c3 * temp * np.random.rand(dim) * (food - female[i])
else:
if random.random() < model_threshold:
# 战斗模式
for i in range(male_number):
fm = math.exp(-female_best / (male_fitness[i] + 1e-30))
direction = random.choice(vec_flag)
new_male[i] = male[i] + direction * c3 * fm * np.random.rand(dim) * (
quantity * male_best_solution - male[i])
for i in range(female_number):
ff = math.exp(-male_best / (female_fitness[i] + 1e-30))
direction = random.choice(vec_flag)
new_female[i] = female[i] + direction * c3 * ff * np.random.rand(dim) * (
quantity * female_best_solution - female[i])
else:
# 交配模式
for i in range(male_number):
mm = math.exp(-female_fitness[i % female_number] / (male_fitness[i] + 1e-30))
direction = random.choice(vec_flag)
new_male[i] = male[i] + direction * c3 * np.random.rand() * mm * (
quantity * female[i % female_number] - male[i])
for i in range(female_number):
mf = math.exp(-male_fitness[i % male_number] / (female_fitness[i] + 1e-30))
direction = random.choice(vec_flag)
new_female[i] = female[i] + direction * c3 * np.random.rand() * mf * (
quantity * male[i % male_number] - female[i])
if random.random() < 0.5:
worst_male = np.argmax(male_fitness)
new_male[worst_male] = lower_bound + np.random.rand(dim) * (
upper_bound - lower_bound)
worst_female = np.argmax(female_fitness)
new_female[worst_female] = lower_bound + np.random.rand(dim) * (
upper_bound - lower_bound)
# 边界处理
new_male = np.clip(new_male, lower_bound, upper_bound)
new_female = np.clip(new_female, lower_bound, upper_bound)
# 更新雄性种群
for i in range(male_number):
new_fit = chung_reynolds(new_male[i],p_train,t_trian,p_valid,t_valid)
if new_fit < male_fitness[i]:
male[i] = new_male[i]
male_fitness[i] = new_fit
# 更新雌性种群
for i in range(female_number):
new_fit = chung_reynolds(new_female[i],p_train,t_trian,p_valid,t_valid)
if new_fit < female_fitness[i]:
female[i] = new_female[i]
female_fitness[i] = new_fit
# 更新全局最优
current_male_best = np.min(male_fitness)
current_female_best = np.min(female_fitness)
gene_best[t] = min(current_male_best, current_female_best)
print('第',t,'轮优化完成，最优适应度值为:',gene_best[t])
if current_male_best < current_female_best:
if current_male_best < gy_best:
gy_best = current_male_best
food[:] = male[np.argmin(male_fitness)]
else:
if current_female_best < gy_best:
gy_best = current_female_best
food[:] = female[np.argmin(female_fitness)]
return food, gy_best, gene_best
# 修正后的邻域生成函数
def neighbor(x, solution_bound):
perturbation = np.random.uniform(-1, 1, size=x.shape)
x_new = x + perturbation
# 确保扰动后仍在边界内
x_new = np.clip(x_new, solution_bound[0], solution_bound[1])
return x_new
'''tent映射'''
def tent(dim):
x=np.zeros([dim])
x[0]=np.random.rand()#初始点
a=0.7 #参数a的值
for i in range(dim-1):
if x[i]<a:
x[i+1]=x[i]/a
else:
x[i+1]=(1-x[i]) / (1 - a)
return x
''' Tent种群初始化函数 '''
def initialnew(pop, dim, ub, lb):
X = np.zeros([pop, dim])
for i in range(pop):
tentvalue=tent(dim)
for j in range(dim):
X[i, j] = tentvalue[j]*(ub[j] - lb[j]) + lb[j]
return X
'''改进的蛇算法优化'''
def SASO(search_agents_no, max_iter, dim, lower_bound,upper_bound,p_train,t_trian,p_valid,t_valid):
# 参数初始化
T = 1000
T_min = 0.01
alpha = 0.9
vec_flag = [1, -1]
food_threshold = 0.25
temp_threshold = 0.6
c1 = 0.5
model_threshold = 0.6
c2 = 0.5
c3 = 2
lower_bound, upper_bound =np.array(lower_bound),np.array(upper_bound)
# 种群初始化（使用numpy数组代替矩阵）
X = initialnew(search_agents_no, dim, upper_bound, lower_bound) #初始化种群
# 适应度计算
fitness = np.array([chung_reynolds(x,p_train,t_trian,p_valid,t_valid) for x in X])
g_best_idx = np.argmin(fitness)
food = X[g_best_idx].copy()
gy_best = fitness[g_best_idx]
# 性别分组
male_number = search_agents_no // 2
female_number = search_agents_no - male_number
male = X[:male_number].copy()
female = X[male_number:].copy()
# 性别相关参数
male_fitness = fitness[:male_number]
female_fitness = fitness[male_number:]
male_best_idx = np.argmin(male_fitness)
male_best = male_fitness[male_best_idx]
male_best_solution = male[male_best_idx].copy()
female_best_idx = np.argmin(female_fitness)
female_best = female_fitness[female_best_idx]
female_best_solution = female[female_best_idx].copy()
gene_best = np.zeros(max_iter)
for t in range(max_iter):
temp = math.exp(-t / max_iter)
quantity = c1 * math.exp((t - max_iter) / max_iter)
quantity = min(quantity, 1)
new_male = male.copy()
new_female = female.copy()
if quantity < food_threshold:
# 无食物模式
for i in range(male_number):
leader = random.choice(male)
am = math.exp(-(chung_reynolds(leader,p_train,t_trian,p_valid,t_valid) / (male_fitness[i] + 1e-30)))
direction = random.choice(vec_flag)
new_male[i] = leader + direction * c2 * am * np.random.rand(dim) * (
upper_bound - lower_bound)
for i in range(female_number):
leader = random.choice(female)
am = math.exp(-(chung_reynolds(leader,p_train,t_trian,p_valid,t_valid) / (female_fitness[i] + 1e-30)))
direction = random.choice(vec_flag)
new_female[i] = leader + direction * c2 * am * np.random.rand(dim) * (
upper_bound - lower_bound)
else:
if temp > temp_threshold:
# 热态探索
for i in range(male_number):
direction = random.choice(vec_flag)
new_male[i] = food + direction * c3 * temp * np.random.rand(dim) * (food - male[i])
for i in range(female_number):
direction = random.choice(vec_flag)
new_female[i] = food + direction * c3 * temp * np.random.rand(dim) * (food - female[i])
else:
if random.random() < model_threshold:
# 战斗模式
for i in range(male_number):
fm = math.exp(-female_best / (male_fitness[i] + 1e-30))
direction = random.choice(vec_flag)
new_male[i] = male[i] + direction * c3 * fm * np.random.rand(dim) * (
quantity * male_best_solution - male[i])
for i in range(female_number):
ff = math.exp(-male_best / (female_fitness[i] + 1e-30))
direction = random.choice(vec_flag)
new_female[i] = female[i] + direction * c3 * ff * np.random.rand(dim) * (
quantity * female_best_solution - female[i])
else:
# 交配模式
for i in range(male_number):
mm = math.exp(-female_fitness[i % female_number] / (male_fitness[i] + 1e-30))
direction = random.choice(vec_flag)
new_male[i] = male[i] + direction * c3 * np.random.rand() * mm * (
quantity * female[i % female_number] - male[i])
for i in range(female_number):
mf = math.exp(-male_fitness[i % male_number] / (female_fitness[i] + 1e-30))
direction = random.choice(vec_flag)
new_female[i] = female[i] + direction * c3 * np.random.rand() * mf * (
quantity * male[i % male_number] - female[i])
if random.random() < 0.5:
worst_male = np.argmax(male_fitness)
new_male[worst_male] = lower_bound + np.random.rand(dim) * (
upper_bound - lower_bound)
worst_female = np.argmax(female_fitness)
new_female[worst_female] = lower_bound + np.random.rand(dim) * (
upper_bound - lower_bound)
# 边界处理
new_male = np.clip(new_male, lower_bound, upper_bound)
new_female = np.clip(new_female, lower_bound, upper_bound)
# 更新雄性种群
for i in range(male_number):
new_fit = chung_reynolds(new_male[i],p_train,t_trian,p_valid,t_valid)
if new_fit < male_fitness[i]:
male[i] = new_male[i]
male_fitness[i] = new_fit
# 更新雌性种群
for i in range(female_number):
new_fit = chung_reynolds(new_female[i],p_train,t_trian,p_valid,t_valid)
if new_fit < female_fitness[i]:
female[i] = new_female[i]
female_fitness[i] = new_fit
# 更新全局最优
current_male_best = np.min(male_fitness)
current_female_best = np.min(female_fitness)
gene_best[t] = min(current_male_best, current_female_best)
print('第',t,'轮优化完成，最优适应度值为:',gene_best[t])
if current_male_best < current_female_best:
if current_male_best < gy_best:
gy_best = current_male_best
food[:] = male[np.argmin(male_fitness)]
else:
if current_female_best < gy_best:
gy_best = current_female_best
food[:] = female[np.argmin(female_fitness)]
# 模拟退火
if T > T_min:
food_new = neighbor(food, [lower_bound,upper_bound])
f_new = chung_reynolds(food_new,p_train,t_trian,p_valid,t_valid)
delta = f_new - gy_best
if delta < 0 or np.random.rand() < math.exp(-delta / T):
food[:] = food_new
gy_best = f_new
T *= alpha
return food, gy_best, gene_best
调用上述代码的程序写在同一个文件中了，只需要改下面的mode='选择的方法'即可。
下面程序中有两个地方设定了False，分别是优化前面和训练前面，因为优化的时间很长，所以我把结果保存了下来，后面只需要调用保存的参数即可，如果要重新优化，改成optimization_again = True。同理如果要重新训练，则改train_again = True训练即可。代码如下：
# coding: utf-8
# In[1]: 导入必要的库函数
import os
import numpy as np
import matplotlib.pyplot as plt
from sklearn.metrics import confusion_matrix
import seaborn as sns
from model import ProposedNet2D,SE_AttentionLayer
import tensorflow as tf
from tensorflow.keras.losses import SparseCategoricalCrossentropy
from tensorflow.keras.optimizers import Adam,SGD
from util import read_directory
from optim import HUATU,Optim
import random
from tensorflow.keras.utils import plot_model
import os
# os.environ["CUDA_DEVICE_ORDER"] = "PCI_BUS_ID"
# os.environ["CUDA_VISIBLE_DEVICES"] = "0"
# In[] 加载数据
height=64
width=64
# 时频图---2D-CNN输入
x_train1,y_train=read_directory('image/train',height,width,normal=1)
x_valid1,y_valid=read_directory('image/test',height,width,normal=1)
# In[] 采用优化算法优化CNN超参数
mode = 'SO' # SA PSO SASO SO 四种算法 用哪种方法就改成哪种方法的名字 比如用SASIO就改成 mode = 'SASO'
optimization_again = False # 为 False 就直接调用之前的优化的参数 True就重新优化
if optimization_again:
best_solution, best_fitness, fitness_history = Optim(mode,x_train1,y_train,x_valid1,y_valid)
np.savez('result/'+mode+'para.npz',best_solution=best_solution,best_fitness=best_fitness,fitness_history=fitness_history)
else:
data = np.load('result/'+mode+'para.npz')
best_solution, best_fitness, fitness_history = data['best_solution'].reshape(-1,),data['best_fitness'],data['fitness_history'].reshape(-1,)
HUATU(mode,fitness_history)
# In[] 利用优化得到的参数进行搭建模型
# 构建网络结构是 输入层-2d卷积-池化-Bn-SE注意力-2d卷积-池化-Bn-SE注意力-2d卷积-池化-Bn-SE注意力-全连接-输出层
# 参数设置，有点多 分别
# 是第1卷积层的卷积核数量k1_num、卷积核大小k1_size*k1_size
# 是第2卷积层的卷积核数量k2_num、卷积核大小k2_size*k2_size
# 是第3卷积层的卷积核数量k3_num、卷积核大小k3_size*k3_size
# 是第4卷积层的卷积核数量k4_num、卷积核大小k4_size*k4_size
# 是第5卷积层的卷积核数量k5_num、卷积核大小k4_size*k5_size
# 全连接层的节点数fc1,fc2
# 训练次数 epoches，batchsize,学习率 learningrate
tf.random.set_seed(0)
np.random.seed(0)
para = [best_solution[i] if i==0 else int(best_solution[i]) for i in range(len(best_solution))]
print(mode,'优化得到的超参数:',para)
learningrate=para[0]#学习率
epochs = para[1]#迭代次数
batch_size = para[2] #BATCH_SIZE
k1_num,k1_size,k2_num,k2_size,k3_num,k3_size,k4_num,k4_size,k5_num,k5_size,fc1,fc2=para[3:]
model= ProposedNet2D(k1_num,5,k2_num,k2_size,k3_num,k3_size,k4_num,k4_size,k5_num,k5_size,fc1,fc2)
model.summary()
checkpoint = tf.keras.callbacks.ModelCheckpoint(filepath='model/'+mode+'best', monitor='val_accuracy', verbose=0, save_best_only=True, mode = 'max')
model.compile(optimizer=Adam(learning_rate=learningrate),
loss=SparseCategoricalCrossentropy(),
metrics=['accuracy'] )
# exit()
# In[] 训练
train_again = False # 为 False 的时候就直接加载训练好的模型进行测试 True就重新训练
# 训练模型
if train_again:
history = model.fit(x_train1,y_train, epochs=epochs,
validation_data=(x_valid1, y_valid),
batch_size=batch_size, verbose=1,callbacks= [checkpoint])
model=tf.keras.models.load_model('model/'+mode+'best',custom_objects = {'SE_AttentionLayer': SE_AttentionLayer})
# 画loss曲线
plt.figure()
plt.ylabel('MSE')
plt.xlabel('Epoch')
plt.plot(history.history['loss'], label='training')
# plt.plot(history.history['val_loss'], label='validation')
plt.title('loss curve')
plt.legend()
plt.savefig('result/'+mode+'model2d_loss_curve.jpg')
else: # 加载模型
model=tf.keras.models.load_model('model/'+mode+'best',custom_objects = {'SE_AttentionLayer': SE_AttentionLayer})
#
test_pred = model.predict(x_valid1)
pred_labels = np.argmax(test_pred,1)
acc=np.sum(pred_labels==y_valid)/len(y_valid)
print('分类精度为:%f '%(acc*100),'%')
C1= confusion_matrix(y_valid, pred_labels)
print('混淆矩阵:')
print(C1)
xtick=[str(i) for i in range(10)]
ytick=[str(i) for i in range(10)]
plt.figure()
# C2=C1/C1.sum(axis=1)
sns.heatmap(C1,fmt='g', cmap='Blues',annot=True,cbar=False,xticklabels=xtick, yticklabels=ytick)
plt.xlabel('Predict label')
plt.ylabel('True label')
plt.title('Confusion_matrix')
plt.savefig('result/'+mode+'Confusion_matrix_2d.jpg')
plt.show()
每一种方法都优化完成后，可以运行下面的代码进行适应度曲线的对比分析，这里要说明一下，每个人优化的结果都不一定一样，这是因为随机数存在的原因，其次因为数据简单，所以都可能优化到100%，所以可以简单的改一改保存的结果再画图，你懂的。：
import numpy as np
import matplotlib.pyplot as plt
data = np.load('result/SApara.npz')
best_solution0, best_fitness0, fitness_history0 = data['best_solution'].reshape(-1,),data['best_fitness'],data['fitness_history'].reshape(-1,)
data = np.load('result/PSOpara.npz')
best_solution1, best_fitness1, fitness_history1 = data['best_solution'].reshape(-1,),data['best_fitness'],data['fitness_history'].reshape(-1,)
data = np.load('result/SOpara.npz')
best_solution2, best_fitness2, fitness_history2 = data['best_solution'].reshape(-1,),data['best_fitness'],data['fitness_history'].reshape(-1,)
data = np.load('result/SASOpara.npz')
best_solution3, best_fitness3, fitness_history3 = data['best_solution'].reshape(-1,),data['best_fitness'],data['fitness_history'].reshape(-1,)
# Plot optimization process
plt.figure(figsize=(10, 6))
plt.plot(fitness_history0, linewidth=2,label='SA')
plt.plot(fitness_history1, linewidth=2,label='PSO')
plt.plot(fitness_history2, linewidth=2,label='SO')
plt.plot(fitness_history3, linewidth=2,label='SASO')
plt.legend()
plt.title("Optimization Process")
plt.xlabel("Iteration")
plt.ylabel("Best Fitness")
plt.grid(True)
plt.savefig('result/网络超参数寻优适应度函数对比.jpg')
plt.show()
print("SA最佳解决方案：", [best_solution0[i] if i==0 else int(best_solution0[i]) for i in range(len(best_solution0))])
print("SA最佳适应度：", best_fitness0)
print("PSO最佳解决方案：", [best_solution1[i] if i==0 else int(best_solution1[i]) for i in range(len(best_solution1))])
print("PSO最佳适应度：", best_fitness1)
print("SO最佳解决方案：", [best_solution2[i] if i==0 else int(best_solution2[i]) for i in range(len(best_solution2))])
print("SO最佳适应度：", best_fitness2)
print("SASO最佳解决方案：", [best_solution3[i] if i==0 else int(best_solution3[i]) for i in range(len(best_solution3))])
print("SASO最佳适应度：", best_fitness3)
3 其他对比
除了上面之外一般还有验证本文方法在不同噪声等级下的有效性分析，就是通过给数据里面添加同等长度的随机数。
还有因为我们用的模型是VGG13+SE的网络结构，还可以验证与不加SE注意力的对比。
4 结语
此博客是轴承故障【从数据处理到模型优化全流程】系列的最后一篇，主要是从多个角度对比所提方法的有效性。我只做了两项，更多的对比等待读者去发掘。
获取更多内容请点击【专栏】，您的点赞与收藏是我持续更新【Python神经网络1000个案例分析】的动力。