find_best_threshold.py

import os
import numpy as np
from sklearn.metrics import (
    f1_score,
    precision_score,
    recall_score,
    roc_curve,
    auc,
    precision_recall_curve,
    det_curve,
)
import matplotlib.pyplot as plt


action = "test"
# unibo_0: "smallnet_mlp_full_as_is_REPLICATE_AGAIN" (epoch=9-step=3160.ckpt)
# unibo_1: "smallnet_mlp_full_as_is_REPLICATE_AGAIN" (epoch=15-step=5056.ckpt)
experiment_name = (
    "smallnet_mlp_full_as_is_REPLICATE_AGAIN_TEST_EPOCH15_RETRY-32-true_TEST"
)
which_folder = f"experiments_data/experiment_{experiment_name}/window_0"
which_predictions = None


def main():
    prediction_type = "validate" if action in {"fit", "validate"} else "test"
    if which_predictions is not None:
        predictions_file = os.path.join(which_folder, which_predictions)
    else:
        predictions_file = os.path.join(
            which_folder, f"predictions_{prediction_type}_all.npz"
        )
    preds_all = np.load(predictions_file)
    labels = preds_all["label"]
    preds = preds_all["prediction"]

    print("Finding best threshold...")
    thr = find_best_threshold(
        y_true=labels,
        y_probs=preds,
        metric=f1_score,
    )

    # AUROC
    fpr, tpr, thresholds = roc_curve(labels, preds)
    roc_auc = auc(fpr, tpr)
    fnr = 1 - tpr
    eer_threshold = thresholds[np.nanargmin(np.absolute((fnr - fpr)))]

    print("ROC AUC:", roc_auc)
    print("Best threshold (F1):", thr)
    print("EER Threshold:", eer_threshold)
    print(
        "EER:",
        fpr[np.nanargmin(np.absolute((fnr - fpr)))],
        fnr[np.nanargmin(np.absolute((fnr - fpr)))],
    )
    # thr = eer_threshold
    print("F1:", f1_score(labels, preds >= thr))
    print(
        "Precision:",
        precision_score(labels, preds >= thr),
    )
    print(
        "Recall:",
        recall_score(labels, preds >= thr),
    )
    print(
        "F1 score:",
        f1_score(labels, preds >= thr),
    )
    print(
        "Accuracy:",
        np.mean(labels == (preds >= thr).astype(int)),
    )

    # Calculate APCER and BPCER
    # Labels: 0 = real/bona fide, 1 = attack
    # Predictions >= threshold = classified as attack (1)
    # Predictions < threshold = classified as real (0)

    predicted_labels = (preds >= thr).astype(int)

    # APCER: Attack Presentation Classification Error Rate
    # Proportion of attack presentations incorrectly classified as bona fide
    attack_samples = labels == 1
    if np.sum(attack_samples) > 0:
        apcer = np.sum((predicted_labels == 0) & attack_samples) / np.sum(
            attack_samples
        )
    else:
        apcer = 0.0

    # BPCER: Bona Fide Presentation Classification Error Rate
    # Proportion of bona fide presentations incorrectly classified as attacks
    bona_fide_samples = labels == 0
    if np.sum(bona_fide_samples) > 0:
        bpcer = np.sum((predicted_labels == 1) & bona_fide_samples) / np.sum(
            bona_fide_samples
        )
    else:
        bpcer = 0.0

    print("APCER:", apcer)
    print("BPCER:", bpcer)
    print("ACER (Average):", (apcer + bpcer) / 2)

    # Save Precision-Recall curve
    precision, recall, pr_thresholds = precision_recall_curve(labels, preds)
    plt.figure(figsize=(8, 6))
    plt.plot(recall, precision, linewidth=2, label=f"PR Curve")

    # Find the point corresponding to the best threshold
    best_pred = (preds >= thr).astype(int)
    best_precision = precision_score(labels, best_pred)
    best_recall = recall_score(labels, best_pred)
    plt.plot(
        best_recall,
        best_precision,
        "ro",
        markersize=8,
        label=f"Best Threshold ({thr:.3f})",
    )

    plt.xlabel("Recall")
    plt.ylabel("Precision")
    plt.title("Precision-Recall Curve")
    plt.legend()
    plt.grid(True, alpha=0.3)
    plt.tight_layout()
    pr_curve_path = "precision_recall_curve.png"
    plt.savefig(pr_curve_path, dpi=300, bbox_inches="tight")
    plt.close()
    print(f"Precision-Recall curve saved to: {pr_curve_path}")

    # Save DET curve
    fpr_det, fnr, det_thresholds = det_curve(labels, preds)
    plt.figure(figsize=(8, 6))
    plt.plot(fpr_det, fnr, linewidth=2, label="DET Curve")

    # Find the point corresponding to the best threshold for DET curve
    best_pred = (preds >= thr).astype(int)
    best_fpr = np.sum((best_pred == 1) & (labels == 0)) / np.sum(
        labels == 0
    )  # False Positive Rate
    best_fnr = np.sum((best_pred == 0) & (labels == 1)) / np.sum(
        labels == 1
    )  # False Negative Rate
    plt.plot(
        best_fpr, best_fnr, "ro", markersize=8, label=f"Best Threshold ({thr:.3f})"
    )

    plt.xlabel("False Positive Rate")
    plt.ylabel("False Negative Rate")
    plt.title("Detection Error Tradeoff (DET) Curve")
    plt.legend()
    plt.grid(True, alpha=0.3)
    plt.tight_layout()
    det_curve_path = "det_curve.png"
    plt.savefig(det_curve_path, dpi=300, bbox_inches="tight")
    plt.close()
    print(f"DET curve saved to: {det_curve_path}")


def find_best_threshold(y_true, y_probs, metric=f1_score):
    thresholds = np.linspace(0, 1, 1000)
    scores = []

    for t in thresholds:
        y_pred = (y_probs >= t).astype(int)
        score = metric(y_true, y_pred)
        scores.append(score)

    best_index = np.argmax(scores)
    return thresholds[best_index]


if __name__ == "__main__":
    main()