Photo by Ryan Jacobson on Unsplash
SVM（核技巧）
决策树
随机森林
非线性支持向量机SVM（核技巧）¶
当数据线性可分时，我们知道如何使用SVM来进行分类(了解线性SVM),而当线性不可分时，SVM通过某种非线性映射（核函数kernel）将数据映射到高维空间，在这个高维空间中数据是线性可分的，在高维空间中构造最优分类超平面。
可以通过这个视频先有个直观感受：SVM 分割超平面
这里我们直接引入线性可分时的超平面方程： \begin{equation}\label{eq:1} f(x)=\sum_{i=1}^{n} \alpha_{i} y_{i} x_{i}^{T} x+b \end{equation}
其中，$y_i \in \{+1,-1\}$表示样本 $x_i$的标签,$\alpha_i$表示各项对应的拉格朗日乘子，$Support Vector$对应的$\alpha_i \ne 0$,非$Support Vector$对应的 $\alpha_i = 0$,法向量$w=\sum_{i=1}^{n} a_{i} y_{i} x_{i}$。
将公式$\ref{eq:1}$表达为内积形式：
\begin{equation}\label{eq:2} f(\boldsymbol{x})=\sum_{i=1}^{n} \alpha_{i} y_{i}\left\langle\boldsymbol{x}_{i}, \boldsymbol{x}\right\rangle+ b \end{equation}
当线性不可分时，我们希望存在一个映射$\phi$能将数据映射到另外的特征空间，在新的特征空间中使用线性分类器：
\begin{equation}\label{eq:3} f(\boldsymbol{x})=\sum_{i=1}^{n} \alpha_{i} y_{i}\left\langle\emptyset\left(\boldsymbol{x}_{i}\right), \emptyset(\boldsymbol{x})\right\rangle+ b \end{equation}
那么，核函数是指这样的一个函数$K$，对于所有的$x,z \in \mathrm{X}$，满足$K(\mathrm{x}, \mathrm{z})=<\phi(\mathrm{x}), \phi(\mathrm{z})>$,使用这样的核函数有什么好处呢？简而言之，如果没有核函数，我们就要分两步，第一，将$x,z$显示的映射到新的特征空间后得到$\phi(\mathrm{x}), \phi(\mathrm{z})$；第二，计算$\phi(\mathrm{x}), \phi(\mathrm{z})$的内积。有了核函数，没有显示映射这一步，直接求最后的内积结果。
这种将内积替换成核函数的方式就叫做核技巧$(Kernel Trick)$。
举个例子，现在二维平面上存在若干点，点集$A$服从$\left\{x, y | x^{\wedge} 2+y^{\wedge} 2=1\right\}$，点集$B$服从$\left\{x, y | x^{\wedge} 2+y^{\wedge} 2=9\right\}$。
很显然，红蓝点集是线性不可分的，现在做一下映射变化$(x, y) \rightarrow\left(x, y, x^{\wedge} 2+y^{\wedge} 2\right)$，变换后在三维空间中点的分布为：
现在再看两个点集就是线性可分的，用一个平面去分割。
再看一个更一般的例子：
现在红蓝点集服从两个椭圆分布，一个理想的分界线应该是一个椭圆，假如二维平面的两个坐标轴是$X_{1}$和$X_{2}$，那么这个分解椭圆的一般方程是：
\begin{equation}\label{eq:4} a_{1} X_{1}+a_{2} X_{1}^{2}+a_{3} X_{2}+a_{4} X_{2}^{2}+a_{5} X_{1} X_{2}+a_{6}=0 \end{equation}
现在，我们换一种形式，令：
\begin{equation}\label{eq:5} \mathrm{Z}_{1}=\mathrm{X}_{1}, \mathrm{Z}_{2}=X_{1}^{2}, \mathrm{Z}_{3}=X_{2}, \mathrm{Z}_{4}=X_{2}^{2}, \mathrm{Z}_{5}=X_{1} X_{2} \end{equation}
那么公式$\ref{eq:4}$就变成了如下形式：
\begin{equation}\label{eq:6} \sum_{i=1}^{5} a_{i} Z_{i}+a_{6}=0 \end{equation}
由公式$\ref{eq:4}$到公式$\ref{eq:6}$,我们相当于做了一个由二维到五维空间的特征映射，二维空间的线性不可分就变成了五维空间的线性可分。
现在，我们再假设一个这样的映射函数：
\begin{equation}\label{eq:7} \emptyset\left(\left(x_{1}, x_{2}\right)\right)=\left(\sqrt{2} x_{1}, x_{1}^{2}, \sqrt{2} x_{2}, x_{2}^{2}, \sqrt{2} x_{1} x_{2}, 1\right) \end{equation}
这是一个二维映射到五维的例子，对于两个向量$a_1=(x_1, x_2)$,$a_2=(y_1, y_2)$我们要计算映射到五维空间后两个向量的内积：
\begin{equation}\label{eq:8} \left\langle\phi\left(\left(x_{1}, x_{2}\right)\right), \phi\left(\left(y_{1}, y_{2}\right)\right)\right\rangle=2 x_{1} y_{1}+x_{1}^{2} y_{1}^{2}+2 x_{2} y_{2}+x_{2}^{2} y_{2}^{2}+2 x_{1} x_{2} y_{1} y_{2}+1 \end{equation}
另外，我们不进行映射计算，直接计算下面的公式：
\begin{equation}\label{eq:9} \left(\left\langle a_{1}, a_{2}\right\rangle+ 1\right)^{2}=2 x_{1} y_{1}+x_{1}^{2} y_{1}^{2}+2 x_{2} y_{2}+x_{2}^{2} y_{2}^{2}+2 x_{1} x_{2} y_{1} y_{2}+1 \end{equation}
你会发现，两个公式的计算结果是相同的，区别是什么呢？
重点来了：
一个是根据映射函数，映射到高维空间中，然后再根据内积的公式进行计算，计算量大；
另一个则直接在原来的低维空间中进行计算，而不需要显式地写出映射后的结果，计算量小。
其实，聪明的你早就看出来了，这个$\left(\left\langle a_{1}, a_{2}\right\rangle+ 1\right)^{2}$不就是一个核函数嘛，没错，我们通过这个核函数$k\left(a_{1}, a_{2}\right)=\left(\left\langle a_{1}, a_{2}\right\rangle+ 1\right)^{2}$的低维运算得到了先映射再内积的高维运算的结果，这就是核函数的神奇之处，它有效减少了我们的计算量。在这个例子中，我们对一个2维空间做映射，映射后的新空间是原来空间所有一阶和二阶的组合，得到了5维的新空间；如果原空间是3维的，那么我们会得到19维的新空间，这个数目是呈爆炸性增长的。如果我们使用$\phi(·)$做映射计算，难度非常大，而且如果遇到无穷维的情况，就根本无从计算了。所以使用核函数进行计算是非常有必要的。
那么问题来了，核函数从哪来？事实上，构造核函数一般来讲是一件很困难的事，通常人们会从一些常用的核函数中进行选择，根据问题和数据的不同，选择不同的参数，得到不同的核函数。径向基核函数是SVM中常用的一个核函数，径向基核函数的高斯版本的公式如下：
\begin{equation}\label{eq:10} \kappa\left(x_{1}, x_{2}\right)=\exp \left\{-\frac{\left\|x_{1}-x_{2}\right\|^{2}}{2 \sigma^{2}}\right\} \end{equation}
参数$\sigma$与$LASSO$算法中的惩罚因子$\alpha$类似，如果$\sigma$选得很大的话，高次特征上的权重实际上衰减得非常快，如果$\sigma$选得很小，则可以将任意的数据映射为线性可分，但随之而来的可能是非常严重的过拟合问题。总的来说，通过调控参数$\sigma$，高斯核实际上具有相当高的灵活性，也是使用最广泛的核函数之一。
决策树(Decision Tree)¶
决策树是一种基本的分类与回归算法，是一种监督学习(Supervised Learning)算法，顾名思义，决策树呈树形结构，在分类问题中，表示基于特征对实例进行分类的过程，树内部节点表示一个特征或者属性，叶节点表示一个类。
图片来源：https://cuijiahua.com/blog/2017/11/ml_2_decision_tree_1.html
可以将决策树看成是一个$if-then$规则的集合，将决策树转换成$if-then$规则的过程是这样的：由决策树的根节点到叶节点的每一条路径构建一条规则；路径上内部节点特征对应着条件，而叶节点的类对应着结论。决策树的路径或其对应的$if-then$规则集合有一个重要的性质：互斥且完备，也就是说每一个实例都被一条路径或者一条规则覆盖，而且只被一条路径或者一条规则覆盖。
决策树学习通常包括三个步骤：
特征选择
决策树生成
决策树修剪
特征选择¶
特征选择在于选取对训练数据具有分类能力的特征，一般遵从两个准侧：信息增益和信息增益比。
在内容展开之前，需要先了解两个概念——熵和条件熵。在热力学中，熵的物理意义是体系混乱程度或者不稳定性的度量，在信息论与概率统计中，熵是表示随机变量不确定性的度量。
假设$X$是一个取有限个值的离散随机变量，其概率分布为：
\begin{equation}\label{eq:11} P\left(X=x_{i}\right)=p_{i}, i=1,2, \cdots, n \end{equation}
那么，随机变量$X$的熵定义为：
\begin{equation}\label{eq:12} H(X)=-\sum_{i=1}^{n} p_{i} \log \left(p_{i}\right) \end{equation}
熵越大，不确定性越大。
条件熵是表示在已知随机变量$X$的条件下随机变量$Y$的不确定性。
\begin{equation}\label{eq:13} H(Y | X)=\sum_{i=1}^{n} p_{i} H\left(Y | X=x_{i}\right)\\ p_{i}=P\left(X=x_{i}\right) \end{equation}
当熵和条件熵中的概率是由数据估计（特别是极大似然估计）得到时，所对应的熵和条件熵分别称为经验熵和经验条件熵。
信息增益¶
信息增益：特征$A$对训练数据集$D$的信息增益$g(D, A)$定义为集合$D$的经验熵$H(D)$与特征$A$给定条件下$D$的经验条件熵$H(D | A)$之差，即
\begin{equation}\label{eq:14} g(D, A)=H(D)-H(D | A) \end{equation}
信息增益表示得知特征$X$的信息而使得类$Y$的信息的不确定性减少的程度。一般地，熵$H(Y)$与条件熵$H(Y | X)$之差称为互信息(mutual information)，决策树学习中的信息增益等价于训练数据集中类与特征的互信息。
信息增益大的特征具有更强的分类能力，所以根据信息增益准则，对训练数据集$D$，计算其每个特征的信息增益，并比它们的大小，从而选择信息增益最大的特征。
假设训练数据集为$D$，样本容量为$|D|$,有$k$个类别$C_{k}$,$\left|C_{k}\right|$为类别$C_{k}$的样本个数。某一特征$A$有$n$个不同的取值$a_{1}, a_{2}, \cdots, a_{n}$。根据特征$A$的取值可将数据集$D$划分为$n$个子集$D_{1}, D_{2}, \cdots, D_{n}$,$\left|D_{i}\right|$为$D_{i}$的样本个数。并记子集$D_{i}$中属于类$C_{k}$的样本的集合为$D_{i k}$，$\left|D_{i k}\right|$为$D_{i k}$的样本个数。
信息增益算法
- 输入：训练数据集$D$和特征$A$；
- 输出：特征$A$对训练数据集$D$的信息增益$g(D, A)$.
(1) 计算数据集$D$的经验熵$H(D)$
- \begin{equation}\label{eq:15} H(D)=-\sum_{k=1}^{K} \frac{\left|C_{k}\right|}{|D|} \log \frac{\left|C_{k}\right|}{|D|} \end{equation}
(2) 计算特征$A$数据集$D$的经验条件熵$H(D | A)$
- \begin{equation}\label{eq:16} H(D | A)=\sum_{i=1}^{n} \frac{\left|D_{i}\right|}{|D|} H\left(D_{i}\right)=\sum_{i=1}^{n} \frac{\left|D_{i}\right|}{|D|} \sum_{k=1}^{K} \frac{\left|D_{i k}\right|}{D_{i}} \log \frac{\left|D_{i k}\right|}{D_{i}} \end{equation}
- (3) 计算信息增益
- \begin{equation}\label{eq:17} g(D, A)=H(D)-H(D | A) \end{equation}
信息增益比¶
以信息增益作为特征选择准则，会存在偏向于选择取值较多的特征的问题。可以采用信息增益比对这一问题进行校正。特征$A$对训练数据集$D$的信息增益比定义为其信息增益与训练集$D$关于特征$A$的值的熵之比，即
\begin{equation}\label{eq:18} g_{R}(D | A)=\frac{g(D, A)}{H_{A}(D)}\\ H_{A}(D)=-\sum_{i=1}^{n} \frac{\left|D_{i}\right|}{|D|} \log \frac{\left|D_{i}\right|}{|D|} \end{equation}
决策树的生成¶
ID3算法¶
- 从根节点开始，对结点计算所有可能特征的信息增益，选择信息增益最大的特征作为结点的特征，并由该特征的不同取值构建子节点；
- 对子节点递归地调用以上方法，构建决策树；
- 直到所有特征的信息增益均很小或者没有特征可选时为止。
C4.5算法¶
C4.5算法对ID3算法进行了改进，它在生成决策树的过程中，使用信息增益比来进行特征选择。
CART算法¶
分类与回归树（classification and regression tree，CART）是应用广泛的决策树学习方法，既可用于分类也可用于回归。CART假设决策树是一个二叉树，递归地二分每个特征，将特征空间划分为有限个单元，并在这些单元上确定预测的概率分布。
递归地构建二叉决策树过程中，对于回归树，采用的是平方误差最小化准则；对于分类树，采用基尼指数最小化准则。
关于CART算法的具体内容此处不再展开，感兴趣的朋友可以查阅参考文献。
决策树修剪¶
决策树的剪枝往往通过极小化决策树整体的损失函数(loss function)或代价函数(cost function)来实现。
设树$T$的叶结点个数为$|T|$，$t$是树$T$的叶结点，该叶结点有$N_t$个样本点，其中$k$类的样本点有$N_{ik}$个，$H_t(T)$为叶结点$t$上的经验熵，$\alpha \ge 0$为参数，则损失函数可定义为： \begin{equation}\label{eq:19} C_{\alpha}(T)=\sum_{t=1}^{|T|} N_{t} H_{t}(T)+\alpha|T| \end{equation}
其中经验熵为 \begin{equation}\label{eq:20} H_{t}(T)=-\sum_{k} \frac{N_{t k}}{N_{t}} \log \frac{N_{t k}}{N_{t}} \end{equation}
损失函数的第一部分表示训练数据的预测误差，第二部分表示模型复杂度。$\alpha$的大小反映了对模型训练集拟合度和模型复杂度的折衷考虑。剪枝的过程就是当$\alpha$确定时，选择损失函数最小的模型。
剪枝过程为：
- 计算每个节点的经验熵；
- 递归地从树的叶节点向上回缩，如果将某一个父节点的所有叶节点合并，能够使得其损失函数减小，则进行剪枝，将父节点变成新的叶节点；
- 返回2，直到不能继续合并。
随机森林（Random Forest）¶
我们知道，如果我们在生成决策树时不去限制其最大深度，就会发生过拟合，这是因为此时决策树有无限的灵活性，它可以一直生长，直到每一个样本都在一个叶子节点上，完美地将它们分类。反过来，如果我们限制其深度，如果设置不当，也会带来巨大偏差。
所以，为了避免限制树的深度，减少方差(好)，增加偏差(坏)，我们可以将许多决策树合并到一个称为随机森林的集成模型中。
随机森林是很多决策树的组合，这个模型不只是简单地平均树的预测值(我们可以称之为“森林”)，它之所以叫随机森林，是因为它使用了两个关键概念：
构建树的时候对训练数据集随机采样。
分割结点的时候随机选择特征集合的子集。
训练数据随机采样¶
训练时，随机森林中的每一棵树从数据集中随机采样，这种采样是有放回的，也就意味着一些样本可能在一棵树中被使用多次。这个想法是通过在不同的样本上训练每棵树，虽然每棵树可能对特定的一组训练数据有很高的方差(拟合的很好)，但总体来说，整个森林的方差会更低，但不会以增加偏差为代价。
测试时，通过平均每棵树的预测来做出最终预测。这种在不同训练子集上训练每一个学习器，然后平均预测输出的过程就是bagging,是bootstrap aggregating的简写。
特征随机采样¶
随机森林的另一个主要概念是，在每个决策树中分解每个节点时，只考虑所有特征的子集。通常将此值设置为$sqrt(n \_features)$用于分类，也就是说如果每个树中的每个节点上有16个特性，则只考虑4个随机特性来拆分节点。
所以，随机森林就是将数百或数千棵决策树组合在一起，使用稍有不同的训练子集去训练每棵决策树，分裂结点时考虑有限数量的特征，最终预测通过平均每棵树的预测得出。
Random Forest in Practice¶
接下来，我们使用$Scikit-Learn$构建一个随机森林模型来解决一个实际问题。
Dataset¶
我们要解决的问题是对一个人是否健康进行分类，特征包括个体的社会经济和生活方式等方面，健康状况不佳的标签为0，健康状况良好的标签为1。数据集点此下载。
import pandas as pd
RSEED = 50
df = pd.read_csv("./2015.csv/2015.csv").sample(100000, random_state=RSEED)
df.head()
df = df.select_dtypes("number")
df["_RFHLTH"] = df["_RFHLTH"].replace({2: 0})
df = df.loc[df["_RFHLTH"].isin([0, 1])].copy()
df = df.rename(columns={"_RFHLTH": "label"})
df["label"].value_counts()
标签不平衡意味着准确性不是最好的度量标准。
我们删除一些建模用不到的属性。
# Remove columns with missing values
df = df.drop(
columns=[
"POORHLTH",
"PHYSHLTH",
"GENHLTH",
"PAINACT2",
"QLMENTL2",
"QLSTRES2",
"QLHLTH2",
"HLTHPLN1",
"MENTHLTH",
]
)
划分数据集合，划分比例{训练集：测试集 = 7 : 3 }
from sklearn.model_selection import train_test_split
import numpy as np
# Extract the labels
labels = np.array(df.pop("label"))
# 30% exapmles in test data
train, test, train_labels, test_labels = train_test_split(
df, labels, stratify=labels, test_size=0.3, random_state=RSEED
)
train = train.fillna(train.mean())
test = test.fillna(test.mean())
# Features for feature importances
features = list(train.columns)
train.shape
test.shape
创建模型并训练。
from sklearn.ensemble import RandomForestClassifier
# Create the model with 100 trees
model = RandomForestClassifier(
n_estimators=100, random_state=RSEED, max_features="sqrt", n_jobs=-1, verbose=1
)
# Fit on training data
model.fit(train, train_labels)
import numpy as np
n_nodes = []
max_depths = []
for ind_tree in model.estimators_:
n_nodes.append(ind_tree.tree_.node_count)
max_depths.append(ind_tree.tree_.max_depth)
print("Average number of nodes %d" % (int(np.mean(n_nodes))))
print("Average maximum depth %d" % (int(np.mean(max_depths))))
我们定义两个函数用来计算一些度量标准，绘制ROC曲线以及混淆矩阵。
import matplotlib.pyplot as plt
from sklearn.metrics import precision_score, recall_score, roc_auc_score, roc_curve
def evaluate_mdoel(predictions, probs, train_predictions, train_probs):
"""Compare machine learning model to baseline performance.
Computes statistics and shows ROC curve."""
baseline = {}
baseline["recall"] = recall_score(test_labels, [1 for _ in range(len(test_labels))])
baseline["precision"] = precision_score(
test_labels, [1 for _ in range(len(test_labels))]
)
baseline["roc"] = 0.5
results = {}
results["recall"] = recall_score(test_labels, predictions)
results["precision"] = precision_score(test_labels, predictions)
results["roc"] = roc_auc_score(test_labels, probs)
train_results = {}
train_results["recall"] = recall_score(train_labels, train_predictions)
train_results["precision"] = precision_score(train_labels, train_predictions)
train_results["roc"] = roc_auc_score(train_labels, train_probs)
for metric in ["recall", "precision", "roc"]:
print(
f"{metric.capitalize()} Baseline:{round(baseline[metric],2)} Test:{round(results[metric],2)} Train:{round(train_results[metric],2)}"
)
# calculate false positive rate and true positive rate
base_fpr, base_tpr, _ = roc_curve(test_labels, [1 for _ in range(len(test_labels))])
model_fpr, model_tpr, _ = roc_curve(test_labels, probs)
plt.figure(figsize=(8, 6))
plt.rcParams["font.size"] = 16
# plot both curves
plt.plot(base_fpr, base_tpr, "b", label="baseline")
plt.plot(model_fpr, model_tpr, "r", label="model")
plt.legend()
plt.xlabel("False Positive Rate")
plt.ylabel("True Positive Rate")
plt.title("ROC Curves")
混淆矩阵
from sklearn.metrics import confusion_matrix
import itertools
def plot_confusion_matrix(
cm, classes, normalize=False, title="Confusion matrix", cmap=plt.cm.Oranges
):
"""
This function prints and plots the confusion matrix.
Normalization can be applied by setting `normalize=True`.
Source: http://scikit-learn.org/stable/auto_examples/model_selection/plot_confusion_matrix.html
"""
if normalize:
cm = cm.astype("float") / cm.sum(axis=1)[:, np.newaxis]
print("Normalized confusion matrix")
else:
print("Confusion matrix, without normalization")
print(cm)
plt.figure(figsize=(10, 10))
plt.imshow(cm, interpolation="nearest", cmap=cmap)
plt.title(title, size=24)
plt.colorbar(aspect=4)
tick_marks = np.arange(len(classes))
plt.xticks(tick_marks, classes, rotation=45, size=14)
plt.yticks(tick_marks, classes, size=14)
fmt = ".2f" if normalize else "d"
thresh = cm.max() / 2.0
# Labeling the plot
for i, j in itertools.product(range(cm.shape[0]), range(cm.shape[1])):
plt.text(
j,
i,
format(cm[i, j], fmt),
fontsize=20,
horizontalalignment="center",
color="white" if cm[i, j] > thresh else "black",
)
plt.grid(None)
plt.tight_layout()
plt.ylabel("True label", size=18)
plt.xlabel("Predicted label", size=18)
随机森林的预测结果
train_rf_predictions = model.predict(train)
train_rf_probs = model.predict_proba(train)[:, 1]
rf_predictions = model.predict(test)
rf_probs = model.predict_proba(test)[:, 1]
evaluate_mdoel(rf_predictions, rf_probs, train_rf_predictions, train_rf_probs)
cm = confusion_matrix(test_labels,rf_predictions)
plot_confusion_matrix(cm,classes=['Pool Health','Good Health'],title='Health Confusion Matrix')
与单决策树相比，该模型假阳性成分(FP)较少，假阴性成分(FN)较多。总的来说，随机森林比单个决策树做得好得多。
最后，我们看一下特征的重要性。
fi = pd.DataFrame(
{"feature": features, "importance": model.feature_importances_}
).sort_values("importance", ascending=False)
fi.head()
Last But Not Least¶
机器学习中常见且常用的分类或回归模型还有GBDT、XGBoost等，感兴趣的同学可以查阅资料，等有机会回头再来补充....