I started to improve the solution by transforming the thresholds = np.arange(0,1,0.1)
into a smarter, dichotomous way of finding the maximum
Then I realized, after 2 hours of work, that getting all the accuracies were far more cheaper than just finding the maximum !! (Yes it is totally counter-intuitive).
I wrote a lot of comments here below to explain my code. Feel free to delete all these to make the code more readable.
import numpy as np
# Definition : we predict True if y_score > threshold
def ROC_curve_data(y_true, y_score):
y_true = np.asarray(y_true, dtype=np.bool_)
y_score = np.asarray(y_score, dtype=np.float_)
assert(y_score.size == y_true.size)
order = np.argsort(y_score) # Just ordering stuffs
y_true = y_true[order]
# The thresholds to consider are just the values of score, and 0 (accept everything)
thresholds = np.insert(y_score[order],0,0)
TP = [sum(y_true)] # Number of True Positives (For Threshold = 0 => We accept everything => TP[0] = # of postive in true y)
FP = [sum(~y_true)] # Number of True Positives (For Threshold = 0 => We accept everything => TP[0] = # of postive in true y)
TN = [0] # Number of True Negatives (For Threshold = 0 => We accept everything => we don't have negatives !)
FN = [0] # Number of True Negatives (For Threshold = 0 => We accept everything => we don't have negatives !)
for i in range(1, thresholds.size) : # "-1" because the last threshold
# At this step, we stop predicting y_score[i-1] as True, but as False.... what y_true value say about it ?
# if y_true was True, that step was a mistake !
TP.append(TP[-1] - int(y_true[i-1]))
FN.append(FN[-1] + int(y_true[i-1]))
# if y_true was False, that step was good !
FP.append(FP[-1] - int(~y_true[i-1]))
TN.append(TN[-1] + int(~y_true[i-1]))
TP = np.asarray(TP, dtype=np.int_)
FP = np.asarray(FP, dtype=np.int_)
TN = np.asarray(TN, dtype=np.int_)
FN = np.asarray(FN, dtype=np.int_)
accuracy = (TP + TN) / (TP + FP + TN + FN)
sensitivity = TP / (TP + FN)
specificity = TN / (FP + TN)
return((thresholds, TP, FP, TN, FN))
The all process is just a single loop, and the algorithm is just trivial.
In fact, the stupidly simple function is 10 times faster than the solution proposed before me (commpute the accuracies for thresholds = np.arange(0,1,0.1)
) and 30 times faster than my previous smart-ass-dychotomous-algorithm...
You can then easily compute ANY KPI you want, for example :
def max_accuracy(thresholds, TP, FP, TN, FN) :
accuracy = (TP + TN) / (TP + FP + TN + FN)
return(max(accuracy))
def max_min_sensitivity_specificity(thresholds, TP, FP, TN, FN) :
sensitivity = TP / (TP + FN)
specificity = TN / (FP + TN)
return(max(np.minimum(sensitivity, specificity)))
If you want to test it :
y_score = np.random.uniform(size = 100)
y_true = [np.random.binomial(1, p) for p in y_score]
data = ROC_curve_data(y_true, y_score)
%matplotlib inline # Because I personnaly use Jupyter, you can remove it otherwise
import matplotlib.pyplot as plt
plt.step(data[0], data[1])
plt.step(data[0], data[2])
plt.step(data[0], data[3])
plt.step(data[0], data[4])
plt.show()
print("Max accuracy is", max_accuracy(*data))
print("Max of Min(Sensitivity, Specificity) is", max_min_sensitivity_specificity(*data))
Enjoy ;)
accuracy = np.array(accuracy)
should have been altered toaccuracy = np.array(accuracies)
or similar:)