利用遗传算法（GA）与模拟退火算法（SA）求目标函数最小值

标签：模拟退火 parents GA fitness np 遗传算法 offspring best population

要求实现一个演化计算的算法，求测试函数的最小值。

要求：群体规模NP=100；最大迭代次数不超过3000代。或者，总的计算次数小于100*3000。算法需独立运行30次，并记录进化的过程。

一、遗传算法原理

遗传算法（Genetic Algorithm, GA）是模拟达尔文生物进化论的自然选择和遗传学机理的生物进化过程的计算模型，是一种通过模拟自然进化过程搜索最优解的方法。

其主要特点是直接对结构对象进行操作，不存在求导和函数连续性的限定；具有内在的隐并行性和更好的全局寻优能力；采用概率化的寻优方法，不需要确定的规则就能自动获取和指导优化的搜索空间，自适应地调整搜索方向。

遗传算法以一种群体中的所有个体为对象，并利用随机化技术指导对一个被编码的参数空间进行高效搜索。其中，选择、交叉和变异构成了遗传算法的遗传操作；参数编码、初始群体的设定、适应度函数的设计、遗传操作设计、控制参数设定五个要素组成了遗传算法的核心内容。

算法原理

初始化种群：

随机生成一定数量的初始解（个体），构成种群。

$P(0) = \{ \mathbf{x}_1, \mathbf{x}_2, \ldots, \mathbf{x}_{NP} \}$

适应度评估：
$f(\mathbf{x}_i) = \text{objective\_function}(\mathbf{x}_i)$

选择操作：

根据适应度值选择优秀的个体作为父代

$p_i = \frac{f(\mathbf{x}_i)}{\sum_{j=1}^{NP} f(\mathbf{x}_j)}$

交叉操作：

通过交叉操作生成新的子代，模拟自然界中的基因重组。常用的交叉方法有单点交叉、多点交叉、均匀交叉等。

$\begin{cases} \text{offspring}[k, 0:\text{crossover\_point}] = \text{parent1}[0:\text{crossover\_point}] \\ \text{offspring}[k, \text{crossover\_point}:] = \text{parent2}[\text{crossover\_point}:] \end{cases}$

变异操作：

对交叉生成的子代进行变异操作，模拟自然界中的基因突变。变异操作有助于维持种群的多样性，防止早熟收敛。

种群更新：

用新生成的子代替换种群中的部分个体，形成新一代种群。

$P(t+1) = \{ \text{offspring} \} \cup \{ \text{elite\_individuals} \}$

终止条件：

达到全局最优。

$\text{if } t \geq \text{max\_generations} \text{ or } \text{convergence\_criterion\_met}$

伪代码

初始化种群 P(0)
评估 P(0) 的适应度
t = 0
while 未达到终止条件 do
    从 P(t) 中选择父代
    执行交叉和变异操作生成子代
    评估子代的适应度
    从 P(t) 和子代中选择个体形成 P(t+1)
    t = t + 1
end while
输出找到的最优解

二、目标函数实现

构造目标函数

对于给定的输入向量x_{i}，函数计算每个元素xi的负值乘以它的平方根的正弦值之和。

def objective_function(x):
    return -sum(xi * np.sin(np.sqrt(abs(xi))) for xi in x)

初始化种群

def initialize_population(pop_size, dim, bounds):
    return np.random.uniform(bounds[0], bounds[1], (pop_size, dim))

initialize_population 函数生成一个随机种群，种群的个体数为 pop_size，每个个体的维度为 dim，个体值在 bounds 范围内均匀分布。

评估种群

def evaluate_population(population):
    return np.array([objective_function(ind) for ind in population])

evaluate_population 函数计算种群中每个个体的适应度值，返回一个包含适应度值的数组。

选择父代

def select_parents(population, fitness, num_parents):
    parents = np.empty((num_parents, population.shape[1]))
    for i in range(num_parents):
        min_fitness_idx = np.argmin(fitness)
        parents[i, :] = population[min_fitness_idx, :]
        fitness[min_fitness_idx] = np.inf  # 排除已选父代
    return parents

select_parents 函数基于适应度值选择 num_parents 个父代。适应度值越高的个体被选中的概率越大。被选中的个体在后续迭代中被排除

交叉操作

def crossover(parents, offspring_size):
    offspring = np.empty(offspring_size)
    crossover_point = np.uint8(offspring_size[1] / 2)
    for k in range(offspring_size[0]):
        parent1_idx = k % parents.shape[0]
        parent2_idx = (k + 1) % parents.shape[0]
        offspring[k, 0:crossover_point] = parents[parent1_idx, 0:crossover_point]
        offspring[k, crossover_point:] = parents[parent2_idx, crossover_point:]
    return offspring

crossover 函数通过交叉操作生成新的子代。每对父代的前半部分和后半部分分别组合生成新的个体。

变异操作

def mutate(offspring, mutation_rate, bounds):
    for idx in range(offspring.shape[0]):
        for _ in range(np.uint8(offspring.shape[1] * mutation_rate)):
            i = np.random.randint(0, offspring.shape[1])
            random_value = np.random.uniform(bounds[0], bounds[1])
            offspring[idx, i] = random_value
    return offspring

mutate 函数对交叉产生的子代进行变异操作。变异是指随机改变子代中某些基因的值，以增加种群的多样性。变异率由 mutation_rate 参数控制。

主循环

# 参数设置
pop_size = 800  # 种群大小
dim = 30   # 维度；等于n=30
bounds = (-500, 500)  # 值范围
num_generations = 500  # 迭代次数
num_parents_mating = 30  # 父代数量
mutation_rate = 0.1  # 变异率
runs = 30  # 运行次数

best_solutions = []
best_fitness_values = []
evolution_records = []
population_evolution_records = []

for run in range(runs):
    print(f"Starting run {run + 1}/{runs}...")
    # 初始化种群
    population = initialize_population(pop_size, dim, bounds)
    best_fitness_record = []
    population_fitness_record = []

    # 进化过程
    for generation in tqdm(range(num_generations), desc=f"Run {run + 1} Progress"):
        fitness = evaluate_population(population)
        parents = select_parents(population, fitness, num_parents_mating)
        offspring_crossover = crossover(parents, (pop_size - parents.shape[0], dim))
        offspring_mutation = mutate(offspring_crossover, mutation_rate, bounds)

        # 精英保留策略
        elite_individual = population[np.argmin(fitness)]
        population[parents.shape[0]:, :] = offspring_mutation
        population[0, :] = elite_individual

        # 记录当前代的最佳适应度值
        fitness = evaluate_population(population)
        best_fitness_record.append(np.min(fitness))
        population_fitness_record.append(fitness.copy())

    # 当前运行的最佳解决方案
    best_match_idx = np.argmin(fitness)
    best_solution = population[best_match_idx]
    best_fitness = np.min(fitness)

    best_solutions.append(best_solution)
    best_fitness_values.append(best_fitness)
    evolution_records.append(best_fitness_record)
    population_evolution_records.append(population_fitness_record)
    print(f"Run {run + 1} completed. Best solution fitness = {best_fitness}")

# 输出结果
for i in range(runs):
    print(f"Run {i + 1}: Best solution fitness = {best_fitness_values[i]}")

# 绘制最佳解决方案适应度值的进化过程
plt.figure(figsize=(10, 6))
for i in range(runs):
    plt.plot(evolution_records[i], label=f'Run {i + 1}')
plt.title('Evolution of Best Solution Fitness over Generations')
plt.xlabel('Generation')
plt.ylabel('Best Fitness Value')
plt.legend()
plt.show()

# 绘制种群适应度分布
for run in range(runs):
    fig, ax = plt.subplots(figsize=(10, 6))
    population_fitness_record = population_evolution_records[run]
    for gen in range(0, num_generations, num_generations // 10):
        fitness_values = population_fitness_record[gen]
        ax.hist(fitness_values, bins=50, alpha=0.5, label=f'Generation {gen}')
    ax.set_title(f'Population Fitness Distribution for Run {run + 1}')
    ax.set_xlabel('Fitness Value')
    ax.set_ylabel('Frequency')
    plt.savefig()
    ax.legend()
    plt.show()

三、测试函数结果分析

因为提前知道最小值为-12569；所以需要手动进行调参。

先观测数据

生成30维度样本并计算目标函数值和直方图分布，观测数据分布。

# -*- coding=utf-8 -*-

import numpy as np
import matplotlib.pyplot as plt

def objective_function(x):
    return -sum((xi * np.sin(np.sqrt(abs(xi)))) for xi in x)

# 生成高维随机样本
num_samples = 50000
dim_high = 30  
bounds = (-500, 500)

# 计算目标函数值
for i in range(50):
    samples_high = np.random.uniform(bounds[0], bounds[1], (num_samples, dim_high))
    objective_values_high = np.array([objective_function(sample) for sample in samples_high])

    print("随机样本的目标函数最小值:", np.min(objective_values_high))
    # 绘制目标函数值的分布直方图
    plt.hist(objective_values_high, bins=50, alpha=0.7)
    plt.title('High-Dimensional Objective Function Value Distribution')
    plt.xlabel('Objective Function Value')
    plt.ylabel('Frequency')
    plt.show()

好吧，这个啥都看不出来。。

探究父代数量是如何影响目标函数值

先看看这个参数有什么影响

# 不同父代数量设置
parent_mating_settings = [10, 30, 50]

results = {}

for num_parents_mating in parent_mating_settings:
    best_solutions = []
    best_fitness_values = []
    evolution_records = []

    for run in range(runs):
        print(f"Starting run {run + 1}/{runs} with num_parents_mating = {num_parents_mating}...")
        # 初始化种群
        population = initialize_population(pop_size, dim, bounds)
        best_fitness_record = []

        # 进化过程
        for generation in tqdm(range(num_generations), desc=f"Run {run + 1} Progress"):
            fitness = evaluate_population(population)
            parents = select_parents(population, fitness, num_parents_mating)
            offspring_crossover = crossover(parents, (pop_size - parents.shape[0], dim))
            offspring_mutation = mutate(offspring_crossover, mutation_rate, bounds)

            # 精英保留策略
            elite_individual = population[np.argmin(fitness)]
            population[parents.shape[0]:, :] = offspring_mutation
            population[0, :] = elite_individual

            # 记录当前代的最佳适应度值
            fitness = evaluate_population(population)
            best_fitness_record.append(np.min(fitness))

        # 当前运行的最佳解决方案
        best_match_idx = np.argmin(fitness)
        best_solution = population[best_match_idx]
        best_fitness = np.min(fitness)

        best_solutions.append(best_solution)
        best_fitness_values.append(best_fitness)
        evolution_records.append(best_fitness_record)
        print(f"Run {run + 1} completed. Best solution fitness = {best_fitness}")

    results[num_parents_mating] = {
        'best_solutions': best_solutions,
        'best_fitness_values': best_fitness_values,
        'evolution_records': evolution_records
    }

整体来看父代数更多确实效果要好些，但是好像，30跟50也没什么区别

探究种群大小是如何影响目标函数值

另外一个参数

dim = 30  # 维度
bounds = (-500, 500)  # 值范围
num_generations = 300  # 迭代次数
num_parents_mating = 30  # 父代数量
mutation_rate = 0.1  # 变异率
runs = 3  # 每组实验的运行次数

# 不同种群数量设置
population_sizes = [600, 1000, 1400]
results = {}

for pop_size in population_sizes:
    best_solutions = []
    best_fitness_values = []
    evolution_records = []

    for run in range(runs):
        print(f"Starting run {run + 1}/{runs} with population size = {pop_size}...")
        # 初始化种群
        population = initialize_population(pop_size, dim, bounds)
        best_fitness_record = []

        # 进化过程
        for generation in tqdm(range(num_generations), desc=f"Run {run + 1} Progress"):
            fitness = evaluate_population(population)
            parents = select_parents(population, fitness, num_parents_mating)
            offspring_crossover = crossover(parents, (pop_size - parents.shape[0], dim))
            offspring_mutation = mutate(offspring_crossover, mutation_rate, bounds)

            # 精英保留策略
            elite_individual = population[np.argmin(fitness)]
            population[parents.shape[0]:, :] = offspring_mutation
            population[0, :] = elite_individual

            # 记录当前代的最佳适应度值
            fitness = evaluate_population(population)
            best_fitness_record.append(np.min(fitness))

        # 当前运行的最佳解决方案
        best_match_idx = np.argmin(fitness)
        best_solution = population[best_match_idx]
        best_fitness = np.min(fitness)

        best_solutions.append(best_solution)
        best_fitness_values.append(best_fitness)
        evolution_records.append(best_fitness_record)
        print(f"Run {run + 1} completed. Best solution fitness = {best_fitness}")

    results[pop_size] = {
        'best_solutions': best_solutions,
        'best_fitness_values': best_fitness_values,
        'evolution_records': evolution_records
    }

Pop Size越高整体也是取值越优，选1000吧那就。

————————————————————

结果还是不好啊。

数学真难，痛苦面具

另外一个变异率，就不调了，感觉GA效果好差啊。试一下模拟退火算法吧。

————————————————————

四、模拟退火算法

用模拟退火算法

原理和局限性

模拟退火算法（Simulated Annealing, SA，1983）是一种用于全局优化问题的概率性算法。模拟物理退火过程中的热力学原理。其基本思想是通过在解空间内随机采样，并在一定概率下接受较差的解，以避免陷入局部最优解，从而逐步找到全局最优解。

(wiki百科)

模拟退火算法的主要步骤如下：

初始解和温度设定：从一个初始解 X_0开始，并设定一个初始温度 T_0。
邻域搜索：在当前解的邻域内随机生成一个新的解 x′。
解的接受准则：
- 如果新的解 x′优于当前解 x（即目标函数值更小），则接受 x′ 作为当前解。
- 如果新的解 x′不优于当前解 x，则以一定的概率接受 x′，该概率通常为 $exp(\frac{f(x)-f(x')))}{T})$ ，其中 f(x) 和 f(x′)分别是当前解和新解的目标函数值，T 为当前温度。
温度降低：按照某种退火策略逐步降低温度 T，例如线性降温或指数降温。
迭代过程：重复上述步骤，直到达到停止条件（如温度降低到某个阈值或达到最大迭代次数）。

模拟退火算法的优劣势

优势

全局优化能力强：由于在一定概率下接受较差解，模拟退火算法可以有效避免陷入局部最优，具有较强的全局搜索能力。
简单易实现：算法简单易于实现，不需要对问题结构有太多假设。
适用范围广：可以应用于各种复杂的优化问题，包括离散和连续问题。

劣势

计算效率较低：模拟退火算法的收敛速度较慢，尤其是在解决大规模优化问题时，计算时间可能较长。
参数敏感性：算法的性能依赖于初始温度、退火策略等参数的设置，不同问题需要不同的参数调优。
结果的不确定性：由于算法具有随机性，每次运行的结果可能不同，需要多次运行取最优结果

代码实现

调用scipy的库，真是简洁

# -*- coding=utf-8 -*-

import numpy as np
from scipy.optimize import dual_annealing
import matplotlib.pyplot as plt

history = []

def test_function(x):
    val = -np.sum(x * np.sin(np.sqrt(np.abs(x))))
    history.append(val)
    return val

bounds = [(-500, 500)] * 30

result = dual_annealing(test_function, bounds)

plt.figure(figsize=(12, 6))
plt.plot(history, label="Function Value")
plt.xlabel("Function Calls")
plt.ylabel("Function Value")
plt.title("Evolution of the Function Value During Optimization")
plt.legend()
plt.grid(True)
plt.show()

print(result)

结果

min为-12569.48；X_i=420.97

好使

标签：模拟退火,parents,GA,fitness,np,遗传算法,offspring,best,population
From： https://blog.csdn.net/weixin_46636042/article/details/139757808