XGBoost · LightGBM · CatBoost — production gradient boosting¶
The industrial descendants of Friedman's GBM (Chen & Guestrin, 2016)¶
The previous notebook built gradient boosting from scratch and matched scikit-learn's GradientBoosting. That algorithm is correct but slow and lightly regularised. The production boosters keep Friedman's additive, gradient-driven loop and rebuild everything around it for speed, regularisation, and scale — they are the default winning model on tabular data and a staple of systematic-trading and credit pipelines. Three dominate, each with a distinct idea:
- XGBoost (Chen & Guestrin, 2016) — second-order (Newton) boosting: it uses both the gradient and the curvature (Hessian) of the loss to choose splits and leaf values, adds explicit L1/L2 regularisation on the leaf weights, and is sparsity-aware. The model that popularised the family.
- LightGBM (Microsoft) — histogram split-finding (bin the features) with leaf-wise (best-first) growth and GOSS gradient sampling. Built for speed on large data.
- CatBoost (Yandex) — ordered boosting (a permutation trick that removes the subtle target leakage ordinary boosting suffers) and native categorical handling via ordered target statistics; grows symmetric (oblivious) trees.
We benchmark all three head-to-head on the same Taiwan credit-default and California-housing data, add the production essentials — early stopping and SHAP interpretability — and place them on the running out-of-sample scoreboard. ROC-AUC is defined in the CART notebook (0.5 = chance, 1 = perfect).
import numpy as np, pandas as pd, matplotlib.pyplot as plt, time, warnings
warnings.filterwarnings("ignore")
from sklearn.model_selection import train_test_split
from sklearn.metrics import roc_auc_score, mean_squared_error
import xgboost as xgb, lightgbm as lgb, catboost as cb, shap
BLUE="#2b6cb0"; RED="#c53030"; GREEN="#2f855a"; ORANGE="#dd6b20"; GREY="#a0aec0"; PURP="#6b46c1"
d=pd.read_csv("credit_default.csv"); feat=[c for c in d.columns if c!="default"]
X=d[feat].to_numpy(float); y=d["default"].to_numpy(int)
Xtr,Xte,ytr,yte=train_test_split(X,y,test_size=0.3,random_state=0,stratify=y)
h=pd.read_csv("cali_housing.csv"); hf=[c for c in h.columns if c!="MedHouseVal"]
Xh=h[hf].to_numpy(float); yh=h["MedHouseVal"].to_numpy(float); Xhtr,Xhte,yhtr,yhte=train_test_split(Xh,yh,test_size=0.3,random_state=0)
print(f"xgboost {xgb.__version__} | lightgbm {lgb.__version__} | catboost {cb.__version__} | shap {shap.__version__}")
print(f"credit default: {len(ytr):,}/{len(yte):,}; California: {len(yhtr):,}/{len(yhte):,}")
xgboost 3.3.0 | lightgbm 4.7.0 | catboost 1.2.10 | shap 0.52.0 credit default: 21,000/9,000; California: 14,448/6,192
1. Head-to-head benchmark — accuracy and speed¶
The same hyperparameters across all three (300 trees, learning rate 0.05, depth 4) on both tasks, timing the fit. They should land within a whisker of each other on accuracy — comfortably above the single tree and the parametric baseline — differing mostly in speed, which is where they were engineered to compete.
def bench(make, Xtr, ytr, Xte, score):
t=time.time(); m=make().fit(Xtr,ytr); ft=time.time()-t; return score(m), ft
clf={"XGBoost":lambda:xgb.XGBClassifier(n_estimators=300,learning_rate=0.05,max_depth=4,eval_metric="auc",verbosity=0,n_jobs=-1),
"LightGBM":lambda:lgb.LGBMClassifier(n_estimators=300,learning_rate=0.05,max_depth=4,verbose=-1,n_jobs=-1),
"CatBoost":lambda:cb.CatBoostClassifier(n_estimators=300,learning_rate=0.05,depth=4,verbose=0)}
reg={"XGBoost":lambda:xgb.XGBRegressor(n_estimators=300,learning_rate=0.05,max_depth=4,verbosity=0,n_jobs=-1),
"LightGBM":lambda:lgb.LGBMRegressor(n_estimators=300,learning_rate=0.05,max_depth=4,verbose=-1,n_jobs=-1),
"CatBoost":lambda:cb.CatBoostRegressor(n_estimators=300,learning_rate=0.05,depth=4,verbose=0)}
rows=[]
for nm in clf:
a,ta=bench(clf[nm],Xtr,ytr,Xte,lambda m:roc_auc_score(yte,m.predict_proba(Xte)[:,1]))
r,tr=bench(reg[nm],Xhtr,yhtr,Xhte,lambda m:mean_squared_error(yhte,m.predict(Xhte))**.5)
rows.append([nm,a,ta,r,tr])
tab=pd.DataFrame(rows,columns=["library","AUC","AUC_fit_s","RMSE","RMSE_fit_s"]).set_index("library")
print(tab.round(3).to_string())
fig,ax=plt.subplots(1,3,figsize=(16,4)); nm=list(tab.index); cols=[BLUE,GREEN,ORANGE]
ax[0].bar(nm,tab["AUC"],color=cols); ax[0].set_ylim(0.75,0.785); ax[0].set_title("Credit default — test AUC"); [ax[0].text(i,v+0.0005,f"{v:.3f}",ha="center",fontsize=8) for i,v in enumerate(tab["AUC"])]
ax[1].bar(nm,tab["RMSE"],color=cols); ax[1].set_title("California — test RMSE"); [ax[1].text(i,v+0.002,f"{v:.3f}",ha="center",fontsize=8) for i,v in enumerate(tab["RMSE"])]
ax[2].bar(nm,tab["AUC_fit_s"],color=cols); ax[2].set_title("Fit time — classification (s)"); ax[2].set_ylabel("seconds")
for a in ax: plt.setp(a.get_xticklabels(),rotation=10)
plt.tight_layout(); plt.show()
print(f"On CLASSIFICATION they are near-identical: AUC {tab['AUC'].min():.3f}-{tab['AUC'].max():.3f}, a spread of {tab['AUC'].max()-tab['AUC'].min():.3f}.")
print(f"On REGRESSION they are not quite: RMSE {tab['RMSE'].min():.3f}-{tab['RMSE'].max():.3f}, with CatBoost {100*(tab['RMSE'].max()/tab['RMSE'].min()-1):.0f}% behind the other two.")
print(f"Fit time separates them far more than accuracy does: {tab['AUC_fit_s'].max()/tab['AUC_fit_s'].min():.1f}x between fastest and slowest on the same task.")
print("So the libraries compete on ENGINEERING -- speed, scale, and the regularisation, categorical and")
print("early-stopping machinery below -- rather than on a headline accuracy gap on clean tabular data.")
AUC AUC_fit_s RMSE RMSE_fit_s library XGBoost 0.774 0.215 0.494 0.105 LightGBM 0.773 0.142 0.494 0.097 CatBoost 0.776 1.204 0.521 0.449
On CLASSIFICATION they are near-identical: AUC 0.773-0.776, a spread of 0.003. On REGRESSION they are not quite: RMSE 0.494-0.521, with CatBoost 6% behind the other two. Fit time separates them far more than accuracy does: 8.5x between fastest and slowest on the same task. So the libraries compete on ENGINEERING -- speed, scale, and the regularisation, categorical and early-stopping machinery below -- rather than on a headline accuracy gap on clean tabular data.
2. Early stopping — the production essential¶
Boosting overfits if you add too many trees (the previous notebook's peak-then-decline). In practice you never guess the count: you hold out a validation set, watch its metric each round, and stop when it stops improving. All three support it in one argument; here XGBoost, plotting the train-vs-validation AUC and the chosen stopping point.
Xt2,Xv,yt2,yv=train_test_split(Xtr,ytr,test_size=0.25,random_state=1,stratify=ytr)
m=xgb.XGBClassifier(n_estimators=1000,learning_rate=0.05,max_depth=4,eval_metric="auc",
early_stopping_rounds=30,verbosity=0,n_jobs=-1).fit(Xt2,yt2,eval_set=[(Xt2,yt2),(Xv,yv)],verbose=False)
ev=m.evals_result(); tr_auc=ev["validation_0"]["auc"]; va_auc=ev["validation_1"]["auc"]
fig,ax=plt.subplots(figsize=(7.5,4))
ax.plot(tr_auc,color=BLUE,label="train AUC"); ax.plot(va_auc,color=RED,label="validation AUC")
ax.axvline(m.best_iteration,color=GREEN,ls="--",label=f"early stop @ {m.best_iteration} trees")
ax.set_xlabel("boosting round"); ax.set_ylabel("AUC"); ax.set_title("Early stopping halts at the validation peak"); ax.legend(frameon=False)
plt.tight_layout(); plt.show()
print(f"Training AUC keeps rising toward 1.0; validation AUC peaks near round {m.best_iteration} then flattens/declines -- so we stop there.")
print(f"Held-out test AUC at the chosen stopping point: {roc_auc_score(yte,m.predict_proba(Xte)[:,1]):.4f}. No manual tree-count tuning needed.")
Training AUC keeps rising toward 1.0; validation AUC peaks near round 98 then flattens/declines -- so we stop there. Held-out test AUC at the chosen stopping point: 0.7724. No manual tree-count tuning needed.
3. Interpretability — SHAP¶
A boosted ensemble of hundreds of trees is opaque; SHAP (SHapley Additive exPlanations) makes it legible. Rooted in cooperative game theory, a SHAP value splits each prediction into per-feature contributions that sum to the prediction, and for trees they are computed exactly and fast (TreeSHAP). The beeswarm shows every feature's impact across clients (colour = feature value); the dependence plot shows how the top feature moves risk.
mx=xgb.XGBClassifier(n_estimators=300,learning_rate=0.05,max_depth=4,verbosity=0,n_jobs=-1).fit(Xtr,ytr)
expl=shap.TreeExplainer(mx); Xs=pd.DataFrame(Xte[:2000],columns=feat); sv=expl.shap_values(Xs)
shap.summary_plot(sv,Xs,show=False,max_display=10,plot_size=(9,5)); plt.title("SHAP: which features drive default risk, and which way"); plt.tight_layout(); plt.show()
fig,ax=plt.subplots(figsize=(7.5,4.4))
shap.dependence_plot("PAY_1",sv,Xs,interaction_index=None,ax=ax,show=False)
ax.set_title("SHAP dependence: recent delinquency (PAY_1) pushes risk up"); plt.tight_layout(); plt.show()
print("SHAP confirms the whole tree family's story quantitatively: recent repayment status (PAY_1) dominates, and being in")
print("arrears pushes the predicted log-odds of default sharply UP. Unlike Gini importance, SHAP is per-prediction, signed,")
print("and consistent -- the modern standard for explaining boosted models to risk committees and regulators.")
SHAP confirms the whole tree family's story quantitatively: recent repayment status (PAY_1) dominates, and being in arrears pushes the predicted log-odds of default sharply UP. Unlike Gini importance, SHAP is per-prediction, signed, and consistent -- the modern standard for explaining boosted models to risk committees and regulators.
Reading the rest of the beeswarm¶
The plot above ranks ten features; the paragraph under it explains one. That is the usual way SHAP
gets reported and it wastes most of the figure, so the cell below reads off the top of the ranking
with a direction for each — the correlation between a feature's value and its own SHAP value —
and then takes the clearest of the "lowers risk" features, PAY_AMT2, back to the raw data.
# Read the beeswarm rather than leaving it decorative: rank by mean |SHAP|, and for each of
# the top features report the DIRECTION -- the correlation between the feature's value and its
# own SHAP value. Positive means high values push risk up; negative means they push it down.
_imp = np.abs(sv).mean(0)
_order = np.argsort(-_imp)[:8]
print(f"{'feature':12s}{'mean |SHAP|':>13s}{'direction':>11s} reading")
for _j in _order:
_xj = Xs.iloc[:, _j].to_numpy(); _sj = sv[:, _j]
_r = np.corrcoef(_xj, _sj)[0, 1] if _xj.std() > 0 else np.nan
_way = "raise" if _r > 0 else "lower"
print(f"{feat[_j]:12s}{_imp[_j]:13.4f}{_r:+11.2f} high values {_way} predicted risk")
# The payment-amount features are the ones the beeswarm shows and the text above skips.
print("\nPAY_AMT2 is the amount actually paid in the second-most-recent month, and it is the")
print("clearest of the 'lowers risk' features. Its behaviour is worth stating on the raw data,")
print("because the SHAP direction alone does not say where the signal sits:\n")
_v = d["PAY_AMT2"]
_q = pd.qcut(_v.rank(method="first"), 5, labels=["Q1 (lowest)", "Q2", "Q3", "Q4", "Q5 (highest)"])
_qtab = d.groupby(_q, observed=True).agg(median_paid=("PAY_AMT2", "median"),
default_rate=("default", "mean"))
print(_qtab.round(3).to_string())
_zero = d[_v == 0]; _paid = d[_v > 0]
print(f"\npaid exactly zero: {len(_zero):,} clients, default rate {_zero['default'].mean():.3f}")
print(f"paid anything : {len(_paid):,} clients, default rate {_paid['default'].mean():.3f}")
print(f"(base rate {d['default'].mean():.3f}; correlation of PAY_AMT2 with BILL_AMT2 "
f"{np.corrcoef(_v, d['BILL_AMT2'])[0,1]:.3f}, with LIMIT_BAL {np.corrcoef(_v, d['LIMIT_BAL'])[0,1]:.3f})")
print("\nThree things follow. The relationship is monotone -- more paid, less default -- so the")
print("negative SHAP direction is the right reading. But the signal is concentrated at ZERO:")
print(f"{100*(_v==0).mean():.0f}% of clients paid nothing that month and they default at "
f"{_zero['default'].mean():.1%} against {_paid['default'].mean():.1%}")
print("for everyone else, so the red points at negative SHAP are largely the mirror of a large")
print("positive push for the non-payers. And it is not merely proxying for account size: the")
print("correlation with the bill and with the credit limit is weak, so the amount paid carries")
print("its own information. Because SHAP is conditional on the rest of the model, this is what")
print("PAY_AMT2 adds ON TOP of the delinquency flags -- which is why it sits in the second tier")
print("rather than at the top.")
feature mean |SHAP| direction reading
PAY_1 0.5343 +0.74 high values raise predicted risk
LIMIT_BAL 0.2061 -0.83 high values lower predicted risk
BILL_AMT1 0.1672 +0.07 high values raise predicted risk
PAY_AMT3 0.1239 -0.29 high values lower predicted risk
PAY_AMT1 0.1022 -0.49 high values lower predicted risk
PAY_AMT2 0.0977 -0.64 high values lower predicted risk
PAY_2 0.0860 +0.79 high values raise predicted risk
PAY_3 0.0715 +0.75 high values raise predicted risk
PAY_AMT2 is the amount actually paid in the second-most-recent month, and it is the
clearest of the 'lowers risk' features. Its behaviour is worth stating on the raw data,
because the SHAP direction alone does not say where the signal sits:
median_paid default_rate
PAY_AMT2
Q1 (lowest) 0.0 0.327
Q2 1165.0 0.238
Q3 2009.0 0.217
Q4 4045.5 0.192
Q5 (highest) 10405.5 0.132
paid exactly zero: 5,396 clients, default rate 0.333
paid anything : 24,604 clients, default rate 0.197
(base rate 0.221; correlation of PAY_AMT2 with BILL_AMT2 0.101, with LIMIT_BAL 0.178)
Three things follow. The relationship is monotone -- more paid, less default -- so the
negative SHAP direction is the right reading. But the signal is concentrated at ZERO:
18% of clients paid nothing that month and they default at 33.3% against 19.7%
for everyone else, so the red points at negative SHAP are largely the mirror of a large
positive push for the non-payers. And it is not merely proxying for account size: the
correlation with the bill and with the credit limit is weak, so the amount paid carries
its own information. Because SHAP is conditional on the rest of the model, this is what
PAY_AMT2 adds ON TOP of the delinquency flags -- which is why it sits in the second tier
rather than at the top.
4. Calibration — proportions vs predictions¶
AUC measures ranking; it says nothing about whether a predicted 0.30 actually defaults 30% of the time. For credit and pricing that calibration is what a decision uses, so we check it directly. Left — credit: sort clients into ten equal-count bins by predicted probability and plot the observed default proportion in each bin against the mean predicted probability; points on the 45° line are perfectly calibrated. The Brier score (mean squared error of the probabilities) summarises it in one number. Right — California: the regression analogue — predicted value vs actual, with decile-binned means over the scatter. Together these are the "proportions vs predictions" view of both data sets.
from sklearn.calibration import calibration_curve
from sklearn.metrics import brier_score_loss
from sklearn.linear_model import LogisticRegression
from sklearn.preprocessing import StandardScaler
from sklearn.pipeline import make_pipeline
pmods={"XGBoost":xgb.XGBClassifier(n_estimators=300,learning_rate=0.05,max_depth=4,verbosity=0,n_jobs=-1).fit(Xtr,ytr),
"CatBoost":cb.CatBoostClassifier(n_estimators=300,learning_rate=0.05,depth=4,verbose=0).fit(Xtr,ytr),
"logistic":make_pipeline(StandardScaler(),LogisticRegression(max_iter=2000)).fit(Xtr,ytr)}
fig,ax=plt.subplots(1,2,figsize=(13.5,5))
pmax=max(m.predict_proba(Xte)[:,1].max() for m in pmods.values())
ax[0].plot([0,pmax],[0,pmax],"k--",lw=1,label="perfect calibration")
for (nm,m),c in zip(pmods.items(),[BLUE,GREEN,ORANGE]):
p=m.predict_proba(Xte)[:,1]; pt,pp=calibration_curve(yte,p,n_bins=10,strategy="quantile")
ax[0].plot(pp,pt,"o-",color=c,lw=2,ms=5,label=f"{nm} (Brier {brier_score_loss(yte,p):.3f})")
ax[0].set_xlabel("predicted P(default) [decile bins]"); ax[0].set_ylabel("observed default proportion")
ax[0].set_title("Credit default — reliability curve"); ax[0].legend(frameon=False,fontsize=9)
# California: predicted vs actual
mh=xgb.XGBRegressor(n_estimators=300,learning_rate=0.05,max_depth=4,verbosity=0,n_jobs=-1).fit(Xhtr,yhtr); ph=mh.predict(Xhte)
ax[1].scatter(ph,yhte,s=5,alpha=.12,color=BLUE)
b=np.quantile(ph,np.linspace(0,1,11)); idx=np.clip(np.digitize(ph,b[1:-1]),0,9)
mp=[ph[idx==k].mean() for k in range(10)]; ma=[yhte[idx==k].mean() for k in range(10)]
lim=[min(ph.min(),yhte.min()),max(ph.max(),yhte.max())]
ax[1].plot(lim,lim,"k--",lw=1,label="perfect"); ax[1].plot(mp,ma,"o-",color=RED,lw=2,ms=6,label="decile means")
ax[1].set_xlabel("predicted value ($100k)"); ax[1].set_ylabel("actual value ($100k)")
ax[1].set_title(f"California — predicted vs actual (RMSE {mean_squared_error(yhte,ph)**.5:.3f})"); ax[1].legend(frameon=False)
plt.tight_layout(); plt.show()
print("Credit: the boosters hug the 45-degree line (low Brier) -- a predicted 30% really does default ~30% of the time; the")
print("linear logistic is slightly less calibrated in the high-risk deciles. California: decile means track the diagonal, with")
print("the visible flattening at the top from the data set's known $500k price CAP, which no model can predict past.")
Credit: the boosters hug the 45-degree line (low Brier) -- a predicted 30% really does default ~30% of the time; the linear logistic is slightly less calibrated in the high-risk deciles. California: decile means track the diagonal, with the visible flattening at the top from the data set's known $500k price CAP, which no model can predict past.
5. The out-of-sample scoreboard — the whole subsection¶
Placing the production boosters beside every earlier method on the same held-out test sets. The ensembles cluster at the top; the modern boosters match the random forest and the from-scratch GBM, winning on engineering (speed, early stopping, regularisation, categoricals) rather than a large accuracy gap on clean tabular data.
from sklearn.linear_model import LogisticRegression, LinearRegression
from sklearn.preprocessing import StandardScaler
from sklearn.pipeline import make_pipeline
from sklearn.tree import DecisionTreeClassifier, DecisionTreeRegressor
from sklearn.ensemble import (RandomForestClassifier, RandomForestRegressor,
GradientBoostingClassifier, GradientBoostingRegressor)
cl={}; rg={}
cl["logistic / linear"]=roc_auc_score(yte,make_pipeline(StandardScaler(),LogisticRegression(max_iter=2000)).fit(Xtr,ytr).predict_proba(Xte)[:,1])
rg["logistic / linear"]=mean_squared_error(yhte,LinearRegression().fit(Xhtr,yhtr).predict(Xhte))**.5
cl["single tree"]=roc_auc_score(yte,DecisionTreeClassifier(min_samples_leaf=50,random_state=0).fit(Xtr,ytr).predict_proba(Xte)[:,1])
rg["single tree"]=mean_squared_error(yhte,DecisionTreeRegressor(max_depth=8,random_state=0).fit(Xhtr,yhtr).predict(Xhte))**.5
cl["random forest"]=roc_auc_score(yte,RandomForestClassifier(n_estimators=400,min_samples_leaf=5,random_state=0,n_jobs=-1).fit(Xtr,ytr).predict_proba(Xte)[:,1])
rg["random forest"]=mean_squared_error(yhte,RandomForestRegressor(n_estimators=300,min_samples_leaf=3,random_state=0,n_jobs=-1).fit(Xhtr,yhtr).predict(Xhte))**.5
cl["sklearn GBM"]=roc_auc_score(yte,GradientBoostingClassifier(n_estimators=300,learning_rate=0.05,max_depth=3,random_state=0).fit(Xtr,ytr).predict_proba(Xte)[:,1])
rg["sklearn GBM"]=mean_squared_error(yhte,GradientBoostingRegressor(n_estimators=300,learning_rate=0.05,max_depth=3,random_state=0).fit(Xhtr,yhtr).predict(Xhte))**.5
for nm in ["XGBoost","LightGBM","CatBoost"]:
cl[nm]=tab.loc[nm,"AUC"]; rg[nm]=tab.loc[nm,"RMSE"]
sc=pd.DataFrame({"classification AUC (higher=better)":cl,"regression RMSE (lower=better)":rg}).round(3)
print(sc.to_string())
fig,ax=plt.subplots(1,2,figsize=(14,4.6)); names=list(cl)
colr=[GREY,GREY,GREEN,BLUE,BLUE,BLUE,BLUE]; colr=[ORANGE,GREY,GREEN,PURP,BLUE,GREEN,ORANGE]
ax[0].barh(names,[cl[k] for k in names],color=colr); ax[0].set_xlim(0.5,0.8); ax[0].set_title("Credit default — test AUC (higher better)"); ax[0].invert_yaxis()
ax[1].barh(names,[rg[k] for k in names],color=colr); ax[1].set_title("California — test RMSE (lower better)"); ax[1].invert_yaxis()
plt.tight_layout(); plt.show()
ens = sc.iloc[2:]
print("The single tree and the linear baseline trail; every ensemble -- forest, from-scratch GBM, and the three production")
print(f"boosters -- sits at the top, within {ens.iloc[:,0].max()-ens.iloc[:,0].min():.3f} AUC and {ens.iloc[:,1].max()-ens.iloc[:,1].min():.3f} RMSE of each other. On clean tabular data the choice is about")
print("SPEED, TOOLING, and UNCERTAINTY (BART), not a headline accuracy race -- exactly the judgement a quant workflow makes.")
classification AUC (higher=better) regression RMSE (lower=better) logistic / linear 0.715 0.737 single tree 0.737 0.666 random forest 0.775 0.523 sklearn GBM 0.773 0.522 XGBoost 0.774 0.494 LightGBM 0.773 0.494 CatBoost 0.776 0.521
The single tree and the linear baseline trail; every ensemble -- forest, from-scratch GBM, and the three production boosters -- sits at the top, within 0.003 AUC and 0.029 RMSE of each other. On clean tabular data the choice is about SPEED, TOOLING, and UNCERTAINTY (BART), not a headline accuracy race -- exactly the judgement a quant workflow makes.
6. Summary¶
XGBoost, LightGBM, and CatBoost are Friedman's gradient boosting, re-engineered for production. They kept the additive-gradient loop and added second-order (Newton) boosting and explicit regularisation (XGBoost), histogram + leaf-wise speed (LightGBM), and ordered boosting + native categoricals (CatBoost). On our two clean tabular tasks all three land within a hair of each other and of the random forest (AUC ≈ 0.77, RMSE ≈ 0.51) — the honest lesson being that on tidy data the accuracy differences are small, and the real dividends are speed, early stopping, regularisation, and categorical handling. We added the two production must-haves: early stopping on a validation set (no manual tree-count tuning) and SHAP (signed, per-prediction, exact TreeSHAP explanations) — which confirmed the whole family's finding that recent delinquency (PAY_1) drives default.
This closes the Tree Ensembles subsection: from a single CART built by hand, through the variance-cutting Random Forest, the Bayesian BART with its credible intervals, Friedman's Gradient Boosting from scratch, to these industrial boosters — the models that win most tabular competitions and underpin credit and systematic-trading pipelines. The R companion fits the same models with the xgboost package. Which to reach for: LightGBM/XGBoost for speed and scale, CatBoost when categoricals abound, BART when calibrated uncertainty matters, and a single tree when you must explain the whole model on one page. Choosing among them well — not just fitting them — is the skill the subsection has been building toward.