1.4 Machine Learning Classifiers/Models

本节我们简单介绍逻辑斯蒂回归 (Logistic Regression, LR), 支持向量机(support vector machine, SVM), 决策树(Decision Tree, DT) 这3个分类器的原理和Python/R的实现。

1) Logistic Regression

逻辑回归(Logistic Regression,LR)。在Kaggle竞赛的统计中,LR算法以63.5%的出产率,是处理军事和安全领域中出场率最高的算法。在实际场景中,逻辑回归同样应用广泛,大到国家各项经济政策的制定,小到计算广告CTR,都能看到LR算的身影。

除了应用广泛外,LR的建模过程还体现了数据建模中很重要的思想:对问题划分层次,并利用非线性变换和线性模型的组合,将未知的复杂问题分解为已知的简单问题。因此,我们可以说:理解好逻辑回归的细节,就掌握了数据建模的精髓。

逻辑斯蒂回归虽然名字叫做"回归",实际上解决的是分类的问题。在统计学中,逻辑斯蒂回归是一种二值响应的广义线性模型,在统计学中拟合广义线性模型都可以称为是"回归",这里的"回归"只是沿用了广义线性模型习惯叫法,并不是说它解决的是回归问题。

1a) 数学原理

逻辑回归假设数据服从伯努利分布,通过极大似然函数的方法,运用梯度下降来求解参数,来达到将数据二分类的目的。

逻辑回归是一个非线性模型,但是是其背后是以线性回归为理论支撑的。

提出一个与线性模型y=θTXby=\theta^T \centerdot X_b长相类似但不同的新公式:假设特征X所对应的y值是在指数上变化,那么就可以将结果y值取对数,作为其线性模型逼近的目标。也就是所谓的“对数线性回归”:ln(y)=θTXbln(y) = \theta^T \centerdot X_b

在“对数线性回归”的公式中,可以改写为 y=eθTXby = e^{\theta^T \centerdot X_b}。实际上是在求输入空间X到输出空间Y的非线性函数映射。对数函数的作用是将线性回归模型的预测值与真实标记联系起来

因此可以得到一个一般意义上的单调可微的“联系函数”g(a)=ln(a)g(a)=ln(a)。其本质就是给原来线性变换加上一个非线性变换(或者说映射),使得模拟的函数有非线性的属性,但本质上调参还是线性的,主体是内部线性的调参

那么对于解决分类问题的逻辑回归来说,我们需要找到一个“联系函数”,将线性回归模型的预测值与真实标记联系起来

将“概率”转换为“分类”的工具是“阶梯函数”

p^=f(x)y^={0p^0.51p^>0.5\hat {p} = f(x) \qquad \hat y = \begin{cases} 0& \hat p \leq 0.5 \\ 1& \hat p > 0.5 \end{cases}

但是这个阶梯函数不连续,不能作为“联系函数”g,因此使用对数几率函数来在一定程度上近似阶梯函数,将线性回归模型的预测值转化为分类所对应的概率。

σ(t)=11+et\sigma(t) = \frac{1} {1+e^{-t}}

如果令y为正例,1-y为负例,所谓的“几率”就是二者的比值y1y\frac {y} {1-y}。几率反映了样本x为正例的相对可能性。

“对数几率”就是对几率取对数lny1yln\frac {y} {1-y},对数几率实际上就是之前提到的sigmoid函数,将线性模型转化为分类。

如果令 y=11+eθTXby=\frac{1} {1+e^{-\theta^T \centerdot X_b}}1y=eθTXb1+eθTXb1-y = \frac{e^{-\theta^T \centerdot X_b}} {1+e^{-\theta^T \centerdot X_b}}。带入到对数几率中lny1y=θTXbln\frac {y} {1-y}=\theta^T \centerdot X_b

可以看出,sigmoid实际上就是用线性回归模型的预测结果取逼近真实值的对数几率,因此逻辑回归也被称为“对数几率回归”。

已经知道逻辑回归的模型:

p^=σ(θTxb)=11+eθTXby^={1p^0.50p^0.5\hat p = \sigma(\theta^T \centerdot x_b) = \frac{1} {1+e^{-\theta^T \centerdot X_b}} \qquad \hat y = \begin{cases} 1& \hat p \geq 0.5 \\ 0& \hat p \leq 0.5 \end{cases}

那么,如何求出未知参数θ\theta呢?

首先回顾一下线性回归。在线性回归中,做法如下:

由于已知 θTx_b\theta^T \centerdot x\_b 是估计值,于是用估计值与真值的差来度量模型的好坏。使用MSE(差值的平方和再平均)作为损失函数。 然后就可以通过导数求极值的方法,找到令损失函数最小的θ\theta了。

那么在逻辑回归中,解决思路也大致类似。

逻辑回归和线性回归最大的区别就是:逻辑回归解决的是分类问题,得到的y要么是1,要么是0。而我们估计出来的p是概率,通过概率决定估计出来的p到底是1还是0。因此,也可以将损失函数分成两类:

  • 如果给定样本的真实类别y=1,则估计出来的概率p越小,损失函数越大(估计错误)

  • 如果给定样本的真实类别y=0,则估计出来的概率p越大,损失函数越大(估计错误)

那么将用什么样的函数表示这两种情况呢,可以使用如下函数:

J={log(p^)ify=1log(1p^)ify=0J = \begin{cases} -log(\hat p) & if \quad y = 1 \\ -log(1 - \hat p) & if \quad y = 0 \end{cases}

由于模型是个二分类问题,分类结果y非0即1,因此我们可以使用一个巧妙的方法,通过控制系数的方式,将上面的两个式子合并成一个:

J(p^,y)=log(p^)ylog(1p^)1yJ(\hat p,y) = -log(\hat p)^{y}-log(1-\hat p)^{1-y}

以上是对于单个样本的误差值,那么求整个集合内的损失可以取平均值:

J(θ)=1mi=1my(i)log(p^(i))+(1y(i))log(1p^(i))J(\theta) = - \frac {1} {m} \sum^m_{i=1} y^{(i)} log(\hat p^{(i)}) + (1-y^{(i)})log(1-\hat p^{(i)})

然后,我们将 p^\hat p 替换成sigmoid函数,得到逻辑回归的损失函数如下

J(θ)=1mi=1my(i)log(σ(θTXb(i)))+(1y(i))log(1σ(θTXb(i)))J(\theta) = - \frac {1} {m} \sum^m_{i=1} y^{(i)} log(\sigma(\theta^T \centerdot X_b^{(i)})) + (1-y^{(i)})log(1-\sigma(\theta^T \centerdot X_b^{(i)}))

1b) 模型优缺点

优点

  • 实现简单,广泛的应用于工业问题上;

  • 分类时计算量非常小,速度很快,存储资源低;

  • 便利的观测样本概率分数;

  • 对逻辑回归而言,多重共线性并不是问题,它可以结合L2正则化来解决该问题;

  • 计算代价不高,易于理解和实现.

缺点

  • 当特征空间很大时,逻辑回归的性能不是很好;

  • 容易欠拟合,一般准确度不太高;

  • 不能很好地处理大量多类特征或变量;

  • 只能处理两分类问题(在此基础上衍生出来的softmax可以用于多分类),且必须线性可分;

  • 对于非线性特征,需要进行转换.

1c) 模型实现(Python)

我们基于python的sklearn.linear_model.LogisticRegression函数实现Iris数据集的二分类。

import numpy as np
import matplotlib.pyplot as plt
from sklearn import datasets
from sklearn.model_selection import train_test_split
from sklearn.linear_model import LogisticRegression
from sklearn.metrics import accuracy_score

# 加载数据
iris = datasets.load_iris()
X = iris.data
y = iris.target
X = X[y<2,:2]
y = y[y<2]

# 数据切分
X_train, X_test, y_train, y_test = train_test_split(X, y, random_state=666,test_size=0.25)
log_reg = LogisticRegression()
log_reg.fit(X_train, y_train)
# 查看测试数据集分类准确度
y_predict=log_reg.predict(X_test)
accuracy_score(y_test, y_predict)
"""
输出:1.0
"""

# 查看逻辑回归得到的概率
log_reg.predict_proba(X_test).T
"""
输出:
array([[0.09391473, 0.07585788, 0.79176992, 0.93705398,0.89573287,0.95218244, 0.86675827, 0.01504783, 0.02499892, 0.29394692,
0.8294444 , 0.9617504 , 0.74466265, 0.89573287, 0.11837917,0.35916592, 0.24085869, 0.68831525, 0.73479882, 0.79387392,0.91851748, 0.68279272, 0.02499892, 0.01448875, 0.89810141],
[0.90608527, 0.92414212, 0.20823008, 0.06294602, 0.10426713,0.04781756, 0.13324173, 0.98495217, 0.97500108, 0.70605308,0.1705556 , 0.0382496 , 0.25533735, 0.10426713, 0.88162083,0.64083408, 0.75914131, 0.31168475, 0.26520118, 0.20612608,0.08148252, 0.31720728, 0.97500108, 0.98551125, 0.10189859]])
"""

# 得到逻辑回归分类结果
log_reg.predict(X_test)
"""
输出:
array([1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 1, 1, 1, 0, 0, 0])
"""

1d) 模型实现(R)

data(iris)
set.seed(1)
######读取数据#######
head(iris)
# 去除类别为setosa的样本
df <- iris[iris$Species != 'setosa', ]
# 总样本数为输入矩阵的行数
n_samples <- nrow(df)
# 计算训练集样本数(80%)
n_train <- floor(n_samples * 0.8)  # Calculate the size of training sets

######切分数据集#######
# 产生所有样本序号的重新排列
indices <- sample(1:n_samples)
# 从以上随机排列中选取n_train个作为训练集样本序号
indices <- indices[1:n_train]
# 选取训练集样本
train_data <- df[indices,]
# 选取测试集样本
test_data <- df[-indices,]

######训练LR模型#######
fit.full <- glm(Species~Sepal.Length+Sepal.Width+Petal.Length+Petal.Width,
                data=train_data,family = binomial())
# 模型
summary(fit.full)
# 输出,查看模型
"""
Call:
glm(formula = Species ~ Sepal.Length + Sepal.Width + Petal.Length + 
    Petal.Width, family = binomial(), data = train_data)

Deviance Residuals: 
     Min        1Q    Median        3Q       Max  
-1.93995  -0.01458  -0.00080   0.00531   1.76567  

Coefficients:
             Estimate Std. Error z value Pr(>|z|)  
(Intercept)   -37.064     26.992  -1.373   0.1697  
Sepal.Length   -2.561      2.432  -1.053   0.2922  
Sepal.Width    -6.070      4.427  -1.371   0.1704  
Petal.Length    8.611      4.726   1.822   0.0685 .
Petal.Width    16.637      9.648   1.724   0.0846 .
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for binomial family taken to be 1)

    Null deviance: 109.097  on 79  degrees of freedom
Residual deviance:  11.635  on 75  degrees of freedom
AIC: 21.635

Number of Fisher Scoring iterations: 10
"""
######预测结果#######
pre <- predict.glm(fit.full, type = 'response', newdata = test_data) 
pre
# 结果输出,查看逻辑回归得到的概率
"""
          58           65           83           98          102          104          106          116 
4.932333e-09 7.691400e-08 3.891353e-07 6.861242e-06 9.992686e-01 9.993687e-01 1.000000e+00 9.999837e-01 
         118          122          124          130          140          141          142          143 
9.999995e-01 9.991476e-01 9.278001e-01 9.453108e-01 9.996240e-01 9.999997e-01 9.998213e-01 9.992686e-01 
         144          146          148          150 
9.999997e-01 9.999753e-01 9.978303e-01 9.700878e-01 
"""

值得注意的是,可以根据summary(fit.full)观察每个特征对应参数的P value(输出中的Pr(>|z|) ),可以进行特征筛选,仅使用显著的特征进行再次建模。本次示例中仅介绍了如何使用R实现逻辑回归的简单过程,如想了解更多可以help(glm),或者看这里

2) SVM

SVM学习的基本想法是求解能够正确划分训练数据集并且几何间隔最大的分离超平面。

2a) 决策边界

一个简单的二分类问题,在二维的特征平面中,所有的数据点分为了两类:蓝色圆形和黄色三角。我们的目标是找到了一条决策边界,将数据分类。但实际上我们可以找到多条决策边界。

分类中的“不适定问题”会影响模型的泛化性。比如在下面的模型中,被黄色箭头标出的点被决策边界划为蓝色圆点,但实际上它和黄色三角更近一些。也就说决策边界的选择,不仅要考虑已经存在的数据上的是否分类正确,还要考虑是否能够更好地划分未出现的测试数据:

逻辑回归算法如何求解决策边界的呢?

首先定义一个概率函数sigmoid函数:

σ(t)=11+et\sigma(t) = \frac{1} {1+e^{-t}}

然后根据概率函数进行建模:

P(Y=1)=11+eθTXbP(Y=1) = \frac {1} {1+e^{-\theta^T \centerdot X_b}}

建立损失函数:

J(θ)=1mi=1my(i)log(σ(θTXb(i)))+(1y(i))log(1σ(θTXb(i)))J(\theta) = - \frac {1} {m} \sum^m_{i=1} y^{(i)} log(\sigma(\theta^T \centerdot X_b^{(i)})) + (1-y^{(i)})log(1-\sigma(\theta^T \centerdot X_b^{(i)}))

最小化损失函数,从而求出一条决策边界(线性):

x2=θ0θ1x1θ2x_2 = \frac {-\theta_0 - \theta_1x_1} {\theta_2}

也就说,逻辑回归的损失函数完全是由训练数据集确定的。

SVM算法如何求解决策边界的呢?

SVM要找到一条泛化性比较好的决策边界,就是这条直线要离两个分类都尽可能的远,我们认为这样的决策边界就是好的。

如上图所示:在上面的两个类别中,离决策边界最近的点到决策边界的距离尽可能地远。

那也就是说,我们可以忽略其他大部分的数据点,只关注这几个特殊的点即可。

2b) SVM的数学推导

在搜索最优决策边界时,SVM要求在两个边界中离决策边界最近的点到决策边界的距离尽可能地远。为了定量描述这个过程我们引入支持向量最大间隔

支撑向量就是当我们将最优决策边界向上&下平移,在遇到第一个点时停下来,支撑着两条平移边界的点。

支撑向量到决策边界的距离是d;这两条平移后的直线的间隔(2d)被称为最大间隔Margin

所谓的支撑向量机,最初就是一个线性分类器,只不过这个线性分类器不仅能把样本分对,可以最大化Margin。这样我们就将SVM转换为了一个最优化问题,下面的工作就是求出Margin的数学表达式

回忆解析几何的知识,点(x,y)(x,y)到直线Ax+By+C=0Ax+By+C=0的距离:

d=Ax+By+CA2+B2d = \frac {|Ax+By+C|} {\sqrt {A^2+B^2}}

扩展到n维空间:点xxwTx+b=0w^Tx + b = 0(其中ww是n维向量,bb是截距)的距离为:

d=wTx+bw,w=w12+w22+...+w32d = \frac {w^Tx + b} {||w||}, \quad ||w||=\sqrt {w_1^2+w_2^2+...+w_3^2}

然后我们去找决策边界的表达式:

求出在满足下面的条件下,ww的值是多少。

{wTx(i)+b1y(i)=1wTx(i)+b1y(i)=1\begin{cases} w^Tx^{(i)}+b \ge 1 \quad \forall y^{(i)}=1 \\ w^Tx^{(i)}+b \le -1 \quad \forall y^{(i)}=-1 \\ \end{cases}

我们将上面的两个式子进行合并,即可以将y(i)y^{(i)}提到前面去。这样,支撑向量机的所有数据点都要满足下面的式子:

y(i)(wTx(i)+b)1y^{(i)}(w^Tx^{(i)}+b) \ge 1

对于任意支撑向量点xx到决策边界的距离为dd,我们要最大化Margin,将前面的式子带入后得到max2wTx+bwmax\frac {2\|w^Tx+b\|} {\|\|w\|\|},也就是max2wmax\frac {2} {\|\|w\|\|}

即相当于最小化min12wTwmin \frac {1}{2} w^Tw

OK,现在我们已经得到了SVM的最优化问题:

minΦ(x)=12wTws.t.y(i)(wTx(i)+b)1min \quad \Phi(x) = \frac {1}{2} w^Tw \\ s.t. \quad y^{(i)}(w^Tx^{(i)}+b) \ge 1

即最优化的目标函数为min12wTwmin\frac {1}{2} w^Tw,还要由限定条件(s.t.subjecttos.t. subject to,受限制于):所有数据满足y(i)(wTx(i)+b)1y^{(i)}(w^Tx^{(i)}+b) \ge 1

SVM的最优化问题是有限制条件的,对于有约束条件的最优化问题,用拉格朗日乘数法来解决,得到(aia_i是拉格朗日系数):

Lp=12w2i=1laiyi(wTxi+b)+i=1laiL_p = \frac{1}{2}||w||^2-\sum_{i=1}^l a_iy_i(w^T \centerdot x_i+b)+\sum_{i=1}^la_i

此时,我们要求LpL_p基于w,bw,b的极小值。

2c) Soft Margin SVM

线性可分问题中,对于样本点来说,存在一根直线可以将样本点划分,我们称之为Hard Margin SVM;但是(同样线性不可分),有时候会出现不那么完美,样本点会有一些噪声或者异常点,并不能完全分开。即没有一条直线可以将样本分成两类。那么就提出了Soft Margin SVM

Soft Margin SVM的思想也很朴素,就是在Hard Soft的基础上,将原来的约束条件放宽一些,增加容错性。

Hard Soft中的约束条件为:

Φ(x)=12wTws.t.y(i)(wTx(i)+b)1\Phi(x) = \frac {1}{2} w^Tw \\ s.t. \quad y^{(i)}(w^Tx^{(i)}+b) \ge 1

对于限制条件y(i)(wTx(i)+b)1y^{(i)}(w^Tx^{(i)}+b)\ge 1,形象地说,Margin区域里必须是任何数据点都没有,所有的数据点都必须在(y(i)wTx(i)+b)=1(y^{(i)}w^Tx^{(i)}+b)=1(y(i)wTx(i)+b)=1(y^{(i)}w^Tx^{(i)}+b)=-1两条直线的外侧。

如果有些数据点不能满足这个要求,就对条件加以宽松,在margin区域外给他一个宽松量ζi\zeta_i(大于等于0):

Φ(x)=12wTws.t.y(i)(wTx(i)+b)1ζi\Phi(x) = \frac {1}{2} w^Tw \\ s.t. \quad y^{(i)}(w^Tx^{(i)}+b) \ge 1-\zeta_i

但是我们很容易地想到:容错空间ζi\zeta_i也不能无限制的放大。在最小化的同时加上所有点的容错空间的和,就可以在最小化的同时又可以容忍一点程度的错误。并且通过参数CC来平衡重要程度。

因此Soft Margin SVM最优化问题对应的数学表达式为:

Φ(x)=12wTw+Ci=1mζis.t.y(i)(wTx(i)+b)1ζiζi0\Phi(x) = \frac {1}{2} w^Tw + C\sum_{i=1}^{m}\zeta_i\\ s.t. \quad y^{(i)}(w^Tx^{(i)}+b) \ge 1-\zeta_i \\ \zeta_i \ge 0

在这里i=1mζi\sum_{i=1}^{m}\zeta_i,其实相当于在Soft Margin SVM中加入了L1正则,避免因为训练出的模型往极端的方式发展,让模型的拥有一定的容错能力,泛化性有所提升。相应的,也有L2正则,即把i=1mζi\sum_{i=1}^{m}\zeta_i换成i=1mζi2\sum_{i=1}^{m}\zeta_i^2

系数CC越大,相应的容错空间越小。如果取正无穷,意味着逼迫着容错空间ζi\zeta_i趋近于零,也就变成了Hard Margin SVM

Soft Margin的目标函数长什么样呢?还是用拉格朗日的方法:

LP=12w2+Ci=1lζii=1lai[yi(wxi+b)1+ζi]i=1lμiζiL_P = \frac{1}{2}||w^2||+C\sum_{i=1}^l\zeta_i - \sum_{i=1}^l a_i[y_i(w \centerdot x_i + b)-1+\zeta_i]-\sum_{i=1}^l \mu_i\zeta_i

由于优化函数增加了惩罚项,且增加了一些约束条件,因此目标函数变得很复杂了。

2d) 模型优缺点

优点

  • 可以解决高维问题,即大型特征空间;

  • 解决小样本下机器学习问题;

  • 能够处理非线性特征的相互作用;

  • 无局部极小值问题;(相对于神经网络等算法)

  • 无需依赖整个数据;

  • 泛化能力比较强.

缺点:

  • 当观测样本很多时,效率并不是很高;

  • 对非线性问题没有通用解决方案,有时候很难找到一个合适的核函数;

  • 对于核函数的高维映射解释力不强,尤其是径向基函数;

  • 常规SVM只支持二分类;

  • 对缺失数据敏感.

2e) 模型实现(Python)

我们基于python的sklearn.svm.LinearSVC函数实现Iris数据集的二分类。值得注意的是,LinearSVC实现了线性分类支持向量机,它是给根据liblinear实现的,可以用于二类分类,也可以用于多类分类。另外python的sklearn.svm.SVC是基于libsvm库实现的支持向量机,支持linear, poly, rbf, sigmoid, precomputed多个核函数,通过kernel参数进行指定。

import numpy as np
import matplotlib.pyplot as plt
from sklearn import datasets
from sklearn.preprocessing import StandardScaler
from sklearn.model_selection import train_test_split
from sklearn.svm import LinearSVC

# 加载数据
iris = datasets.load_iris()
X = iris.data
y = iris.target
X = X[y<2,:2]
y = y[y<2]
# 数据归一化
standardScaler = StandardScaler()
standardScaler.fit(X)
X_std = standardScaler.transform(X)
# 数据切分
X_train, X_test, y_train, y_test = train_test_split(X, y, random_state=666,test_size=0.25)
svc = LinearSVC(C=1)
svc.fit(X_train,y_train)
# 查看测试数据集分类准确度
svc.score(X_test, y_test)
"""
输出:1.0
"""
# 得到svm分类结果
svc.predict(X_test)
"""
输出:
array([1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 0,1, 1, 0])
"""

2f) 模型实现(R)

这里我们利用e1071这个包进行SVM计算,这个包是基于libsvm库实现的,支持多个核函数,通过参数kernel进行指定。这里我们使用linear核函数,实现线性分类支持向量机

data(iris)
set.seed(1)
library('e1071')
######读取数据#######
head(iris)
# 去除类别为setosa的样本
df <- iris[iris$Species != 'setosa', ]
all_data <- df[, 1:(ncol(df) - 1)]
all_classes <- factor(df$Species)
# 总样本数为输入矩阵的行数
n_samples <- nrow(all_data)
# 计算训练集样本数(80%)
n_train <- floor(n_samples * 0.8)  # Calculate the size of training sets

######切分数据集#######
# 产生所有样本序号的重新排列
indices <- sample(1:n_samples)
# 从以上随机排列中选取n_train个作为训练集样本序号
indices <- indices[1:n_train]
# 选取训练集样本
train_data <- all_data[indices,]
train_classes <- all_classes[indices]
# 选取测试集样本
test_data <- all_data[-indices,]
test_classes <- all_classes[-indices]

######训练SVM模型#######
model<-svm(train_data, train_classes,type='C',kernel='linear',probability=T)
# 模型
summary(model)
# 输出,查看模型
"""
Call:
svm.default(x = train_data, y = train_classes, type = "C", kernel = "linear", probability = T)


Parameters:
   SVM-Type:  C-classification 
 SVM-Kernel:  linear 
       cost:  1 

Number of Support Vectors:  15

 ( 7 8 )


Number of Classes:  2 

Levels: 
 versicolor virginica

"""
######预测结果#######
pred<-predict(model,test_data,probability=T)
pred
# 结果输出,查看SVM预测得到标签和概率
"""
        58         65         83         98        102        104        106        116        118        122 
versicolor versicolor versicolor versicolor  virginica  virginica  virginica  virginica  virginica  virginica 
       124        130        140        141        142        143        144        146        148        150 
 virginica  virginica  virginica  virginica  virginica  virginica  virginica  virginica  virginica  virginica 
attr(,"probabilities")
     versicolor   virginica
58  0.998263719 0.001736281
65  0.996828184 0.003171816
83  0.993788977 0.006211023
98  0.985421225 0.014578775
102 0.048124518 0.951875482
104 0.056467147 0.943532853
106 0.002234467 0.997765533
116 0.017955780 0.982044220
118 0.009627553 0.990372447
122 0.051792276 0.948207724
124 0.230168648 0.769831352
130 0.278480368 0.721519632
140 0.055826303 0.944173697
141 0.003998612 0.996001388
142 0.041970244 0.958029756
143 0.048124518 0.951875482
144 0.004313567 0.995686433
146 0.019069236 0.980930764
148 0.090848574 0.909151426
150 0.197592888 0.802407112
Levels: versicolor virginica
"""

3) Decision Tree

决策树是一个非常有意思的模型,它的建模思路是尽可能模拟人做决策的过程。因此决策树几乎没有任何抽象,完全通过生成决策规则来解决分类和回归问题。因为它的运行机制能很直接地被翻译成人类语言,即使对建模领域完全不了解的非技术人员也能很好地理解它。因此在学术上被归为白盒模型(white box model)。

3a) DT的基本思想

决策树是一种常见的机器学习算法,类似流程图的结构,其中每个内部节点表示一个测试功能,即类似做出决策的过程(动作),比如下面这个“相亲决策树”:

决策树的思想还是非常直观的。用决策树分类:从根节点开始,对实例的某一特征进行测试,根据测试结果将实例分配到其子节点,此时每个子节点对应着该特征的一个取值,如此递归的对实例进行测试并分配,直到到达叶节点,最后将实例分到叶节点的类中。

3b) DT的数学原理

假设给定训练数据集D=(x1,y1),(x2,y2),..(xN,yN)D={(x_1,y_1),(x_2,y_2),……..(x_N,y_N)},其中xi=(xi(1),xi(2),.,xi(n))Txi=(x_i(1) ,x_i(2),…….,x_i(n))^T为输入实例(特征向量),n为特征个数,yi12..Kyi∈{1,2,…..K}为类标记(label),i=123,Ni=1,2,3,……,N,N为样本容量。

学习目标:根据给定的训练数据集构建一个决策模型,使它能够对实例进行正确的分类。

决策树学习本质上是从训练数据集中归纳出一组分类规则。与训练数据集不相矛盾的决策树(即能对训练数据进行正确分类的决策树)可能是0个或多个。我们需要找到一个与训练数据矛盾较小的决策树,同时具有很好的泛化能力

与其他模型相同,决策树学习用损失函数表示这一目标。决策树学习的损失函数通常是正则化的极大似然函数。决策树学习的策略是以损失函数为目标函数的最小化

当损失函数确定以后,学习问题就变为在损失函数意义下选择最优决策树的问题。因为从所有可能的决策树中选取最优决策树是NP完全问题,所以现实中决策树学习算法通常采用启发式方法,近似求解这一最优化问题。

3c) DT的构造步骤

决策树学习的算法通常是一个递归地选择最优特征,并根据该特征对训练数据进行分割,使得对各个子数据集有一个最好的分类的过程。这一过程对应着对特征空间的划分,也对应着决策树的构建。决策树通常有三个步骤

  • 特征选择

  • 决策树的生成

    (1) 开始构建根节点,将所有训练数据都放在根节点,选择一个最优特征,按着这一特征将训练数据集分割成子集,使得各个子集有一个在当前条件下最好的分类。 (2) 如果这些子集已经能够被基本正确分类,那么构建叶节点,并将这些子集分到所对应的叶子节点去。 (3) 如果还有子集不能够被正确的分类,那么就对这些子集选择新的最优特征,继续对其进行分割,构建相应的节点。 (4) 如此递归进行,直至所有训练数据子集被基本正确的分类,或者没有合适的特征为止。每个子集都被分到叶节点上,即都有了明确的类,这样就生成了一颗决策树。

  • 决策树的修剪

    生成的决策树可能发生过拟合现象。我们需要对已生成的树自下而上进行剪枝,将树变得更简单,从而使其具有更好的泛化能力。具体地,就是去掉过于细分的叶结点,使其回退到父结点,甚至更高的结点,然后将父结点或更高的结点改为新的叶结点,从而使得模型有较好的泛化能力。

根据不同的特征重要性的衡量方式,常用的决策树有ID3C4.5

3d) 模型优缺点

优点:

  • 决策树易于理解和解释,可以可视化分析,容易提取出规则;

  • 可以同时处理标称型和数值型数据;

  • 比较适合处理有缺失属性的样本;

  • 能够处理不相关的特征;

  • 在相对短的时间内能够对大型数据源做出可行且效果良好的结果.

缺点:

  • 容易发生过拟合(随机森林可以很大程度上减少过拟合);

  • 容易忽略数据集中属性的相互关联;

  • 对于那些各类别样本数量不一致的数据,在决策树中,进行属性划分时,不同的判定准则会带来不同的属性选择倾向.

3e) 模型实现(Python)

我们基于python的sklearn.tree.DecisionTreeClassifier函数实现Iris数据集的二分类。值得注意的是决策树算法有很多种,例如常用的CART、C4.5算法,sklearn.tree.DecisionTreeClassifier函数使用CART算法的优化版本,scikit-learn实现暂时不支持分类变量。

import numpy as np
import matplotlib.pyplot as plt
from sklearn import datasets
from sklearn.model_selection import train_test_split
from sklearn.tree import DecisionTreeClassifier
from sklearn.metrics import accuracy_score

# 加载数据
iris = datasets.load_iris()
X = iris.data
y = iris.target
X = X[y<2,:2]
y = y[y<2]

# 数据切分
X_train, X_test, y_train, y_test = train_test_split(X, y, random_state=666,test_size=0.25)
# 创建决策树对象,最大深度max_depth为2层,criterion评判标准为entropy(熵)
dt_clt = DecisionTreeClassifier(max_depth=2,criterion='entropy')
dt_clt.fit(X_train,y_train)
# 查看测试数据集分类准确度
y_predict=dt_clt.predict(X_test)
accuracy_score(y_test, y_predict)
"""
输出:1.0
"""

# 查看决策树得到的概率
dt_clt.predict_proba(X_test).T
"""
输出:
array([[0., 0., 0.96774194, 0.96774194, 0.96774194,0.96774194, 0.96774194, 0., 0., 0.,0.96774194, 0.96774194, 0.96774194, 0.96774194,0.0.2, 0., 0.96774194, 1., 0.96774194,0.96774194, 0.96774194, 0., 0., 0.96774194],
[1.,1., 0.03225806, 0.03225806, 0.03225806,0.03225806,0.03225806,1.,1.,1.,0.03225806,0.03225806, 0.03225806, 0.03225806, 1.,0.8,1., 0.03225806, 0., 0.03225806,0.03225806, 0.03225806, 1., 1., 0.03225806]])
"""

# 得到决策树分类结果
dt_clt.predict(X_test)
"""
输出:
array([1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 0,1, 1, 0])
"""

3f) 模型实现(R)

针对不同的决策树算法,例如常用的CART、C4.5算法,R语言中有不同的包可以实现。这里我们以基于CART算法的决策树为例。

data(iris)
set.seed(1)
library(rpart)
######读取数据#######
head(iris)
# 去除类别为setosa的样本
df <- iris[iris$Species != 'setosa', ]
# 总样本数为输入矩阵的行数
n_samples <- nrow(all_data)
# 计算训练集样本数(80%)
n_train <- floor(n_samples * 0.8)  # Calculate the size of training sets

######切分数据集#######
# 产生所有样本序号的重新排列
indices <- sample(1:n_samples)
# 从以上随机排列中选取n_train个作为训练集样本序号
indices <- indices[1:n_train]
# 选取训练集样本
train_data <- df[indices,]
# 选取测试集样本
test_data <- df[-indices,]

######训练决策树模型#######
# 初步训练DT模型
dtree<-rpart(Species~.,data=train_data,method="class", minsplit=2,parms=list(split="gini"))  #使用CART算法,若split="information",则为ID3算法。
# 模型
printcp(dtree)
# 输出,查看模型
"""
Classification tree:
rpart(formula = Species ~ ., data = train_data, method = "class", 
    parms = list(split = "gini"), minsplit = 2)

Variables actually used in tree construction:
[1] Petal.Length Petal.Width  Sepal.Length

Root node error: 34/80 = 0.425

n= 80 

        CP nsplit rel error  xerror     xstd
1 0.852941      0  1.000000 1.00000 0.130045
2 0.044118      1  0.147059 0.26471 0.083124
3 0.029412      3  0.058824 0.26471 0.083124
4 0.014706      4  0.029412 0.20588 0.074334
5 0.010000      6  0.000000 0.23529 0.078920
"""

# prune()剪枝,提高模型的泛化能力
dtree$cptable  #剪枝前的复杂度 ,nsplit分枝数,rel error训练集对应的误差 ,xerror 10折交叉验证误差,xstd交叉验证误差的标准差
"""
          CP nsplit  rel error    xerror       xstd
1 0.85294118      0 1.00000000 1.0000000 0.13004524
2 0.04411765      1 0.14705882 0.2647059 0.08312402
3 0.02941176      3 0.05882353 0.2647059 0.08312402
4 0.01470588      4 0.02941176 0.2058824 0.07433384
5 0.01000000      6 0.00000000 0.2352941 0.07892005
"""
tree<-prune(dtree,cp=0.02)
tree$cptable #剪枝后的复杂度
"""
          CP nsplit  rel error    xerror       xstd
1 0.85294118      0 1.00000000 1.0000000 0.13004524
2 0.04411765      1 0.14705882 0.2647059 0.08312402
3 0.02941176      3 0.05882353 0.2647059 0.08312402
4 0.02000000      4 0.02941176 0.2058824 0.07433384
"""

######预测结果#######
prob <- predict(tree,newdata=test_data,type="prob")
prob
# 结果输出,查看决策树得到的概率
"""
    setosa versicolor virginica
58       0 0.00000000 1.0000000
65       0 1.00000000 0.0000000
83       0 1.00000000 0.0000000
98       0 1.00000000 0.0000000
102      0 0.03225806 0.9677419
104      0 0.03225806 0.9677419
106      0 0.03225806 0.9677419
116      0 0.03225806 0.9677419
118      0 0.03225806 0.9677419
122      0 0.03225806 0.9677419
124      0 0.03225806 0.9677419
130      0 1.00000000 0.0000000
140      0 0.03225806 0.9677419
141      0 0.03225806 0.9677419
142      0 0.03225806 0.9677419
143      0 0.03225806 0.9677419
144      0 0.03225806 0.9677419
146      0 0.03225806 0.9677419
148      0 0.03225806 0.9677419
150      0 0.03225806 0.9677419
"""
pr_label <- predict(tree,newdata=test_data,type="class")
pr_label
# 结果输出,查看决策树得到的标签
"""
        58         65         83         98        102        104        106        116        118        122 
 virginica versicolor versicolor versicolor  virginica  virginica  virginica  virginica  virginica  virginica 
       124        130        140        141        142        143        144        146        148        150 
 virginica versicolor  virginica  virginica  virginica  virginica  virginica  virginica  virginica  virginica 
Levels: setosa versicolor virginica
"""

4) More Classifiers/Models with Python

方法含义python function

KNN

在k-NN分类中,输出是一个分类族群。一个对象的分类是由其邻居的“多数表决”确定的,k个最近邻居中最常见的分类决定了赋予该对象的类别。

from sklearn.neighbors.KNeighborsClassifier

RF

随机森林(RandomForest), 指的是利用多棵树对样本进行训练并预测的一种分类器。对于每棵树,它们使用的训练集是从总的训练集中有放回采样出来的。

from sklearn.ensemble.RandomForestClassifier

NB

朴素贝叶斯的思想基础是这样的:对于给出的待分类项,求解在此项出现的条件下各个类别出现的概率,哪个最大,就认为此待分类项属于哪个类别。

sklearn.naive_bayes.GaussianNB

BPNN

BP神经网络是一种多层网络算法,其核心是反向传播误差,即:使用梯度下降法(或其他算法),通过反向传播来不断调整网络的权值和阈值,使网络的误差平方和最小。

from sklearn.neural_network.MLPClassifier

5) References

Take a break

Bayes 贝叶斯

不论是学习概率统计还是机器学习的过程中,贝叶斯 (Bayes) 总是是绕不过去的一道坎。

据现有资料显示,贝叶斯全名为托马斯·贝叶斯(Thomas Bayes,1701?-1761),是一位与牛顿同时代的牧师,是一位业余数学家,平时就思考些有关上帝的事情,当然,统计学家都认为概率这个东西就是上帝在掷骰子。当时贝叶斯发现了古典统计学当中的一些缺点,从而提出了自己的“贝叶斯统计学”,也就是贝叶斯定理。托马斯·贝叶斯的洞察力非常简明。一个假设是真实的概率取决于两个标准:1.根据当前的知识(“先验”),判断它的合理程度;2.评估它与新的证据的契合程度。

但贝叶斯统计当中由于引入了一个主观因素(先验概率),一点都不被当时的人认可。直到20世纪中期,也就是快200年后了,统计学家在古典统计学中遇到了瓶颈,伴随着计算机技术的发展,当统计学家使用贝叶斯统计理论时发现能解决很多之前不能解决的问题,从而贝叶斯统计学一下子火了起来。

运用先前的经验和记忆中积累的知识、和意识中提炼出的新证据,我们对日常事物的概率进行分配和生活进行管理。举一个生活中的简单事件:

接听手机。通常在工作时你将它放在办公桌上,而在家里时把它放在在充电器上。现在你在家里的阳台上浇花,听到屋内电话声响起。新的数据会告诉你它处于室内任何地方,但你仍会直接走向充电器。因为你将的先前对手机位置的认知(通常在办公桌上或家中的充电器上)与新的证据(房屋的某处)相结合,从而确定了它的位置。如果手机不在充电器上,那么你会唤起先前你在某些放置过手机的位置的认知来缩小搜索范围。 你会忽略房子里大部分的地方,如冰箱、袜子抽屉等等,因为这些地方在你先前所积累的认知中被认定为极不可能的位置,你会在最终找到电话之前思考最可能的地方。

在这个找电话的过程中,你便正在使用贝叶斯定理。

贝叶斯之谜:

1763年12月23日,托马斯·贝叶斯( Thomas Bayes,1701?-1761 )的遗产受赠者R.Price牧师在英国皇家学会宣读了贝叶斯的遗作《论机会学说中一个问题的求解》,其中给出了贝叶斯定理,这一天现在被当作贝叶斯定理的诞生日。虽然贝叶斯定理在今天已成为概率论统计最经典的内容之一,但贝叶斯本人却笼罩在谜团之中。

现有资料表明,贝叶斯是一位神职人员,长期担任英国坦布里奇韦尔斯地方教堂的牧师,他从事数学研究的目的是为了证明上帝的存在。他在1742年当选为英国皇家学会会士,但没有记录表明他此前发表过任何科学或数学论文。他的提名是由皇家学会的重量级人物签署的,但为什么提名以及他为何能当选,至今任是个谜。贝叶斯的研究工作和他本人在他生活的时代很少有人关注,贝叶斯定理出现后很快就被遗忘了,后来大数学家拉普拉斯使它重新被科学界所熟悉,但直到二十世纪随着统计学的广泛应用才备受瞩目。贝叶斯的出生年份至今也没有清楚确定,甚至关于如今广泛流传的他的画像是不是贝叶斯本人,也仍然存在争议。

—— 摘自《机器学习》(周志华)

Last updated