混淆矩阵和精确性可以帮助我们了解概率类模型的分类结果。然而,我们选择概率类模型进行分类,大多数时候都不是为了单单追求效果,而是希望看到预测的相关概率。这种概率给出预测的可信度,所以对于概率类模型,我们希望能够由其他的模型评估指标来帮助我们判断,模型在"概率预测"
这项工作上,完成得如何。本文介绍概率类模型
独有的评估指标。本文字数8216,建议收藏。
概率预测的准确程度被称为"校准程度",是衡量算法预测出的概率和真实结果的差异的一种方式。一种比较常用的指标叫做布里尔分数,它被计算为是概率预测相对于测试样本的均方误差,表示为:
其中
是样本数量,
为概率类模型预测出的概率,
是样本所对应的真实结果,只能取到0或者1,如果事件发生则为1,如果不发生则为0。
这个指标衡量了概率距离真实标签结果的差异。布里尔分数的范围是从0到1,分数越高则预测结果越差劲,校准程度越差,因此布里尔分数越接近0越好。
sklearn.metrics.brier_score_loss(y_true,
y_prob,
*,
sample_weight=None,
pos_label=None)
y_true : array, shape (n_samples,) 真实标签 y_prob : array, shape (n_samples,) 预测出的概率值 sample_weight : array-like of shape = [n_samples], optional 样本权重 pos_label : int or str, default=None 默认情况下如果label={0,1}或者label={-1,1}两种情况,参数pos_label可以省略,使用默认值。若label={1,2},不是标准型,则pos_label必须等于2。一般会认为较小值为负样本label,较大值为正样本label
布里尔分数可以用于任何可以使用predict_proba
接口调用概率的模型。
import numpy as np
import matplotlib.pyplot as plt
from sklearn.naive_bayes import GaussianNB
from sklearn.datasets import load_digits
from sklearn.model_selection import train_test_split
from sklearn.metrics import brier_score_loss
from sklearn.svm import SVC
from sklearn.linear_model import LogisticRegression as LR
digits = load_digits()
X, y = digits.data, digits.target
Xtrain,Xtest,Ytrain,Ytest = train_test_split(X,y,test_size=0.3,random_state=420)
gnb = GaussianNB().fit(Xtrain,Ytrain)
# 注意,第一个参数是真实标签,第二个参数是预测出的概率值
# 在二分类情况下,接口predict_proba会返回两列,
# 但SVC的接口decision_function却只会返回一列
# 要随时注意,使用了怎样的概率分类器,以辨别查找置信度的接口,以及这些接口的结构
brier_score_loss(Ytest, gnb.predict_proba(Xtest)[:,1], pos_label=1)
>>> 0.032619662406118764
lr = LR(C=1., solver='lbfgs',max_iter=3000,multi_class="auto").fit(Xtrain,Ytrain)
svc = SVC(kernel = "linear",gamma=1).fit(Xtrain,Ytrain)
brier_score_loss(Ytest,lr.predict_proba(Xtest)[:,1],pos_label=1)
>>> 0.011386864786268255
#由于SVC的置信度并不是概率,为了可比性,
# 需要将SVC的置信度“距离”归一化,压缩到[0,1]之间
svc_prob = (svc.decision_function(Xtest) -
svc.decision_function(Xtest).min())/(svc.decision_function(Xtest).max() -
svc.decision_function(Xtest).min())
brier_score_loss(Ytest,svc_prob[:,1],pos_label=1)
>>> 0.23818950248917947
将每个分类器每个标签类别下的布里尔分数可视化:
import pandas as pd
name = ["Bayes","Logistic","SVC"]
color = ["red","black","orange"]
df = pd.DataFrame(index=range(10),columns=name)
for i in range(10):
df.loc[i,name[0]] = brier_score_loss(Ytest,gnb.predict_proba(Xtest)[:,i],pos_label=i)
df.loc[i,name[1]] = brier_score_loss(Ytest,lr.predict_proba(Xtest)[:,i],pos_label=i)
df.loc[i,name[2]] = brier_score_loss(Ytest,svc_prob[:,i],pos_label=i)
plt.figure(figsize=(16,6))
for i in range(df.shape[1]):
plt.plot(range(10),df.iloc[:,i],c=color[i],label=name[i])
plt.legend(fontsize=20,loc='upper left')
plt.show()
可以观察到,逻辑回归的布里尔分数有着压倒性优势,SVC的效果明显弱于贝叶斯和逻辑回归(SVC是强行利用sigmoid函数来压缩概率,因此SVC产出的概率结果并不那么可靠)。贝叶斯位于逻辑回归和SVC之间,效果也不错,但比起逻辑回归,还是不够精确和稳定。
对于一个给定的概率分类器,在预测概率为条件的情况下,真实概率发生的可能性的负对数。这里可以参考逻辑回归。
由于是损失,因此对数似然函数的取值越小,则证明概率估计越准确,模型越理想。对数损失只能用于评估分类型模型。
对于一个样本,如果样本的真实标签
在{0,1}中取值,并且这个样本在类别1下的概率估计为
,则这个样本所对应的对数损失是:
sklearn.metrics.log_loss(y_true,
y_pred,
*,
eps=1e-15,
normalize=True,
sample_weight=None,
labels=None)
y_true : array-like or label indicator matrix 样本的真实标签 y_pred : array-like of float, shape = (n_samples, n_classes) or (n_samples,) 预测的概率,注意不要误以为预测标签,要以接口predict_proba来调用 eps : float 精度 normalize : bool, optional (default=True) 如果为true,则返回每个样本的平均损失。否则,返回每个样本损失的总和。 sample_weight : array-like of shape (n_samples,), default=None 样本权重 labels : array-like, optional (default=None) 如果未提供,将从y_true推断标签。如果标签为None且y_pred具有形状(n_samples),则假定标签为二进制,并从y_true推论得出。
>>> from sklearn.metrics import log_loss
# 高斯贝叶斯
>>> log_loss(Ytest,gnb.predict_proba(Xtest))
2.4725653911460683
# 逻辑回归
>>> log_loss(Ytest,lr.predict_proba(Xtest))
0.12738186898609435
# 支持向量机
>>> log_loss(Ytest,svc_prob)
1.625556312147472
用log_loss
得出的结论和布里尔分数得出的结论不一致:当使用布里尔分数作为评判标准的时候,SVC的估计效果是最差的,逻辑回归和贝叶斯的结果相接近。而使用对数似然的时候,虽然依然是逻辑回归最强大,但贝叶斯却没有SVC的效果好。为什么会有这样的不同呢?
因为逻辑回归和SVC都是以最优化为目的来求解模型,然后进行分类的算法。而朴素贝叶斯中,却没有最优化的过程。对数似然函数直接指向模型最优化的方向,甚至就是逻辑回归的损失函数本身,因此在逻辑回归和SVC上表现得更好。
那什么时候使用对数似然,什么时候使用布里尔分数?
在现实应用中,对数似然函数是概率类模型评估的黄金指标,是评估概率类模型的优先选择。
但是它也有一些缺点。
是产出概率的算法有自己的调节方式,就是调节概率的校准程度。校准程度越高,模型对概率的预测越准确,算法在做判断时就越有自信,模型就会更稳定。
可靠性曲线(reliability curve),又叫做概率校准曲线(probability calibration curve),可靠性图(reliability diagrams),这是一条以预测概率为横坐标,真实标签为纵坐标的曲线。
希望预测概率和真实值越接近越好,最好两者相等,因此一个模型/算法的概率校准曲线越靠近对角线越好。校准曲线因此也是模型评估指标之一。和布里尔分数相似,概率校准曲线是对于标签的某一类来说的,因此一类标签就会有一条曲线,或者我们可以使用一个多类标签下的平均来表示一整个模型的概率校准曲线。
sklearn.calibration.calibration_curve(y_true,
y_prob,
*,
normalize=False,
n_bins=5,
strategy='uniform')
y_true : array, shape (n_samples,) 真实标签 y_prob : array, shape (n_samples,) 预测返回的,正类别下的概率值或置信度 normalize : bool, optional, default=False 布尔值,默认False。是否将y_prob中输入的内容归一化到[0, 1]之间,当y_prob并不是真正的概率的时候可以使用。如果该参数为True,则会将y_prob中最小的值归一化为0,最大值归一化为1。 n_bins : int 整数值,表示分箱的个数。如果箱数很大,则需要更多的数据。
返回
trueproba 可靠性曲线的纵坐标,结构为(n_bins, ),是每个箱子中少数类(Y=1)的占比 predproba 可靠性曲线的横坐标,结构为(n_bins, ),是每个箱子中概率的均值
import numpy as np
import matplotlib.pyplot as plt
from sklearn.datasets import make_classification as mc
from sklearn.naive_bayes import GaussianNB
from sklearn.svm import SVC
from sklearn.linear_model import LogisticRegression as LR
from sklearn.metrics import brier_score_loss
from sklearn.model_selection import train_test_split
from sklearn.calibration import calibration_curve
X, y = mc(n_samples=100000,n_features=20 #总共20个特征
,n_classes=2 #标签为2分类
,n_informative=2 #其中两个代表较多信息
,n_redundant=10 #10个都是冗余特征
,random_state=42)
#样本量足够大,因此使用1%的样本作为训练集
Xtrain, Xtest, Ytrain, Ytest = train_test_split(X, y, test_size=0.99, random_state=42)
gnb = GaussianNB()
gnb.fit(Xtrain,Ytrain)
gnb_score = gnb.score(Xtest, Ytest)
y_pred = gnb.predict(Xtest)
prob = gnb.predict_proba(Xtest)[:,1] #我们的预测概率 - 横坐标
#Ytest - 我们的真实标签 - 横坐标
#从类calibiration_curve中获取横坐标和纵坐标
trueproba, predproba = calibration_curve(Ytest, prob ,
n_bins=10 #输入希望分箱的个数
)
fig = plt.figure(figsize=(10,6))
ax1 = plt.subplot()
ax1.plot([0, 1], [0, 1], "k:", label="Perfectly calibrated")
ax1.plot(predproba, trueproba,"s-",label="%s (%1.3f)" % ("Bayes",gnb_score))
ax1.set_ylabel("True probability for class 1", fontsize=20)
ax1.set_xlabel("Mean predcited probability", fontsize=20)
ax1.set_ylim([-0.05, 1.05])
ax1.legend(fontsize=15)
plt.show()
fig, axes = plt.subplots(1,3,figsize=(18,5))
for ind,i in enumerate([3,10,100]):
ax = axes[ind]
ax.plot([0, 1], [0, 1], "k:", label="Perfectly calibrated")
trueproba, predproba = calibration_curve(Ytest, prob,n_bins=i)
ax.plot(predproba, trueproba,"s-",label="n_bins = {}".format(i))
ax1.set_ylabel("True probability for class 1")
ax1.set_xlabel("Mean predcited probability")
ax1.set_ylim([-0.05, 1.05])
ax.legend()
plt.show()
很明显可以看出,
n_bins
越大,箱子越多,概率校准曲线就越精确,但是太过精确的曲线不够平滑,无法和我们希望的完美概率密度曲线相比较。
n_bins
越小,箱子越少,概率校准曲线就越粗糙,虽然靠近完美概率密度曲线,但是无法真实地展现模型概率预测的结果。
因此我们需要取一个既不是太大,也不是太小的箱子个数,让概率校准曲线既不是太精确,也不是太粗糙,而是一条相对平滑,又可以反应出模型对概率预测的趋势的曲线。通常来说,建议先试试看箱子数等于10的情况。箱子的数目越大,所需要的样本量也越多,否则曲线就会太过精确。
X, y = mc(n_samples=100000,n_features=20 #总共20个特征
,n_classes=2 #标签为2分类
,n_informative=2 #其中两个代表较多信息
,n_redundant=10 #10个都是冗余特征
,random_state=42)
#样本量足够大,因此使用1%的样本作为训练集
Xtrain, Xtest, Ytrain, Ytest = train_test_split(X, y,test_size=0.99,random_state=42)
gnb = GaussianNB()
lr = LR(C=1., solver='lbfgs',max_iter=3000,multi_class="auto")
svc = SVC(kernel = "linear",gamma=1)
names = ["GaussianBayes","Logistic","SVC"]
clfs=[gnb,lr,svc]
fig, ax1 = plt.subplots(figsize=(10,5))
ax1.plot([0, 1], [0, 1], "k:", label="Perfectly calibrated")
for clf, name in zip(clfs,names):
print(name)
clf.fit(Xtrain,Ytrain)
y_pred = clf.predict(Xtest)
#hasattr(obj,name):查看一个类obj中是否存在名字为name的接口,存在则返回True
if hasattr(clf, "predict_proba"):
prob = clf.predict_proba(Xtest)[:,1]
else: # use decision function
prob= clf.decision_function(Xtest)
prob= (prob - prob.min()) / (prob.max() - prob.min())
#返回布里尔分数
clf_score = brier_score_loss(Ytest, prob, pos_label=y.max())
trueproba, predproba = calibration_curve(Ytest, prob,n_bins=10)
ax1.plot(predproba, trueproba,"s-",label="%s (%1.3f)" % (name, clf_score))
ax1.set_ylabel("True probability for class 1")
ax1.set_xlabel("Mean predcited probability")
ax1.set_ylim([-0.05, 1.05])
ax1.legend()
ax1.set_title('Calibration plots (reliability curve)')
plt.show()
sigmoid函数
的形状。它的概率校准曲线效果是典型的置信度不足的分类器(under-confident classifier)
的表现:大量的样本点集中在决策边界的附近,它们的置信度靠近0.5左右,即便决策边界能够将样本点判断正确,模型本身对这个结果也不是非常确信的。相对的,离决策边界很远的点的置信度就会很高,因为它很大可能性上不会被判断错误。支持向量机在面对混合度较高的数据的时候,有着天生的置信度不足的缺点。Sigmoid函数
相反的形状,说明数据集中的特征不是相互条件独立的。这与我们设定的数据集相关,其中 的10个冗余特征之间不是完全独立的。可以通过绘制直方图来查看模型的预测概率的分布。直方图是以样本的预测概率分箱后的结果为横坐标,每个箱中的样本数量为纵坐标的一个图像。
fig, ax2 = plt.subplots(figsize=(10,6))
for clf, name in zip(clfs,names):
clf.fit(Xtrain,Ytrain)
y_pred = clf.predict(Xtest)
if hasattr(clf, "predict_proba"):
prob = clf.predict_proba(Xtest)[:,1]
else: # use decision function
prob = clf.decision_function(Xtest)
prob = (prob - prob.min()) / (prob.max() - prob.min())
ax2.hist(prob
,bins=10
,label=name
,histtype="step" #设置直方图为透明
,lw=2 #设置直方图每个柱子描边的粗细
)
ax2.set_ylabel("Distribution of probability", fontsize=20)
ax2.set_xlabel("Mean predicted probability", fontsize=20)
ax2.set_xlim([-0.05, 1.05])
ax2.set_xticks([0,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1])
ax2.legend(loc=9, fontsize=15)
plt.show()
以上的样本都在0和1的附近,可以说是置信度最高的算法,但是贝叶斯的布里尔分数却不如逻辑回归,这证明贝叶斯中在0和1附近的样本中有一部分是被分错的。
等近似回归有两种回归可以使用,一种是基于Platt的Sigmoid模型
的参数校准方法,一种是基于等渗回归(isotonic calibration)的非参数的校准方法。概率校准应该发生在测试集上,必须是模型未曾见过的数据。
在数学上,使用这两种方式来对概率进行校准的原理十分复杂,而此过程我们在sklearn中无法进行干涉,大家不必过于去深究。
class sklearn.calibration.CalibratedClassifierCV (base_estimator=None,
method='sigmoid',
cv='warn')
这是一个带交叉验证的概率校准类,它使用交叉验证生成器,对交叉验证中的每一份数据,它都在训练样本上进行模型参数估计,在测试样本上进行概率校准,然后为我们返回最佳的一组参数估计和校准结果。每一份数据的预测概率会被求解平均。注意,类CalibratedClassifierCV没有接口decision_function,要查看这个类下校准过后的模型生成的概率,必须调用predict_proba接口。
base_estimator=None, 需要校准其输出决策功能的分类器,必须存在"predict_proba"或"decision_function"接口。如果参数"cv = prefit",分类器必须已经拟合数据完毕。 method='sigmoid' 进行概率校准的方法,可输入"sigmoid"或者"isotonic"。
当校准的样本量太少(比如,小于等于1000个测试样本)的时候,不建议使用等渗回归,因为它倾向于过拟合。样本量过少时请使用sigmoids,即Platt校准。 cv='warn' 整数,确定交叉验证的策略。可能输入是:
sklearn.model_selection.StratifiedKFold
进行折数分割。如果y是连续型变量,则使用sklearn.model_selection.KFold
进行分割。# 封装可视化函数
def plot_calib(models,names,Xtrain,Xtest,Ytrain,Ytest,n_bins=10):
import matplotlib.pyplot as plt
from sklearn.metrics import brier_score_loss
from sklearn.calibration import calibration_curve
fig, (ax1, ax2) = plt.subplots(1, 2,figsize=(20,6))
ax1.plot([0, 1], [0, 1], "k:", label="Perfectly calibrated")
for clf, name in zip(models,names):
clf.fit(Xtrain,Ytrain)
y_pred = clf.predict(Xtest)
#hasattr(obj,name):查看一个类obj中是否存在名字为name的接口,存在则返回True
if hasattr(clf, "predict_proba"):
prob = clf.predict_proba(Xtest)[:,1]
else: # use decision function
prob = clf.decision_function(Xtest)
prob = (prob - prob.min()) / (prob.max() - prob.min())
#返回布里尔分数
clf_score = brier_score_loss(Ytest, prob, pos_label=y.max())
trueproba, predproba = calibration_curve(Ytest, prob,n_bins=n_bins)
ax1.plot(predproba, trueproba,"s-",label="%s (%1.3f)" % (name, clf_score))
ax2.hist(prob, range=(0, 1), bins=n_bins, label=name, histtype="step", lw=2)
ax2.set_ylabel("Distribution of probability", fontsize=20)
ax2.set_xlabel("Mean predicted probability", fontsize=20)
ax2.set_xlim([-0.05, 1.05])
ax2.legend(loc=9)
ax2.set_title("Distribution of probablity", fontsize=20)
ax1.set_ylabel("True probability for class 1", fontsize=20)
ax1.set_xlabel("Mean predcited probability", fontsize=20)
ax1.set_ylim([-0.05, 1.05])
ax1.legend(fontsize=15)
ax1.set_title('Calibration plots(reliability curve)', fontsize=20)
plt.show()
# 模型准备
from sklearn.calibration import CalibratedClassifierCV
names = ["GaussianBayes","Logistic","SVC"]
models = [GaussianNB()
,LR(C=1., solver='lbfgs',max_iter=3000,multi_class="auto")
#定义两种校准方式
,SVC(kernel = 'rbf')]
# 数据准备
X, y = mc(n_samples=100000,n_features=20 #总共20个特征
,n_classes=2 #标签为2分类
,n_informative=2 #其中两个代表较多信息
,n_redundant=10 #10个都是冗余特征
,random_state=42)
Xtrain,Xtest,Ytrain,Ytest = train_test_split(X,y,test_size = 0.9,random_state = 0)
# 函数执行
plot_calib(models,names,Xtrain,Xtest,Ytrain,Ytest)
gnb = GaussianNB().fit(Xtrain,Ytrain)
gnb.score(Xtest,Ytest)
>>> 0.8669222222222223
brier_score_loss(Ytest,gnb.predict_proba(Xtest)[:,1],pos_label = 1)
>>> 0.11514876553209971
gnbisotonic = CalibratedClassifierCV(gnb, cv=2, method='isotonic').fit(Xtrain,Ytrain)
gnbisotonic.score(Xtest,Ytest)
>>> 0.8665111111111111
brier_score_loss(Ytest,gnbisotonic.predict_proba(Xtest)[:,1],pos_label = 1)
>>> 0.09576140551554904
可以看出,校准概率后,布里尔分数明显变小了,但整体的准确率却略有下降,这证明算法在校准之后,尽管对概率的预测更准确了,但模型的判断力略有降低。
思考: 布里尔分数衡量模型概率预测的准确率,布里尔分数越低,代表模型的概率越接近真实概率,当进行概率校准后,本来标签是1的样本的概率应该会更接近1,而标签本来是0的样本应该会更接近0,没有理由布里尔分数提升了,模型的判断准确率居然下降了。但从结果来看,模型的准确率和概率预测的正确性并不是完全一致的,为什么会这样呢? 对于不同的概率类模型,原因是不同的。
当然,可能还有更多更深层的原因,比如概率校准过程中的数学细节如何影响了我们的校准,class calibration_curve中是如何分箱,如何通过真实标签和预测值来生成校准曲线使用的横纵坐标的,这些过程中也可能有着让布里尔分数和准确率向两个方向移动的过程。 当两者相悖的时候,请务必以准确率为标准。但是这不代表说布里尔分数和概率校准曲线就无效了。概率类模型几乎没有参数可以调整,除了换模型之外,鲜有更好的方式帮助我们提升模型的表现,概率校准是难得的可以帮助我们针对概率提升模型的方法。
GaussianNB_CCC = [GaussianNB()
,CalibratedClassifierCV(GaussianNB(),method = 'sigmoid',cv=5)
,CalibratedClassifierCV(GaussianNB(),method = 'isotonic',cv=5)
]
LR_CCC = [LR(solver = 'lbfgs',max_iter = 3000,multi_class = 'auto')
,CalibratedClassifierCV(LR(solver = 'lbfgs',max_iter = 3000,multi_class = 'auto'),method = 'sigmoid',cv=5)
,CalibratedClassifierCV(LR(solver = 'lbfgs',max_iter = 3000,multi_class = 'auto'),method = 'isotonic',cv=5)
]
SVC_CCC = [SVC(kernel = 'rbf',gamma='scale')
,CalibratedClassifierCV(SVC(kernel = 'rbf',gamma='scale'),method = 'sigmoid',cv = 5)
,CalibratedClassifierCV(SVC(kernel = 'rbf',gamma='scale'),method = 'isotonic',cv = 5)
]
estimators_CCC = [GaussianNB_CCC,LR_CCC,SVC_CCC]
GaussianNB_names =['GaussianNB','GaussianNB+sigmoid','GaussianNB+isotonic']
LR_names = ['LogisticRegression','LogisticRegression+sigmoid','LogisticRegression+isotonic']
SVC_names = ['SVC','SVC+sigmoid','SVC+isotonic']
names = [GaussianNB_names,LR_names,SVC_names]
for i in range(len(estimators_CCC)):
plot_calib(estimators_CCC[i],names[i],Xtrain,Xtest,Ytrain,Ytest)
从校正结果来看:
在现实中,可以根据自己的需求来调整模型。
概率拟合
,并使用概率校准的方式来调节模型。准确率和Recall
,可以考虑使用天生就非常准确的概率类模型逻辑回归,也可以考虑使用除了概率校准之外还有很多其他参数可调的支持向量机分类器。