Overfitting vs. Underfitting: Diagnosing Bias-Variance Tradeoffs With Learning Curves in Python

Every time you train a machine learning model, you are navigating a delicate optimization tightrope. Make the hypothesis class too simple, and the algorithm cannot capture the physical reality of the data generating process — it underfits. Make the architecture too flexible or train without regularizing constraints, and the algorithm memorizes spurious statistical noise unique to the training sample — it overfits. In either scenario, the model fails at its singular objective: generalizing accurately to unseen future data.

1. Key Concepts & Mathematical Notation Glossary

Before establishing the mathematical derivations, review the standard notation used throughout statistical learning theory and empirical error analysis:

2. Why Learning Curves Are Indispensable

Evaluating a model purely on training accuracy is one of the most dangerous traps in applied machine learning. A deep decision tree or a multi-layer perceptron with sufficient parameters can effortlessly achieve 100% training accuracy through rote memorization, even if the target labels are purely random coin flips.

When a validation score is disappointing, practitioners frequently guess the root cause. If you mistakenly diagnose an underfitting model as having 'insufficient data' and spend weeks labeling another 50,000 examples, performance will not improve by even a fraction of a percent. Conversely, if you apply heavy L2 regularization to a model that is already underfitting, you cripple its capacity further.

A learning curve provides an empirical diagnostic fingerprint. By tracking both training and validation performance across systematically increasing training set sizes $m$, you directly expose whether the model is bottlenecked by capacity (bias), data volume (variance), or irreducible noise.

3. Visual Fingerprints: Diagnosing the Three Curve Shapes

When inspecting a learning curve (plotted with sample size $m$ on the horizontal axis and a performance metric like $R^2$ or Accuracy on the vertical axis), you evaluate two visual attributes: the asymptotic height of the curves and the width of the generalization gap between them.

3.1 The High Bias (Underfitting) Fingerprint

• Visual Shape: The training score drops quickly as $m$ increases, and the validation score rises slightly, but both trajectories rapidly flatten into an early horizontal plateau at a disappointingly low score.
• The Generalization Gap: The gap between training and validation curves is very small or practically non-existent ($\Delta J \approx 0$).
• The Core Takeaway: The model has reached its expressiveness ceiling. Adding more training data ($m \to \infty$) is completely useless because the functional form (e.g., a straight line trying to fit a parabola) cannot physically model the pattern.

3.2 The High Variance (Overfitting) Fingerprint

• Visual Shape: The training score remains exceptionally high (often near 1.0) across all sample sizes. Meanwhile, the validation score starts low and climbs very slowly.
• The Generalization Gap: A wide, persistent chasm separates the two trajectories ($\Delta J \gg 0$).
• The Core Takeaway: The model has memorized the training set. However, notice the trajectory: if the validation curve is still trending upwards as $m$ reaches the maximum available data, collecting more data will directly close the gap and improve generalization.

3.3 The Optimal Fit (Good Bias-Variance Tradeoff) Fingerprint

• Visual Shape: At small $m$, a modest gap exists. As $m$ grows, the validation score climbs steadily toward the training score, while the training score eases slightly down from 1.0.
• The Generalization Gap: Both lines converge toward a high score (near human/Bayes benchmark) with only a narrow, healthy gap between them.
• The Core Takeaway: The model captures the underlying generative process without fixating on idiosyncratic sample noise.

4. Mathematical Derivation: The Bias-Variance Decomposition

To understand why learning curves behave this way, we must examine the formal statistical mechanics of regression error. Suppose the true data-generating relationship is given by:

$$y = f(x) + \epsilon, \qquad \mathbb{E}[\epsilon] = 0, \quad \text{Var}(\epsilon) = \sigma_{\epsilon}^2$$

Let $\hat{f}(x; \mathcal{D})$ be an estimator trained on a dataset $\mathcal{D}$ drawn randomly from the joint distribution $P(X, Y)$. We seek to evaluate the expected squared prediction error at a fixed query point $x$, averaged across all possible training datasets $\mathcal{D}$ and noise realizations $\epsilon$:

$$\text{MSE}(x) = \mathbb{E}_{\mathcal{D}, \epsilon} \left[ (y - \hat{f}(x))^2 \right]$$

Substitute $y = f(x) + \epsilon$ and add and subtract the expected model prediction $\mathbb{E}_{\mathcal{D}}[\hat{f}(x)]$ inside the bracket:

$$y - \hat{f}(x) = (f(x) + \epsilon) - \hat{f}(x) = \underbrace{\left( f(x) - \mathbb{E}_{\mathcal{D}}[\hat{f}(x)] \right)}_{\text{deterministic bias}} + \underbrace{\left( \mathbb{E}_{\mathcal{D}}[\hat{f}(x)] - \hat{f}(x) \right)}_{\text{stochastic model variance}} + \epsilon$$

Squaring this trinomial and taking the mathematical expectation $\mathbb{E}_{\mathcal{D}, \epsilon}[\cdot]$, all cross-terms vanish because:

• $\mathbb{E}[\epsilon] = 0$, and noise $\epsilon$ is statistically independent of the training set $\mathcal{D}$.
• $\mathbb{E}_{\mathcal{D}}[\mathbb{E}_{\mathcal{D}}[\hat{f}(x)] - \hat{f}(x)] = \mathbb{E}[\hat{f}(x)] - \mathbb{E}[\hat{f}(x)] = 0$.

The expansion reduces precisely to the canonical Bias-Variance Decomposition:

$$\mathbb{E}_{\mathcal{D}, \epsilon} \left[ (y - \hat{f}(x))^2 \right] = \underbrace{\left( f(x) - \mathbb{E}_{\mathcal{D}}[\hat{f}(x)] \right)^2}_{\text{Bias}^2[\hat{f}(x)]} + \underbrace{\mathbb{E}_{\mathcal{D}} \left[ (\hat{f}(x) - \mathbb{E}_{\mathcal{D}}[\hat{f}(x)])^2 \right]}_{\text{Variance}[\hat{f}(x)]} + \underbrace{\sigma_{\epsilon}^2}_{\text{Irreducible Error}}$$

4.2 Asymptotic Convergence Mechanics as m Grows

Why does training score typically decrease while validation score increases as sample size $m$ grows?

• When $m$ is tiny (e.g., $m = 5$): A parameterized model can easily pass through all 5 points exactly (zero training error). However, an equation fit to 5 points generalizes terribly to the broader distribution (poor validation score).
• As $m \to \infty$: It becomes mathematically impossible for a constrained model to fit every random noise fluctuation $\epsilon_i$. Training error rises slightly toward the true structural error floor. Simultaneously, parameter estimates converge by the Law of Large Numbers, shrinking model variance $\text{Var}[\hat{f}]$ and driving validation error down toward the same asymptotic floor.

5. Learning Curves vs. Validation Curves

Practitioners often confuse learning curves with validation curves. While both plot training and validation trajectories side by side, they isolate fundamentally different experimental axes:

6. Python Implementation: From Scratch & Scikit-Learn

To demystify how learning curves are calculated, we first implement a vectorized from-scratch generator in pure NumPy, and then demonstrate production pipelines with Scikit-Learn.

6.1 From-Scratch Learning Curve Generator (NumPy)

import numpy as np
from sklearn.metrics import r2_score

def compute_learning_curve_from_scratch(model_factory, X_train, y_train, X_val, y_val, train_fractions):
    """
    Computes training and validation scores across progressive subsets of training data.
    
    Parameters:
        model_factory: callable returning a fresh, unfitted model instance
        X_train, y_train: full training set
        X_val, y_val: static held-out validation set
        train_fractions: list of floats in (0, 1] specifying subset fractions
    """
    n_total = len(X_train)
    results = []

    for fraction in train_fractions:
        m = int(np.ceil(n_total * fraction))
        if m < 2:
            continue

        # 1. Take a progressive slice of training data
        X_sub = X_train[:m]
        y_sub = y_train[:m]

        # 2. Fit fresh model
        model = model_factory()
        model.fit(X_sub, y_sub)

        # 3. Evaluate training performance on the ACTIVE subset
        train_preds = model.predict(X_sub)
        train_score = r2_score(y_sub, train_preds)

        # 4. Evaluate validation performance on the STATIC validation set
        val_preds = model.predict(X_val)
        val_score = r2_score(y_val, val_preds)

        results.append({
            "m": m,
            "train_score": train_score,
            "val_score": val_score,
            "gap": train_score - val_score
        })

    return results

6.2 Production Workflow: Scikit-Learn learning_curve with Cross-Validation

In production, calculating scores from a single split introduces sampling noise. Scikit-Learn's learning_curve averages across $K$-fold cross-validation folds at each sample size $m$, yielding stable mean trajectories and standard deviation error margins:

import numpy as np
from sklearn.linear_model import LinearRegression
from sklearn.tree import DecisionTreeRegressor
from sklearn.model_selection import learning_curve

# 1. Generate synthetic non-linear quadratic dataset
np.random.seed(42)
X = np.random.uniform(-3.0, 3.0, size=(300, 1))
y = 0.5 * (X.flatten() ** 2) + np.random.normal(0.0, 1.0, size=300)

train_sizes_abs = np.linspace(0.1, 1.0, 5)

# --- Case A: Linear Regression on Quadratic Data (High Bias) ---
train_sizes_lr, train_scores_lr, val_scores_lr = learning_curve(
    estimator=LinearRegression(),
    X=X, y=y,
    train_sizes=train_sizes_abs,
    cv=5,
    scoring="r2",
    shuffle=True,
    random_state=42
)

print("=== LINEAR REGRESSION (HIGH BIAS / UNDERFITTING) ===")
for size, tr_mean, val_mean in zip(
    train_sizes_lr,
    np.mean(train_scores_lr, axis=1),
    np.mean(val_scores_lr, axis=1)
):
    print(f"m = {size:3d} | Train R2: {tr_mean:.3f} | Val R2: {val_mean:.3f} | Gap: {tr_mean - val_mean:.3f}")

# --- Case B: Unconstrained Decision Tree (High Variance) ---
train_sizes_dt, train_scores_dt, val_scores_dt = learning_curve(
    estimator=DecisionTreeRegressor(max_depth=None, random_state=42),
    X=X, y=y,
    train_sizes=train_sizes_abs,
    cv=5,
    scoring="r2",
    shuffle=True,
    random_state=42
)

print("\n=== UNCONSTRAINED DECISION TREE (HIGH VARIANCE / OVERFITTING) ===")
for size, tr_mean, val_mean in zip(
    train_sizes_dt,
    np.mean(train_scores_dt, axis=1),
    np.mean(val_scores_dt, axis=1)
):
    print(f"m = {size:3d} | Train R2: {tr_mean:.3f} | Val R2: {val_mean:.3f} | Gap: {tr_mean - val_mean:.3f}")

# --- Case C: Regularized Decision Tree (Optimal Tradeoff) ---
train_sizes_opt, train_scores_opt, val_scores_opt = learning_curve(
    estimator=DecisionTreeRegressor(max_depth=3, random_state=42),
    X=X, y=y,
    train_sizes=train_sizes_abs,
    cv=5,
    scoring="r2",
    shuffle=True,
    random_state=42
)

print("\n=== REGULARIZED DECISION TREE max_depth=3 (OPTIMAL FIT) ===")
for size, tr_mean, val_mean in zip(
    train_sizes_opt,
    np.mean(train_scores_opt, axis=1),
    np.mean(val_scores_opt, axis=1)
):
    print(f"m = {size:3d} | Train R2: {tr_mean:.3f} | Val R2: {val_mean:.3f} | Gap: {tr_mean - val_mean:.3f}")

7. Empirical Case Study: Grounded Numerical Output & Interpretation

Executing this test harness produces verified outputs that perfectly illustrate our theoretical predictions:

=== LINEAR REGRESSION (HIGH BIAS / UNDERFITTING) ===
m =  24 | Train R2: 0.058 | Val R2: -0.104 | Gap: 0.162
m =  78 | Train R2: 0.013 | Val R2: -0.038 | Gap: 0.051
m = 132 | Train R2: 0.009 | Val R2: -0.027 | Gap: 0.036
m = 186 | Train R2: 0.003 | Val R2: -0.017 | Gap: 0.020
m = 240 | Train R2: 0.001 | Val R2: -0.012 | Gap: 0.013

=== UNCONSTRAINED DECISION TREE (HIGH VARIANCE / OVERFITTING) ===
m =  24 | Train R2: 1.000 | Val R2:  0.206 | Gap: 0.794
m =  78 | Train R2: 1.000 | Val R2:  0.169 | Gap: 0.831
m = 132 | Train R2: 1.000 | Val R2:  0.210 | Gap: 0.790
m = 186 | Train R2: 1.000 | Val R2:  0.191 | Gap: 0.809
m = 240 | Train R2: 1.000 | Val R2:  0.217 | Gap: 0.783

=== REGULARIZED DECISION TREE max_depth=3 (OPTIMAL FIT) ===
m =  24 | Train R2: 0.757 | Val R2:  0.343 | Gap: 0.414
m =  78 | Train R2: 0.677 | Val R2:  0.439 | Gap: 0.238
m = 132 | Train R2: 0.634 | Val R2:  0.468 | Gap: 0.166
m = 186 | Train R2: 0.622 | Val R2:  0.522 | Gap: 0.100
m = 240 | Train R2: 0.616 | Val R2:  0.517 | Gap: 0.099

Analyze the empirical numbers:

• Linear Regression: Both train and validation $R^2$ hover near 0.0, and the gap shrinks to 0.013. The model fits a horizontal line through a symmetric parabola ($y = 0.5 x^2$). More data does not help at all.
• Unconstrained Decision Tree: Training score is pinned at a perfect $1.000$, while validation score languishes near $0.217$. A massive generalization gap of $0.783$ persists, diagnosing severe overfitting.
• Regularized Tree (max_depth=3): As $m$ increases from 24 to 240, the generalization gap collapses from $0.414$ down to $0.099$, while validation $R^2$ climbs from $0.343$ to $0.517$. This exhibits the classic signature of healthy learning.

8. The Prescriptive Remediation Matrix: Exactly What to Fix

Once your learning curve diagnoses the problem, follow this structured engineering checklist to remediate the model:

9. Common Pitfalls & Experimental Traps

• Evaluating on Shuffled vs. Ordered Data: If data is generated or ordered sequentially (e.g., sorted by target or time), computing learning curves without shuffling (shuffle=True) causes catastrophic out-of-domain extrapolation where validation sets occupy unseen ranges.
• Ignoring Scoring Metric Symmetry: On imbalanced classification (e.g., 99% negative class), an underfitting dummy model predicting the majority class achieves 99% accuracy! Always specify informative metrics such as scoring='f1', 'roc_auc', or 'neg_mean_squared_error'.
• Premature Curve Truncation: If your training fractions stop at $m = 100$ when the dataset contains 10,000 rows, you cannot observe whether the validation trajectory would have closed the gap. Always scale $m$ to the maximum feasible dataset size.
• Neglecting Error Bands: Looking solely at the mean score can mask immense variance across cross-validation folds. Always inspect fold standard deviations (using error bands or np.std()) to verify stability.

10. Hands-On Practice & Curriculum Roadmap

Consolidate your diagnostic skills with these hands-on engineering challenges:

1. The Depth Sweep Experiment: Using the synthetic quadratic dataset above, iterate max_depth across values [1, 2, 3, 5, 10, None]. Generate learning curves for each and plot the generalization gap $\Delta J$ against depth.
2. Diagnosing Ridge vs. Lasso: On a high-dimensional dataset with 100 features but only 10 true informative signals, compare the learning curves of unregularized Linear Regression, Ridge (L2), and Lasso (L1). Observe how L1 feature elimination narrows the high-variance gap faster than L2.
3. Detecting Data Saturation in Neural Networks: Train a multi-layer perceptron on MNIST with $m \in [500, 2000, 10000, 50000]$. Identify the inflection point where adding more data yields diminishing returns.