LearnML Algorithms & EvaluationLogistic Regression and Metrics

📊MediumEvaluation & Benchmarks

Logistic Regression and Metrics

Route access-change requests with logistic regression from scratch: derive sigmoid and log loss, fit NumPy weights, select a cost-aware threshold on validation data, audit ranking and calibration, then compare with scikit-learn.

22 min read

Learning path

Step 32 of 177 in the full curriculum

Linear Regression from Scratch Decision Trees, Forests, and Boosting

Personalize this lesson

Adapt explanations and teaching visuals to your background and preferred voice.

Logistic regression turns features into a binary-routing score, then challenges you to justify treating that score as a probability.^{[1]Reference 1The Elements of Statistical Learning.https://hastie.su.domains/ElemStatLearn/} In the previous lesson, access-change request REQ-10234 needed a numeric estimate: would retrieved evidence fit the response-latency budget? Now it needs a decision: may the workflow auto-process the permission change under policy line P-7, or must a human review it first?

You'll fit a model by hand on eight historical access-change requests, derive its gradient, implement each line in NumPy, and check it against scikit-learn. A separate eight-request validation slice then carries the decisions that matter: threshold cost, ranking, and calibration. Training output explains the mechanism. Held-out labels are where you begin to trust an operating policy.

Logistic-regression probability field over ambiguity and policy-support features, with eight labeled training requests, decision contours at thresholds 0.20, 0.50, and 0.80, and held-out score, cost, AUC, and calibration audits. — Training requests shape a probability field and its threshold contours. Held-out requests then test the routing cost, ranking quality, and probability meaning without leaking validation labels into the fit.

You need the linear-regression chapter (design matrix, gradient descent loop, NumPy broadcasting) and comfort with basic probability. Each derivative gets worked with real numbers on paper before the code appears.

From latency estimates to access-review routing

In the previous chapter you predicted a number: "how many milliseconds will this evidence cost?" Now the target is binary: 1 means an access request needs human review before any permission change; 0 means policy evidence supports automatic handling. An incorrect automatic approval can cost $120 in permission rollback and investigation. An unnecessary review costs $18 of queue time.

Linear regression on a 0/1 label would happily predict 1.3 or −0.2. Those numbers are meaningless as probabilities. We need something that squeezes any real number into (0, 1).

This is also how you evaluate LLM-era classifiers. A prompt-injection guardrail, a toxicity filter, or an agent action router emits a score; a threshold turns that score into an intervention. Precision, recall, cost, and calibration belong in its model card and on-call review.

The sigmoid turns any score into a probability

The sigmoid (or logistic) function does that:

\sigma(z) = \frac{1}{1 + e^{-z}}

It has four properties:

As z → +∞, σ(z) → 1
As z → −∞, σ(z) → 0
At z = 0, σ(0) = 0.5
Its derivative is σ'(z) = σ(z) · (1 − σ(z)), which is positive for finite z and peaks at 0.25 when z=0. This tidy form is why the gradient of log-loss later becomes so simple.

Compute five concrete values by hand before reaching for code:

z	e^{-z}	1 + e^{-z}	σ(z)
−4.0	54.6	55.6	0.018
−1.0	2.718	3.718	0.269
0.0	1.0	2.0	0.500
1.0	0.368	1.368	0.731
4.0	0.018	1.018	0.982

sigmoid-probabilities.py

import numpy as np

def sigmoid(z: np.ndarray) -> np.ndarray:
    return 1.0 / (1.0 + np.exp(-z))

scores = np.array([-4.0, -1.0, 0.0, 1.0, 4.0])
print("scores:", scores)
print("probabilities:", sigmoid(scores).round(3))

Output

scores: [-4. -1.  0.  1.  4.]
probabilities: [0.018 0.269 0.5   0.731 0.982]

Smooth sigmoid and derivative curves from scores minus six to six, showing maximum slope 0.25 at zero and comparing an equal score move from minus one to one with a move from four to six, where probability changes by 0.462 versus 0.016. — The derivative peaks at `0.25` when `z = 0`, so a two-point score move near the center changes probability by `0.462`. The same score move in the positive tail changes it by only `0.016`.

The curve saturates quickly. A score of +4 produces p = 0.982; whether that number represents an honest frequency still needs validation.

The model itself is still linear inside:

z = w_0 + w_1 x_1 + w_2 x_2

Then p = sigma(z). The positive class is needs human review. This one line is the whole logistic regression model.

The score z is also the logit, or log-odds:^{[1]Reference 1The Elements of Statistical Learning.https://hastie.su.domains/ElemStatLearn/}

z = \log \left(\frac{p}{1 - p}\right)

At z = 0, the odds are 1:1 and p = 0.5. Adding one to z multiplies the odds by e; it doesn't guarantee that the displayed probability is calibrated on production traffic.

Our running example: eight access-change requests

We continue the REQ-10234 access-control workflow. Each historical request has two scaled features produced from retrieved evidence:

x1 = ambiguity score: conflicting audit logs, missing requester context, or unclear scope description
x2 = auto-policy support score: strength of evidence that policy line P-7 permits automatic handling

The label y = 1 when the completed case required human review before a permission change, and y = 0 when automatic handling was acceptable.

This tiny dataset stays with the entire chapter:

Request	`x1` ambiguity	`x2` auto-policy support	`y` needs review
R1	1.0	1	0
R2	1.5	2	0
R3	2.0	3	1
R4	3.5	6	1
R5	4.0	7	1
R6	4.5	5	1
R7	2.8	4	0
R8	3.2	8	0

We kept this training fixture small so you can verify each number. R3 needed review despite only moderate ambiguity. R8 had substantial ambiguity, but strong policy support made automatic processing acceptable. The resulting coefficient signs should now tell a coherent story: ambiguity raises review risk; supporting policy evidence lowers it.

Add the intercept column of ones and build the design matrix X (8 × 3) and the label vector y.

Log-loss (binary cross-entropy) is the right loss for probabilities

If the model outputs p = 0.93 for a request with y = 0, it made a confident mistake that deserves a large penalty. If it outputs p = 0.07 for that request, the loss is small. Negative log-likelihood encodes that asymmetry.

For one example the loss is the negative log-likelihood (also called binary cross-entropy).^{[2]Reference 2Pattern Recognition and Machine Learning.https://www.microsoft.com/en-us/research/publication/pattern-recognition-machine-learning/}^{[3]Reference 3Machine Learning: A Probabilistic Perspective.https://probml.github.io/pml-book/book0.html}

L = − \bigl[ y \log p + (1 − y) \log (1 − p) \bigr]

When y=1 the second term vanishes and we pay −log p (large when p is near 0). When y=0 we pay −log(1−p).

Why not squared error? Reusing MSE from linear regression makes the objective non-convex after composition with a sigmoid. For unregularized linear logistic regression, log loss is convex in the weights and yields the clean gradient below.^{[1]Reference 1The Elements of Statistical Learning.https://hastie.su.domains/ElemStatLearn/} Convexity makes optimization behavior easier to reason about, but it never replaces evaluation on unseen data.

Compute by hand for R4 (x1=3.5, x2=6, y=1) with a starting guess w = [0, 0, 0]:

z = 0 + 0·3.5 + 0·6 = 0
p = σ(0) = 0.5
L = − [1 · log(0.5) + 0] = − log(0.5) ≈ 0.693

Now suppose the model is slightly better and predicts p = 0.82 for the same row:

L = − log(0.82) ≈ 0.198
The loss dropped dramatically because the model became more confident in the correct direction, exactly the behavior we want.

For code, compute the same loss from logits rather than taking log(sigmoid(z)). The expression log(1 + exp(z)) - y*z is the same loss, but evaluating exp(z) directly can overflow. NumPy's logaddexp(0, z) - y*z evaluates that softplus term stably even when z is very large in either direction.

stable-log-loss-from-logits.py

import numpy as np

def log_loss_from_logits(y: np.ndarray, z: np.ndarray) -> float:
    return float(np.mean(np.logaddexp(0.0, z) - y * z))

print("correct positive:", round(log_loss_from_logits(np.array([1.0]), np.array([2.0])), 3))
print("confident wrong:", round(log_loss_from_logits(np.array([1.0]), np.array([-4.0])), 3))
extreme = log_loss_from_logits(np.array([0.0, 1.0]), np.array([1000.0, -1000.0]))
print("extreme logits finite:", np.isfinite(extreme))

Output

correct positive: 0.127
confident wrong: 4.018
extreme logits finite: True

The log-loss gradient collapses to (p − y) · x

Derive this once, with real numbers, on R4. The same pattern will reappear in neural-network classifiers later.

We have L(p, y) and p = σ(z), z = w·x (including w₀).

By chain rule: ∂L/∂wⱼ = (∂L/∂p) · (∂p/∂z) · (∂z/∂wⱼ)
First, ∂L/∂p = − (y/p − (1−y)/(1−p) )
Second, ∂p/∂z = p · (1 − p) (the sigmoid derivative)
Third, ∂z/∂wⱼ = xⱼ because increasing weight wⱼ changes the score in proportion to feature xⱼ.

The algebra collapses more neatly than it looks. Look at the two label cases:

label	loss	derivative with respect to score z	meaning
y = 1	−log(p)	p − 1	prediction is too low unless p is already near 1
y = 0	−log(1−p)	p	prediction is too high unless p is already near 0

Both rows are the same formula: ∂L/∂z = p − y. Then the weight gradient is:

\frac{\partial L}{\partial w_j} = (p - y) \cdot x_j

This holds for each weight, including the intercept (x₀ = 1).^{[2]Reference 2Pattern Recognition and Machine Learning.https://www.microsoft.com/en-us/research/publication/pattern-recognition-machine-learning/}^{[3]Reference 3Machine Learning: A Probabilistic Perspective.https://probml.github.io/pml-book/book0.html}

Plug in numbers for the starting point on R4 (p=0.5, y=1, x=[1, 3.5, 6]):

grad_w0 = (0.5 − 1) · 1 = −0.5
grad_w1 = (0.5 − 1) · 3.5 = −1.75
grad_w2 = (0.5 − 1) · 6 = −3.0

The gradient points uphill (direction of increasing loss). We subtract it (scaled by the learning rate) to descend.

That sign is the entire learning signal. If y=1 and p=0.5, p-y is negative, so subtracting the gradient increases weights attached to positive features and raises the future score for similar requests. If y=0 and p=0.8, p-y is positive, so subtracting the gradient lowers those weights and reduces the score.

Do the same arithmetic for each row, average the gradients, and you have one step of gradient descent for logistic regression. The code is close to the linear-regression GD loop you already wrote; the prediction and gradient formula change.

check-gradient-with-finite-differences.py

import numpy as np

def sigmoid(z: float) -> float:
    return float(1.0 / (1.0 + np.exp(-z)))

def loss(w: np.ndarray, x: np.ndarray, y: float) -> float:
    z = x @ w
    return float(np.logaddexp(0.0, z) - y * z)

x = np.array([1.0, 3.5, 6.0])
y = 1.0
w = np.zeros(3)
analytic = (sigmoid(x @ w) - y) * x
eps = 1e-6
numeric = np.array([
    (loss(w + eps * np.eye(3)[j], x, y) - loss(w - eps * np.eye(3)[j], x, y)) / (2 * eps)
    for j in range(3)
])
print("analytic:", analytic.round(3))
print("numeric :", numeric.round(3))
print("match:", np.allclose(analytic, numeric, atol=1e-6))

Output

analytic: [-0.5  -1.75 -3.  ]
numeric : [-0.5  -1.75 -3.  ]
match: True

Logistic-regression gradient heatmap for eight access-change requests at zero weights, showing each row contribution to intercept, ambiguity, and policy gradients, their average gradient, and the first update to weights 0, 0.069, and 0.075. — At zero weights, every request starts at `p = 0.5`. Multiplying each signed error by its feature vector creates the row heatmap; column averages give `g = [0, -0.344, -0.375]`, so the first `lr = 0.2` update becomes `w = [0, 0.069, 0.075]`.

Implementing logistic regression from scratch in NumPy

This complete, self-contained implementation can be copied into a file and run.

logistic-regression-numpy.py

import numpy as np

def sigmoid(z: np.ndarray) -> np.ndarray:
    return 1.0 / (1.0 + np.exp(-np.clip(z, -50, 50)))

def log_loss_from_logits(y_true: np.ndarray, z: np.ndarray) -> float:
    return float(np.mean(np.logaddexp(0.0, z) - y_true * z))

def fit_logistic(X: np.ndarray, y: np.ndarray, lr: float = 0.2, epochs: int = 3000, verbose: bool = True) -> np.ndarray:
    n, d = X.shape
    w = np.zeros(d)
    for epoch in range(epochs):
        z = X @ w
        p = sigmoid(z)
        grad = (1 / n) * (X.T @ (p - y))
        w -= lr * grad
        if verbose and epoch in (0, 200, 500, 1000, 2999):
            loss = log_loss_from_logits(y, X @ w)
            print(f"after_epoch={epoch + 1:4d}  loss={loss:.4f}  w={np.round(w, 3)}")
    return w

def predict_proba(X: np.ndarray, w: np.ndarray) -> np.ndarray:
    return sigmoid(X @ w)

def predict_with_threshold(p: np.ndarray, threshold: float) -> np.ndarray:
    return (p >= threshold).astype(int)

def confusion_matrix_at_threshold(y_true: np.ndarray, p: np.ndarray, threshold: float) -> dict[str, int]:
    pred = predict_with_threshold(p, threshold)
    tp = int(((pred == 1) & (y_true == 1)).sum())
    fp = int(((pred == 1) & (y_true == 0)).sum())
    fn = int(((pred == 0) & (y_true == 1)).sum())
    tn = int(((pred == 0) & (y_true == 0)).sum())
    return {"tp": tp, "fp": fp, "fn": fn, "tn": tn}

def precision_recall_f1(cm: dict[str, int]) -> tuple[float, float, float]:
    tp, fp, fn = cm["tp"], cm["fp"], cm["fn"]
    precision = tp / (tp + fp) if tp + fp else 0.0
    recall = tp / (tp + fn) if tp + fn else 0.0
    f1 = 2 * precision * recall / (precision + recall) if precision + recall else 0.0
    return precision, recall, f1

def roc_points(y_true: np.ndarray, p: np.ndarray, thresholds: list[float]) -> list[tuple[float, float, float]]:
    points = []
    for threshold in thresholds:
        cm = confusion_matrix_at_threshold(y_true, p, threshold)
        fpr = cm["fp"] / (cm["fp"] + cm["tn"]) if cm["fp"] + cm["tn"] else 0.0
        tpr = cm["tp"] / (cm["tp"] + cm["fn"]) if cm["tp"] + cm["fn"] else 0.0
        points.append((threshold, fpr, tpr))
    return points

def expected_calibration_error(y_true: np.ndarray, p: np.ndarray, n_bins: int = 4) -> float:
    edges = np.linspace(0.0, 1.0, n_bins + 1)
    ece = 0.0
    for lo, hi in zip(edges[:-1], edges[1:]):
        mask = (p >= lo) & ((p < hi) if hi < 1.0 else (p <= hi))
        if not mask.any():
            continue
        avg_pred = p[mask].mean()
        frac_pos = y_true[mask].mean()
        ece += mask.mean() * abs(avg_pred - frac_pos)
    return float(ece)

# Historical access-change requests: [ambiguity, auto-policy support]
X_raw = np.array([
    [1.0, 1.0], [1.5, 2.0], [2.0, 3.0], [3.5, 6.0],
    [4.0, 7.0], [4.5, 5.0], [2.8, 4.0], [3.2, 8.0],
])
y = np.array([0, 0, 1, 1, 1, 1, 0, 0])

X = np.c_[np.ones((len(X_raw), 1)), X_raw]   # add intercept column

w = fit_logistic(X, y, lr=0.2, epochs=3000)
p = predict_proba(X, w)
cm = confusion_matrix_at_threshold(y, p, threshold=0.5)
precision, recall, f1 = precision_recall_f1(cm)
ece = expected_calibration_error(y, p, n_bins=4)

print("\nFinal weights (w0, w1, w2):", w.round(3))
print("Probabilities:", p.round(3))
print("Confusion @ 0.5:", cm)
print("Precision / recall / F1:", tuple(round(v, 3) for v in (precision, recall, f1)))
print("Training ECE (4 bins):", round(ece, 3))
print("weights_match_reference", np.allclose(w, [-4.894, 3.399, -0.959], atol=0.01))
print("confusion_matches_reference", cm == {"tp": 3, "fp": 1, "fn": 1, "tn": 3})
print("f1_matches_reference", np.isclose(f1, 0.75))
print("ece_in_expected_band", 0.26 < ece < 0.28)

Output

after_epoch=   1  loss=0.6830  w=[0.    0.069 0.075]
after_epoch= 201  loss=0.4648  w=[-2.205  1.615 -0.462]
after_epoch= 501  loss=0.4206  w=[-3.524  2.452 -0.686]
after_epoch=1001  loss=0.4092  w=[-4.356  3.023 -0.85 ]
after_epoch=3000  loss=0.4072  w=[-4.894  3.399 -0.959]

Final weights (w0, w1, w2): [-4.894  3.399 -0.959]
Probabilities: [0.079 0.153 0.274 0.777 0.88  0.996 0.687 0.156]
Confusion @ 0.5: {'tp': 3, 'fp': 1, 'fn': 1, 'tn': 3}
Precision / recall / F1: (0.75, 0.75, 0.75)
Training ECE (4 bins): 0.268
weights_match_reference True
confusion_matches_reference True
f1_matches_reference True
ece_in_expected_band True

The negative intercept starts low when both signals are absent. Ambiguity's positive coefficient raises review risk, while the negative policy-support coefficient lowers it, matching the feature definition. Each progress line reports loss and weights after the same update, so snapshots stay comparable. These eight rows teach model mechanics; they never justify a shipping threshold by themselves.

The code builds a working logistic regression model whose scoring, gradients, and threshold behavior are visible.

Select a threshold on held-out requests

The fitted model emits a score in (0, 1). A business action needs a threshold t: route to review when p >= t. The training confusion matrix at t = 0.5 is a useful code check (TP=3, FP=1, FN=1, TN=3), but selecting t on those same rows would reward overfitting.

Use a validation slice with cases the optimizer never saw. Here a false negative means automatic processing when review was required and costs $120. A false positive means unnecessary human review and costs $18.

choose-threshold-on-validation.py

import numpy as np

def sigmoid(z: np.ndarray) -> np.ndarray:
    return 1.0 / (1.0 + np.exp(-z))

def metrics(y: np.ndarray, p: np.ndarray, threshold: float) -> tuple[float, float, float, int]:
    pred = (p >= threshold).astype(int)
    tp = int(((pred == 1) & (y == 1)).sum())
    fp = int(((pred == 1) & (y == 0)).sum())
    fn = int(((pred == 0) & (y == 1)).sum())
    precision = tp / (tp + fp) if tp + fp else 0.0
    recall = tp / (tp + fn) if tp + fn else 0.0
    f1 = 2 * precision * recall / (precision + recall) if precision + recall else 0.0
    cost = 120 * fn + 18 * fp
    return precision, recall, f1, cost

w = np.array([-4.894, 3.399, -0.959])
X_val_raw = np.array([
    [1.2, 2.0], [2.1, 4.0], [2.4, 3.0], [3.0, 5.0],
    [3.6, 4.0], [3.2, 7.0], [4.2, 5.0], [2.8, 6.0],
])
y_val = np.array([0, 0, 1, 0, 1, 0, 1, 1])
X_val = np.c_[np.ones(len(X_val_raw)), X_val_raw]
p_val = sigmoid(X_val @ w)

for t in (0.20, 0.50, 0.80):
    precision, recall, f1, cost = metrics(y_val, p_val, t)
    print(f"t={t:.2f} precision={precision:.3f} recall={recall:.3f} f1={f1:.3f} cost=${cost}")

Output

t=0.20 precision=0.667 recall=1.000 f1=0.800 cost=$36
t=0.50 precision=0.750 recall=0.750 f1=0.750 cost=$138
t=0.80 precision=1.000 recall=0.500 f1=0.667 cost=$240

On this validation slice, t = 0.20 costs least because avoiding missed reviews outweighs two extra queue items. In production, estimate these costs with stakeholders and re-evaluate thresholds as request mix shifts.

Precision, recall, F1, and class imbalance

For any threshold, fill the 2x2 confusion matrix:

Eight validation scores split at threshold 0.20 into TN 2, FP 2, FN 0, and TP 4, yielding precision 0.667, recall 1.0, and 36 dollars of weighted error cost. — At `t = 0.20`, all four required reviews are caught. Two routine requests enter the queue, producing `$36` of validation cost while avoiding the `$120` false-negative penalty.

Definitions:

Precision = TP / (TP + FP): among requests routed to review, what fraction required it?
Recall = TP / (TP + FN): among requests requiring review, what fraction reached it?
F1 = 2 * precision * recall / (precision + recall): one balance of those rates.

Accuracy alone can be nearly useless when positive cases are rare:

accuracy-fails-on-rare-review-cases.py

import numpy as np

y_true = np.array([1] * 4 + [0] * 96)
predict_auto_process_everything = np.zeros(100, dtype=int)
accuracy = (predict_auto_process_everything == y_true).mean()
recall = ((predict_auto_process_everything == 1) & (y_true == 1)).sum() / (y_true == 1).sum()
print("accuracy:", round(float(accuracy), 3))
print("review recall:", round(float(recall), 3))

Output

accuracy: 0.96
review recall: 0.0

ROC and precision-recall curves measure ranking

To draw an ROC (Receiver Operating Characteristic) curve, vary the threshold from high to low and plot recall (true-positive rate) on the y-axis against false-positive rate on the x-axis. AUC answers a ranking question: how often does a random positive score rank above a random negative score? Ties get half credit. AUC never chooses a business threshold for you.

auc-from-positive-negative-pairs.py

import numpy as np

p_val = np.array([0.061, 0.169, 0.595, 0.624, 0.971, 0.325, 0.990, 0.244])
y_val = np.array([0, 0, 1, 0, 1, 0, 1, 1])
positive = p_val[y_val == 1]
negative = p_val[y_val == 0]
pair_scores = [(pos > neg) + 0.5 * (pos == neg) for pos in positive for neg in negative]
auc = float(np.mean(pair_scores))
print("positive-negative pairs:", len(pair_scores))
print("validation AUC:", round(auc, 3))

Output

positive-negative pairs: 16
validation AUC: 0.812

Validation ROC staircase with threshold checkpoints 0.80, 0.50, and 0.20 alongside a four-by-four positive-versus-negative score matrix with 13 ranking wins out of 16 pairs, producing AUC 0.8125. — The staircase sweeps every validation score. The matrix gives the equivalent ranking view: positive requests outrank negatives in `13` of `16` pairs, so AUC is `13 / 16 = 0.8125`.

For rare positives, ROC can look generous because many true negatives keep false-positive rate small. Plot a precision-recall curve as well: it exposes how much review load accompanies higher recall.^{[1]Reference 1The Elements of Statistical Learning.https://hastie.su.domains/ElemStatLearn/}

Validation precision-recall curve beside review queues at thresholds 0.80, 0.50, and 0.20, showing two true positives and zero false positives, three and one, then four and two. — Lowering the threshold grows the review queue from `2 TP / 0 FP` to `4 TP / 2 FP`. The final required catch comes with one additional false alarm, leaving precision at `0.667`.

Average precision (AP) summarizes the precision-recall curve with a stepwise weighted average of precision at each positive-label recall increase. People sometimes label AP as PR-AUC, but trapezoidal area under the same curve can produce a different number. Name the statistic you compute. Unlike ROC-AUC, AP focuses on the positive class and is usually more informative when review-required cases are rare. On the same eight validation scores:

average-precision-on-validation.py

import numpy as np

p_val = np.array([0.061, 0.169, 0.595, 0.624, 0.971, 0.325, 0.990, 0.244])
y_val = np.array([0, 0, 1, 0, 1, 0, 1, 1])

# AP = average of precision@k for each true positive, in descending score order
order = np.argsort(-p_val)
y_sorted = y_val[order]
tp_cum = np.cumsum(y_sorted)
fp_cum = np.cumsum(1 - y_sorted)
precision_at = tp_cum / (tp_cum + fp_cum)
ap = float(precision_at[y_sorted == 1].mean())
print("validation average precision:", round(ap, 3))
print("validation ROC-AUC:   ", 0.812)

Output

validation average precision: 0.854
validation ROC-AUC:    0.812

On this balanced slice AP and ROC-AUC are close. When positives become rare, report both: ROC-AUC can stay high while AP falls with the positive prevalence and the ranking among the top of the queue.

Calibration: do scores represent frequencies?

A model can rank access-change requests well while producing over-confident or under-confident numbers. If requests scored near 0.80 require review only half the time, exposing 80% to an operator misstates risk.

A reliability diagram compares average score with observed positive frequency per bin. Expected Calibration Error (ECE) summarizes binned gaps:

\mathrm{ECE} = \sum_b \frac{|B_b|}{n}\left|\operatorname{avg}_{i\in B_b}(p_i) - \operatorname{avg}_{i\in B_b}(y_i)\right|

ECE depends on sample size, prevalence, and bin choices. It has no universal pass threshold. State the binning scheme, measure on held-out labels, and use a separate calibration split when fitting Platt scaling, isotonic regression, or temperature scaling. Temperature scaling is a common post-hoc method for neural classifiers.^{[4]Reference 4On Calibration of Modern Neural Networkshttps://arxiv.org/abs/1706.04599}

calibration-on-validation-labels.py

import numpy as np

p_val = np.array([0.061, 0.169, 0.595, 0.624, 0.971, 0.325, 0.990, 0.244])
y_val = np.array([0, 0, 1, 0, 1, 0, 1, 1])
edges = np.linspace(0.0, 1.0, 5)
ece = 0.0
for lo, hi in zip(edges[:-1], edges[1:]):
    mask = (p_val >= lo) & ((p_val < hi) if hi < 1.0 else (p_val <= hi))
    if mask.any():
        avg_p = p_val[mask].mean()
        frac_pos = y_val[mask].mean()
        ece += mask.mean() * abs(avg_p - frac_pos)
        print(f"[{lo:.2f}, {hi:.2f}) n={mask.sum()} avg_p={avg_p:.3f} frac_pos={frac_pos:.3f}")
print("validation ECE:", round(float(ece), 3))

Output

[0.00, 0.25) n=3 avg_p=0.158 frac_pos=0.333
[0.25, 0.50) n=1 avg_p=0.325 frac_pos=0.000
[0.50, 0.75) n=2 avg_p=0.609 frac_pos=0.500
[0.75, 1.00) n=2 avg_p=0.980 frac_pos=1.000
validation ECE: 0.139

Eight validation requests are enough for arithmetic, not a release gate. A single disputed label moves both metrics noticeably:

one-label-changes-small-sample-conclusions.py

import numpy as np

def f1_at_half(y: np.ndarray, p: np.ndarray) -> float:
    pred = p >= 0.5
    tp = ((pred == 1) & (y == 1)).sum()
    fp = ((pred == 1) & (y == 0)).sum()
    fn = ((pred == 0) & (y == 1)).sum()
    return float(2 * tp / (2 * tp + fp + fn))

def ece(y: np.ndarray, p: np.ndarray) -> float:
    edges = np.linspace(0.0, 1.0, 5)
    total = 0.0
    for lo, hi in zip(edges[:-1], edges[1:]):
        mask = (p >= lo) & ((p < hi) if hi < 1.0 else (p <= hi))
        if mask.any():
            total += mask.mean() * abs(p[mask].mean() - y[mask].mean())
    return float(total)

p = np.array([0.061, 0.169, 0.595, 0.624, 0.971, 0.325, 0.990, 0.244])
base = np.array([0, 0, 1, 0, 1, 0, 1, 1])
one_flip = base.copy()
one_flip[3] = 1
for name, labels in (("base", base), ("one flip", one_flip)):
    print(f"{name:8s} f1={f1_at_half(labels, p):.3f} ece={ece(labels, p):.3f}")

Output

base     f1=0.750 ece=0.139
one flip f1=0.889 ece=0.209

First audit action selection:

Validation threshold comparison with stacked error-cost bars: threshold 0.20 costs 36 dollars from two false positives, 0.50 costs 138 dollars from one false positive and one false negative, and 0.80 costs 240 dollars from two false negatives. — Threshold `0.20` wins this validation audit because two `$18` reviews cost less than one `$120` missed review. It's a candidate policy, not a permanent constant.

Then audit probability meaning:

Four validation calibration bins with paired average-score and observed-positive-rate bars, plus weighted ECE contributions 0.066, 0.041, 0.027, and 0.005 summing to 0.139. — Each bin contributes its score-frequency gap weighted by its share of the eight validation requests: `0.066 + 0.041 + 0.027 + 0.005 ≈ 0.139`. The arithmetic is exact for this slice; the evidence is still too small for a release claim.

Verify the NumPy fit with scikit-learn

A from-scratch implementation earns trust by matching a maintained implementation under equivalent settings.^{[5]Reference 5Scikit-learn: Machine Learning in Python.https://jmlr.org/papers/v12/pedregosa11a.html} Scikit-learn applies regularization by default; set C=np.inf here to compare against our unregularized gradient-descent fit.

compare-with-scikit-learn.py

import numpy as np
from sklearn.linear_model import LogisticRegression
from sklearn.metrics import f1_score, roc_auc_score

X_raw = np.array([
    [1.0, 1.0], [1.5, 2.0], [2.0, 3.0], [3.5, 6.0],
    [4.0, 7.0], [4.5, 5.0], [2.8, 4.0], [3.2, 8.0],
])
y = np.array([0, 0, 1, 1, 1, 1, 0, 0])
X = np.c_[np.ones((len(X_raw), 1)), X_raw]
w_ours = np.array([-4.894, 3.399, -0.959])
X_val_raw = np.array([
    [1.2, 2.0], [2.1, 4.0], [2.4, 3.0], [3.0, 5.0],
    [3.6, 4.0], [3.2, 7.0], [4.2, 5.0], [2.8, 6.0],
])
y_val = np.array([0, 0, 1, 0, 1, 0, 1, 1])
X_val = np.c_[np.ones(len(X_val_raw)), X_val_raw]

clf = LogisticRegression(C=np.inf, solver="lbfgs", max_iter=5000, fit_intercept=False)
clf.fit(X, y)
w_sk = clf.coef_[0]
print("sklearn weights:", np.round(w_sk, 3))
print("our weights   :", w_ours)
print("weights close :", np.allclose(w_sk, w_ours, atol=0.05))
p_val = clf.predict_proba(X_val)[:, 1]
print("validation F1 @ 0.20:", round(f1_score(y_val, (p_val >= 0.20).astype(int)), 3))
print("validation AUC:", round(roc_auc_score(y_val, p_val), 3))

Output

sklearn weights: [-4.914  3.414 -0.964]
our weights   : [-4.894  3.399 -0.959]
weights close : True
validation F1 @ 0.20: 0.8
validation AUC: 0.812

The weight match validates the mechanics. The validation metrics validate only this small held-out exercise; larger and temporally later data is still required before deployment.

Symptom, cause, fix table (the debugging checklist)

Symptom	Likely cause	Fix
All predicted probabilities hover around 0.4-0.6	Features are weak, data is tiny, or GD is under-trained	Add informative evidence features, increase epochs, tune learning rate, or evaluate a stronger model
F1 at threshold 0.5 is mediocre but recall at 0.2 is excellent	Class imbalance + wrong operating point	Read the ROC or cost-sensitive curve; pick threshold by expected business cost, not by accuracy
Good AUC but large reliability gaps on substantial held-out data	Model ranks well but scores misstate frequency	Withhold raw percentages; calibrate on a reserved split and re-audit
Your NumPy weights differ from sklearn by >0.3	Forgot the intercept column, used different regularization, or feature scaling mismatch	Compare design matrices, disable sklearn regularization, and verify the same X
Loss decreases then suddenly increases	Learning rate too large or probability-form loss overflowed	Lower `lr`, compute loss with `logaddexp`, and print loss during fitting

Try it yourself

These exercises turn the chapter from reading into mastery. Use training rows for fitting and validation rows for policy selection.

Start with w = [0,0,0] and compute the gradient vector on the full batch by hand for the initial step (lr = 0.5). Verify it matches the initial printed line of your NumPy run.
Change the learning rate to 2.0 and watch the loss explode or oscillate. Explain why in one sentence.
Sweep thresholds [0.2, 0.35, 0.5, 0.65, 0.8] on validation and print precision, recall, F1, and cost. Which threshold minimizes the stated cost?
Bin validation probabilities into four equal-width bins, compute ECE by hand, and compare with code output. Repeat with eight bins and explain why the number changes.
Flip one validation label and measure F1 and ECE again. Write the argument for gathering a larger, later-in-time validation set before publishing probability claims.

If any answer feels hand-wavy, print the arrays, loss numbers, and concrete cells of the confusion matrix. Concrete numbers build classification skill.

Eight validation requests remain far too few for a release claim. Next, decision trees and ensembles attack the same access-review routing task with nonlinear rules; the later cross-validation lesson deepens the evidence required before shipping a model.

Next Step

Continue to Decision Trees, Forests, and Boosting

There you'll model the same access-review decision with nonlinear rules and inspect why high-scoring rankings still need calibration audits.

PreviousLinear Regression from Scratch

Share this article

X Facebook LinkedIn Bluesky Reddit Hacker News Email

References

The Elements of Statistical Learning.

Hastie, T., Tibshirani, R., Friedman, J. · 2009

https://hastie.su.domains/ElemStatLearn/

Pattern Recognition and Machine Learning.

Bishop, C. M. · 2006

https://www.microsoft.com/en-us/research/publication/pattern-recognition-machine-learning/

Machine Learning: A Probabilistic Perspective.

Murphy, K. P. · 2012

https://probml.github.io/pml-book/book0.html

On Calibration of Modern Neural Networks

Guo, C., Pleiss, G., Sun, Y., & Weinberger, K. Q. · 2017

https://arxiv.org/abs/1706.04599

Scikit-learn: Machine Learning in Python.

Pedregosa, F., et al. · 2011 · JMLR

https://jmlr.org/papers/v12/pedregosa11a.html

Discussion

Questions and insights from fellow learners.

Discussion loads when you reach this section.

Logistic Regression and Metrics

From latency estimates to access-review routing

The sigmoid turns any score into a probability

Compute $\sigma(2.0)$ by hand. What decision does it suggest for an access-change request if the human-review threshold is `0.5`?

Our running example: eight access-change requests

Log-loss (binary cross-entropy) is the right loss for probabilities

The log-loss gradient collapses to (p − y) · x

An access-change request has $y = 0$ but the model predicts $p = 0.8$ . Is $p-y$ positive or negative, and what does gradient descent do to the score for similar requests?

You fit the model and the final weight vector is $w = [-4.9, 3.4, -0.96]$ . For `REQ-10234`, suppose $x_1 = 3.0$ (ambiguity) and $x_2 = 5.0$ (auto-policy support). Compute $z$ and its human-review score. Would you route it to review at threshold `0.5`?

Implementing logistic regression from scratch in NumPy

Select a threshold on held-out requests

At threshold `0.50`, validation has `FP=1` and `FN=1`. At threshold `0.20`, it has `FP=2` and `FN=0`. If an automatic-processing mistake costs `$120` and an unnecessary review costs `$18`, which threshold wins on validation?

Precision, recall, F1, and class imbalance

ROC and precision-recall curves measure ranking

Calibration: do scores represent frequencies?

Validation AUC is `0.812` and four-bin ECE is `0.139`, both measured on only eight requests. A product manager wants to display "This access request has an 87% chance of manual review." What is your response?

Verify the NumPy fit with scikit-learn

Symptom, cause, fix table (the debugging checklist)

Try it yourself

Mastery Check

Discussion

From latency estimates to access-review routing

The sigmoid turns any score into a probability

Compute σ(2.0)\sigma(2.0)σ(2.0) by hand. What decision does it suggest for an access-change request if the human-review threshold is 0.5?

Compute σ(2.0)\sigma(2.0)σ(2.0) by hand. What decision does it suggest for an access-change request if the human-review threshold is 0.5?

Our running example: eight access-change requests

Log-loss (binary cross-entropy) is the right loss for probabilities

The log-loss gradient collapses to (p − y) · x

An access-change request has y=0y = 0y=0 but the model predicts p=0.8p = 0.8p=0.8. Is p−yp-yp−y positive or negative, and what does gradient descent do to the score for similar requests?

An access-change request has y=0y = 0y=0 but the model predicts p=0.8p = 0.8p=0.8. Is p−yp-yp−y positive or negative, and what does gradient descent do to the score for similar requests?

Implementing logistic regression from scratch in NumPy

Select a threshold on held-out requests

At threshold 0.50, validation has FP=1 and FN=1. At threshold 0.20, it has FP=2 and FN=0. If an automatic-processing mistake costs $120 and an unnecessary review costs $18, which threshold wins on validation?

At threshold 0.50, validation has FP=1 and FN=1. At threshold 0.20, it has FP=2 and FN=0. If an automatic-processing mistake costs $120 and an unnecessary review costs $18, which threshold wins on validation?

Precision, recall, F1, and class imbalance

ROC and precision-recall curves measure ranking

Calibration: do scores represent frequencies?

Validation AUC is 0.812 and four-bin ECE is 0.139, both measured on only eight requests. A product manager wants to display "This access request has an 87% chance of manual review." What is your response?

Validation AUC is 0.812 and four-bin ECE is 0.139, both measured on only eight requests. A product manager wants to display "This access request has an 87% chance of manual review." What is your response?

Verify the NumPy fit with scikit-learn

Symptom, cause, fix table (the debugging checklist)

Try it yourself

Mastery Check

Discussion

Compute $\sigma(2.0)$ by hand. What decision does it suggest for an access-change request if the human-review threshold is `0.5`?

An access-change request has $y = 0$ but the model predicts $p = 0.8$ . Is $p-y$ positive or negative, and what does gradient descent do to the score for similar requests?

At threshold `0.50`, validation has `FP=1` and `FN=1`. At threshold `0.20`, it has `FP=2` and `FN=0`. If an automatic-processing mistake costs `$120` and an unnecessary review costs `$18`, which threshold wins on validation?

Validation AUC is `0.812` and four-bin ECE is `0.139`, both measured on only eight requests. A product manager wants to display "This access request has an 87% chance of manual review." What is your response?