<jupyter_start><jupyter_text>Precision-Recall and Receiver Operating Characteristic CurvesIn this notebook, we will learn about constructing and interpreting precision-recall (PR) and receiver operating characteristic (ROC) curves. These curves are typically used to judge the performances of probabilistic classifiers beyond the confusion matrix and other performance measures that are derived from the confusion matrix. Note that the confusion matrix primarily judges categorical decisions, while PR and ROC curves are assessing probabilistic decisions. If a classifier does indeed make only categorical decisions (e.g. True or False, 1 or 0, yes or no, etc..), it may be sufficient to stick with a confusion matrix (although it is possible to turn categorical decisions to probabilistic ones, see [isotonic regression](https://en.wikipedia.org/wiki/Isotonic_regression)). However, oftentimes, a classifier makes probabilistic predictions (e.g. 80% True and 20% False, 75% 1 and 25% 0, 33% yes and 67% no, etc...) and the tragedy of judging classifier performance starts by choosing an arbitrary threshold (typically, > 50%) to categorize the prediction. For example, if a probabilistic classifier predicts True at 80%, one might say that the prediction is categorically True (since 80% > 50%); if a prediction is True at 49%, one might say that the prediction is categorically False (since 49% < 50%). The 50% cut-off is arbitrary, and one may choose any cut-off or threshold to impact the confusion matrix. With PR and ROC curves, the idea is to vary the cut-off from $[0, 1]$ and measure the performance over all these thresholds. In fact, after varying the cut-off, one may use PR and ROC curves to choose an optimal threshold to balance for trade-offs between precision vs recall or true positive rate (TPR) vs false positive rate (FPR). To understand PR and ROC curves, first understand the confusion matrix. Assume that there are only 2 classes to predict; True and False. Then a 2 x 2 matrix may be created such that the rows represent the predictions and the columns represent the truth. ```+-------+------+-------+|       | True | False |+-------+------+-------+| True  | tp   | fp    |+-------+------+-------+| False | fn   | tn    |+-------+------+-------+```Each element in the confusion matrix stands for something.* tp is the number of true positive (model predicts the observation as true when the observation is really true)* fp is the number of false positive (model predicts the observation as true when the observation is really false)* fn is the number of false negative (model predicts the observation as false when the observation is really true)* tn is the number of true negative (model predicts the observation as false when the observation is really false)From these simple counts, a wonderful list of performance measures may be created (refer to the [Wikipedia site](https://en.wikipedia.org/wiki/Confusion_matrix)). Here, we only focus on precision, recall, TPR and FPR, which are defined as follows.* $\text{precision}=\frac{\text{TP}}{\text{TP} + \text{FP}}$* $\text{recall}=\frac{\text{TP}}{\text{TP} + \text{FN}}$* $\text{TPR}=\frac{\text{TP}}{\text{TP} + \text{FN}}$* $\text{FPR}=\frac{\text{FP}}{\text{FP} + \text{TN}}$Note that recall and TPR are the same thing and also referred to by many other names (e.g. sensitivity).If a classifier is a probabilistic one, we may choose an arbitrary threshold, above which we will categorically say the prediction is True, to induce a confusion matrix and hence precision, recall, TPR and FPR. Changing the threshold will change the confusion matrix and hence, precision, recall, TPR and FPR (typically). We can `cheat` and choose one threshold to produce a set of desired precision, recall, TPR and FPR and report out that our classifier performs really well. However, we know better that we should look at these performance measures over all possible thresholds and report the aggregate (typically average) performance. The PR and ROC curves are visualizations of this operation, and the area under these curves (integrating over these curves) are the expected performance across all thresholds. The PR curve plots precision vs recall and we may observe the trade-off between the two; the ROC curve plots TPR vs FPR and we also observe the trade-off between the two. Note that precision, recall, TPR and FPR are scalar values that all lie in the range of $[0, 1]$. As such, the PR and ROC curves' domains in the x and y axis are also in the range $[0, 1]$. Meaning, the whole area is a square and a square with side length of 1 has an area of 1. This geometric understanding is important because when we draw the curve and integrate, the integration value will always be $[0, 1]$ where a lower value indicates `bad` performance and a higher value indicates `good` performance. The integration of the curve is often called the area under the curve (AUC), and the AUC for PR is denoted as AUC-PR and for ROC is AUC-ROC. Simulate dataLet's start out by simulating some data. We have 2 classes, 0 and 1, and we sample from 2 different multivariate gaussians.* $\mu_0 = [1.5, 2.5, 3.3]$, $\Sigma_0 = \begin{bmatrix} 1 & 0 & 0 \\ 0 & 1 & 0 \\ 0 & 0 & 1 \end{bmatrix}$* $\mu_1 = [1.7, 2.9, 3.0]$, $\Sigma_1 = \begin{bmatrix} 1.0 & 0.5 & 0.2 \\ 0.5 & 1.0 & 2.0 \\ 0.2 & 2.0 & 1.0 \end{bmatrix}$We purposefully make the means between the two models very similar to give the classifier problems deciding. We take 10,000 samples from each (total 20,000 samples) for training, `T`, and 25 samples from each (total 50 samples) for validation, `V`.<jupyter_code>%matplotlib inline
import matplotlib.pyplot as plt
import seaborn as sns
import numpy as np
import pandas as pd
from collections import namedtuple
from scipy.stats import multivariate_normal
import warnings

warnings.filterwarnings('ignore')
np.random.seed(37)
plt.style.use('ggplot')

MVN = namedtuple('MVN', 'mean cov')
DATA = namedtuple('DATA', 'X, y')

mvn0 = MVN(np.array([1.5, 2.5, 3.3]), np.array([[1.0, 0.0, 0.0], [0.0, 1.0, 0.0], [0.0, 0.0, 1.0]]))
mvn1 = MVN(np.array([1.7, 2.9, 3.0]), np.array([[1.0, 0.5, 0.2], [0.5, 1.0, 2.0], [0.2, 2.0, 1.0]]))

N = 10000
X0 = np.array([multivariate_normal.rvs(mvn0.mean, mvn0.cov) for _ in range(N)])
X1 = np.array([multivariate_normal.rvs(mvn1.mean, mvn1.cov) for _ in range(N)])
y0 = np.full((N, 1), 0, dtype=np.int32)
y1 = np.full((N, 1), 1, dtype=np.int32)

X = np.vstack([X0, X1])
y = np.vstack([y0, y1])

T = DATA(X, y)

N = 25
X0 = np.array([multivariate_normal.rvs(mvn0.mean, mvn0.cov) for _ in range(N)])
X1 = np.array([multivariate_normal.rvs(mvn1.mean, mvn1.cov) for _ in range(N)])
y0 = np.full((N, 1), 0, dtype=np.int32)
y1 = np.full((N, 1), 1, dtype=np.int32)

X = np.vstack([X0, X1])
y = np.vstack([y0, y1])

V = DATA(X, y)<jupyter_output><empty_output><jupyter_text>Learn a classifierAfter we have our training `T` and validation `V` data, we plug in the data into a random forest classifier. We use the `predict_proba` function to retrieve the prediction probabilities of the validation example being 1 (or True). Below, `outcomes` is a list of tuples, where each tuple is a pair composed of the true label `y_t` and predicted label `y_p`.<jupyter_code>from sklearn.ensemble import RandomForestClassifier

rf_model = RandomForestClassifier(random_state=37)
rf_model.fit(T.X, T.y)

y_pred = rf_model.predict_proba(V.X)[:,1]

outcomes = sorted([(y_t, y_p) for y_t, y_p in zip(V.y[:,0], y_pred)], 
                  key=lambda tup: (-tup[1], -tup[0]))<jupyter_output><empty_output><jupyter_text>Construct the precision-recall curveNow with `y_t, y_p`, we create count tp, fp and fn and compute precision and recall. Note that `get_pr()` computes the precision and recall for one threshold `t`, and that `get_prs` returns a list of precisions, recalls, and thresholds. To make a precision-recall curve, we use the `step` function. We make two precision-recall curves; one where the precision and recall values are computed from Scikit and another one where these values are computed using our own code.  This side-by-side comparison is just to show that we know how to compute and graph these precision-recall curves ourselves.<jupyter_code>from sklearn.metrics import precision_recall_curve

def get_pr(outcomes, t):
    tp = 0
    fp = 0
    fn = 0

    for y_t, y_p in outcomes:
        clazz_t = int(y_t)
        clazz_p = 1 if y_p >= t else 0
        tp = tp + (1 if clazz_t == 1 and clazz_p == 1 else 0)
        fp = fp + (1 if clazz_t == 0 and clazz_p == 1 else 0)
        fn = fn + (1 if clazz_t == 1 and clazz_p == 0 else 0)

    pre = tp / (tp + fp)
    rec = tp / (tp + fn)
    return pre, rec

def get_prs(outcomes, thresholds=np.linspace(0.2, 1, 9)):
    pres = []
    recs = []
    thrs = []
    for t in thresholds:
        pre, rec = get_pr(outcomes, t)
        pres.append(pre)
        recs.append(rec)
        thrs.append(t)
    pres.append(1.0)
    recs.append(0.0)
    return np.array(pres), np.array(recs), np.array(thrs)

pre_m, rec_m, _ = get_prs(outcomes, thresholds=np.linspace(0.2, 1, 100))
pre_s, rec_s, _ = precision_recall_curve(V.y, y_pred)
baseline = np.sum(V.y) / len(V.y)

fig, ax = plt.subplots(1, 2, figsize=(15, 5))

ax[0].step(rec_s, pre_s, color='b', alpha=0.5, where='post', label='PR curve')
ax[0].set_xlabel('recall')
ax[0].set_ylabel('precision')
ax[0].set_title('Scikit Precision-Recall Curve')
ax[0].plot((0, 1), (baseline, baseline), 'r--', alpha=0.3, label='baseline')
ax[0].legend()

ax[1].step(rec_m, pre_m, color='b', alpha=0.5, where='post', label='PR curve')
ax[1].set_xlabel('recall')
ax[1].set_ylabel('precision')
ax[1].set_title('Manual Precision-Recall Curve')
ax[1].plot((0, 1), (baseline, baseline), 'r--', alpha=0.3, label='baseline')
ax[1].legend()

plt.tight_layout()<jupyter_output><empty_output><jupyter_text>Integrate the precision-recall curve (using trapezoid method)Now we integrate over these curves.* `apr` is the average precision value computed from Scikit's API* `apr_s` is the average precision value computed from our code (precision and recall from Scikit)* `apr_m` is the average precision value computed from our code (precision and recall from our code)Note how they are off? Our code uses trapezoid integration which is optimistic.<jupyter_code>from sklearn.metrics import average_precision_score

def get_apr(pre, rec):
    x = np.flip(rec)
    y = np.flip(pre)
    return np.trapz(y, x)

apr = average_precision_score(V.y, y_pred)
apr_s = get_apr(pre_s, rec_s)
apr_m = get_apr(pre_m, rec_m)

print('apr = {:.5f}, apr_s = {:.5f}, apr_m = {:.5f}'.format(apr, apr_s, apr_m))<jupyter_output>apr = 0.70137, apr_s = 0.71357, apr_m = 0.71357<jupyter_text>Integrate the precision-recall curve (using weighted average)Here, we use the suggested conservative integration approach and now all average percision agree.<jupyter_code>def get_apr(pre, rec):
    x = np.flip(rec)
    y = np.flip(pre)
    
    total = 0
    for i in range(len(x)):
        r_c = x[i]
        r_p = x[i if i - 1 < 0 else i - 1]
        p_i = y[i]
        a = (r_c - r_p) * p_i
        total = total + a
    return total

apr_s = get_apr(pre_s, rec_s)
apr_m = get_apr(pre_m, rec_m)

print('apr = {:.5f}, apr_s = {:.5f}, apr_m = {:.5f}'.format(apr, apr_s, apr_m))<jupyter_output>apr = 0.70137, apr_s = 0.70137, apr_m = 0.70137<jupyter_text>Construct the receiver operating characterics (ROC) curveWe turn our attention to the ROC curve now. The algorithms for generating the data for PR and ROC curves are nearly identical and differ only in what we return. In PR, we return the precision and recall, but in ROC, we return TPR and FPR.Note that here, we pass in thresholds in a descending order and there is a 2 prepended to the thresholds in $[0, 1]$. The 2 prepended to the thresholds is easy to explain; it's simply there to level off the curve (we actually did a similar thing when constructing the PR curve, although we manipulated the TPR and FPR instead of the thresholds). The descending order of the thresholds makes the FPR vector sort in an increasing sequence.<jupyter_code>from sklearn.metrics import roc_curve

def get_fpr_tpr(outcomes, t):
    tp = 0
    tn = 0
    fp = 0
    fn = 0

    for y_t, y_p in outcomes:
        clazz_t = int(y_t)
        clazz_p = 1 if y_p >= t else 0
        tp = tp + (1 if clazz_t == 1 and clazz_p == 1 else 0)
        fp = fp + (1 if clazz_t == 0 and clazz_p == 1 else 0)
        fn = fn + (1 if clazz_t == 1 and clazz_p == 0 else 0)
        tn = tn + (1 if clazz_t == 0 and clazz_p == 0 else 0)

    fpr = fp / (fp + tn)
    tpr = tp / (tp + fn)
    return fpr, tpr

def get_all_fpr_tpr(outcomes, thresholds=np.flip(np.append(np.linspace(0, 1, 11), [2]))):
    fprs = []
    tprs = []
    thrs = []
    for t in thresholds:
        fpr, tpr = get_fpr_tpr(outcomes, t)
        fprs.append(fpr)
        tprs.append(tpr)
        thrs.append(t)
    return np.array(fprs), np.array(tprs), np.array(thrs)

fpr_m, tpr_m, t_m = get_all_fpr_tpr(outcomes, thresholds=np.flip(np.append(np.linspace(0, 1, 100), [2])))
fpr_s, tpr_s, t_s = roc_curve(V.y, y_pred)

fig, ax = plt.subplots(1, 2, figsize=(15, 5))

ax[0].step(fpr_s, tpr_s, color='b', alpha=0.5, where='post', label='ROC curve')
ax[0].set_xlabel('fpr')
ax[0].set_ylabel('tpr')
ax[0].set_title('Scikit ROC')
ax[0].plot((0, 1), (0, 1), 'r--', alpha=0.2, label='baseline')
ax[0].legend()

ax[1].step(fpr_m, tpr_m, color='b', alpha=0.5, where='post', label='ROC curve')
ax[1].set_xlabel('fpr')
ax[1].set_ylabel('tpr')
ax[1].set_title('Manual ROC')
ax[1].plot((0, 1), (0, 1), 'r--', alpha=0.2, label='baseline')
ax[1].legend()

plt.tight_layout()<jupyter_output><empty_output><jupyter_text>Integrate the ROC curve (using trapezoid method)Here we use trapezoid integration over the ROC curves and observe that our manual approach agrees with Scikit's.<jupyter_code>from sklearn.metrics import roc_auc_score

def get_auc(tpr, fpr):
    return np.trapz(tpr, fpr)

auc = roc_auc_score(V.y, y_pred)
auc_s = get_auc(tpr_s, fpr_s)
auc_m = get_auc(tpr_m, fpr_m)

print('auc = {:.5f}, auc_s = {:.5f}, auc_m = {:.5f}'.format(auc, auc_s, auc_m))<jupyter_output>auc = 0.71200, auc_s = 0.71200, auc_m = 0.71200
{'fraction_non_alphanumeric': 0.12318076355199326, 'fraction_numerical': 0.028756239893130843, 'mean_word_length': 4.0084507042253525, 'pattern_counts': {'":': 0, '<': 28, '<?xml version=': 0, '>': 30, 'https://': 2, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '7', 'dir': 'jupyter-structured-clean-dedup', 'id': '283279', 'license': 'CC-BY-4.0', 'n_tokens_mistral': 5141, 'n_tokens_neox': 4737, 'n_words': 1830, 'path': 'sphinx/datascience/source/precision-recall-roc.ipynb', 'repo_name': 'oneoffcoder/books', 'size': '14223'}