本文详细介绍了遗传规划算法的原理以及如何使用它进行因子挖掘。具体到实操流程上，我们使用到了Python中专门实现遗传规划的gplearn包进行代码实现。在对gplearn包进行修改优化后，我们以2016至2022年作为输入数据在商品期货上挖掘了5个具有较好选期能力的截面alpha因子，并回测了其在全样本上的表现。回测结果表明，遗传规划挖掘出的因子在样本内外均具备一定的有效性。
一、因子挖掘方法论
在我们之前的期货多因子系列报告中，各类因子往往是基于一定的经济学逻辑或是历史经验人工构造出来的。如经典的截面动量因子，就是我们认为资产之间的相对强弱关系会延续（动量效应），那么依据“强者恒强，弱智恒弱”的逻辑我们就可以构造出动量因子。这样的构建因子的方式属于“先有逻辑，后有公式”的“演绎法”。在因子构建过程中，还有一种方式，即我们可以先通过机器学习的方法将大量的历史数据与收益率序列相拟合，生成海量因子，再从中归纳总结出背后的逻辑。这种方式就属于“先有公式，后有逻辑”的“归纳法”。在“归纳法”法中，遗传规划作为一种仿生启发式算法，在我们对目标因子知之甚少的情况下可以帮助我们挖掘因子。它的优势在于通过启发式的搜索，可以挖掘出人工难以构造、复杂的因子。本节我们将重点介绍遗传规划算法以及如何使用Python的gplearn包进行因子挖掘。
（一）遗传算法与遗传规划简介
遗传算法（Genetic Algorithm）最初由美国密歇根大学的J.Holland提出，是一种通过模拟自然界生物进化过程（“物竞天择、适者生存”）搜索最优解的算法，其本质上也是一种监督学习算法。对于一个最优化的问题，它借鉴了生物学中的现象（遗传、突变、杂交等）使一定数量的初始解按照适应度的方向进化为更优的解。进化过程从完全随机生成的个体种群开始，一代代进化。每一代中，会基于个体的适应度筛选出较优个体并在个体中发生变异、进化进而生成新的种群，新种群则继续进行迭代直至生成最优（适应度最高）种群。遗传规划/遗传编程（Genetic Progamming）是遗传算法中的一个分支，相较于遗传算法，它最大的不同在于将进化对象编码成树而不是简单的线性串。得益于此，遗传规划可以表达更为复杂的语义。下图简单展示了在遗传规划算法中可以作为进化个体的公式树（）。
我们沿着公式树从下往上即可还原整个公式，其中红色的节点表征函数，其他的粉红色叶子则表征变量和常数。在因子挖掘过程中，公式树也可以视作我们的因子表达式，所以公式树的进化过程实际上就是挖掘因子。
（二）适应度函数
适应度函数用来计算每一代种群中个体的适应度，适应度则衡量了个体与最终的目标个体的相符程度。每一代进化过程中，只保留适应度较高的个体作为后续进化的种群，适应度较低的个体则被淘汰。从因子挖掘的角度来说，我们可以使用用来评判因子有效程度的指标做为个体适应度，如使用因子的RankIC、RankICIR可以挖掘与未来收益率相关程度较高的线性因子；互信息（Mutual information）可以用来挖掘与识别非线性的因子；或是直接使用因子回测后的夏普率（sharpe ratio）作为适应度。
1.杂交变异
杂交变异是遗传规化中的主要变异方式，是将两个个体中的遗传物质进行混合交换。在公式树进化过程中，杂交变异需要选择两个种群，其中一个种群中的最高适应度个体作为“父代”，从“父代”中随机挑选一个子树以待替换；并从另一个种群中随机挑选一个个体作为“捐赠者”，也从“捐赠者”中随机挑选一个子树用以替换“父代”中的子树。最终生成的公式树用以后续进化。2.子树变异
子树变异是一种更为激进的变异方式，称它更激进是因为子树变异如同对公式树进行“移植手术”般，直接将公式树中的子树替换成最初生成的随机子树。这样做可以将已经淘汰的子树重新引入进化过程中以增强遗传物质的多样性。3.Hoist变异
Hoist（抬升）变异也是一种较为激进的变异方式。主要是为了防止公式树过于“臃肿庞大”而对公式树中的随机子树进行“截断”，并保留子树中的一小部分用以替换原子树。4.点变异
点变异是遗传规划中最常见的变异方式，与子树变异类似都是为了增添遗传物质的多样性而将过去淘汰的节点重新引入进化过程中。不同于子树变异激进的变异方式，点变异较为温和，仅仅只替换公式数中的随机某个节点而不是整个子树。以上即为公式树的主要变异方式，它们四者的概率之和应小于1。从因子的可解释性角度来说，我们应该尽量减少因子挖掘过程中子树变异、Hoist变异与点变异的概率，而提升杂交变异概率。这么做的目的也是为了防止因子表达式过于难以解释。
1.参数说明与设置
本文主要通过Python中的专门实现遗传规划的gplearn包进行商品期货的截面因子挖掘，有关gplearn的详细介绍，感兴趣的读者可以查阅官方文档（https://gplearn.readthedocs.io/en/stable/）。
gplearn中的特征转化器SymbolicTransformer，可以将初代种群的原始特征通过公式树转换为适应度最优的新特征，因此可以用于因子的挖掘。其中SymbolicTransformer模型的主要参数说明与设置如下：
2.gplearn的改进
直接使用gplearn挖掘截面选期因子，是无法实现或者效果往往是不尽人意的，其原因有三点：lgplearn原始模型输入变量只能是2维数据，即NumPy中的2维数组。对应到因子挖掘过程中，输入变量只能是包含时间序列以及各项特征的单一资产，无法扩充到截面上多个资产。lgplearn所提供的原函数集（加、减、乘、除、开方、取对数、取绝对值等）较为简单且无法考虑时间序列上的运算（滚动窗口运算）。lgplearn自带适应度函数（IC、RankIC）较为简单，且出来的挖掘因子往往样本内的效果较好，样本外效果较差，即过拟合程度较高。针对第一点，我们将输入的训练变量更改为3维数组（3D array），增添截面这一维度。修改过后的3维数组的shape为（x,y,z），其中x表征OHLC、成交量、持仓量、仓单数等变量，y表征时间序列，z表征截面上的品种个数。我们通过修改原包中的大量函数，用以适配3维数组的运算。针对第二点，我们对gplearn源代码进行修改后引入了一系列可以进行时间序列运算的函数以及Ta-lib包中的函数并将之汇总扩充到原函数集中，用以提高因子挖掘的能力，具体算子函数表请见附录。
最后我们针对第三点，自定义了适应度函数，用以提升因子挖掘效果。在初步尝试中，我们直接使用单因子回测后策略的夏普比率作为个体额适应度指标，但多次测试下来，因子也存在样本外失效较快的情况。对于此，我们进行简单的交叉验证，即将样本内的训练集进一步划分为80%训练集+20%验证集，并分别计算训练集与验证集的夏普率。得出两者值后，再根据以下条件筛选因子：训练集夏普率绝对值大于1（取绝对值为了考虑因子方向）；以上便是我们对于gplearn的改造与优化，下一节我们将展示具体的因子挖掘流程与回测结果。
（一）截面选期因子挖掘流程
本节我们将遵循以下步骤进行因子挖掘：
1.样本选择与回测细节
商品期货池：我们从国内期货市场筛选出历史流动性较好的41个品种具体如下表所示。交易价格：我们使用各个商品期货品种中主力合约的复权收盘价计算收益。
回测区间：我们使用2016/1/1-2022/1/1作为样本内训练集，2022/1/2-2023/3/10作为样本外测试集，总回测区间为2016/1/1-2023/3/10。输入变量（初代种群）：如上文图表7所提，各个品种开盘价、收盘价、最高价、最低价、成交量、成交额、成交均价、持仓量、当日收盘价收益率、仓单数、仓单变化数11个原始因子。目标变量：各个品种2个交易日后收益率，我们回测时在T日收盘计算因子、T+1日做入、T+2日产生收益。单因子回测方法：五组分层回测法，默认使用因子排序前20%品种与排序后20%品种构建多空组合并计算收益（等权）。调仓周期：每日调仓。
2.使用gplearn进行因子挖掘
将我们的训练集与测试集输入模型，使用SymbolicTransformer.fit()方法即可开始训练模型。在训练过程中，首先通过我们图表7中的算子函数集随机生成大量公式（因子表达式）。通过自定义适应度函数计算每个个体的适应度值，并按照图表1中的遗传规划流程循环进化。在上一节自定义适应度函数中已经提到，我们会将样本内数据进一步拆分并进行简单交叉验证以对抗过拟合。对挖掘出来的因子进行单因子分层回测，考察因子在样本内与样本外的表现。尝试对回测效果较好的因子内在含义进行归因解释。
（二）截面选期因子挖掘结果及整体回测结果
经过多次挖掘，我们总结出了以下5个在全样本上都表现较好的截面选期Alpha因子。下面是Alpha因子在样本内、样本外、以及全样本上的单因子回测结果。
（三）截面选期因子回测表现及简单归因
Alpha1的因子表达式为ts_midpoint(ts_pct_change(ts_inverse_cv(ts_ AROONOSC(volume,avgprice,10),10),5),486)，该因子属于量价相关性类因子。它描述了短期量价相关性趋势的波动的变化率的长期走势。该因子为正向因子，则表明量价相关性趋势的波动的变化程度越高越好。下面是该因子在全样本上的具体回测结果，其中样本线以左表明样本内，样本线以右表明样本外。
回测下来，该因子分层的单调性虽不明显，但是在时序上整体的RankIC累计值与净值的走势都比较流畅，稳定性强，回测基本上小于4%。2.Alpha2因子
Alpha2的因子表达式为ts_sum(ts_maxmin(ts_maxmin(warehouse,126),126),63)，该因子属于基本面仓单类因子的衍生物。回测结果表明，该因子分层的单调性强，RankIC累计值与净值走势也呈现出强稳定性。样本外选期能力虽有所衰减，但仍然具备一定的有效性。3.Alpha3因子
Alpha3的因子表达式为ts_ema(ts_sum(ts_rsi(ts_kama(openinterest, 486), 63), 243), 63)，它描述了持仓量的库夫曼移动均线在过去一段时间内的相对强弱程度之和的平均水平。该因子为负向因子，表明品种长期持仓量的均线持续走高时，倾向于做空该品种。Alpha3因子在整个时序上分层的单调性较为一般，从RankIC累计值走势来看，样本内该因子预测方向也发生过较长时间的反转。但其在样本外的表现较为稳定，回测幅度较小。
4.Alpha4因子
Alpha4的因子表达式为ts_ema(ts_inverse_cv（ts_corr(volume,avgprice,21),42）,105)。它反映了回看期21日（一个月）的成交均价与成交量的相关性变化程度的指数移动平均。根据因子回测结果，该因子为负向因子，表明品种过去一个月的量价相关性在过去一段时间的变化程度越低越好。从回测结果来看，Alpha1因子不仅在分层回测中表现出来较好的单调性（分层明显）且整个时序上保持了较好的负向选期能力（RankIC累计值稳定下行）。在样本外，Alpha1因子也保持了较强的alpha能力且整体的回撤水平较低。5.Alpha5因子
Alpha5的因子表达式为ts_dema(ts_median(ts_cov(open, volume,21),21),15)，它同样可以归属为量价相关性类因子。它描述了开盘价与成交量的相关性中枢的移动平均水平。回测结果表明，该因子可以很好的区分截面上品种的强弱关系，且RankIC累计值在全样本上能维持稳定的下行走势。Alpha5虽在样本外发生了最大回撤，但回撤时间不长且近期仍维持了较好的走势。三、总结与思考
本文详细介绍了遗传规划算法的原理以及如何使用它进行因子挖掘。具体到实操流程上，我们使用到了Python中专门实现遗传规划的gplearn包进行代码实现。在因子挖掘的过程中，我们对gplearn包中的源代码进行了改进与优化：扩充了原始用于进化的算子函数集，使因子在进化过程中可以进行时间序列上运算。自定义了适应度函数，以对抗过拟合程度并提高因子的可解释性。在完成了以上修改后，我们利用了遗传规划算法，以2016至2022年作为输入数据在商品期货上挖掘了5个具有较好选期能力的截面alpha因子，并回测了5个因子在全样本上的表现：Alpha1：年化收益:9.32%，年化波动:5.42%，夏普:1.72，最大回撤:4.15%，卡玛比率:2.25；Alpha2：年化收益:11.07%，年化波动:7.00%，夏普:1.58，最大回撤:6.05%，卡玛比率:1.83；Alpha3：年化收益:10.24%，年化波动:7.74%，夏普:1.32，最大回撤:9.86%，卡玛比率: 1.05；Alpha4：年化收益:9.35 %，年化波动:7.41%，夏普:1.26，最大回撤:8.16%，卡玛比率: 1.15；Alpha5：年化收益:8.40%，年化波动:7.11%，夏普:1.18，最大回撤:9.00%，卡玛比率:0.93；回测下来，由算法挖掘出来的因子在样本内外均具备一定的有效性，表明利用遗传规划可以帮助我们归纳并总结出具有一定alpha能力的因子。但这种方法仍然存在着许多缺点。第一，因子的挖掘过程中具有很强的随机性，不同的输入特征变量、模型的超参数、算子函数、适应度函数都会影响最终的因子挖掘结果。第二，很容易产生过拟合的问题，即挖掘出来的因子往往在样本内表现较好而样本外失效较快，即使本文使用了简单的交叉验证以对抗过拟合，但也无法保证样本外因子表现持续较好。第三，遗传规划挖掘的因子复杂程度普遍较高，绝大部分并不具有明确的经济意义。因此在后续的研究中，可以从两方面进行改进以提升我们因子挖掘的效果，首先对输入特征变量进行特征工程使其转变成有一定的经济含义的因子，再进行挖掘，缩小挖掘的随机性并提高因子的可解释性。其次是可以在挖掘过程中加入时序交叉验证的步骤以更好地对抗因子过拟合的问题。
四、附录
我们在遗传规划挖掘因子的过程中使用到了以下算子函数，其中X1，X2为表征不同特征变量（open、close、…）的2维数组，每行对应每一个交易日，每列对应每一个品种，d为时序中的回看天数。