Python遗传规划实战：用gplearn和DEAP解决符号回归问题（附完整代码）

Python遗传规划实战：用gplearn和DEAP解决符号回归问题

在机器学习领域，符号回归是一种独特的技术，它不依赖于预设的模型结构，而是通过进化计算自动发现数据背后的数学表达式。本文将深入探讨如何使用Python生态中的两大工具——gplearn和DEAP库，从零开始构建符号回归解决方案。

1. 符号回归与遗传规划基础

符号回归（Symbolic Regression）是一种特殊的回归分析方法，与传统线性回归或神经网络不同，它不需要预先定义模型形式。遗传规划（Genetic Programming, GP）则是实现符号回归的核心算法，通过模拟自然进化过程来自动发现最优数学模型。

遗传规划的工作流程包括：

• 种群初始化：随机生成一组数学表达式作为初始种群
• 适应度评估：计算每个表达式对数据的拟合程度
• 选择操作：根据适应度选择优秀个体进入下一代
• 遗传操作：通过交叉和变异产生新的表达式
• 迭代优化：重复上述过程直到满足终止条件

与传统机器学习方法相比，符号回归的优势在于：

• 可解释性强：生成的数学表达式直观易懂
• 无需特征工程：自动发现变量间的关系
• 灵活性高：适应各种复杂非线性关系

2. gplearn库实战应用

gplearn是一个专为符号回归设计的Python库，其API设计与scikit-learn保持高度一致，使得熟悉scikit-learn的用户能够快速上手。

2.1 环境配置与数据准备

首先安装gplearn库：

pip install gplearn

准备示例数据，我们创建一个具有明确数学关系的合成数据集：

import numpy as np
from sklearn.model_selection import train_test_split

# 生成合成数据
np.random.seed(42)
X = np.random.uniform(-10, 10, size=(1000, 2))
y = X[:, 0]**2 + 3*X[:, 1] - 5  # 真实关系式

# 划分训练测试集
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2)

2.2 基础模型训练

使用gplearn的SymbolicRegressor进行基础训练：

from gplearn.genetic import SymbolicRegressor

# 初始化符号回归器
est = SymbolicRegressor(population_size=1000,
                        generations=20,
                        stopping_criteria=0.01,
                        random_state=42)

# 训练模型
est.fit(X_train, y_train)

# 评估模型
print("训练R²:", est.score(X_train, y_train))
print("测试R²:", est.score(X_test, y_test))
print("最佳表达式:", est._program)

2.3 高级配置与优化

gplearn提供了丰富的参数用于控制进化过程：

种群参数配置

est = SymbolicRegressor(
    population_size=5000,      # 种群规模
    generations=50,            # 进化代数
    tournament_size=20,        # 锦标赛选择规模
    init_method='half and half', # 初始化方法
    init_depth=(2,6),          # 初始树深度范围
    function_set=('add', 'sub', 'mul', 'div', 'sqrt', 'log')  # 函数集
)

遗传操作概率调整

est.set_params(
    p_crossover=0.7,          # 交叉概率
    p_subtree_mutation=0.1,   # 子树变异概率
    p_hoist_mutation=0.05,    # 提升变异概率
    p_point_mutation=0.1,     # 点变异概率
    parsimony_coefficient=0.01 # 简洁性系数
)

自定义适应度函数

from gplearn.fitness import make_fitness

def _mape(y, y_pred, w):
    """自定义MAPE指标"""
    diffs = np.abs((y - y_pred) / np.maximum(1e-6, np.