Random Forest & Bagging

Why this session exists

A single tree can change dramatically when a few training rows change. Bagging asks whether averaging many such fits can stabilize the prediction. A random forest adds another source of diversity: each split considers a random subset of features rather than letting every tree repeatedly choose the same dominant predictor.

The useful result is not “more trees never make mistakes.” We need to understand the covariance that averaging cannot remove, what out-of-bag evaluation actually excludes, and why feature importance remains a property of a model and dataset. We will build a small forest with an inspectable per-node splitter, compute OOB probabilities, and compare with scikit-learn without claiming identical trees or a guaranteed benchmark score.

Prerequisites

Review decision-tree impurity and recursion, bias and variance, sampling with replacement, and binary probability averaging. The standalone program below supports finite numeric inputs and binary labels 0/1. Its leaves store the positive-class fraction, so even a bootstrap sample containing only one class has an unambiguous probability representation.

(a) Intuition · 5 min

Imagine a committee whose members examine overlapping, randomly selected evidence. If every member makes the same error, voting adds no information. If their prediction errors vary, averaging can reduce the variable component. This analogy is about error dependence, not an assertion that a majority of any set of experts is automatically better than the best one.

Averaging preserves the mean of identically distributed member predictions, but a bagged base learner is fitted on a bootstrap sample rather than the original full dataset. That can change its bias relative to a single full-data fit. Feature restrictions can change both bias and variance too. “Forest only changes variance” is a helpful first intuition, not a universal decomposition theorem.

The three ideas to own before we build

• Bootstrap rows: draw $n$ row indices with replacement from $n$ training rows. The sample has $n$ draws, not $0.632n$ draws; repeated indices account for the difference.
• Sample features per split: draw a fresh candidate subset at each node. Sampling a single subset once per tree is a random-subspace variant, not the same algorithm.
• Aggregate consistently: regression commonly averages leaf predictions. Scikit-learn's classifier averages class probabilities then takes the largest average; that can differ from majority voting over each tree's hard label.

A brief history

The ensemble guide distinguishes pasting, bagging, random subspaces and random patches by whether rows and/or features are sampled and whether replacement is used. Its forest section cites Breiman's Random Forests (2001), while the Extra Trees section describes additional randomization of thresholds.

That taxonomy is the history we need for implementation: resampling rows, restricting candidate features and randomizing thresholds are separate mechanisms. A former per-tree feature-subset example should not be called a faithful per-split forest merely because both contain randomness. There is no need to invent a competition timeline or a company's internal model-replacement story to motivate the distinction.

The bagging pipeline

A learner can parallelize trees because one tree does not depend on another tree's predictions. Their statistical errors are not therefore independent: they share an underlying training sample and correlated features. Keep computational independence distinct from the covariance assumption used in the variance formula.

Why bootstrap covers ~63 %

For one fixed row, the chance of avoiding it on a single draw is $1-1/n$. Across $n$ independent draws with replacement, exclusion probability is $(1-1/n)^n$. By linearity of expectation,

As $n\to\infty$, exclusion approaches $e^{-1}\approx0.3679$ and inclusion approaches 0.6321. For finite $n$ the exact expression differs, and any one bootstrap has a random distinct count. At $n=1$, the only row is always included and there is no OOB prediction.

For $B$ independent bootstrap index draws, the number of trees excluding a particular row has expected value $B(1-1/n)^n$. At small tree counts, some rows may have no OOB voters. Report the covered fraction instead of treating zero-filled prediction slots as valid estimates.

Bagging vs Random Forest — one tweak makes a big difference

The sqrt(d) candidate count is the scikit-learn classifier default. The 1.7 regressor default max_features=1.0 considers all features, not the older rule of $d/3$. Defaults are versioned choices, not derivations of an optimum. Reducing candidate features may lower inter-tree correlation but can deprive trees of useful splits. Depth, leaf size, bootstrap sample size and feature count all deserve attention; no theorem makes max_features the only real hyperparameter.

Extra Trees do not simply choose a completely arbitrary final split. They randomize candidate thresholds and then compare candidate quality. Additional randomness can reduce correlation while increasing approximation error. Any speed or accuracy advantage depends on data and implementation, not a fixed multiplier.

The math of variance reduction

Let $Z_b$ be a tree's prediction at a fixed input, with variance $\sigma^2$ and equal pairwise correlation $\rho$. Then

There are $B$ variance terms and $B(B-1)$ ordered covariance terms. This is a covariance identity, not an invocation of the central limit theorem. Independence implies $\rho=0$, giving $\sigma^2/B$. For $\rho=1$, the variance is unchanged. With $\sigma^2=1$, $\rho=0.2$ and $B=100$, it is 0.208, not 0.01.

The equal-correlation/equal-variance model is illustrative. The variance is of numeric predictions over repeated training/randomization, not finite-test classification error after thresholding. Adding a tree can flip one finite-test label the wrong way, so accuracy need not improve monotonically. As the number of randomized trees grows, Monte Carlo variation can stabilize while distribution bias, leakage and a poor representation remain. Select an adequate count under memory and latency limits, not a “cannot overfit” slogan.

Out-of-bag evaluation — useful, not a free replacement for every split

For each training row, average predictions only from trees whose bootstrap indices excluded that row. Compare those averages with labels for covered rows. This reuses fitting work efficiently but is not literally $k$-fold CV: different rows have overlapping, differently sized sets of OOB trees.

OOB isolation applies to a tree's row fit, not automatically to every upstream transformation. If a supervised encoder or feature selection step used all training labels before bootstrapping, OOB rows influenced that preprocessing. OOB also does not enforce temporal or group independence. A group/time holdout and final untouched test set remain necessary when those match the deployment question. Repeatedly tuning on OOB results can overfit that estimate just as repeated validation search can.

Feature importance — the free bonus (with caveats)

MDI sums sample-weighted training impurity reductions attributed to splits on a feature and normalizes across features. It is cheap, but high-cardinality features have more opportunities for lucky reductions. The importance is conditional on the fitted model, and correlated inputs can share or substitute for one another's apparent contribution.

Permutation importance measures a chosen score drop after shuffling a feature on evaluation data. It can be negative in a finite sample. It answers “how much does this model rely on this column under this perturbation?” rather than “what is the causal effect of this variable?” A weak model can assign low importance to genuinely useful information it failed to learn.

If two columns duplicate each other, shuffling only one leaves the other available. Both individual importances may look small. n_repeats repeats shuffles to measure variability; it does not group correlated features. To measure joint reliance, explicitly permute the whole group using the same row permutation, preserving relationships within that group while disrupting its relationship with other columns and the target. Even that can create out-of-distribution combinations. SHAP likewise requires a background/reference and assumptions about dependence; it is not an automatic ground truth.

(c) Hands-on · 25 min

This standalone small forest reuses the previous lesson's Gini/midpoint logic, adds a fresh feature draw inside recursion, and records in-bag/OOB indices. The brute-force splitter is deliberately bounded to 120 generated rows and 24 trees. Scikit-learn supplies a separate larger comparison on the same data, not the hidden implementation of the scratch forest.

import numpy as np
from sklearn.datasets import make_classification
from sklearn.model_selection import train_test_split
from sklearn.ensemble import RandomForestClassifier
from sklearn.tree import DecisionTreeClassifier
from sklearn.inspection import permutation_importance
 
 
def gini(y):
    p = y.mean()
    return 2*p*(1-p)
 
 
def grow(X, y, rng, max_features, depth=0, max_depth=6, min_leaf=2):
    node = {'p': float(y.mean()), 'n': len(y)}
    if depth >= max_depth or len(y) < 2*min_leaf or np.unique(y).size == 1:
        return node
    # Crucial: this executes AGAIN at every recursive node.
    features = rng.choice(X.shape[1], max_features, replace=False)
    node['candidates'] = features.tolist()
    best = None
    for j in features:
        vals = np.unique(X[:, j])
        for lo, hi in zip(vals[:-1], vals[1:]):
            t = lo/2 + hi/2  # avoids overflow from lo + hi
            if t <= lo:
                t = hi  # strict < preserves the intended partition
            mask = X[:, j] < t
            nl = int(mask.sum())
            if min(nl, len(y)-nl) < min_leaf:
                continue
            gain = gini(y)-(nl*gini(y[mask])+(len(y)-nl)*gini(y[~mask]))/len(y)
            if best is None or gain > best[0]:
                best = (gain, int(j), float(t), mask)
    if best is None or best[0] <= 1e-12:
        return node
    _, j, t, mask = best
    node.update(j=j, t=t,
        left=grow(X[mask], y[mask], rng, max_features, depth+1, max_depth, min_leaf),
        right=grow(X[~mask], y[~mask], rng, max_features, depth+1, max_depth, min_leaf))
    return node
 
 
def tree_prob(node, X):
    out = []
    for x in X:
        current = node
        while 'j' in current:
            current = current['left'] if x[current['j']] < current['t'] else current['right']
        out.append(current['p'])
    return np.asarray(out)
 
 
class Forest:
    def __init__(self, trees=24, max_features=2, seed=42):
        self.count, self.max_features, self.seed = trees, max_features, seed
 
    def fit(self, X, y):
        X, y = np.asarray(X, float), np.asarray(y)
        if X.ndim != 2 or len(X) == 0 or y.shape != (len(X),):
            raise ValueError('nonempty numeric matrix and matching labels required')
        if not np.isfinite(X).all() or not np.isin(y, [0, 1]).all():
            raise ValueError('finite numeric features and binary 0/1 labels required')
        if self.count < 1 or not 1 <= self.max_features <= X.shape[1]:
            raise ValueError('invalid tree or feature count')
        rng = np.random.default_rng(self.seed)
        self.trees, self.inbag, self.oob = [], [], []  # repeated fit resets state
        totals, counts = np.zeros(len(y)), np.zeros(len(y), dtype=int)
        for _ in range(self.count):
            idx = rng.integers(0, len(y), len(y))
            excluded = np.flatnonzero(np.bincount(idx, minlength=len(y)) == 0)
            tree = grow(X[idx], y[idx], rng, self.max_features)
            self.trees.append(tree); self.inbag.append(idx); self.oob.append(excluded)
            totals[excluded] += tree_prob(tree, X[excluded])
            counts[excluded] += 1
        self.oob_counts_ = counts
        self.oob_probability_ = np.divide(totals, counts,
            out=np.full(len(y), np.nan), where=counts>0)
        return self
 
    def predict_proba(self, X):
        positive = np.mean([tree_prob(t, X) for t in self.trees], axis=0)
        return np.column_stack([1-positive, positive])
 
    def predict(self, X):
        return self.predict_proba(X).argmax(axis=1)  # ties go to class 0
 
 
# A rounded midpoint must not erase a valid split between adjacent floats.
adjacent_X = np.array([1., np.nextafter(1., np.inf)]).reshape(-1, 1)
adjacent_y = np.array([0, 1])
adjacent_tree = grow(adjacent_X, adjacent_y, np.random.default_rng(0), 1, min_leaf=1)
assert np.allclose(tree_prob(adjacent_tree, adjacent_X), adjacent_y)
 
X, y = make_classification(n_samples=120, n_features=6, n_informative=3,
    n_redundant=1, flip_y=.04, random_state=42)
Xd, Xt, yd, yt = train_test_split(X, y, test_size=.25, stratify=y, random_state=43)
f = Forest().fit(Xd, yd)
covered = f.oob_counts_ > 0
print('OOB coverage', covered.mean(), 'OOB accuracy',
      np.mean((f.oob_probability_[covered] > .5) == yd[covered]))
print('scratch test', np.mean(f.predict(Xt)==yt))
for inside, outside in zip(f.inbag, f.oob):
    assert not np.intersect1d(inside, outside).size
assert np.allclose(f.predict_proba(Xt).sum(axis=1), 1)
pred = f.predict_proba(Xt).copy()
f.fit(Xd, yd)
assert len(f.trees) == 24 and np.allclose(pred, f.predict_proba(Xt))
 
# One-class bootstraps and absent OOB predictions are handled explicitly.
tiny = Forest(trees=2, max_features=1).fit(np.array([[0.]]), np.array([1]))
assert np.allclose(tiny.predict_proba([[2.]]), [[0., 1.]])
assert tiny.oob_counts_[0] == 0 and np.isnan(tiny.oob_probability_[0])
 
base = DecisionTreeClassifier(random_state=42).fit(Xd, yd)
reference = RandomForestClassifier(n_estimators=64, max_features=2,
    min_samples_leaf=2, max_depth=6, oob_score=True, n_jobs=1, random_state=42).fit(Xd, yd)
print('single tree / sklearn forest test', base.score(Xt, yt), reference.score(Xt, yt))
print('sklearn OOB', reference.oob_score_)
 
# Tree-count comparison uses OOB development scores, never selects from test accuracy.
for b in [16, 32, 64]:
    r = RandomForestClassifier(n_estimators=b, oob_score=True, n_jobs=1, random_state=42).fit(Xd, yd)
    print('trees', b, 'OOB', r.oob_score_)
 
# This held-out importance is descriptive; using it to select features needs a new evaluation.
p = permutation_importance(reference, Xt, yt, n_repeats=5, n_jobs=1, random_state=42)
print('MDI', reference.feature_importances_)
print('permutation mean / SD', p.importances_mean, p.importances_std)
 
# Exact bootstrap and covariance arithmetic checks.
assert np.isclose((1-1/10000)**10000, np.exp(-1), atol=2e-5)
assert np.isclose(.2+(1-.2)/100, .208)

What each block does

grow samples candidate features on every recursive call, then applies the familiar split loop only to those candidates. If none of the sampled features produces positive gain, it stops; scikit-learn's exact behavior when candidate features cannot yield a valid partition is more nuanced. This teaching implementation is therefore not expected to reproduce the library tree-for-tree. The midpoint calculation avoids summing two large endpoints directly. If rounding makes the midpoint equal to the lower endpoint, it uses the upper endpoint instead: under strict <, this still separates the two observed values. The adjacent-float fixture checks this numerical boundary independently of bootstrap randomness.

Leaves always represent class-1 probability with a scalar, avoiding misaligned probability columns when a bootstrap happens to omit a class. A general multiclass implementation must map each tree's classes_ columns back to a global class order. It cannot assume every tree has the same columns or that class labels equal positional indices.

fit clears prior state, otherwise fitting twice would silently double the forest. OOB predictions are computed while the original fitting rows are available; no later API accepts an unrelated row order and pretends its indices match. NaN marks uncovered rows, and coverage is printed before accuracy. Hard predictions use a documented tie rule, with OOB > .5 matching class-zero ties.

The single-tree and sklearn comparisons are demonstrations on one generated dataset, not estimates of an industrial ranking. The count sweep uses OOB rather than repeatedly consulting the final test. The importance values describe this fixed fitted model on test; using them to redesign features would make that set developmental and require another final evaluation.

(d) Production reality · 15 min

Forest or booster under real constraints?

A forest is often an effective baseline when threshold interactions matter and independent tree training is convenient. It can still require encoding, missing-value policy, group-aware splitting, class-imbalance handling and calibration. Hundreds of deep trees can create a large serving artifact. A booster may obtain a useful representation with smaller trees and sequential corrections, but tree count, depth and learning rate interact. Neither family is always smaller, faster or more accurate at a given quality target.

On a dataset with millions of rows, first measure fit time, peak memory, prediction latency and validation quality. Candidate interventions include reducing depth or leaf count, increasing minimum leaf size, limiting bootstrap sample size, reducing tree count once prediction variability stabilizes, and using a bounded number of worker threads. Extra Trees or histogram-based alternatives are different algorithms/configurations, not drop-in promises of identical statistical properties. Subsampling must retain important rare classes and temporal/group structure.

For class imbalance, start by declaring the positive class and reporting relevant metrics, not by assuming bagging solves it. A tree may receive few or no positive examples in a bootstrap; class weighting, sampling strategy and thresholds change different parts of the system. Probability averaging does not guarantee calibration.

Correlated-feature practice with an answer

Duplicate a useful column in both development and evaluation arrays, optionally with small noise. Refit, then compare MDI and individual permutation importances for the pair. Do not expect an exact 50/50 split: competing thresholds, feature sampling, ties and noise affect allocation.

For joint permutation, construct one row permutation order, then assign the two-column block with Xp[:, group] = Xeval[order][:, group]. Measure the original score minus the perturbed score over several independently drawn orders. Reusing the same order for the group preserves the pair's internal relationship. This estimates reliance on the group under this perturbation, not a causal effect; correlated proxies outside the group can still mask reliance.

Where this shows up in the rest of the plan

Gradient boosting replaces independent fits with sequential loss corrections. Evaluation metrics explains why stable average probabilities do not directly specify an alerting policy. Feature engineering carries the OOB/CV boundary back into preprocessing. The forest is both a modeling candidate and a way to learn covariance-aware aggregation.

(e) Recall + stretch · 10 min

1. What is bagging? Answer: bootstrap aggregation—fit on with-replacement row samples and combine predictions. Variance reduction depends on covariance, and bias relative to a full-data learner can change.
2. Why about 63% distinct rows? Answer: inclusion probability is $1-(1-1/n)^n$, tending to $1-e^{-1}$. It is an expectation, not an exact sample size.
3. What makes this a per-split forest? Answer: a new rng.choice occurs inside grow, not just once in the tree-building loop.
4. Is OOB k-fold CV? Answer: no; each row uses its own excluded-tree subset. It is useful but does not solve preprocessing, group/time leakage or repeated-search bias.
5. Do more trees guarantee better finite-test accuracy? Answer: no. The covariance identity describes numeric prediction variance under stated assumptions; thresholded accuracy can fluctuate.
6. What is wrong with default importance as truth? Answer: MDI reflects training split use and can favor high-cardinality variables. Permutation reflects a model's reliance under a perturbation, with correlated-feature and distribution limits.
7. Do trees require standard scaling? Answer: usually not; numeric ordering governs threshold partitions. Encoding, missingness and finite-precision qualifications still matter.
8. When prefer a forest? Answer: when its measured validation/serving trade-off and parallel fitting suit the problem. Small data or a need for explanations alone does not guarantee superiority to boosting.
9. Extra Trees? Answer: additional threshold randomization, and no bootstrap by sklearn default. It trades search effort and diversity against fit quality; measure rather than promise a win.
10. Stretch: slow 10-million-row fit. Answer: profile, test leaf/depth/sample/tree budgets, cap parallelism to available memory, and compare alternatives under equal evaluation. A fixed downsample size or library switch is not a guarantee of quality or equivalent behavior.

Explain-out-loud test

Derive the covariance formula including its cross terms. Then explain why bootstrapping changes rows while feature sampling changes local split opportunities. Finally identify two ways an apparently good OOB score can leak information. A complete answer must mention the upstream fit boundary, not only the bootstrap indices.

What comes next

Continue to Gradient Boosting — XGBoost, LightGBM, where each tree has a different target because it corrects the current ensemble rather than independently predicting the original labels.

Sources and execution boundary

Checked 2026-10-07: scikit-learn 1.7 ensembles, RandomForestClassifier API, and permutation importance. The lab is a bounded local CPU verification of the scratch mechanism, OOB bookkeeping and API use. No large-scale throughput, production ranking, fairness or calibration benchmark is implied.