Gradient Boosting: Regression and Classification

AdaBoost reweights training samples so the next stump focuses on hard examples. Gradient Boosting takes a different path: each new tree directly predicts the errors of the current ensemble. The framework generalizes to any differentiable loss — MSE for regression, log loss for classification.

Anchor: 6-sample house prices (regression trace) and 8-sample loan defaults (classification trace).

import numpy as np

# Regression anchor
X_reg = np.array([650, 850, 1100, 1400, 1600, 1900])
y_reg = np.array([180, 220, 280, 340, 370, 430])

# Classification anchor: [income_$k, credit_score]
X_clf = np.array([[25,580],[32,610],[45,650],[60,680],[70,710],[80,730],[90,750],[110,780]])
y_clf = np.array([1, 1, 1, 0, 0, 0, 0, 0])  # 1=default

AdaBoost vs Gradient Boosting

Gradient Boosting is the more general framework: AdaBoost is a special case of GB with exponential loss and stumps.

Gradient Boosting Regression — 3-Round Trace

Initial Prediction

$F_0(x)=\bar{y}=(180+220+280+340+370+430)/6=303.3$

Start with the mean. MSE loss is minimized by the mean, so this is the optimal constant prediction.

Round 1: Fit a Tree to Residuals

Pseudo-residuals: $r_i^{(1)}=y_i-F_0(x_i)=y_i-303.3$

Train a regression stump on $(X_{\text{reg}},r_1)$. Best split at sq_ft ≤ 1250:

• Left leaf (samples 1,2,3): $\bar{r}=(-123.3-83.3-23.3)/3=-76.6$
• Right leaf (samples 4,5,6): $\bar{r}=(36.7+66.7+126.7)/3=76.7$

Update with learning rate $\nu =0.1$:

$F_1(x)=F_0(x)+\nu \cdot h_1(x)$

• sq_ft ≤ 1250: $F_1=303.3+0.1\times (-76.6)=303.3-7.66=295.6$
• sq_ft > 1250: $F_1=303.3+0.1\times 76.7=303.3+7.67=311.0$

Round 2: Fit Tree to New Residuals

$r_i^{(2)}=y_i-F_1(x_i)$:

• Sample 1: $180-295.6=-115.6$ (was −123.3 — shrinking ✓)
• Sample 4: $340-311.0=+29.0$ (was +36.7 — shrinking ✓)

Same split threshold wins again (sq_ft ≤ 1250). New leaf means: left = −75.7, right = +69.3.

$F_2:\text{left}=295.6+0.1\times (-75.7)=288.0,\text{right}=311.0+0.1\times 69.3=317.9$

Round 3 and Convergence

Residuals keep shrinking each round by 10% (ν=0.1). After T rounds:

$F_T(x)=F_0(x)+\nu {\sum }_{t=1}^T{h}_t(x)$

For x_new = sq_ft=1250 (boundary → left branch):

Why Small Learning Rate Works

With $\nu =1.0$: each tree fully corrects the residual in one step → fast convergence but memorizes training data quickly.

With $\nu =0.1$: each step corrects 10% of the residual → smooth path through the loss surface → better generalization at convergence.

Rule of thumb: $\nu \le 0.1$, compensate with higher n_estimators (500–1000+). The tradeoff is identical to gradient descent step size.

Gradient Boosting for Classification — Log-Odds View

For binary classification, GB minimizes cross-entropy. Predictions live in log-odds space; the pseudo-residuals are probability errors.

Initial Prediction

$P(\text{default})=3/8=0.375$ (3 defaults in 8 samples).

$F_0=\log \left(\frac{P(y=1)}{P(y=0)}\right)=\log \left(\frac{3}{5}\right)=-0.511$

Initial probability for all samples: ${\hat{p}}_0=\sigma (-0.511)=\frac{1}{1+e^{0.511}}=0.375$.

Round 1: Probability Residuals

$r_i^{(1)}=y_i-{\hat{p}}_0=y_i-0.375$

Best split on residuals: income ≤ 55k (same clean boundary as before).

Leaf value formula for log-loss (second-order approximation):

${\gamma }_{\text{leaf}}=\frac{{\sum }_i{r}_i}{{\sum }_i{\hat{p}}_i(1-{\hat{p}}_i)}$

• Left (samples 1,2,3): $\frac{3\times 0.625}{3\times 0.375\times 0.625}=\frac{1.875}{0.703}=2.667$
• Right (samples 4–8): $\frac{5\times (-0.375)}{5\times 0.375\times 0.625}=\frac{-1.875}{1.172}=-1.600$

Update ($\nu =0.1$):

• Left: $F_1=-0.511+0.1\times 2.667=-0.244\to \hat{p}=\sigma (-0.244)=0.439$
• Right: $F_1=-0.511+0.1\times (-1.600)=-0.671\to \hat{p}=\sigma (-0.671)=0.338$

Defaults (y=1) increase from 0.375 → 0.439 ✓. Non-defaults (y=0) decrease from 0.375 → 0.338 ✓. Each round nudges probabilities in the right direction.

sklearn Implementation

from sklearn.ensemble import GradientBoostingClassifier, GradientBoostingRegressor
from sklearn.datasets import fetch_california_housing, load_breast_cancer
from sklearn.model_selection import train_test_split
from sklearn.metrics import mean_squared_error, r2_score, roc_auc_score
import numpy as np

# Regression: California Housing
ch = fetch_california_housing()
X_r, y_r = ch.data, ch.target
Xr_tr, Xr_te, yr_tr, yr_te = train_test_split(X_r, y_r, test_size=0.2, random_state=42)

gb_reg = GradientBoostingRegressor(
    n_estimators=200, learning_rate=0.1, max_depth=3,
    subsample=0.8, random_state=42
)
gb_reg.fit(Xr_tr, yr_tr)
y_pred_r = gb_reg.predict(Xr_te)
print(f"GB Regressor: RMSE={np.sqrt(mean_squared_error(yr_te, y_pred_r)):.4f}, R²={r2_score(yr_te, y_pred_r):.4f}")

# Classification: Breast Cancer
bc = load_breast_cancer()
X_c, y_c = bc.data, bc.target
Xc_tr, Xc_te, yc_tr, yc_te = train_test_split(X_c, y_c, test_size=0.2, random_state=42, stratify=y_c)

gb_clf = GradientBoostingClassifier(
    n_estimators=200, learning_rate=0.1, max_depth=3,
    subsample=0.8, random_state=42
)
gb_clf.fit(Xc_tr, yc_tr)
y_prob_c = gb_clf.predict_proba(Xc_te)[:, 1]
print(f"GB Classifier: Test={gb_clf.score(Xc_te, yc_te):.4f}, AUC={roc_auc_score(yc_te, y_prob_c):.4f}")

GB Regressor (RMSE=0.451) beats Random Forest (RMSE=0.503) on California Housing — GB's sequential residual fitting extracts more signal given enough rounds.

Hyperparameter Sweep: n_estimators × learning_rate

configs = [(50, 0.5), (100, 0.2), (200, 0.1), (500, 0.05), (1000, 0.01)]
print(f"{'n':>6} | {'lr':>5} | {'RMSE':>8} | {'R²':>8}")
for n, lr in configs:
    gb = GradientBoostingRegressor(n_estimators=n, learning_rate=lr,
                                    max_depth=3, random_state=42)
    gb.fit(Xr_tr, yr_tr)
    pred = gb.predict(Xr_te)
    rmse = np.sqrt(mean_squared_error(yr_te, pred))
    r2   = r2_score(yr_te, pred)
    print(f"{n:>6} | {lr:>5.2f} | {rmse:>8.4f} | {r2:>8.4f}")

n=200/lr=0.1 and n=500/lr=0.05 give nearly identical results — the product n×lr governs the effective step budget. n=1000/lr=0.01 degrades because 1000 rounds at 0.01 is equivalent to 100 rounds at 0.1, not enough budget.

max_depth Sweep

print(f"{'depth':>8} | {'RMSE':>8}")
for depth in [1, 2, 3, 5, 7]:
    gb = GradientBoostingRegressor(n_estimators=200, learning_rate=0.1,
                                    max_depth=depth, random_state=42)
    gb.fit(Xr_tr, yr_tr)
    rmse = np.sqrt(mean_squared_error(yr_te, gb.predict(Xr_te)))
    print(f"{depth:>8} | {rmse:>8.4f}")

GB with max_depth=3 (8 leaves): trees capture pairwise interactions. Unlike AdaBoost (which needs depth=1), GB benefits from slightly deeper trees — but depth=5+ starts overfitting even with learning rate regularization.

Stochastic Gradient Boosting: subsample

print(f"{'subsample':>10} | {'RMSE':>8}")
for ss in [0.5, 0.6, 0.8, 1.0]:
    gb = GradientBoostingRegressor(n_estimators=200, learning_rate=0.1,
                                    max_depth=3, subsample=ss, random_state=42)
    gb.fit(Xr_tr, yr_tr)
    rmse = np.sqrt(mean_squared_error(yr_te, gb.predict(Xr_te)))
    print(f"{ss:>10} | {rmse:>8.4f}")

subsample=0.8: train each tree on a random 80% of the training data. This injects noise — different from bootstrap (with replacement) but similar effect. The randomness decorrelates consecutive trees, acting as regularization. subsample=1.0 (full dataset) is slightly worse because consecutive trees see the same data and can overfit the same patterns.

Feature Importance

importances = gb_reg.feature_importances_
for name, imp in sorted(zip(ch.feature_names, importances), key=lambda x: -x[1]):
    print(f"  {name:20s}: {imp:.4f}")

MedInc (median income) accounts for 38% of feature importance — the dominant predictor of California housing prices. Geography (Latitude + Longitude = 34%) is the second strongest signal.

Staged Prediction: RMSE Over Rounds

from sklearn.metrics import mean_squared_error
import numpy as np

train_rmse = []
test_rmse  = []

for y_pred_tr, y_pred_te in zip(
    gb_reg.staged_predict(Xr_tr),
    gb_reg.staged_predict(Xr_te)
):
    train_rmse.append(np.sqrt(mean_squared_error(yr_tr, y_pred_tr)))
    test_rmse.append(np.sqrt(mean_squared_error(yr_te, y_pred_te)))

print(f"Round   1: Train={train_rmse[0]:.4f}, Test={test_rmse[0]:.4f}")
print(f"Round  50: Train={train_rmse[49]:.4f}, Test={test_rmse[49]:.4f}")
print(f"Round 100: Train={train_rmse[99]:.4f}, Test={test_rmse[99]:.4f}")
print(f"Round 200: Train={train_rmse[199]:.4f}, Test={test_rmse[199]:.4f}")

Train RMSE decreases monotonically (every new tree reduces training error). Test RMSE plateaus around round 120–150 — adding more trees beyond the plateau risks overfitting. The shaded region marks the optimal n_estimators for this dataset.

GB vs AdaBoost vs Random Forest

Test Your Understanding

1. The leaf value formula for classification is $\gamma =\frac{\sum r_i}{\sum {\hat{p}}_i(1-{\hat{p}}_i)}$. For the left leaf (3 samples, all defaults, $r=0.625$, $\hat{p}=0.375$): verify the computation. Why does the denominator use $\hat{p}(1-\hat{p})$ instead of simply $n_{\text{leaf}}$ (as in regression)? What does this term represent in the second-order Taylor expansion of cross-entropy loss?

2. The n_estimators × learning_rate table shows n=1000/lr=0.01 gives RMSE=0.462 — worse than n=200/lr=0.1 (RMSE=0.451). The effective step budget is 1000×0.01=10 vs 200×0.1=20. Why is n=1000/lr=0.01 with budget 10 worse than n=200/lr=0.1 with budget 20, even though 1000 trees > 200 trees?

3. subsample=0.8 outperforms subsample=1.0. Each tree in Stochastic GB sees a random 80% subset — no replacement (unlike bootstrap). At round $t$, two consecutive trees share 80%×80%=64% of the data in expectation. How does this compare to Random Forest's bootstrap overlap (~63%), and what different regularization effect does each achieve?

4. GB with max_depth=3 gives RMSE=0.4512, while max_depth=1 (stumps, like AdaBoost) gives RMSE=0.5612. But AdaBoost with max_depth=1 achieves test accuracy comparable to GB on classification tasks. Why does GB benefit more from max_depth=3 than AdaBoost does — even though both are boosting methods?

5. The staged prediction curve shows train RMSE monotonically decreasing while test RMSE plateaus. In theory, once training error reaches a minimum, adding more trees cannot decrease test error — but it also shouldn't increase it (the new trees only add to the existing sum). What breaks this reasoning and causes test error to eventually increase with too many rounds?