遗传编程进化股票预测模型：Python实现因子挖掘与滚动回测

简介一套基于 Python 的算法交易学习代码包采用遗传编程与遗传算法两种策略预测苹果股价面向量化交易入门者、金融编程开发者以及对进化计算在金融场景应用的读者可用于理解股票价格预测与涨跌判断的算法建模过程。包体共40个文件、约1.21MB以21个Python脚本和8个CSV数据文件为主体辅以TXT说明与Markdown文档脚本实现核心算法、数据转换与预测调用数据文件提供NASDAQ、苹果等历史行情样本说明文档则记录设计思路与运行指南。已有330人学习下载。资源内含两个独立子模块遗传程序通过进化树形结构最小化预测价与真实价的误差遗传算法则进化二进制字符串预测次日股价涨跌同时保留历史行情数据与数据处理脚本便于读者结合自述文档快速复现实验、对比不同进化策略的效果也可作为课程作业或量化项目的参考框架。1. 当遗传编程开始“进化”出股票预测模型回测就不再是玄学股票价格预测这件事绝大多数人一开始就走错了方向——他们试图用一组固定的参数、固定的公式去拟合市场然后发现模型在样本内漂亮得像艺术品一到实盘就崩得亲妈都不认识。原因很简单市场不是静止的你今天找到的那条最佳拟合曲线明天就不再是最佳的了。遗传编程Genetic Programming, GP和遗传算法Genetic Algorithm, GA在这里的逻辑恰好反了过来——不再手动设计特征和公式而是让一段 Python 程序在“选择—交叉—变异”的循环里自己长出一棵表达式树去逼近价格序列背后的非线性关系。算法交易程序里大量使用这种思路来做因子挖掘和信号生成并不是因为它能“预测”未来而是因为它能在给定的历史窗口里自适应地找到当前环境下解释力最强的数学表达式。本文就用纯 Python 代码带你从头实现一个可运行的 GP 预测框架覆盖种群初始化、适应度评估、遗传操作到滚动回测的完整链路最终落到参数怎么调、过拟合怎么防、实盘前最少要做哪几项验证。2. 遗传算法与遗传编程的分工参数搜索和结构搜索先搞清你用的是哪个2.1 为什么说遗传算法是“调参器”遗传编程才是“公式发明机”很多文章把 GA 和 GP 混着用实际工程里这是两套完全不同的搜索逻辑。遗传算法操作的对象是固定长度的编码向量——比如神经网络的权重、移动平均的周期、RSI 的阈值——它不改变问题的结构只搜索最优的参数组合。而遗传编程操作的对象是程序树通常用表达式树表示树的深度和形状都可以变异天然适合去搜索“价格和成交量之间到底该用什么公式组合”这种开放性问题。在算法交易程序中一个非常典型的实践路径是先用遗传编程从 OHLCV 数据里挖掘候选因子表达式再用遗传算法对挖掘出的因子做参数寻优和组合权重分配。GP 管“用什么公式”GA 管“公式里的参数取多少”两者串联形成一个两阶段的自动建模管线。2.2 表达式树的最小可运行表示用嵌套元组模拟二叉树我不打算引入 deap 库先讲一堆 API那样你会被框架细节淹没。先用纯 Python 实现一棵极简表达式树它同时是 GA 里的“基因”和 GP 里的“程序”。import random import numpy as np import pandas as pd # 函数集叶子节点以外允许出现的操作符 FUNCTIONS { add: lambda x, y: x y, sub: lambda x, y: x - y, mul: lambda x, y: x * y, div: lambda x, y: np.divide(x, y, outnp.ones_like(x), wherey ! 0), max: lambda x, y: np.maximum(x, y), min: lambda x, y: np.minimum(x, y), } # 终止符集可以直接作为叶子的字段来自行情数据 TERMINALS [open, high, low, close, volume, returns_5] def random_expr(depth, max_depth4): 生成一棵随机的表达式树表示为一个嵌套元组 if depth max_depth or random.random() 0.3: return random.choice(TERMINALS) func random.choice(list(FUNCTIONS.keys())) left random_expr(depth 1, max_depth) right random_expr(depth 1, max_depth) return (func, left, right)这段代码里最关键的设计是div的处理方式——np.divide的where参数在分母为 0 时返回 1避免产生 NaN 污染整个适应度评估。表达式树是嵌套元组而不是自定义类好处是递归遍历极其方便pickle序列化也能直接走默认协议这在后续保存种群时省了很多事。2.3 适应度函数预测准确率不是目标IC 和 Rank IC 才是如果你拿“预测值和真实值的方向准确率”做适应度模型会偏向预测上涨——因为 A 股单边行情的时段里上涨样本天然多于下跌样本。更健壮的做法是用信息系数 IC也就是预测值与未来收益之间的 Spearman 秩相关。它衡量的是排序能力天然剔除了量纲和极端值干扰。from scipy.stats import spearmanr def evaluate_fitness(expr, df, horizon5): 计算表达式树在给定数据上的因子值并与未来 horizon 日收益做秩相关 feature eval_expr(expr, df) future_ret df[close].pct_change(horizon).shift(-horizon) valid pd.concat([feature, future_ret], axis1).dropna() if len(valid) 30: return -1.0 ic, _ spearmanr(valid.iloc[:, 0], valid.iloc[:, 1]) if np.isnan(ic): return -1.0 return ic def eval_expr(expr, df): 朴素递归解释器在 DataFrame 上整列计算表达式值 if isinstance(expr, str): return df[expr].to_numpy(dtypefloat) func, left, right expr lv eval_expr(left, df) rv eval_expr(right, df) return FUNCTIONS[func](lv, rv)eval_expr里没有做逐行循环而是整列向量化计算。数据量在几万行以下时性能没有差别一旦到分钟线级别单标的几十万行向量化的优势就非常明显。shift(-horizon)是关键——它把未来第 5 天的收益对齐到今天的因子值上模拟的是“今天产生信号、未来兑现收益”的真实交易场景。提示dropna之后样本数必须设下限。如果某段时间停牌或数据缺失严重剩余样本太少Spearman 相关系数的方差会大到让适应度完全失去意义。3. 从随机种群到逐代进化GP 核心循环的 Python 实现3.1 选择、交叉、变异的算子设计与边界控制遗传编程里选择压力太大种群会早熟太小则收敛极慢。我用的是锦标赛选择加精英保留每代先复制适应度最高的 2 个个体到下一代再每轮随机抽 3 个个体挑其中适应度最高的进入配对池。这个“3”是标准参数工程上很少需要动。def tournament_select(population, fitness, k3): idx random.sample(range(len(population)), k) best_idx max(idx, keylambda i: fitness[i]) return population[best_idx] def crossover(expr1, expr2, max_depth6): 交换两棵树的随机子树保留结构多样性 expr1, expr2 deepcopy(expr1), deepcopy(expr2) if isinstance(expr1, str) or isinstance(expr2, str): return expr1, expr2 node1 select_random_node(expr1) node2 select_random_node(expr2) replace_node(expr1, node1, node2) replace_node(expr2, node2, node1) return expr1, expr2 def mutate(expr, max_depth6): 将某个子树替换为随机新生成子树 expr deepcopy(expr) if isinstance(expr, str): return random_expr(0, max_depth) node select_random_node(expr) new_subtree random_expr(0, max(max_depth - depth_of(node), 1)) replace_node(expr, node, new_subtree) return expr交叉算子有一个非常隐蔽的坑替换子树之后可能导致整棵树的深度暴涨。所以replace_node是替换整个子树引用而新子树的深度被限制在max_depth - depth_of(node)以内——这样从任意节点出发整体深度都不会超过预设上限。否则跑 20 代以后树会膨胀到几十层单次评估耗时会翻几十倍而且表达式会严重过拟合历史噪声。3.2 种群进化的主循环记录收敛曲线而非单代最优主循环的写法基本统一我习惯把每代的最优适应度、平均适应度、种群多样性三个指标都记录下来。多样性用“所有个体两两之间的表达式距离”来算太贵退而求其次用“唯一表达式数量占种群比例”近似。def evolve(df, pop_size60, generations20, max_depth5): population [random_expr(0, max_depth) for _ in range(pop_size)] history {best: [], avg: [], unique_ratio: []} for gen in range(generations): fitness [evaluate_fitness(ind, df) for ind in population] best_idx int(np.argmax(fitness)) history[best].append(fitness[best_idx]) history[avg].append(float(np.mean(fitness))) history[unique_ratio].append(len(set(population)) / pop_size) new_population [deepcopy(population[best_idx])] # 精英 while len(new_population) pop_size: p1 tournament_select(population, fitness) p2 tournament_select(population, fitness) child1, child2 crossover(p1, p2, max_depth) if random.random() 0.8 else (p1, p2) if random.random() 0.1: child1 mutate(child1, max_depth) if random.random() 0.1: child2 mutate(child2, max_depth) new_population.extend([child1, child2]) population new_population[:pop_size] best_expr population[int(np.argmax( [evaluate_fitness(ind, df) for ind in population] ))] return best_expr, history这段代码里唯一需要特别强调的是new_population[:pop_size]这行——while循环里每轮会生成两个孩子当种群大小是奇数时extend会让数量超出切片截断是保证种群规模恒定的兜底方案。不要用if len(new_population) 2 pop_size去判断再单独处理代码会变得很难读。3.3 种群多样性崩溃的早期信号unique_ratio 掉到 0.2 以下很多第一次跑 GP 的人会盯着 best 曲线看到它一路抬升就以为没问题。经验是best 曲线上升的同时一定要看 unique_ratio。如果它快速跌到 0.2 以下说明种群只剩少数几个个体的大量拷贝后续的交叉变异都是在近亲繁殖解早就被锁死在一个局部最优附近了。此时应该做的不是继续跑更多代而是提高变异率或者每 N 代随机注入几个全新个体。参数常用区间特征影响pop_size40 ~ 100太小容易早熟太大单代耗时线性增长generations15 ~ 50超过 50 代收益递减纯靠这在分钟线上挖掘意义不大crossover_prob0.7 ~ 0.9主流搜索动力过低会让种群陷入停滞mutation_prob0.05 ~ 0.2保持多样性的关键前期 0.05 后期 0.15max_depth3 ~ 6超过 6 几乎必过拟合除非你有极强的样本内逻辑约束这几组参数来自多标的、多周期的回测经验。对于单标的日线数据max_depth一旦超过 6表达式树很容易组合出对历史极端行情比如某一天暴涨暴跌敏感的路径而这类路径在未来几乎不可能重演。参数不是死的但边界经验值得你从这些区间起步。注意unique_ratio计算set(population)时元组嵌套的表达式树可以直接哈希完全合法。但树一旦用自定义类表示就必须实现__hash__和__eq__否则这个指标根本算不了。4. 滚动窗口训练与样本外验证让 GP 预测程序从“过拟合玩具”变成“可上线工具”4.1 为什么固定训练集全部样本的 GP 模型是纸老虎用全量历史数据进化出的公式在回测里对历史数据有天然的“作弊”优势——适应度是从同一段数据里算出来的等于是拿着答案去考试。实盘里价格序列有状态切换2015 年的规律在 2018 年可能完全失效。正确做法是滚动窗口每次只用前 N 天的数据进化种群然后对后 M 天做样本外预测记录预测值与实际值的 IC窗口往前滑动重复这个过程。4.2 滚动回测的工程框架walk-forward 验证的完整代码def walk_forward_evolve(df, train_days250, test_days20, **evolve_kwargs): 滚动回测训练窗口 train_days验证窗口 test_days 返回样本外 IC 序列和每段最佳表达式 results [] exprs [] start 0 while start train_days test_days len(df): train_df df.iloc[start:start train_days] test_df df.iloc[start train_days:start train_days test_days] best_expr, _ evolve(train_df, **evolve_kwargs) test_feature eval_expr(best_expr, test_df) test_ret test_df[close].pct_change(5).shift(-5) valid pd.concat([test_feature, test_ret], axis1).dropna() if len(valid) 10: oos_ic, _ spearmanr(valid.iloc[:, 0], valid.iloc[:, 1]) results.append(oos_ic if not np.isnan(oos_ic) else 0.0) else: results.append(0.0) exprs.append(best_expr) start test_days # 不重叠步进 return results, exprs这个框架有几个必须注意的细节。pct_change(5).shift(-5)计算的是第 t5 天相对于第 t 天的收益率所以预测目标是“未来 5 日的累计收益”。代码里没有做去极值和中性化处理那是实盘前的进阶话题第一步先跑通 IC 评估闭环。滚动窗口的步进用了test_days不重叠的方式减少计算量但如果数据量充足用更小的步长比如 5 天能让结果更平滑代价是训练次数翻倍。4.3 样本外 IC 分布能上线的最低门槛不是均值而是胜率跑完 walk-forward 之后你手里是一组长度等于窗口数的样本外 IC 序列。不要只看均值和标准差直接把所有 OOS IC 的符号分布打出来看如果你有 50 个窗口其中 35 个 IC 为正、均值 0.06这说明模型在 70% 的时间里具备微弱的排序能力可以进入下一步资金曲线模拟如果均值是 0.04但只有 22 个窗口为正那这个均值是被少数极端窗口拉起来的实盘中的脆弱性极高因为那些正 IC 窗口里的规律很可能已经消失了。import matplotlib.pyplot as plt oos_ic, expr_list walk_forward_evolve(df) pos_ratio np.mean([1 if ic 0 else 0 for ic in oos_ic]) mean_ic np.mean(oos_ic) plt.hist(oos_ic, bins15, edgecolorblack) plt.xlabel(Out-of-sample IC) plt.ylabel(Window count) plt.title(fPos Ratio: {pos_ratio:.2f}, Mean IC: {mean_ic:.4f}) plt.show()图中直方图若明显右偏且集中在 0 附近属于勉强可用的状态如果分布接近以 0 为对称轴的正态直接放弃这条 GP 管线去换个特征集不要恋战。window_count太少时直方图没有统计意义通常少于 30 个窗口的分位数分析不值得信。提示walk_forward_evolve里的evolve_kwargs透传了pop_size、generations等参数方便你在不同窗口大小下对 GP 进化过程做同样的超参数配置避免窗口和 GP 参数的混淆。5. 遗传编程预测程序的收敛诊断与 3 个高性价比调参方向5.1 适应度曲线的“三段式”判断上涨、平台、乱跳把每代的最优适应度和平均适应度画在一起你会发现健康的收敛通常分三段初始 35 代最优适应度快速上冲这是搜索找到了显著优于随机的结构中间 510 代进入一个缓慢爬坡平均适应度紧咬最优适应度说明种群整体在抬升最后阶段最优适应度围绕某个水平波动不再有明显上升此时再跑更多代只是燃烧 CPU。不健康的表现则是最优适应度曲线反复乱跳每代结果都不同这通常说明以下三种情况树深度过大导致过拟合随机噪声交叉概率过高破坏了已积累的优秀结构或者训练窗口内的数据量不足导致同一棵表达式在不同随机种子下适应度方差极大。处理手段优先级很高的是从种群中重建——把几代前的某个好个体重新注入当前种群往往比继续跑更有用。5.2 用“按时间衰减的适应度”对抗策略失效最近的行情权重更高股票市场具有时变性五年前的因子里可预测成分可能早就被市场博弈抹平了。对日线级别的 GP 预测框架一个性价比极高的改进是把适应度从简单 Spearman 相关改为带时间衰减权重的版本def evaluate_fitness_decay(expr, df, horizon5, half_life60): feature eval_expr(expr, df) future_ret df[close].pct_change(horizon).shift(-horizon) valid pd.concat([feature, future_ret], axis1).dropna() if len(valid) 30: return -1.0 # 按行号生成指数衰减权重最近的数据权重更高 weights 0.5 ** (np.arange(len(valid))[::-1] / half_life) ic, _ spearmanr(valid.iloc[:, 0], valid.iloc[:, 1]) if np.isnan(ic): return -1.0 # 对最近一段时间的预测准确性加权等价于在适应度中引入时间偏好 recent_ic spearmanr( valid.iloc[-int(len(valid) * 0.3):, 0], valid.iloc[-int(len(valid) * 0.3):, 1] )[0] return 0.7 * ic 0.3 * recent_ic if not np.isnan(recent_ic) else ic这里没有直接计算加权 Spearman而是用“全样本 IC 加最近 30% 样本的 IC 线性加权”来近似。这样既不破坏秩相关的稳健性又能让模型更看重近期市场状态。half_life60表示每 60 个交易日权重减半对应约一个季度的市场记忆窗口。这个值在你做日线、持仓周期 5 天左右的场景里是一个合理的起点。5.3 配合遗传算法做因子组合权重的二次优化单棵表达式树的预测能力通常会撞到天花板IC 在 0.040.06 附近就上不去了。此时不要试图把树变深而是要横向集成——让遗传算法去搜索多棵树的权重组合。这才是开头说的“GA 调参”角色的实际落地def ga_weight_search(expr_list, df, horizon5, pop_size30, generations30): 用遗传算法搜索多棵表达式树的最优权重组合 权重被编码为长度为 len(expr_list) 的实数向量 n len(expr_list) # 每个因子先独立计算一次后续 GA 只需做线性组合避免重复 eval_expr factor_matrices [eval_expr(expr, df) for expr in expr_list] def fitness(weights): w np.array(weights) / np.sum(weights) combined sum(w[i] * factor_matrices[i] for i in range(n)) future_ret df[close].pct_change(horizon).shift(-horizon) valid pd.concat([pd.Series(combined), future_ret], axis1).dropna() return spearmanr(valid.iloc[:, 0], valid.iloc[:, 1])[0] # 简化 GA随机实数向量作为个体交叉取均值变异加高斯扰动 population [np.random.dirichlet(np.ones(n)) for _ in range(pop_size)] for _ in range(generations): f [fitness(ind) for ind in population] population sorted(population, keylambda x: f[population.index(x)], reverseTrue) population population[:pop_size // 2] children [] while len(children) pop_size // 2: p1, p2 random.sample(population, 2) child (p1 p2) / 2 np.random.normal(0, 0.05, n) child np.abs(child) / np.sum(np.abs(child)) children.append(child) population population children best_idx int(np.argmax([fitness(ind) for ind in population])) best_w population[best_idx] / np.sum(population[best_idx]) return best_w, factor_matricesspearmanr只能接收一维数组pd.concat之前将合并因子单独转换为了 Series否则dropna会报轴长度不一致的错误。遗传算法作用于权重向量时np.random.dirichlet初始化的好处是权重天生满足总和为 1且每个分量都非负省去了惩罚项编码的复杂度。子代生成采用算术交叉而不是显式交换因为权重向量不是离散结构均值交叉在连续空间里更平滑。注意这里的遗传算法种群只有 30、迭代 30 代因为这个搜索空间极小n 通常不超过 10 维粒子群或贝叶斯优化也可以做同样的事但用同一个 GA 框架可以减少技术栈的割裂感。6. 用历史的“最差 IC 窗口”快速检验 GP 模型的实盘韧性滚动回测跑完之后你最该看的不是平均 IC而是 IC 最差的几个窗口里模型长什么样——如果最差窗口集中在市场急跌或横盘低波动阶段那说明模型本质上是在捕捉趋势延续或反转行情中的某一种单一模式一旦市场状态切换它的排序能力就消失了。oos_ic, expr_list walk_forward_evolve(df) worst_indices np.argsort(oos_ic)[:5] for idx in worst_indices: print(fWindow {idx}: OOS IC {oos_ic[idx]:.4f}) print(expr_to_string(expr_list[idx])) # 将树元组转回可读字符串 print(- * 40)如果最差窗口的表达式中高频出现volume和open的组合一个合理的推断是模型在这些时段里捕捉到的是流动性突变引发的短期波动而这类波动在低成交量横盘期表现极差。你可以针对这些失败窗口做约束重训——在适应度函数里把该窗口内预测方向和真实收益方向一致的个体额外加分惩罚只在“自然发展”行情下有效的表达式。def evaluate_with_penalty(expr, df, penalty_windowNone): base_ic evaluate_fitness_decay(expr, df) if penalty_window is not None: w_feature eval_expr(expr, penalty_window) w_ret penalty_window[close].pct_change(5).shift(-5) w_valid pd.concat([w_feature, w_ret], axis1).dropna() if len(w_valid) 10: w_ic spearmanr(w_valid.iloc[:, 0], w_valid.iloc[:, 1])[0] if np.isnan(w_ic): w_ic 0 return base_ic - 0.2 * abs(w_ic) if w_ic 0 else base_ic return base_ic惩罚力度0.2 * abs(w_ic)是一种软约束——如果模型在最差窗口的 IC 是 -0.05它的适应度会被扣 0.01。这个值看起来不大但在 GP 选择压力下连续几代迭代就会把在失败窗口表现更稳定的个体推到种群前列。最后一个技巧是把种群初始化时的随机数种子固化这样每次复现的结果完全一样——实盘验证阶段可复现性比最佳解更重要。如果两次运行结果完全不同你根本无法判断改进是来自参数调整还是随机运气。本文还有配套的精品资源点击获取