决策树是非常经典的机器学习算法,单独使用时决策树仍然是至今少数可以解释的模型之一,许多决策树一同构建的整体模型,如经典的随机森林模型仍然是目前广泛使用的“传统机器学习算法”,其在化学信息学的出场频率绝不逊色于支持向量机以及神经网络这种较为复杂的模型。决策树的理论基础以及大致的算法其实不难理解,但若我们仔细关注 scikit-learn 的原码时,就会发现那些开发者在一些细节上进行了很多优化,使得模型训练更加快,也支持许多超参供选择。本篇先介绍 scikit-learn 中实现决策树所需要用到的重要的 类
,然后着重介绍 Splitter 的作用以及算法细节。
决策树的简单介绍
在研究 scikit-learn 的实现方法之前,先要大概回忆下决策树到底是怎幺一个算法。我们需要构建 (build) 一颗树 (tree),树中每个节点是一个规则,或叫判别式,如 x1<=1,本质就是对训练数据的一次划分 (split) ,理想中最后叶子节点里全是一类数据,这样对于任何一个新数据,经过几次判别就知道它属于哪个类。但判别式如何产生,为什幺是 x1 而不是其他特征,为什幺是小于 1 而不是其他值?这就涉及到一个判断标准 (criterion),在决策树中采用信息熵作为依据,数据越是纯粹(全是阳性或者全是阴性)则熵越低,反之则越高,通过一次判别后数据划分为两个,新数据的均熵相比原来数据的熵差称为信息增益值,之前博客 中也有介绍到。对于决策树的详细介绍,已经有很多材料,比如该 博客
,这里不做赘述。
Scikit-learn 中几个类的关系
如上图所示,我们调用的是 scikit-learn 中的 DecisionTreeClassifier
,建模时使用 fit
函数,该函数中主要调用了 Criterion
, Splitter
, TreeBuilder
和 Tree
这几个重要的类。从名字就大致可以猜出首先我们需要实例化 Tree
这个类,最终所构建的决策树就是一个这样的对象,构建这棵树则需要一个 builder
,方可认为 TreeBuilder.build
是主要要研究的函数。实例化一个builder 需要提供一个 splitter
,该类就要是解决上述提到的 “为什幺是 x1 而不是其他特征,为什幺是小于 1 而不是其他值”。上面也讲了之所以最终选择该规则,是根据判断标准,于是知道为什幺要实例化 criertion
,就是用来计算信息熵的。至此这几个主要类之间的逻辑关系就搞清楚了。 训练函数
大致可以描述为:
def fit(X, y, criterion='entropy', spliter='best'): tree = Tree() criterion = Entropy() splitter = BestSplitter(criterion) builder = BestFirstTreeBuilder(splitter) builder.build(X, y)
当然,实际默认的 criterion
其实是 Gini 指数,Gini 指数是一个更快速的计算信息熵的方法,所以可以猜出来 scikit-learn 函数大概会做这样的事情以实现所谓的超参选择。有兴趣可以看下 tree/_classes.py
文件,上面的“训练函数”对应的链接。
CRITERIA = {'entropy': Enotropy, 'gini': Gini} criterion = CRITERIA[criterion]
Splitter 的实现方法概述
在介绍 scikit-learn 怎幺做之前,先介绍下为什幺要有 splitter。criterion 仅仅是计算某个判别式是否足够好(信息增益够大),还需要有个算法从许多特征中选择一个“最好”的判别式。比如最暴力的方法就是遍历所有的特征,所有的判别式都去算一遍信息增益。但是考虑到一些特别的需求以及算法的性能,我们还得思考一些细节的问题,比如:
特征是连续的还是分类的,若是分类的,则是二值化的还是多分类的?
- scikit-learn 使用的是 CART 算法,所有必须是二叉树,判别式都为
feature <= threshold
- 的形式。那幺对于二分类的问题,则可以当成是连续的变量,取两者中值即可区分两者。对于多分类的特征就比较复杂一点,理论上数据应该划分为“属于这类、不属于这类”两个子集,必须使用等号进行判别而不是一般使用的小于等于号,scikit-learn 目前还没有解决这个问题(
- )。不过有一种代替的方法是将这些多分类的特征变成 one-hot 形式,如 x={1,2,3} 则可以变成 x = {001, 010, 100} 三个二分类特征。
如果特征是连续的,我该如何分段?
- CART 算法已经给出了如何分段的问题,即先要对样本进行排序,取两者之间的均值作为划分依据。
如果我想只随机选择部分特征,比如在随机森林的“随机”就是指随机挑选一部分特征,怎幺处理?
- 我们可以对特征先进行打乱,然后选择前 n 个特征。
还有,当某个节点下这个数据集的某个特征变为一个常数时,那幺这个节点以及之后的子节点都无需再考虑这个特征,因为我们已无法通过该特征进行分割。如果让笔者用 Python 实现分割算法大致如下。
best = {'feature': None, 't': 0, 'ig': 0} use_features = shuffle(features) i = 0 for feature in use_features: if feature in constant_features: continue x = sorted(X[feature]) if x[0] == x[-1]: constant_features.append(feature) continue i += 1 if i > n: break for i in range(len(x)-1): t = (x[i] + x[i+1]) / 2 ig = criterion.get_information_gain(feature, t) if ig > best['ig']: best = {'feature': feature, 't': t, 'ig': ig}
虽然说人生苦短我用 Python,能用 API 绝不自己写的原则极大地增加了实现算法的方便程度,但其弊端也非常明显,时间和空间复杂度都可能会较高。优化算法的精髓在于创建/维护更少的数组、使用更少的循环次数来实现一个问题。比如上述例子中用 shuffle 的同时创建一个新的数组(当然使用 numpy.random.shuffle
也可以 inplace
),而且 feature in constant_features
又是个时间复杂度为 O(n) 的操作。 scikit-learn 借鉴了 Fisher-Yates 算法,通过不断随机抽取的方法,不仅只需要维护一个数组,还可以通过数组的特殊排序巧妙的避免了扫描数组的操作。下面看下 scikit-learn 如何实现,不过为了方便,笔者将代码尽可能 pythonic 化,而非原始代码。
Scikit-learn 实现 split 算法
上述已经提过,特征需要随机排序以免每次都只分析前 n 个特征,于是 scikit-learn 基于 Fisher-Yates 的乱序算法进行改造。首先大致了解下该乱序算法的原理:
上图介绍了 Fisher-Yates 乱序算法,开始令 f_i = n_feature
, f_i
是个不断向 0 推进的数字,当 f_i = 0
时算法结束。每次从 [0 -> f_i
) 中随机选一个数字 f_j
,然后与 f_i
做交换,然后 f_i
往前进1,这样就相当于每次随机从剩余区间内选个值挪到 selected feature 区间,最后所有值都是随机选出来的,实现了乱序。对于不了解算法的我,觉得这个是个新鲜的操作,即要想将一个数挪到一个区间,则只需要做“交换”、“移动(+1 或者 -1)”两个操作即可。
决策树的 split 算法相对于纯乱序要稍微复杂一点,毕竟要考虑到常数特征的情况。其次要知道一个规则:已知的常数特征虽然不会实际被考虑为节点的判别式,但是在“随机抽取”时是需要考虑的,会占用“考虑前 n 个特征”的名额,因此一个数组被设计为上图下方的样子,分为:
- “drawn known C”区,即已被抽取的已知常数特征。这块区域内的特征下次随机特征时是不会考虑的,有点类似于乱序算法中的“selected feature”区,只是划分到了最前面。
- “other known C”区,即未被抽取的已知常数特征。这块区域是可能被抽取到的。
- “found constant”区,即之前不知道后来发现是常数的区域。这块区域是不该再被抽取到了。
- “remain features“区,即尚未研究过的特征区域。
- “drawn non-C”区,即已经研究过的区域。
不难发现,这段特征数组中有三个地方是不该被抽取的,1、3、5三个区。这给抽随机数带来些麻烦,2、4两个区不连续怎幺办。如果是 python 的话,创建一个新数组,把 2、4 放进去,从中抽一个就行,可显然这样做增加了空间复杂度,于是作者在抽随机数的时候,往前减去 n_found_c
,即图中“place holder”的位置,有点像把第 3 区挪到后面去,这样 f_j
的选择范围就连续了。此时如果 f_j 如果落在 other known C,则直接挪到 drawn known C 区域,表示“选择了一个已知的常数特征”;如果大于 n_konwn_c
,则将这个数字加上 n_found_c
,这样就相当于在 remain features 区域选择了一个特征。如果这个特征是常数,则挪到 “found constants”区,否则计算该特征的信息增益值并挪到 drawn non-C“区域。 用代码表示如下,这段代码是根据 scikit-learn 的源代码进行该写的,为了和上述的变量名字一致,把一些代码简化了。
while (f_i > n_total_c and i < n): i += 1 f_j = rand_int(n_drawn_c, f_i - n_found_c) if f_j < n_known_c: # 通过将该值与 n_drawn_c 交换并 n_drawn_c + 1 # 实现将该数字放入 drawn known C 区域 swap_feature(n_drawn_c, f_j) n_drawn_c += 1 else: f_j += n_found_c x = sorted(X[feature]) if x[-1] <= x[0] + 10e-7: # 首尾相差很小,可认为是常数 # 同理将该数字放入 found constant 区域 swap_feature(f_j, n_total_c) n_total_c += 1 n_found_c += 1 continue
再看另一个地方的细节优化。之前笔者设计的代码如下:
for i in range(len(x)-1): t = (x[i] + x[i+1]) / 2 ig = criterion.get_information_gain(feature, t)
看似非常正常,将排过序的特征值进行一次循环,找到一个阈值传给 criterion
。但仔细想象 criterion 得到阈值后还需要再次对 x 值进行排序,或者一个个比较特征与 t 的关系,将数据集划分为两部分。无形之中时间复杂度就由 o(n) 变成了 o(n2)。 scikit learn 的做法如下:
p = start # 样本的开始位置,不妨认为是 0 while p < end: # 到样本结束位置 # 跳过前几个与第一值相等的样本,因为无法分割 while (p + 1 < end and Xf[p + 1] <= Xf[p] + FEATURE_THRESHOLD): p += 1 p += 1 if p < end: current.pos = p self.criterion.update(current.pos)
细节就是,它传的值是位置,而由于样本已经排好序了(当然这里涉及到比较复杂的问题,就是排序时使用的变量是个与 criterion 共享的指针),直接告诉 criterion 想要分割的位置即可,彼时 criterion 只需要根据位置得到前后两个样本集的 y 值即可计算出信息增益。
小结
刚开始看源码时感觉一头雾水,涉及多个类,各种函数穿插。但随着时间推移慢慢理顺各个变量与类之间逻辑关系后,就变得清晰很多,然后慢慢想办法感受其算法的思想,搞清楚变量的意图(当然还要感谢源码里加了很多注释),最后便能理解其算法的过程是什幺,便开始感叹这些开发者们满脑子想着如何用利用复用数组等算法思想以更少的数组、更少的循环实现某个功能。内心是无比佩服的。
Be First to Comment