八、支持向量机(下)

8.1 二分类SVC的进阶

8.1.1 SVC用于二分类的原理复习

工作原理：支持向量机分类器，是在数据空间中找出一个超平面作为决策边界，利用这个决策边界来对数据进行分类，并使分类误差尽量小的模型;

平行于决策边界的两条虚线是距离决策边界相对距离为1的超平面，他们分别压过两类样本中距离决策边界最近的样本点，这些样本点就被成为支持向量*;

8.1.2 参数C的理解进阶

• 存在软间隔的数据：有一些数据，可能是线性可分，但在线性可分状况下训练准确率不能达到100%，即无法让训练误差为0：

此时，我们得要让我们决策边界能够忍受一小部分训练误差，我们就不能单纯地寻求最大边际了。因为边际越大被分错的样本也就会越多，因此我们 需要找出一个”最大边际“与”被分错的样本数量“之间的平衡 ——>松弛系数ζ及其惩罚项系数C；

而在软间隔情况下我们的边际依然由支持向量决定，但此时此刻的支持向量可能就不再是来自两种标签类别的点了，而是分布在决策边界两边的，同类别的点：

要想法子把下面那个误分类的紫色点弄成分类正确的

此时我们的虚线超平面ωx_i+b=1-ζ_i是由混杂在红色点中间的紫色点来决定的，所以此时此刻，这个紫色点就是我们的支持向量了；
所以软间隔让决定两条虚线超平面的支持向量可能是来自于同一个类别的样本点，而硬间隔的时候两条虚线超平面必须是由来自两个不同类别的支持向量决定的；
而C值则是决定了会决定我们究竟是依赖红色点作为支持向量（只追求最大边界），还是我们要依赖软间隔中，混杂在红色点中的紫色点来作为支持向量（追求最大边界和判断正确的平衡）：如果C值设定比较大，那SVC可能会选择边际较小的，能够更好地分类所有训练点的决策边界，不过模型的训练时间也会更长。如果C的设定值较小，那SVC会尽量最大化边界，尽量将掉落在决策边界另一方的样本点预测正确，决策功能会更简单，但代价是训练的准确度。

观察一下我们对不同数据集分类时，支持向量都有哪些？软间隔如何影响了超平面和支持向量，就一目了然了：

import numpy as np
import matplotlib.pyplot as plt
from matplotlib.colors import ListedColormap
from sklearn import svm
from sklearn.datasets import make_circles, make_moons, make_blobs,make_classification
n_samples = 100
datasets = [
        make_moons(n_samples=n_samples, noise=0.2, random_state=0),
        make_circles(n_samples=n_samples, noise=0.2, factor=0.5, random_state=1),
        make_blobs(n_samples=n_samples, centers=2, random_state=5),
        make_classification(n_samples=n_samples,n_features = 
2,n_informative=2,n_redundant=0, random_state=5)
 ]
Kernel = ["linear"] #四个数据集分别是什么样子呢？
for X,Y in datasets:
        plt.figure(figsize=(5,4))
        plt.scatter(X[:,0],X[:,1],c=Y,s=50,cmap="rainbow")
nrows=len(datasets)
ncols=len(Kernel) + 1
fig, axes = plt.subplots(nrows, ncols,figsize=(10,16))
#第一层循环：在不同的数据集中循环
for ds_cnt, (X,Y) in enumerate(datasets):

        ax = axes[ds_cnt, 0]
        if ds_cnt == 0:
            ax.set_title("Input data")
        ax.scatter(X[:, 0], X[:, 1], c=Y, zorder=10, cmap=plt.cm.Paired,edgecolors='k')
        ax.set_xticks(())
        ax.set_yticks(())

        for est_idx, kernel in enumerate(Kernel):
            ax = axes[ds_cnt, est_idx + 1]

            clf = svm.SVC(kernel=kernel, gamma=2).fit(X, Y)
            score = clf.score(X, Y)

            ax.scatter(X[:, 0], X[:, 1], c=Y
                       ,zorder=10
                       ,cmap=plt.cm.Paired,edgecolors='k')
            ax.scatter(clf.support_vectors_[:, 0], clf.support_vectors_[:, 1], s=100,
                        facecolors='none', zorder=10, edgecolors='white')

            x_min, x_max = X[:, 0].min() - .5, X[:, 0].max() + .5
            y_min, y_max = X[:, 1].min() - .5, X[:, 1].max() + .5

            XX, YY = np.mgrid[x_min:x_max:200j, y_min:y_max:200j]
            Z = clf.decision_function(np.c_[XX.ravel(), YY.ravel()]).reshape(XX.shape)
            ax.pcolormesh(XX, YY, Z > 0, cmap=plt.cm.Paired)
            ax.contour(XX, YY, Z, colors=['k', 'k', 'k'], linestyles=['--', '-', '--'],
                        levels=[-1, 0, 1])
            ax.set_xticks(())
            ax.set_yticks(())

            if ds_cnt == 0:
                ax.set_title(kernel)

            ax.text(0.95, 0.06, ('%.2f' % score).lstrip('0')
                   , size=15
                   , bbox=dict(boxstyle='round', alpha=0.8, facecolor='white')
                   #为分数添加一个白色的格子作为底色
                   , transform=ax.transAxes #确定文字所对应的坐标轴，就是ax子图的坐标轴本身
                   , horizontalalignment='right' #位于坐标轴的什么方向
                   )
            
plt.tight_layout()
plt.show()

数据的原始分布形式：（分别是月亮型、环形、混合以及软间隔

以及绘制超平面后的样子：

上面的图中可以看到，两条超平面内所有的点和虚线超平面外，但属于另一个类别的点全都被当作支持向量了

8.1.3 二分类SVC中的样本不均衡问题：重要参数class_weight

• 样本不均衡问题：样本不均衡是指在一组数据集中，标签的一类天生占有很大的比例，但我们有着捕捉出某种特定的分类的需求的状况；
• 带来的问题：
  (1)分类模型天生会倾向于多数的类，让多数类更容易被判断正确，少数类被牺牲掉;
  (2)模型评估指标会失去意义；

先要让算法意识到数据的标签是不均衡的，通过施加一些惩罚或者改变样本本身，来让模型向着捕获少数类的方向建模；然后，我们要改进我们的模型评估指标，使用更加针对于少数类的指标来优化模型；

• 在逻辑回归中，参数class_weight默认None，此模式表示假设数据集中的所有标签是均衡的，即自动认为标签的比例是1：1，当样本不均衡的时候，我们可以使用形如{"标签的值1"：权重1，"标签的值2"：权重2}的字典来输入真实的样本标签比例，来让算法意识到样本是不平衡的。或者使用”balanced“模式，直接使用n_samples/(n_classes * np.bincount(y))作为权重，可以比较好地修正我们的样本不均衡情况；

但在SVM中，我们的分类判断是基于决策边界的，而最终决定究竟使用怎样的支持向量和决策边界的参数是参数C，所以所有的样本均衡都是通过参数C来调整的：

• SVC的参数：class_weight：可输入字典或者"balanced”，可不填，默认None，对SVC，将类i的参数C设置为class_weight [i] * C。如果希望改善样本不均衡状况，请输入形如{"标签的值1"：权重1，"标签的值2"：权重2}的字典，则参数C将会自动被设为：标签的值1的C：权重1 * C，标签的值2的C：权重2 * C。或者使用balance模式；
• SVC的接口fit的参数：sample_weight：结构为 (n_samples, )，必须对应输入fit中的特征矩阵的每个样本，对应着每个样本的权重（class_weight对应类别权重，sample_weight对应样本权重）。每个样本在fit时的权重，让权重 * 每个样本对应的C值来迫使分类器强调设定的权重更大的样本。通常，较大的权重加在少数类的样本上，以 迫使模型向着少数类的方向建模；

通常来说，这两个参数我们只选取一个来设置。如果我们同时设置了两个参数，则C会同时受到两个参数的影响，即 class_weight中设定的权重 * sample_weight中设定的权重 * C

如何使用参数：首先，我们来自建一组样本不平衡的数据集。我们在这组数据集上 建两个SVC模型，一个设置有class_weight参数，一个不设置class_weight参数 。我们对两个模型分别进行评估并画出他们的决策边界，以此来观察class_weight带来的效果

1. 导入需要的库和模块

import numpy as np
import matplotlib.pyplot as plt
from sklearn import svm
from sklearn.datasets import make_blobs

2. 创建样本不均衡的数据集

class_1 = 500 #类别1有500个样本，10：1
class_2 = 50 #类别2只有50个
centers = [[0.0, 0.0], [2.0, 2.0]] #设定两个类别的中心
clusters_std = [1.5, 0.5] #设定两个类别的方差，通常来说，样本量比较大的类别会更加松散，方差就大
#n_samples可以输入列表，列表中有几类值就说明有几类标签
X, y = make_blobs(n_samples=[class_1, class_2],
                      centers=centers,
                      cluster_std=clusters_std,
                      random_state=0, shuffle=False)

#看看数据集长什么样
plt.scatter(X[:, 0], X[:, 1], c=y, cmap="rainbow",s=10) #其中红色点是少数类，紫色点是多数类

3. 在数据集上分别建模

#不设定class_weight，C默认是1，这里不写也行
clf = svm.SVC(kernel='linear', C=1.0)
clf.fit(X, y) 

#设定class_weight
wclf = svm.SVC(kernel='linear', class_weight={1: 10})
wclf.fit(X, y) 

#给两个模型分别打分看看，这个分数是accuracy准确度
clf.score(X,y)
#0.9418181818181818

wclf.score(X,y)
#0.9127272727272727
#发现做了样本均衡之后准确率下降

4. 绘制两个模型下数据的决策边界
   用Contour我们只需要在我们的样本构成的平面上，把所有到决策边界的距离为0的点相连，就是我们的决策边界。而到决策边界的距离可以使用我们说明过的接口decision_function来调用：

#4. 绘制两个模型下数据的决策边界

#首先要有数据分布
plt.figure(figsize=(6,5))
plt.scatter(X[:, 0], X[:, 1], c=y, cmap="rainbow",s=10)
ax = plt.gca() #获取当前的子图，如果不存在，则创建新的子图

#绘制决策边界的第一步：要有网格
xlim = ax.get_xlim()
ylim = ax.get_ylim()

xx = np.linspace(xlim[0], xlim[1], 30)
yy = np.linspace(ylim[0], ylim[1], 30)
YY, XX = np.meshgrid(yy, xx)
xy = np.vstack([XX.ravel(), YY.ravel()]).T 

#第二步：找出我们的样本点到决策边界的距离
Z_clf = clf.decision_function(xy).reshape(XX.shape) 
a = ax.contour(XX, YY, Z_clf, colors='black', levels=[0], alpha=0.5, linestyles=['-'])

Z_wclf = wclf.decision_function(xy).reshape(XX.shape) 
b = ax.contour(XX, YY, Z_wclf, colors='red', levels=[0], alpha=0.5, linestyles=['-'])

#第三步：画图例
plt.legend([[*a.legend_elements()][0][0], [*b.legend_elements()][0][0]], ["non weighted", "weighted"],
           loc="upper right")#upper right表示要先是在右上角
plt.show()

上面代码的新东西：

a.collections #调用这个等高线对象中画的所有线，返回一个惰性对象
#<a list of 1 PathCollection objects>

#用[*]把它打开试试看
[*a.collections] #返回了一个linecollection对象，其实就是我们等高线里所有的线的列表
#[<matplotlib.collections.PathCollection at 0x245ad285a90>]

#现在我们只有一条线，所以我们可以使用索引0来锁定这个对象
a.collections[0]
#<matplotlib.collections.PathCollection at 0x245ad285a90>

#plt.legend([对象列表],[图例列表],loc)
#只要对象列表和图例列表相对应，就可以显示出图例

灰色是我们做样本平衡之前的决策边界。灰色线上方的点被分为一类，下方的点被分为另一类。红色是我们做样本平衡之后的决策边界，同样是红色线上方一类，红色线下方一类。可以看到，平衡前大约有一半少数类（红色）被分错，多数类（紫色点）几乎都被分类正确了；

不难得出，在成功分类了大多数红色点的同时，把更多的紫色点给误分类了。做了样本平衡后，少数类几乎全部都被分类正确了，但是多数类有许多被分错了。
可以看出，从准确率的角度来看，不做样本平衡的时候准确率反而更高，做了样本平衡准确率反而变低了，这是因为做了样本平衡后，为了要更有效地捕捉出少数类，模型误伤了许多多数类样本，而多数类被分错的样本数量 > 少数类被分类正确的样本数量，使得模型整体的精确性下降；

而在现实中，我们往往都在追求捕捉少数类，因为在很多情况下，将少数类判断错的代价是巨大的。

8.2 SVC的模型评估指标

如果我们的目标是希望尽量捕获少数类，那准确率这个模型评估逐渐失效，所以我们需要新的模型评估指标来帮助我们。其实我们只需要查看模型在少数类上的准确率就好了，只要能够将少数类尽量捕捉出来，就能够达到我们的目的。
新问题：对多数类判断错误后的补救行为往往伴随着很高的成本；
所以在现实中，我们往往在寻找捕获少数类的能力和将多数类判错后需要付出的成本的平衡。为了评估这样的能力，我们将引入新的模型评估指标：混淆矩阵和ROC曲线来帮助我们；

8.2.1 混淆矩阵（Confusion Matrix）

混淆矩阵是二分类问题的多维衡量指标体系，在样本不平衡时极其有用。在混淆矩阵中，我们将少数类认为是正例，多数类认为是负例。在决策树，随机森林这些普通的分类算法里，即是说少数类是1，多数类是0。在SVM里，就是说少数类是1，多数类是-1。不过混淆矩阵十分乱，下面用简化方式来看看：

混淆矩阵中，永远是真实值在前，预测值在后。其实可以很容易看出，11和00的对角线就是全部预测正确的，01和10的对角线就是全部预测错误的

基于混淆矩阵，我们有六个不同的模型评估指标，所有以11和00为分子的指标(判断正确的类别所占比例)都是越接近1越好，所以以01和10为分子的指标（判断错误的类别所占比例）都是越接近0越好。对于所有的指标，我们用橙色表示分母，用绿色表示分子，则我们有：

(1)模型整体效果：准确率

准确率Accuracy就是所有预测正确的所有样本除以总样本，通常来说越接近1越好；

(2)捕捉少数类的艺术：精确度，召回率和F1 score

• 精确度/查准率Precision：表示所有被我们预测为是少数类的样本中，真正的少数类所占的比例（也就是纵向的精确度）。在支持向量机中，精确度可以被形象地表示为决策边界上方的所有点（所有被预测为1的点）中，红色点（真正为1的少数类点）所占的比例。精确度越高，代表我们捕捉正确的红色点越多，对少数类的预测越精确。精确度越低，则代表我们误伤了过多的多数类。精确度是”将多数类判错后所需付出成本“的衡量。

#所有判断正确并确实为1的样本 / 所有被判断为1的样本
#对于没有class_weight，没有做样本平衡的灰色决策边界来说：
(y[y == clf.predict(X)] == 1).sum()/(clf.predict(X) == 1).sum()
#0.7142857142857143

#对于有class_weight，做了样本平衡的红色决策边界来说：
(y[y == wclf.predict(X)] == 1).sum()/(wclf.predict(X) == 1).sum()
#0.5102040816326531

可以看出，做了样本平衡之后，精确度是下降的。因为很明显，样本平衡之后，有更多的多数类紫色点被我们误伤了

在现实的样本不平衡例子中，当每一次将多数类判断错误的成本非常高昂的时候（比如大众召回车辆的例子），我们会追求高精确度。但如果我们的目标是不计一切代价捕获少数类，那我们并不在意精确度；

• 召回率Recall/敏感度(sensitivity)/真正率/查全率：表示所有真实为1的样本中(11)，被我们预测正确的样本所占的比例(11+10)。在支持向量机中，召回率可以被表示为，决策边界上方的所有红色点占全部样本中的红色点的比例。召回率越高，代表我们尽量捕捉出了越多的少数类，召回率越低，代表我们没有捕捉出足够的少数类：

#所有predict为1的点 / 全部为1的点的比例
#对于没有class_weight，没有做样本平衡的灰色决策边界来说：
(y[y == clf.predict(X)] == 1).sum()/(y == 1).sum()
#0.6

#对于有class_weight，做了样本平衡的红色决策边界来说：
(y[y == wclf.predict(X)] == 1).sum()/(y == 1).sum()
#1.0

召回率和精确度是此消彼长的，两者之间的平衡代表了捕捉少数类的需求和尽量不要误伤多数类的需求的平衡。究竟要偏向于哪一方，取决于我们的业务需求：究竟是误伤多数类的成本更高，还是无法捕捉少数类的代价更高。

• F1 measure：为了同时兼顾精确度和召回率，我们创造了两者的调和平均数作为考量两者平衡的综合性指标。两个数之间的调和平均倾向于靠近两个数中比较小的那一个数，因此我们追求尽量高的F1 measure，能够保证我们的精确度和召回率都比较高:

越接近1越好，代表精确度和召回率越高；

• 假负率（False Negative Rate）：从Recall延申出来的另一个评估指标，它等于 1 - Recall，用于衡量所有真实为1的样本中，被我们错误判断为0的(图中灰色线下方的红点），通常用得不多：

(3) 判错多数类的考量：特异度与假正率

• 特异度(Specificity)：表示所有真实为0的样本中，被正确预测为0的样本所占的比例，度衡量了一个模型将多数类判断正确的能力。在支持向量机中，可以形象地表示为，决策边界下方的点占所有紫色点的比例：
  

#所有被正确预测为0的样本 / 所有的0样本
#对于没有class_weight，没有做样本平衡的灰色决策边界来说：
(y[y == clf.predict(X)] == 0).sum()/(y == 0).sum()
#0.976

#对于有class_weight，做了样本平衡的红色决策边界来说：
(y[y == wclf.predict(X)] == 0).sum()/(y == 0).sum()
#0.904

• 假正率（False Positive Rate）：即1 - specificity，就是一个模型将多数类判断错误的能力。在支持向量机中，假正率就是决策边界上方的紫色点（所有被判断错误的多数类）占所有紫色点的比例：
  

(4)sklearn中的混淆矩阵

没有特异度和假正率，需要调用混淆矩阵自己计算。精确度和召回率曲线一般都不用；

8.2.2 ROC曲线以及其相关问题

• 假正率的应用：帮助理解每判断正确一个少数类，就有多少个多数类会被判断错误的变化。recall比假正率增加的快，则说明判断出一个少数类的代价不大；反之亦反；、
• ROC曲线：让我们衡量模型在尽量捕捉少数类的时候，误伤多数类的情况如何变化。这是一条以不同阈值下的假正率FPR为横坐标，不同阈值下的召回率Recall为纵坐标的曲线；

(1)概率(probability)与阈值(threshold)

• 似然规定刚刚：当一个样本所对应的这个标签类别下的似然大于0.5的时候，这个样本就被分为这一类；

1. 自建数据集

class_1_ = 7
class_2_ = 4
centers_ = [[0.0, 0.0], [1,1]]
clusters_std = [0.5, 1]
X_, y_ = make_blobs(n_samples=[class_1_, class_2_],
                      centers=centers_,
                      cluster_std=clusters_std,
                      random_state=0, shuffle=False)
plt.scatter(X_[:, 0], X_[:, 1], c=y_, cmap="rainbow",s=30)
plt.show()

2. 建模，调用概率

from sklearn.linear_model import LogisticRegression as LogiR

clf_lo = LogiR().fit(X_,y_)

prob = clf_lo.predict_proba(X_) 
prob

该数组代表了上图中11个样本在类别0/1的概率

#将样本和概率放到一个DataFrame中
import pandas as pd
prob = pd.DataFrame(prob)

prob.columns = ["0","1"]

prob

3. 使用阈值0.5，大于0.5的样本被预测为1，小于0.5的样本被预测为0

#手动调节阈值，来改变我们的模型效果
for i in range(prob.shape[0]):
        if prob.loc[i,"1"] > 0.5:
            prob.loc[i,"pred"] = 1
        else:
            prob.loc[i,"pred"] = 0
            
prob["y_true"] = y_

#根据1的概率排序，概率越大越前面（ascending表示逆序）
prob = prob.sort_values(by="1",ascending=False)

prob

4. 使用混淆矩阵查看结果

from sklearn.metrics import confusion_matrix as CM, precision_score as P, recall_score as R

#第一个参数是真实值，第二个参数是预测值，第三个参数表示前面那个是少数类
CM(prob.loc[:,"y_true"],prob.loc[:,"pred"],labels=[1,0])
#array([[2, 2],
#       [0, 7]], dtype=int64)
#正对角线表示判断正确，副对角线表示判断错误

#试试看手动计算Precision和Recall?

#因为混淆矩阵在真实值和预测值上有区别，所以只能把填真实值/少数类写在前面
P(prob.loc[:,"y_true"],prob.loc[:,"pred"],labels=[1,0])
#1.0

R(prob.loc[:,"y_true"],prob.loc[:,"pred"],labels=[1,0])
#0.5

5. 假如我们使用0.4作为阈值呢？

for i in range(prob.shape[0]):
        if prob.loc[i,"1"] > 0.4:
            prob.loc[i,"pred"] = 1
        else:
            prob.loc[i,"pred"] = 0
            
prob

CM(prob.loc[:,"y_true"],prob.loc[:,"pred"],labels=[1,0])
#array([[2, 2],
#       [1, 6]], dtype=int64)

P(prob.loc[:,"y_true"],prob.loc[:,"pred"],labels=[1,0])
#0.6666666666666666

R(prob.loc[:,"y_true"],prob.loc[:,"pred"],labels=[1,0])
#0.5

#注意，降低或者升高阈值并不一定能够让模型的效果变好，一切都基于我们要追求怎样的模型效果
#通常来说，降低阈值能够升高recall

• 注意：并不是升高阈值，就一定能够增加或者减少Recall，一切要根据数据的实际分布来进行判断；

(2)SVM实现概率预测：重要参数probability，接口predict_proba以及decision_function

到超平面的距离一定程度上反应了样本归属于某个标签类的可能性。接口decision_function返回的值也因此被我们认为是SVM中的置信度（confidence）；

不难看出，最右边的那个点是红色类的概率要大得多，而压在超平面上的点则不一定是红/蓝色，即越靠近超平面的点，其属于那一类的概率就越小

• 重要参数probability：置信度始终不是概率，它没有边界，可以无限大，大部分时候也不是以百分比或者小数的形式呈现，而SVC的判断过程又不像决策树一样可以求解出一个比例，因此较难设置阈值；

#重要参数probability

#使用最初的X和y，样本不均衡的这个模型
class_1 = 500 #类别1有500个样本
class_2 = 50 #类别2只有50个
centers = [[0.0, 0.0], [2.0, 2.0]] #设定两个类别的中心
clusters_std = [1.5, 0.5] #设定两个类别的方差，通常来说，样本量比较大的类别会更加松散
X, y = make_blobs(n_samples=[class_1, class_2],
                      centers=centers,
                      cluster_std=clusters_std,
                      random_state=0, shuffle=False) 

#看看数据集长什么样
plt.scatter(X[:, 0], X[:, 1], c=y, cmap="rainbow",s=10) 
#其中红色点是少数类，紫色点是多数类

clf_proba = svm.SVC(kernel="linear",C=1.0,probability=True).fit(X,y)

#生成的各类标签下的概率
clf_proba.predict_proba(X)

clf_proba.predict_proba(X).shape
#(550, 2)

clf_proba.decision_function(X)

clf_proba.decision_function(X).shape

值得注意的是，在二分类过程中，decision_function只会生成一列距离，样本的类别由距离的符号来判断，但是predict_proba会生成两个类别分别对应的概率。SVM也可以生成概率，所以我们可以使用和逻辑回归同样的方式来在SVM上设定和调节我们的阈值

如果我们的确需要置信度分数，但不一定非要是概率形式的话，那建议可以将probability设置为False，使用decision_function这个接口而不predict_proba；

(3)绘制SVM的ROC曲线

画一个ROC曲线需要概率、阈值以及混淆矩阵，并从混淆矩阵中不断获取FPR以及recall组成横纵坐标，进而画图

#绘制SVM的ROC曲线

#首先来看看如何从混淆矩阵中获取FPR和Recall
cm = CM(prob.loc[:,"y_true"],prob.loc[:,"pred"],labels=[1,0])

cm
#array([[2, 2],
#       [1, 6]], dtype=int64)

#FPR，被预测错误的0占所有真正的0的比例
cm[1,0]/cm[1,:].sum()
#0.14285714285714285

#Recall
cm[0,0]/cm[0,:].sum()
#0.5

#开始绘图
recall = []
FPR = []

probrange = np.linspace(clf_proba.predict_proba(X)
[:,1].min(),clf_proba.predict_proba(X)[:,1].max(),num=50,endpoint=False)

from sklearn.metrics import confusion_matrix as CM, recall_score as R
import matplotlib.pyplot as plot

for i in probrange:
        y_predict = []
        for j in range(X.shape[0]):
            if clf_proba.predict_proba(X)[j,1] > i:
                y_predict.append(1)
            else:
                y_predict.append(0)
        cm = CM(y,y_predict,labels=[1,0])
        recall.append(cm[0,0]/cm[0,:].sum())
        FPR.append(cm[1,0]/cm[1,:].sum())
        
recall.sort()
FPR.sort()

plt.plot(FPR,recall,c="red")
plt.plot(probrange+0.05,probrange+0.05,c="black",linestyle="--")
plt.show()

我们希望随着Recall的不断提升，FPR增加得越慢越好，这说明我们可以尽量高效地捕捉出少数类，而不会将很多地多数类判断错误——>希望的图像：纵坐标急速上升，横坐标缓慢增长(就是上图所示的情况）；
中间的虚线代表着，当recall增加1%，我们的FPR也增加1%，也就是说，我们每捕捉出一个少数类，就会有一个多数类被判错，这种情况下，模型的效果就不好，这种模型捕获少数类的结果，会让许多多数类被误伤，从而增加我们的成本；
ROC曲线通常都是凸型的。对于一条凸型ROC曲线来说，曲线越靠近左上角越好，越往下越糟糕，曲线如果在虚线的下方，则证明模型完全无法使用。但是它也有可能是一条凹形的ROC曲线。对于一条凹型ROC曲线来说，应该越靠近右下角越好，凹形曲线代表模型的预测结果与真实情况完全相反，但是可以人为逆转（调整混淆矩阵的label）来调整；
最糟糕的就是，无论曲线是凹形还是凸型，曲线位于图像中间，和虚线非常靠近，那我们拿它无能为力；

• AUC面积：代表了ROC曲线下方的面积，这个面积越大，代表ROC曲线越接近左上角，模型就越好；

(4)sklearn中的ROC曲线和AUC面积

• sklearn.metrics.roc_curve：帮助计算ROC曲线的横坐标假正率FPR，纵坐标Recall和对应的阈值的类；
• sklearn.metrics.roc_auc_score：帮助我们计算AUC面积的类；

sklearn.metrics.roc_curve (y_true, y_score, pos_label=None, sample_weight=None, drop_intermediate=True)

y_true: 数组，形状 = [n_samples]，真实标签，一定是第一个输的；
y_score: 数组，形状 = [n_samples]，置信度分数，可以是正类样本的概率值，或置信度分数，或者decision_function返回的距离；
pos_label : 整数或者字符串, 默认None，表示被认为是正类样本的类别；
sample_weight: 形如 [n_samples]的类数组结构，可不填，表示样本的权重；
drop_intermediate: 布尔值，默认True，如果设置为True，表示会舍弃一些ROC曲线上不显示的阈值点，这对于计算一个比较轻量的ROC曲线来说非常有用；

这个类以此返回：FPR，Recall以及阈值。

from sklearn.metrics import roc_curve

FPR, recall, thresholds = roc_curve(y,clf_proba.decision_function(X), pos_label=1)

FPR

recall

thresholds 
#此时的threshold就不是一个概率值，而是距离值中的阈值了，所以它可以大于1，也可以为负

sklearn.metrics.roc_auc_score (y_true, y_score, average=’macro’, sample_weight=None, max_fpr=None)

AUC面积的分数使用以上类来进行计算，输入的参数也比较简单，就是真实标签，和与roc_curve中一致的置信度分数或者概率值

from sklearn.metrics import roc_auc_score as AUC

area = AUC(y,clf_proba.decision_function(X))
area
#0.9696400000000001

#开始画图

plt.figure()
plt.plot(FPR, recall, color='red',
             label='ROC curve (area = %0.2f)' % area)
plt.plot([0, 1], [0, 1], color='black', linestyle='--')
#若是从0到1，那么图就会紧贴着坐标轴
plt.xlim([-0.05, 1.05])
plt.ylim([-0.05, 1.05])
plt.xlabel('False Positive Rate')
plt.ylabel('Recall')
plt.title('Receiver operating characteristic example')
plt.legend(loc="lower right")
plt.show()

(5)利用ROC曲线找出最佳阈值

我们的希望是，模型在捕获少数类的能力变强的时候，尽量不误伤多数类，也就是说，随着recall的变大，FPR的大小越小越好；

• 约登指数：最优点，其实是Recall和FPR差距最大的点；

#利用ROC曲线找出最佳阈值

maxindex = (recall - FPR).tolist().index(max(recall - FPR))
maxindex
#43

#decision_founction生成的置信度来说，最佳阈值是
thresholds[maxindex] 
#-1.0860191749391461


#我们可以在图像上来看看这个点在哪里
plt.scatter(FPR[maxindex],recall[maxindex],c="black",s=30) 

#把上述代码放入这段代码中：
plt.figure()
plt.plot(FPR, recall, color='red',
             label='ROC curve (area = %0.2f)' % area)
plt.plot([0, 1], [0, 1], color='black', linestyle='--')
plt.scatter(FPR[maxindex],recall[maxindex],c="black",s=30)
plt.xlim([-0.05, 1.05])
plt.ylim([-0.05, 1.05])
plt.xlabel('False Positive Rate')
plt.ylabel('Recall')
plt.title('Receiver operating characteristic example')
plt.legend(loc="lower right")
plt.show()

最佳阈值就这样选取出来了，由于现在我们是使用decision_function来画ROC曲线，所以我们选择出来的最佳阈值其实是最佳距离。

8.3 使用SVC时的其他考虑

• decision_function_shape：决定我们究竟使用哪一种分类模式；
  因为支持向量机是天生二分类的模型，所以支持向量机在处理多分类问题的时候，是把多分类问题转换成了二分类问题来解决；

8.4 SVC真实数据案例：预测明天是否会下雨

• 核心目的：通过巧妙的预处理和特征工程来向大家展示，在现实数据集上我们往往如何做数据预处理，或者我们都有哪些预处理的方式和思路；

(1)导库导数据，探索特征

导入需要的库

import pandas as pd
import numpy as np
from sklearn.model_selection import train_test_split

导入数据，探索数据

weather = pd.read_csv(r"E:\jupyter\TCR\8\weatherAUS5000.csv",index_col=0)
weather.head()

各个变量的含义：

#将特征矩阵和标签Y分开
X = weather.iloc[:,:-1] 
Y = weather.iloc[:,-1] 

X.shape
#(5000, 21)

#探索数据类型
X.info()

#探索缺失值
#缺失值所占总值的比例 isnull().sum()/X.shape[0]
X.isnull().mean()

粗略观察可以发现，这个特征矩阵由一部分分类变量和一部分连续变量组成，其中 云层遮蔽程度虽然是以数字表示，但是本质却是分类变量。大多数特征都是采集的自然数据

#探索标签的分类
np.unique(Y)
#array(['No', 'Yes'], dtype=object)
#说明标签是二分类

(2)分集，优先探索标签

分训练集和测试集，并做描述性统计：
现实中需要先分训练/测试集，然后再进行预处理，因为测试集在现实中往往是不可得的。但是这样会使得工作量翻倍：

#分训练/测试集
Xtrain, Xtest, Ytrain, Ytest = train_test_splt(X, Y, test_size=0.3, random_state=420)

#恢复索引
for i in [Xtrain, Xtest, Ytrain, Ytest]:
    i.index = range(i.shape[0])

#查看是否有样本不均衡问题
Ytrain.value_counts()
#No     2704
#Yes     796
#Name: RainTomorrow, dtype: int64

Ytest.value_counts()
#No     1157
#Yes     343
#Name: RainTomorrow, dtype: int64

#有轻微的样本不均衡

#查看样本中两个类别的比例
Ytrain.value_counts()[0]/Ytrain.value_counts()[1]
#3.3969849246231156

#将标签编码
from sklearn.preprocessing import LabelEncoder #标签专用编码类，允许输入一维数据
encorder = LabelEncoder().fit(Ytrain)
#使用训练集进行训练，然后再训练集和测试集上分别进行transform
#相当于让类根据Ytrain的样子来构建模型
Ytrain = pd.DataFrame(encorder.transform(Ytrain))
Ytest = pd.DataFrame(encorder.transform(Ytest))

#如果测试集中出现了训练集中未出现的标签类别，就会报错
#这种时候就要重新建模

(3)探索特征，开始处理特征矩阵

1. 描述性统计与异常值

#用于观察数据是否有偏态/异常等，了解数据

#描述性统计
#即使什么都不写，也会默认自动返回25，50，75的值
#可以大概知道是否有偏态/异常，比如99%分布是10，但是最大值却是1000，那么就是异常/严重偏态
#转置T是为了方便观看
Xtrain.describe([0.01,0.05,0.1,0.25,0.5,0.75,0.9,0.99]).T

很明显可以看出均值的量纲不统一，不过方差的差异就没那么大。然后要再看看各种值可否为负，以及各种数据的合理性；
可以看到，降雨量有着明显的偏态（99%才41，而最大值却有115）

Xtest.describe([0.01,0.05,0.1,0.25,0.5,0.75,0.9,0.99])

对于去kaggle下载数据/使用完整数据的容易出现异常值
处理思路：要先观察异常值出现频率，如果只出现一次，多半是输入错误，直接删除即可
如果出现了多次，那么就是去和业务人员沟通，人为造成的异常值是没有用的，只会影响模型
如果异常值占有总体数据10%及以上，就要把异常值换成非异常但是无干扰的项，比如用0进行替换/当缺失值处理

#先查看原始的数据结构
Xtrain.shape
#(3500, 21)
Xtest.shape
#(1500, 21)

2.处理困难特征：日期

• 特征处理顺序：先处理比较难的特征（算法难理解的那种，比如日期，时间等），将其转化为带有有效信息的标准；然后进行填补缺失值；最后进行编码；
• 日期难处理的原因：(1)较为难处理；(2)其所代表的含义是SVM无法理解的，可能得用其他算法；(3)得思考日期是否影响标签并转化为算法能理解的影响

Xtrainc = Xtrain.copy()

Xtrainc.sort_values(by="Location")

第一列日期的类型为字符串，那么就得把它转为数字来给模型学习

#现有的日期特征，是连续性还是分类型
#因为日期没得取小数，所以它是一年分了365天的分类型变量
#要先判断日期特征中是否有重复
#若不重复就直接编码即可，但是有的话就得考虑其他情况

Xtrain.iloc[:,0].value_counts()

不难看出，有很多日期后面的数字不是1，说明他们出现了不止一次，即有重复

#首先，日期不是独一无二的，日期有重复
#其次，在我们分训练集和测试集之后，日期也不是连续的，而是分散的
#某一年的某一天倾向于会下雨？或者倾向于不会下雨吗？
#不是日期影响了下雨与否，反而更多的是这一天的日照时间，湿度，温度等等这些因素影响了是否会下雨
#光看日期，其实感觉它对我们的判断并无直接影响
#如果我们把它当作连续型变量处理，那算法会人为它是一系列1~3000左右的数字，不会意识到这是日期

Xtrain.iloc[:,0].value_counts().count()
#2141
#如果我们把它当作分类型变量处理，类别太多，有2141类，如果换成数值型，会被直接当成连续型变量，如果做成哑变量，我们特征的维度会爆炸

如果我们的思考简单一些，我们可以直接删除日期这个特征。首先它不是一个直接影响我们标签的特征，并且要处理日期其实是非常困难的

#日期没啥用，直接删了吧

Xtrain = Xtrain.drop(["Date"],axis=1)
Xtest = Xtest.drop(["Date"],axis=1)

一般的算法是不会考虑样本之间（行间）的联系的，只会考虑特征之间（列间）的关系。若要考虑样本间的影响，就要用时间序列分析，它的主要目的是 根据已有的历史数据对未来进行预测。但是，数据中的时间是有重复的，所以这个好像不太行。
我们可以将时间对气候的连续影响，转换为”今天是否下雨“这个特征，巧妙地将样本对应标签之间的联系，转换成是特征与标签之间的联系了。

#将时间对气候的连续影响，转换为”今天是否下雨“这个特征

Xtrain["Rainfall"].head(20)

Xtrain.loc[Xtrain["Rainfall"] >= 1,"RainToday"] = "Yes"
Xtrain.loc[Xtrain["Rainfall"] < 1,"RainToday"] = "No"
Xtrain.loc[Xtrain["Rainfall"] == np.nan,"RainToday"] = np.nan

Xtest.loc[Xtrain["Rainfall"] >= 1,"RainToday"] = "Yes"
Xtest.loc[Xtrain["Rainfall"] < 1,"RainToday"] = "No"
Xtest.loc[Xtrain["Rainfall"] == np.nan,"RainToday"] = np.nan

Xtrain.head()

Xtest.head()

这样一来，就创造了一个特征，今天是否下雨“RainToday”。但是日期所代表的季节特征也会影响下雨，我们可以将月份或者季节提取出来，作为一个特征使用，而舍弃掉具体的日期

#我们可以将月份或者季节提取出来，作为一个特征使用，而舍弃掉具体的日期

#提取月份
int(Xtrain.loc[0,"Date"].split("-")[1])

#apply是对DataFrame某一列进行处理的一个函数
#lambda x匿名函数，表示再DataFrame某一列的每一行执行冒号后的命令
Xtrain["Date"] = Xtrain["Date"].apply(lambda x:int(x.split("-")[1]))
#替换完毕后，我们需要修改列的名称
#rename是比较少有的，可以用来修改单个列名的函数
#我们通常都直接使用 df.columns = 某个列表 这样的形式来一次修改所有的列名
#但rename允许我们只修改某个单独的列，columns={现在的名字:新名字}。用index修改行名
Xtrain = Xtrain.rename(columns={"Date":"Month"})

Xtrain.head()

Xtest["Date"] = Xtest["Date"].apply(lambda x:int(x.split("-")[1]))
Xtest = Xtest.rename(columns={"Date":"Month"})

Xtest.head()

3.处理困难特征：地点
地点的难难在无法让模型理解其气候特点。我们需要让算法意识到，不同的地点因为气候
不同，所以对“明天是否会下雨”有着不同的影响。如果我们能够将地点转换为这个地方的气候的话，我们就可以将不同城市打包到同一个气候中，而同一个气候下反应的降雨情况应该是相似的

#超过25个类别的分类型变量，都会被算法当成是连续性变量
Xtrain.loc[:,"Location"].value_counts().count()
#49

cityll = pd.read_csv(r"E:\jupyter\TCR\8\cityll.csv",index_col=0)
city_climate = pd.read_csv(r"E:\jupyter\TCR\8\Cityclimate.csv")

#每个城市对应的经纬度
cityll.head()
#每个城市对应的气候
city_climate.head()

#处理上面的两张表

#去掉度数符号
cityll["Latitudenum"] = cityll["Latitude"].apply(lambda x:float(x[:-1]))
cityll["Longitudenum"] = cityll["Longitude"].apply(lambda x:float(x[:-1]))

#观察一下所有的经纬度方向都是一致的，全部是南纬，东经，因为澳大利亚在南半球，东半球
#所以经纬度的方向我们可以舍弃了
citylld = cityll.iloc[:,[0,5,6]]

#将city_climate中的气候添加到我们的citylld中
citylld["climate"] = city_climate.iloc[:,-1]

citylld.head()

接下来，我们如果想要计算距离，我们就会需要所有样本数据中的城市

samplecity = pd.read_csv(r"E:\jupyter\TCR\8\samplecity.csv")

#我们对samplecity也执行同样的处理：去掉经纬度中度数的符号，并且舍弃我们的经纬度的方向

samplecity["Latitudenum"] = samplecity["Latitude"].apply(lambda x:float(x[:-1]))
samplecity["Longitudenum"] = samplecity["Longitude"].apply(lambda x:float(x[:-1]))

samplecityd = samplecity.iloc[:,[0,5,6]]

samplecityd.head()

接下来就是算距离了，公式如下所示：

R是地球半径，slat/slon是起始地的纬/经度，elat/elon是结束地点的维/经度

#算两地间的距离

#首先用radians将角度转化为弧度,且每次只能转换一个数字
from math import radians, sin, cos, acos 
citylld.loc[:,"slat"] = citylld.iloc[:,1].apply(lambda x : radians(x))
citylld.loc[:,"slon"] = citylld.iloc[:,2].apply(lambda x : radians(x))
samplecityd.loc[:,"elat"] = samplecityd.iloc[:,1].apply(lambda x : radians(x))
samplecityd.loc[:,"elon"] = samplecityd.iloc[:,2].apply(lambda x : radians(x))

import sys
for i in range(samplecityd.shape[0]):
    slat = citylld.loc[:,"slat"]
    slon = citylld.loc[:,"slon"]
    elat = samplecityd.loc[:,"elat"]
    elon = samplecityd.loc[:,"elon"]
    dist = 6371.01 * np.arccos(np.sin(slat)*np.sin(elat) + 
                              np.cos(slat)*np.cos(elat)*np,cos(slon.values - elon))
    city_index = np.argsort(dist)[0]
    #每次计算后，取距离最近的城市，然后将最近的城市和城市对应的气候都匹配到samplecityd中
    samplecityd.loc[i,"closest_city"] = citylld.loc[city_index,"City"]
    samplecityd.loc[i,"climate"] = citylld.loc[city_index,"climate"]
    
#查看最后的结果，需要检查城市匹配是否基本正确
samplecityd.head()

#查看气候的分布
samplecityd["climate"].value_counts()

#确认无误后，取出样本城市所对应的气候，并保存
locafinal = samplecityd.iloc[:,[0,-1]]
locafinal.head()
locafinal.columns = ["Location","Climate"] #在这里设定locafinal的索引为地点，是为了之后进行map的匹配
locafinal = locafinal.set_index(keys="Location")
locafinal.to_csv(r"C:\work\learnbetter\micro-class\week 8 SVM (2)\samplelocation.csv")
locafinal.head()

可以使用map功能，map能够将特征中的值一一对应到我们设定的字典中，并且用字典中的值来替换样本中原本的值

#将location中的内容替换，并且确保匹配进入的气候字符串中不含有逗号，气候两边不含有空格
#我们使用re这个模块来消除逗号
#re.sub(希望替换的值，希望被替换成的值，要操作的字符串)
#x.strip()是去掉空格的函数
import re
Xtrain["Location"] = Xtrain["Location"].map([locafinal.iloc[:,0]]).apply(lambda x : re.sub(",","",x.strip()))
Xtest["Location"] = Xtrain["Location"].map([locafinal.iloc[:,0]]).apply(lambda x : re.sub(",","",x.strip()))

#修改特征内容之后，我们使用新列名“Climate”来替换之前的列名“Location”
#注意这个命令一旦执行之后，就再没有列"Location"了，使用索引时要特别注意
Xtrain = Xtrain.rename(columns={"Location":"Climate"})
Xtest = Xtest.rename(columns={"Location":"Climate"})

Xtrain.head()
Xtest.head()

4.处理分类型变量：缺失值
对于分类型变量用众数填补，而连续型变量则用均值填补。都是用训练集的众数/均值来进行填补，

#4.处理分类型变量：缺失值

#查看缺失值的缺失情况
Xtrain.isnull().mean()

#首先找出，分类型特征都有哪些
cate = Xtrain.columns[Xtrain.dtypes == "object"].tolist()

#除了特征类型为"object"的特征们，还有虽然用数字表示，但是本质为分类型特征的云层遮蔽程度
cloud = ["Cloud9am","Cloud3pm"]
cate = cate + cloud
cate
#['Date',
# 'Climate',
# 'WindGustDir',
# 'WindDir9am',
# 'WindDir3pm',
# 'Cloud9am',
# 'Cloud3pm']

#对于分类型特征，我们使用众数来进行填补
from sklearn.impute import SimpleImputer

si = SimpleImputer(missing_values=np.nan,strategy="most_frequent")
#注意，我们使用训练集数据来训练我们的填补器，本质是在生成训练集中的众数
si.fit(Xtrain.loc[:,cate])

#然后我们用训练集中的众数来同时填补训练集和测试集
Xtrain.loc[:,cate] = si.transform(Xtrain.loc[:,cate])
Xtest.loc[:,cate] = si.transform(Xtest.loc[:,cate])

Xtrain.head()
Xtest.head()

#查看分类型特征是否依然存在缺失值
Xtrain.loc[:,cate].isnull().mean()
Xtest.loc[:,cate].isnull().mean()

5.处理分类型变量：将分类型变量编码
只有分类型变量需要编码。也需要先用训练集数据来训练模型，然后对模型中的训练/测试集进行transform；
可以先不普通编码，先看看直接拿去用效果怎么用，不好的话再来其他编码

#5.处理分类型变量：将分类型变量编码

#将所有的分类型变量编码为数字，一个类别是一个数字
#只允许二维以上数据进行输入的类
from sklearn.preprocessing import OrdinalEncoder
oe = OrdinalEncoder()

#利用训练集进行fit
oe = oe.fit(Xtrain.loc[:,cate])

#用训练集的编码结果来编码训练和测试特征矩阵
#在这里如果测试特征矩阵报错，就说明测试集中出现了训练集中从未见过的类别
Xtrain.loc[:,cate] = oe.transform(Xtrain.loc[:,cate])
Xtest.loc[:,cate] = oe.transform(Xtest.loc[:,cate])

Xtrain.loc[:,cate].head()
Xtest.loc[:,cate].head()

一定要先填缺失值再来编码

6.处理连续型变量：填补缺失值

#6.处理连续型变量：填补缺失值

#要处理连续性变量，就要先去掉分类型变量
col = Xtrain.columns.tolist()

for i in cate:
    col.remove(i)
    
col

#实例化模型，填补策略为"mean"表示均值
impmean = SimpleImputer(missing_values=np.nan,strategy = "mean") #用训练集来fit模型
impmean = impmean.fit(Xtrain.loc[:,col])
#分别在训练集和测试集上进行均值填补
Xtrain.loc[:,col] = impmean.transform(Xtrain.loc[:,col])
Xtest.loc[:,col] = impmean.transform(Xtest.loc[:,col])

Xtrain.head()
Xtest.head()

7.处理连续型变量：无量纲化
数据的无量纲化是SVM执行前的重要步骤，因此我们需要对数据进行无量纲化。但注意，这个操作我们不对分类型变量进行

#7.处理连续型变量：无量纲化
#也不对分类型变量进行

col.remove("Month")
col

from sklearn.preprocessing import StandardScaler
ss = StandardScaler().fit(Xtrain.loc[:,col])
Xtrain.loc[:,col] = ss.transform(Xtrain.loc[:,col])
Xtest.loc[:,col] = ss.transform(Xtest.loc[:,col])

Xtrain.head()
Xtest.head()

至此，特征工程做完了

(4)建模与模型评估

#建模与模型评估

from time import time
import datetime
from sklearn.svm import SVC
from sklearn.model_selection import cross_val_score
from sklearn.metrics import roc_auc_score, recall_score

Ytrain = Ytrain.iloc[:,0].ravel()
Ytest = Ytrain.iloc[:,0].ravel()

#建模选择自然是我们的支持向量机SVC，首先用核函数的学习曲线来选择核函数
#我们希望同时观察，精确性，recall以及AUC分数

times = time() #因为SVM是计算量很大的模型，所以我们需要时刻监控我们的模型运行时间

for kernel in ["linear","poly","rbf","sigmoid"]:
    clf = SVC(kernel=kernel
             ,gamma="auto"
             ,degree=1
             ,cache_size=5000#代表内存使用量
             ).fit(Xtrain,Ytrain)
    result = clf.predict(Xtest)
    score = clf.score(Xtest,Ytest)#返回准确度
    recall = recall_score(Ytest,result)
    #第一个参数是真实值，第二个参数是置信度
    auc = roc_auc_score(Ytest,clf.decision_function(Xtest))
    print("%s 's testing accuracy %f, recall is %f', auc is %f" % (kernel,score,recall,auc))
    print(datetime.datetime.fromtimestamp(time()-times).strftime("%M:%S:%f"))

归一化后的数据用rbf/sigmoid仍然表现不太好说明是线性的

模型的准确度和auc面积还是勉勉强强，但是每个核函数下的recall都不太高。此时有三种目标：
一，我希望不计一切代价判断出少数类，得到最高的recall；
二，我们希望追求最高的预测准确率，一切目的都是为了让accuracy更高，我们不在意recall或者AUC；
三，我们希望达到recall，ROC和accuracy之间的平衡，不追求任何一个也不牺牲任何一个；

(4)模型调参

1.追求最高Recall

#第一种方式：调整class_weight
times = time()
for kernel in ["linear","poly","rbf","sigmoid"]:
    clf = SVC(kernel = kernel
              ,gamma="auto"
              ,degree = 1
              ,cache_size = 5000
              ,class_weight = "balanced"
             ).fit(Xtrain, Ytrain)
    result = clf.predict(Xtest)
    score = clf.score(Xtest,Ytest)
    recall = recall_score(Ytest, result)
    auc = roc_auc_score(Ytest,clf.decision_function(Xtest))
    print("%s 's testing accuracy %f, recall is %f', auc is %f" % (kernel,score,recall,auc))
    print(datetime.datetime.fromtimestamp(time()-times).strftime("%M:%S:%f"))

在锁定了线性核函数之后，我甚至可以将class_weight调节得更加倾向于少数类，来不计代价提升recall；

times = time()
for kernel in ["linear","poly","rbf","sigmoid"]:
    clf = SVC(kernel = kernel
              ,gamma="auto"
              ,degree = 1
              ,cache_size = 5000
              ,class_weight = {1:10}#注意，这里写的其实是，类别1：10，隐藏了类别0：1这个比例
             ).fit(Xtrain, Ytrain)
    result = clf.predict(Xtest)
    score = clf.score(Xtest,Ytest)
    recall = recall_score(Ytest, result)
    auc = roc_auc_score(Ytest,clf.decision_function(Xtest))
    print("%s 's testing accuracy %f, recall is %f', auc is %f" % (kernel,score,recall,auc))
    print(datetime.datetime.fromtimestamp(time()-times).strftime("%M:%S:%f"))

随着recall地无节制上升，我们的精确度下降得十分厉害，不过看起来AUC面积却还好，稳定保持在0.86左右。如果此时我们的目的就是追求一个比较高的AUC分数和比较好的recall，那我们的模型此时就算是很不错了;

2.追求最高准确率
此时此刻我们不在意我们的Recall了，那我们首先要观察一下，我们的样本不均衡状况。如果我们的样本非常不均衡，但是此时却有很多多数类被判错的话，那我们可以让模型任性地把所有地样本都判断为0，完全不顾少数类：

valuec = pd.Series(Ytest).value_counts()
valuec
#0    1157
#1     343
#dtype: int64

valuec[0]/valuec.sum()
#0.7713333333333333

初步判断，可以认为我们其实已经将大部分的多数类判断正确了，所以才能够得到现在的正确率。为了证明我们的判断，我们可以使用混淆矩阵来计算我们的特异度，如果特异度非常高，则证明多数类上已经很难被操作了

#查看模型的特异度

from sklearn.metrics import confusion_matrix as CM
clf = SVC(kernel = "linear"
          ,gamma="auto"
          ,cache_size = 5000
         ).fit(Xtrain, Ytrain)
result = clf.predict(Xtest)

cm = CM(Ytest,result,labels=(1,0))
cm
#array([[ 149,  194],
#       [  47, 1110]], dtype=int64)

specificity = cm[1,1]/cm[1,:].sum()

specificity #几乎所有的0都被判断正确了，还有不少1也被判断正确了
#0.9593777009507347

可以看到，特异度非常高，此时此刻如果要求模型将所有的类都判断为0，则已经被判断正确的少数类会被误伤，整体的准确率一定会下降。想让模型多判断些少数类就不太行，因为一旦我们让模型更加倾向于少数类，就会有更多的多数类被判错

可以试试看使用class_weight将模型向少数类的方向稍微调整，已查看我们是否有更多的空间来提升我们的准确率。如果在轻微向少数类方向调整过程中，出现了更高的准确率，则说明模型还没有到极限：

irange = np.linspace(0.01,0.05,10)

for i in irange:
    times = time()
    clf = SVC(kernel = "linear"
              ,gamma="auto"
              ,cache_size = 5000
              ,class_weight = {1:1+i}
              ).fit(Xtrain, Ytrain)
    result = clf.predict(Xtest)
    score = clf.score(Xtest,Ytest)
    recall = recall_score(Ytest, result)
    auc = roc_auc_score(Ytest,clf.decision_function(Xtest))
    print("under ratio 1:%f testing accuracy %f, recall is %f', auc is %f" % (1+i,score,recall,auc))
    print(datetime.datetime.fromtimestamp(time()-times).strftime("%M:%S:%f"))

我们的最高准确度是84.53%，超过了我们之前什么都不做的时候得到的84.40%。可见，模型还是有潜力的。

我们可以继续细化我们的学习曲线来进行调整

irange = np.linspace(0.018889,0.027778,10)

for i in irange:
    times = time()
    clf = SVC(kernel = "linear"
              ,gamma="auto"
              ,cache_size = 5000
              ,class_weight = {1:1+i}
              ).fit(Xtrain, Ytrain)
    result = clf.predict(Xtest)
    score = clf.score(Xtest,Ytest)
    recall = recall_score(Ytest, result)
    auc = roc_auc_score(Ytest,clf.decision_function(Xtest))
    print("under ratio 1:%f testing accuracy %f, recall is %f', auc is %f" % (1+i,score,recall,auc))
    print(datetime.datetime.fromtimestamp(time()-times).strftime("%M:%S:%f"))

模型的效果没有太好，并没有再出现比我们的84.53%精确度更高的取值。想继续优化只能换模型了：

from sklearn.linear_model import LogisticRegression as LR

logclf = LR(solver="liblinear").fit(Xtrain, Ytrain)
logclf.score(Xtest,Ytest)

C_range = np.linspace(3,5,10)

for C in C_range:
    logclf = LR(solver="liblinear",C=C).fit(Xtrain, Ytrain)
    print(C,logclf.score(Xtest,Ytest))

尽管我们实现了非常小的提升，但可以看出来，模型的精确度还是没有能够实现质变。也许，要将模型的精确度提升到90%以上，我们需要集成算法：比如，梯度提升树

3.追求平衡
来试试看调节线性核函数的C值能否有效果：

#3.追求平衡

###======【TIME WARNING：10mins】======###
import matplotlib.pyplot as plt
C_range = np.linspace(0.01,20,20)

recallall = []
aucall = []
scoreall = []
for C in C_range:
    times = time()
    clf = SVC(kernel = "linear",C=C,cache_size = 5000
              ,class_weight = "balanced"
             ).fit(Xtrain, Ytrain)
    result = clf.predict(Xtest)
    score = clf.score(Xtest,Ytest)
    recall = recall_score(Ytest, result)
    auc = roc_auc_score(Ytest,clf.decision_function(Xtest))
    recallall.append(recall)
    aucall.append(auc)
    scoreall.append(score)
    print("under C %f, testing accuracy is %f,recall is %f', auc is %f" % (C,score,recall,auc))
    print(datetime.datetime.fromtimestamp(time()-times).strftime("%M:%S:%f"))
    
print(max(aucall),C_range[aucall.index(max(aucall))])
plt.figure()
plt.plot(C_range,recallall,c="red",label="recall")
plt.plot(C_range,aucall,c="black",label="auc")
plt.plot(C_range,scoreall,c="orange",label="accuracy")
plt.legend()
plt.show()

可以看到，一开始随着C值不断变大，recall也不断增加，但是随后就开始下降；而且每次循环的时间也不断上升

可以观察到几个现象：
1.首先，我们注意到，随着C值逐渐增大，模型的运行速度变得越来越慢。所以正常来说，我们应该设定一个较小的C值范围来进行调整；
2.其次，C很小的时候，模型的各项指标都很低，但当C到1以上之后，模型的表现开始逐渐稳定，在C逐渐变大之后，模型的效果并没有显著地提高。可以认为我们设定的C值范围太大了，然而再继续增大或者缩小C值的范围，AUC面积也只能够在0.86上下进行变化了，调节C值不能够让模型的任何指标实现质变；

我们把目前为止最佳的C值带入模型，看看我们的准确率，Recall的具体值：

times = time()
clf = SVC(kernel = "linear",C=0.01,cache_size = 5000
          ,class_weight = "balanced"
         ).fit(Xtrain, Ytrain)
result = clf.predict(Xtest)
score = clf.score(Xtest,Ytest)
recall = recall_score(Ytest, result)
auc = roc_auc_score(Ytest,clf.decision_function(Xtest))
print("testing accuracy %f,recall is %f', auc is %f" % (score,recall,auc))
print(datetime.datetime.fromtimestamp(time()-times).strftime("%M:%S:%f"))
#testing accuracy 0.801333,recall is 0.755102', auc is 0.869490
#00:00:386998

这种情况下模型的准确率，Recall和AUC都没有太差，但是也没有太好，这也许就是模型平衡后的一种结果。那么现在就看看ROC曲线，查看我们是否可以通过调整阈值来对这个模型进行改进：

from sklearn.metrics import roc_curve as ROC
import matplotlib.pyplot as plt

#输入真实值，概率值/置信度
FPR, Recall, thresholds = ROC(Ytest,clf.decision_function(Xtest),pos_label=1)

area = roc_auc_score(Ytest,clf.decision_function(Xtest))
#0.8694900605012965

plt.figure()
plt.plot(FPR, Recall, color='red',
         label='ROC curve (area = %0.2f)' % area)
plt.plot([0, 1], [0, 1], color='black', linestyle='--')
plt.xlim([-0.05, 1.05])
plt.ylim([-0.05, 1.05])
plt.xlabel('False Positive Rate')
plt.ylabel('Recall')
plt.title('Receiver operating characteristic example')
plt.legend(loc="lower right")
plt.show()

#以此模型作为基础，我们来求解最佳阈值
maxindex = (Recall - FPR).tolist().index(max(Recall - FPR))
thresholds[maxindex]
#-0.15553836465583348

要找的是：捕捉少数类的能力recall与多数类判错FPR之间差异最大的点，找到其所对应的索引

基于我们选出的最佳阈值，我们来认为确定y_predict，并确定在这个阈值下的recall和准确度的值

from sklearn.metrics import accuracy_score as AC

times = time()
clf = SVC(kernel = "linear",C=3.1663157894736838,cache_size = 5000
          ,class_weight = "balanced"
         ).fit(Xtrain, Ytrain)

#置信度的dataframe，用于根据阈值来进行分类
prob = pd.DataFrame(clf.decision_function(Xtest))

prob.loc[prob.iloc[:,0] >= thresholds[maxindex],"y_pred"]=1
prob.loc[prob.iloc[:,0] < thresholds[maxindex],"y_pred"]=0

prob.loc[:,"y_pred"].isnull().sum()
#0

#检查模型本身的准确度
score = AC(Ytest,prob.loc[:,"y_pred"].values)
recall = recall_score(Ytest, prob.loc[:,"y_pred"])
print("testing accuracy %f,recall is %f" % (score,recall))
print(datetime.datetime.fromtimestamp(time()-times).strftime("%M:%S:%f"))
#testing accuracy 0.783333,recall is 0.810496
#00:00:002000

反而还不如我们不调整时的效果好。可见，如果我们追求平衡，那SVC本身的结果就已经非常接近最优结果了

8.6 总结

我们了解了SVC类的各种重要参数，属性和接口，其中参数包括软间隔的惩罚系数C，核函数kernel，核函数的相关参数gamma，coef0和degree，解决样本不均衡的参数class_weight，解决多分类问题的参数decision_function_shape，控制概率的参数probability，控制计算内存的参数cache_size，属性主要包括调用支持向量的属性support_vectors_和查看特征重要性的属性coef_。
接口中，我们学习了最核心的decision_function。除此之外，我们介绍了分类模型的模型评估指标：混淆矩阵和ROC曲线，还介绍了部分特征工程和数据预处理的思路