我运行了一个逻辑回归模型，并预测了logit值。我用这个来得到ROC曲线上的点：
from sklearn import metrics
fpr, tpr, thresholds = metrics.roc_curve(Y_test,p)
我知道metrics.roc_auc_score可以给出ROC曲线下面积，但有没有命令可以找到最佳的截断点（阈值）？
ROC曲线和截断点。Python
85
- Shiva Prakash
6
5个回答
82
您可以使用R语言中的epi包 来完成这个任务，但是我在Python中没有找到类似的包或示例。
根据“真正率”高且“假正率”低的逻辑，最优切割点应该在此处。基于这个逻辑，我编写了下面的示例以查找最佳阈值。
Python代码：
import pandas as pd
import statsmodels.api as sm
import pylab as pl
import numpy as np
from sklearn.metrics import roc_curve, auc
# read the data in
df = pd.read_csv("http://www.ats.ucla.edu/stat/data/binary.csv")
# rename the 'rank' column because there is also a DataFrame method called 'rank'
df.columns = ["admit", "gre", "gpa", "prestige"]
# dummify rank
dummy_ranks = pd.get_dummies(df['prestige'], prefix='prestige')
# create a clean data frame for the regression
cols_to_keep = ['admit', 'gre', 'gpa']
data = df[cols_to_keep].join(dummy_ranks.iloc[:, 'prestige_2':])
# manually add the intercept
data['intercept'] = 1.0
train_cols = data.columns[1:]
# fit the model
result = sm.Logit(data['admit'], data[train_cols]).fit()
print result.summary()
# Add prediction to dataframe
data['pred'] = result.predict(data[train_cols])
fpr, tpr, thresholds =roc_curve(data['admit'], data['pred'])
roc_auc = auc(fpr, tpr)
print("Area under the ROC curve : %f" % roc_auc)
####################################
# The optimal cut off would be where tpr is high and fpr is low
# tpr - (1-fpr) is zero or near to zero is the optimal cut off point
####################################
i = np.arange(len(tpr)) # index for df
roc = pd.DataFrame({'fpr' : pd.Series(fpr, index=i),'tpr' : pd.Series(tpr, index = i), '1-fpr' : pd.Series(1-fpr, index = i), 'tf' : pd.Series(tpr - (1-fpr), index = i), 'thresholds' : pd.Series(thresholds, index = i)})
roc.iloc[(roc.tf-0).abs().argsort()[:1]]
# Plot tpr vs 1-fpr
fig, ax = pl.subplots()
pl.plot(roc['tpr'])
pl.plot(roc['1-fpr'], color = 'red')
pl.xlabel('1-False Positive Rate')
pl.ylabel('True Positive Rate')
pl.title('Receiver operating characteristic')
ax.set_xticklabels([])
最佳截断点为0.317628，因此任何高于此值的可以标记为1，否则为0。从输出/图表中可以看出，在TPR交叉1-FPR时，TPR为63%，FPR为36%，在当前示例中，TPR-(1-FPR)最接近零。
输出：
最佳截断点为0.317628，因此任何高于此值的可以标记为1，否则为0。从输出/图表中可以看出，在TPR交叉1-FPR时，TPR为63%，FPR为36%，在当前示例中，TPR-(1-FPR)最接近零。
1-fpr fpr tf thresholds tpr
171 0.637363 0.362637 0.000433 0.317628 0.637795
希望这对你有所帮助。
Edit
为了简化和使代码更易重用，我已经编写了一个函数来查找最佳概率截止点。
Python 代码:
def Find_Optimal_Cutoff(target, predicted):
""" Find the optimal probability cutoff point for a classification model related to event rate
Parameters
----------
target : Matrix with dependent or target data, where rows are observations
predicted : Matrix with predicted data, where rows are observations
Returns
-------
list type, with optimal cutoff value
"""
fpr, tpr, threshold = roc_curve(target, predicted)
i = np.arange(len(tpr))
roc = pd.DataFrame({'tf' : pd.Series(tpr-(1-fpr), index=i), 'threshold' : pd.Series(threshold, index=i)})
roc_t = roc.iloc[(roc.tf-0).abs().argsort()[:1]]
return list(roc_t['threshold'])
# Add prediction probability to dataframe
data['pred_proba'] = result.predict(data[train_cols])
# Find optimal probability threshold
threshold = Find_Optimal_Cutoff(data['admit'], data['pred_proba'])
print threshold
# [0.31762762459360921]
# Find prediction to the dataframe applying threshold
data['pred'] = data['pred_proba'].map(lambda x: 1 if x > threshold else 0)
# Print confusion Matrix
from sklearn.metrics import confusion_matrix
confusion_matrix(data['admit'], data['pred'])
# array([[175, 98],
# [ 46, 81]])
- Manohar Swamynathan
7
5@skmathur，我将其制作为一个可重复使用且简化的函数。希望这有所帮助。 - Manohar Swamynathan
5在“Find_Optimal_Cutoff”函数中，您的Youden指数公式存在问题。
roc_curve返回的是假阳性率（1-特异度）即fpr。您正在减去(1-fpr)。您需要将tpr-(1-fpr)更改为tpr-fpr。 - John BonfardeciEpi包在R中选择最大化（特异性+敏感性）的截断值。因此，它应该是tpr +（1-fpr），而不是代码中给出的tpr-（1-fpr）。 - Ravi
4@JohnBonfardeci 这只是我一个人的感觉吗？我觉得OP的解决方案产生了错误的结果..难道不应该是
pd.Series(tpr-fpr, index=thresholds, name='tf').idxmax() 吗？ - Stefan Falk1我不明白为什么这个回复得到了最多的赞。@JohnBonfardeci和@StefanFalk是正确的，公式是错误的。除此之外，图表的轴描述也是错误的，应该有一个双y轴或者tpr和fpr的图例，而x轴实际上是阈值。 - Zwiebak
显示剩余2条评论
70
根据您的问题中给出的tpr、fpr和阈值，最佳阈值的答案只是：
optimal_idx = np.argmax(tpr - fpr)
optimal_threshold = thresholds[optimal_idx]
- cgnorthcutt
6
1如果我得到了负的最优阈值，那该怎么办？我的输出预测值在[0,1]范围内。 - Allan Ruin
1使用建议的方法：
optimal_idx = np.argmax(tpr - fpr)
optimal_threshold = thresholds[optimal_idx]对我来说不起作用。阈值数组包含负值，但我预期的值应该在0到1之间。 - RafaelCaballero1@rafaelcaballero “不起作用”？你描述的一切都听起来是在正常工作。它不应该在0和1之间。这只是一个得分。 - cgnorthcutt
1那我误解了问题。我以为阈值在0和1之间移动，目标是找到在这个范围内最大化tpr-fpr值的数值。 - RafaelCaballero
@cgnorthcutt 您的代码是正确的。但是TPR = TP /（真正例），FPR = FP /（真负例）。TPR + FPR！= 1。 - Qinsi
@Qinsi - 是的。这就是为什么我在上面说“它不应该在0和1之间。它只是一个分数。”这是一个很好的观点，谢谢提醒。 - cgnorthcutt
19
Youden's J-Score的原生Python实现
def cutoff_youdens_j(fpr,tpr,thresholds):
j_scores = tpr-fpr
j_ordered = sorted(zip(j_scores,thresholds))
return j_ordered[-1][1]
- lee
7
另一个可能的解决方案。
我将创建一些随机数据。
import numpy as np
import pandas as pd
import scipy.stats as sps
from sklearn import linear_model
from sklearn.metrics import roc_curve, RocCurveDisplay, auc
from sklearn.model_selection import train_test_split
import matplotlib.pyplot as plt
import seaborn as sns
# define data distributions
N0 = 300
N1 = 250
dist0 = sps.gamma(a=8, scale=1/10)
x0 = np.linspace(dist0.ppf(0), dist0.ppf(1-1e-5), 100)
y0 = dist0.pdf(x0)
dist1 = sps.gamma(a=15, scale=1/10)
x1 = np.linspace(dist1.ppf(0), dist1.ppf(1-1e-5), 100)
y1 = dist1.pdf(x1)
with plt.style.context("bmh"):
plt.plot(x0, y0, label="NEG")
plt.plot(x1, y1, label="POS")
plt.legend()
plt.title("Gamma distributions")
# create a random dataset
rvs0 = dist0.rvs(N0, random_state=0)
rvs1 = dist1.rvs(N1, random_state=1)
with plt.style.context("bmh"):
plt.hist(rvs0, alpha=.5, label="NEG")
plt.hist(rvs1, alpha=.5, label="POS")
plt.legend()
plt.title("Random dataset")
用观测值（x特征和y目标）初始化数据框架。
df = pd.DataFrame({
"y": np.concatenate(( np.repeat(0, N0) , np.repeat(1, N1) )),
"x": np.concatenate(( rvs0 , rvs1 )),
})
并用箱线图显示它
# plot the data
with plt.style.context("bmh"):
g = sns.catplot(
kind="box",
data=df,
x="y", y="x"
)
ax = g.axes.flat[0]
sns.stripplot(
data=df,
x="y", y="x",
ax=ax, color='k',
alpha=.25
)
plt.show()
# split dataset into train-test
X_train, X_test, y_train, y_test = train_test_split(
df[["x"]], df.y.values, test_size=0.5, random_state=1)
# init and fit Logistic Regression on train set
clf = linear_model.LogisticRegression()
clf.fit(X_train, y_train)
# predict probabilities on x test set
y_proba = clf.predict_proba(X_test)
# compute FPR and TPR from y test set and predicted probabilities
fpr, tpr, thresholds = roc_curve(
y_test, y_proba[:,1], drop_intermediate=False)
# compute ROC AUC
roc_auc = auc(fpr, tpr)
# init a dataframe for results
df_test = pd.DataFrame({
"x": X_test.x.values.flatten(),
"y": y_test,
"proba": y_proba[:,1]
})
# sort it by predicted probabilities
# because thresholds[1:] = y_proba[::-1]
df_test.sort_values(by="proba", inplace=True)
# add reversed TPR and FPR
df_test["tpr"] = tpr[1:][::-1]
df_test["fpr"] = fpr[1:][::-1]
# optional: add thresholds to check
#df_test["thresholds"] = thresholds[1:][::-1]
# add Youden's j index
df_test["youden_j"] = df_test.tpr - df_test.fpr
# define the cut_off and diplay it
cut_off = df_test.sort_values(
by="youden_j", ascending=False, ignore_index=True).iloc[0]
print("CUT-OFF:")
print(cut_off)
# plot everything
with plt.style.context("bmh"):
fig, ax = plt.subplots(1, 3, figsize=(15, 5))
RocCurveDisplay(
fpr=df_test.fpr, tpr=df_test.tpr,
roc_auc=roc_auc).plot(ax=ax[0])
ax[0].set_title("ROC curve")
ax[0].axline(xy1=(0,0), slope=1, color="r", ls=":")
ax[0].plot(cut_off.fpr, cut_off.tpr, 'ko', ms=10)
df_test.plot(
x="youden_j", y="proba", ax=ax[1],
ylabel="Predicted Probabilities", xlabel="Youden j",
title="Youden's index", legend=False
)
ax[1].axvline(cut_off.youden_j, color="k", ls="--")
ax[1].axhline(cut_off.proba, color="k", ls="--")
df_test.plot(
x="x", y="proba", ax=ax[2],
ylabel="Predicted Probabilities", xlabel="X Feature",
title="Cut-Off", legend=False
)
ax[2].axvline(cut_off.x, color="k", ls="--")
ax[2].axhline(cut_off.proba, color="k", ls="--")
plt.show()
我们得到
CUT-OFF:
x 1.065712
y 1.000000
proba 0.378543
tpr 0.852713
fpr 0.143836
youden_j 0.708878
# check results
TP = df_test[(df_test.x>=cut_off.x)&(df_test.y==1)].index.size
FP = df_test[(df_test.x>=cut_off.x)&(df_test.y==0)].index.size
TN = df_test[(df_test.x< cut_off.x)&(df_test.y==0)].index.size
FN = df_test[(df_test.x< cut_off.x)&(df_test.y==1)].index.size
print("True Positive Rate: ", TP / (TP + FN))
print("False Positive Rate:", 1 - TN / (TN + FP))
True Positive Rate: 0.8527131782945736
False Positive Rate: 0.14383561643835618
- Max Pierini
2
2这是一个不错的方法，但它并没有考虑到非唯一概率。roc_curve为每个唯一概率返回一个值，但是对于使用唯一y_proba创建的任何df，
df_test["tpr"] = tpr[1:][::-1]会出现错误。 - Alexander Vocaet我建议在创建
df_test后立即删除所有关于列"proba"的重复项。 - Alexander Vocaet5
虽然我来晚了，但您也可以使用几何平均数来确定最优阈值，如此处所述: threshold tuning for imbalance classification
它的计算方法为:
# calculate the g-mean for each threshold
gmeans = sqrt(tpr * (1-fpr))
# locate the index of the largest g-mean
ix = argmax(gmeans)
print('Best Threshold=%f, G-Mean=%.3f' % (thresholds[ix], gmeans[ix]))
- Ramesh Kumar
2
这里似乎不需要使用
sqrt。Argmax在没有它的情况下也能正常工作。 - Roland Pihlakas它给我们带来了什么最优解？这是否会与Youden指数在此示例中产生矛盾：
tpr=0.5＆fpr=0.5和tpr=0.55＆fpr=0.45？ - GabrielChu网页内容由stack overflow 提供, 点击上面的可以查看英文原文，
原文链接
原文链接
thresholds[np.argmax(tpr - fpr)]和最受欢迎的答案函数threshold = Find_Optimal_Cutoff(data['true'], data['pred'])之间的区别吗？这些阈值在实际计算中非常接近，但是略有不同。 - kevin_theinfinityfundthresholds[np.argmin(abs(tpr-(1-fpr)))]。 - Danthresholds[np.argmin((1 - tpr) ** 2 + fpr ** 2)]。但@cgnorthcutt的解决方案最大化了Youden's J统计量，这似乎是更为被接受的方法。对于你的情况来说，真正的“最佳”取决于假阳性和假阴性的相对成本。 - rcauvin