36 - Random Forest: why a hundred trees err less than one

Ask a person to estimate the weight of an ox at a country fair: they will be wrong, one way or the other. Ask a hundred independent people and take the average: the overestimates and underestimates tend to cancel out, and the average is much closer to the truth than almost any single estimate. It is the "wisdom of the crowd". But it works only on one condition: the errors must be weakly correlated. If a hundred people copy from the same neighbor, the average is wrong just like the neighbor.

The lesson on decision trees ended with bad news: a single tree is unstable, it changes a lot if the data sample changes a little. The Random Forest turns that defect into an advantage: it trains many different trees and averages them. In this lesson we derive why it works (with the formula for the variance of the mean), we build it by hand and we compare it with the library.

Two terms to clarify: ensemble, bias and variance

An ensemble is a model built by combining many simpler models, called base learners. To understand why combining them helps, the error must be decomposed. Imagine playing darts: the bias is how far, on average, your throws are shifted from the center of the target (systematic error: the aim is crooked); the variance is how scattered the throws are among themselves (random error: the hand shakes). A deep tree is like a thrower with straight aim but a shaky hand: on average it hits the center, but each single throw lands wherever it happens to. Averaging many throwers cancels out the tremor.

Why the average of many models is wrong less often

The error of a model decomposes into bias (systematic error, the model is too simple) and variance (how much the prediction changes as the training sample changes). A deep tree has low bias and high variance: the profile of a model that would do well to be averaged.

Take B predictions X₁...X_B, each with variance σ² and pairwise correlation ρ. The variance of their mean is:

Var( (1/B) Σ X_i ) = (1/B²) · [ B·σ² + B(B−1)·ρσ² ]
                   = σ²/B + ((B−1)/B)·ρσ²
                   = ρσ² + (1−ρ)σ²/B

The first step counts the B variances and the B(B−1) covariances (each one worth ρσ²). Two levers can be read from this formula:

• The term (1−ρ)σ²/B vanishes as B grows: more trees, less variance. This is why adding trees does not cause overfitting.
• The term ρσ² does not depend on B: if the trees are highly correlated, averaging helps little. The only way to lower it is to reduce the correlation ρ between the trees.

The two levers correspond to two different ways of being wrong. With few trees you are wrong because of the dispersion in the first term, and you cure it with quantity. With trees that are too similar you are wrong because of the second term, which no quantity cures. It follows that "more trees" is an incomplete recipe: what is new in a Random Forest compared with a simple average is everything that makes the trees different, and there are two ingredients of randomness.

Ingredient 1: bagging (bootstrap aggregating)

Each tree is trained on a bootstrap sample: n examples drawn with replacement from the training set of n examples. Some examples appear several times, others never. How many are left out? The probability that a given example is not drawn in a single draw is 1 − 1/n; in n independent draws it is (1 − 1/n)ⁿ, which for large n tends to 1/e ≈ 0.368. So each tree does not see about 37% of the examples. Let us verify it:

import numpy as np

rng = np.random.default_rng(0)
for n in (10, 100, 1000, 100000):
    fuori = [len(np.setdiff1d(np.arange(n), rng.integers(0, n, size=n))) / n for _ in range(20)]
    print("n=%6d  share left out by the bootstrap: %.3f   (theory (1-1/n)^n = %.3f)" % (n, np.mean(fuori), (1 - 1 / n) ** n))
print("limit 1/e = %.3f" % np.exp(-1))

n=    10  share left out by the bootstrap: 0.400   (theory (1-1/n)^n = 0.349)
n=   100  share left out by the bootstrap: 0.367   (theory (1-1/n)^n = 0.366)
n=  1000  share left out by the bootstrap: 0.371   (theory (1-1/n)^n = 0.368)
n=100000  share left out by the bootstrap: 0.368   (theory (1-1/n)^n = 0.368)
limit 1/e = 0.368

For each size n the listing draws a bootstrap (rng.integers(0, n, size=n): n indices with replacement), counts with np.setdiff1d how many indices between 0 and n−1 were not drawn, and repeats the experiment 20 times to take the average. The comparison with (1−1/n)ⁿ shows the theory: already with n = 100 the excluded fraction is 0.367, very close to the limit 1/e = 0.368.

These excluded examples are a gift: for each tree they form a small ready-made validation set, and they are called out-of-bag (OOB). To evaluate the forest, each example is voted on only by the trees that did not see it in training: an estimate of generalization without a separate validation set.

The flow of bagging: the training set generates several bootstrap samples, one tree per sample, and the final prediction is the mean (or the vote) of all the trees.

Ingredient 2: a random subset of features at each split

Bagging alone is not enough. If there is a very strong feature, all the trees will choose it in the first splits and end up resembling each other (high ρ). The Random Forest adds a rule: at each node, the best split is searched only among m randomly drawn features (usually m = √d for classification). In this way different trees are forced to use different features, and the correlation drops. Each individual tree becomes a little worse, but the average improves: it is the compromise in the formula above.

Implementation by hand

We already have a tree, in scikit-learn, with the max_features parameter that implements exactly the random draw at each node. We build the forest with a loop, the majority vote and the OOB estimate:

from sklearn.datasets import make_classification
from sklearn.model_selection import train_test_split
from sklearn.tree import DecisionTreeClassifier
from sklearn.ensemble import RandomForestClassifier

X, y = make_classification(n_samples=2000, n_features=20, n_informative=5, n_redundant=3,
                           flip_y=0.05, class_sep=0.8, random_state=0)
Xtr, Xte, ytr, yte = train_test_split(X, y, test_size=0.3, stratify=y, random_state=0)

class ForestaDaZero:
    def __init__(self, n_alberi=100, max_features="sqrt", seed=0):
        self.n_alberi, self.max_features, self.seed = n_alberi, max_features, seed

    def fit(self, X, y):
        rng = np.random.default_rng(self.seed)
        n = len(X)
        self.alberi, self.fuori = [], []
        for i in range(self.n_alberi):
            idx = rng.integers(0, n, size=n)                      # bootstrap
            mask_fuori = np.ones(n, dtype=bool); mask_fuori[idx] = False
            t = DecisionTreeClassifier(max_features=self.max_features,
                                       random_state=int(rng.integers(1_000_000)))
            self.alberi.append(t.fit(X[idx], y[idx]))
            self.fuori.append(mask_fuori)
        return self

    def predict_proba(self, X):
        return np.mean([t.predict_proba(X) for t in self.alberi], axis=0)   # mean of the probabilities

    def predict(self, X):
        return self.predict_proba(X).argmax(axis=1)

    def oob_score(self, X, y):
        somma = np.zeros((len(X), 2)); voti = np.zeros(len(X))
        for t, f in zip(self.alberi, self.fuori):
            somma[f] += t.predict_proba(X[f]); voti[f] += 1
        ok = voti > 0                                              # examples left out by at least one tree
        return (somma[ok].argmax(axis=1) == y[ok]).mean()

mia = ForestaDaZero(n_alberi=200).fit(Xtr, ytr)
lib = RandomForestClassifier(n_estimators=200, random_state=0, oob_score=True).fit(Xtr, ytr)
singolo = DecisionTreeClassifier(random_state=0).fit(Xtr, ytr)

print("single tree          test %.3f" % singolo.score(Xte, yte))
print("hand-made forest     test %.3f   OOB %.3f" % ((mia.predict(Xte) == yte).mean(), mia.oob_score(Xtr, ytr)))
print("sklearn forest       test %.3f   OOB %.3f" % (lib.score(Xte, yte), lib.oob_score_))

single tree          test 0.778
hand-made forest     test 0.875   OOB 0.866
sklearn forest       test 0.863   OOB 0.864

The ForestaDaZero class has three roles. fit: for each of the n_alberi trees it draws a bootstrap (idx), marks in mask_fuori which examples were not drawn (the out-of-bag ones), and trains a DecisionTreeClassifier with max_features="sqrt" on the bootstrap examples only; a different random_state for each tree. predict_proba: mean of the probabilities of the individual trees. oob_score: for each example it lets vote only the trees for which it was out-of-bag (somma[f] += ... uses the mask f to sum the probabilities only over the excluded examples of that tree) and compares with the true label.

The single tree scores 0.778 on the test; the hand-made forest 0.875 and the scikit-learn one 0.863 (the small differences between the two depend on the random seeds). Note that the OOB estimate (computed on the training set only) is close to the test accuracy: it works as a free validation. The scikit-learn forest averages the probabilities of the individual trees, like ours, instead of the majority vote of Breiman's original paper: the ensemble guide says so.

How many trees are needed?

for B in (1, 5, 25, 100, 400):
    f = RandomForestClassifier(n_estimators=B, random_state=0).fit(Xtr, ytr)
    print("trees %4d   train %.3f   test %.3f" % (B, f.score(Xtr, ytr), f.score(Xte, yte)))

trees    1   train 0.905   test 0.743
trees    5   train 0.980   test 0.817
trees   25   train 0.999   test 0.857
trees  100   train 1.000   test 0.863
trees  400   train 1.000   test 0.865

The listing trains forests with a growing number of trees and prints training and test. With B = 1 the forest is a single tree on a bootstrap (worse than a tree trained on the whole training set, because it sees only 63% of the distinct examples).

The test goes from 0.743 with a single tree to 0.857 with 25, 0.863 with 100 and 0.865 with 400: the accuracy grows and then levels off: beyond a certain number of trees the gain is negligible and you pay only in time and memory. Adding trees does not cause overfitting, as the formula anticipated; but a training accuracy close to 1.0 means nothing by itself, because every deep tree memorizes its own bootstrap.

Measuring the correlation between the trees

Theory says that feature subsampling serves to decorrelate. Let us measure it: we compare forests with max_features=None (pure bagging, all features at each node) and with sqrt, computing the mean correlation between the probability predictions of pairs of trees on the test set:

def correlazione_media(foresta, X):
    P = np.array([t.predict_proba(X)[:, 1] for t in foresta.estimators_])      # (trees, samples)
    C = np.corrcoef(P)
    return C[np.triu_indices_from(C, k=1)].mean()

for nome, mf in (("pure bagging (all features)", None), ("random forest (sqrt)", "sqrt")):
    f = RandomForestClassifier(n_estimators=100, max_features=mf, random_state=0).fit(Xtr, ytr)
    print("%-32s mean correlation between trees %.3f   test %.3f" % (nome, correlazione_media(f, Xte), f.score(Xte, yte)))

pure bagging (all features)      mean correlation between trees 0.549   test 0.878
random forest (sqrt)             mean correlation between trees 0.440   test 0.863

The correlazione_media function: estimators_ is the list of the trees of the forest. For each one it obtains the probability of class 1 on the test set (a trees × samples matrix); np.corrcoef returns the correlation matrix between trees, and with np.triu_indices_from(C, k=1) the upper triangle without the diagonal is taken (each pair only once) to compute the mean.

The mean correlation between the trees falls from 0.549 to 0.440 when the features are subsampled: it is the ρ term of the formula that goes down, as expected. But look at the accuracy: in this dataset pure bagging scores 0.878 against 0.863 for the forest. With 600 test samples the difference is 9 examples, so it is not enough to declare a winner. The honest message is this: decorrelating is a means, not a guarantee. The advantage of max_features depends on the dataset (it generally grows when a few features dominate), and must be chosen with validation like any hyperparameter.

Feature importance: two methods, one pitfall

The forest offers feature_importances_, computed as the mean impurity reduction brought by each feature. It is fast but has a known defect: it favors features with many distinct values, even useless ones, because they offer more thresholds to try. The alternative is the permutation importance: it randomly shuffles the values of a column in the validation data and measures how much the accuracy drops; if the column is not needed, it does not drop. Let us add a column of pure noise with many distinct values:

from sklearn.inspection import permutation_importance

rng = np.random.default_rng(1)
Xr_tr = np.hstack([Xtr, rng.normal(size=(len(Xtr), 1))])
Xr_te = np.hstack([Xte, rng.normal(size=(len(Xte), 1))])
f = RandomForestClassifier(n_estimators=300, random_state=0).fit(Xr_tr, ytr)

imp_gini = f.feature_importances_
perm = permutation_importance(f, Xr_te, yte, n_repeats=10, random_state=0).importances_mean
print("noise column (last)  impurity importance: %.4f   (mean of the other 20: %.4f)" % (imp_gini[-1], imp_gini[:-1].mean()))
print("noise column (last)  permutation importance:  %.4f   (mean of the 5 highest: %.4f)" % (perm[-1], np.sort(perm[:-1])[-5:].mean()))

noise column (last)  impurity importance: 0.0204   (mean of the other 20: 0.0490)
noise column (last)  permutation importance:  -0.0017   (mean of the 5 highest: 0.0669)

The listing adds a column of continuous Gaussian noise (all distinct values, hence with many possible thresholds) and trains the forest. feature_importances_ is the impurity importance; permutation_importance shuffles one column at a time 10 times on the test set and measures the mean drop in accuracy.

The impurity-based importance attributes to the noise column 0.0204, that is about 2% of the total importance (which sums to 1): a nonzero value for a column that by construction says nothing. The permutation importance on the test set puts it at −0.0017, that is zero within the noise, while the five most important features exceed 0.066 on average. Practical rule: use the permutation one to understand what really matters; remember that, with strongly correlated features, even that one shares the importance among "twin" columns. The scikit-learn guide confirms both points: the impurity importance for trees is strongly biased toward features with many distinct values, and is computed on the training set, while permutation can be computed on unseen data.

Trade-offs

• Pros: an excellent starting point on tabular data, few critical hyperparameters, robust to noise, requires no scaling, free OOB estimate, parallelizable (n_jobs=-1), estimates feature importance.
• Cons: loses the interpretability of the single tree; a heavy model in memory and slower at prediction; does not extrapolate outside the training range; often beaten by boosting (next lesson) when you need the last point of accuracy.
• Hyperparameters that matter: max_features (decorrelation), min_samples_leaf (shallower trees, less noise), n_estimators (more is better until it levels off).

The common thread of the pros and cons is that the forest is an average of many deep trees: averaging makes it robust and almost free of tuning, but removes the readability of the single tree and brings with it the weight of hundreds of trees. Of the three hyperparameters, max_features and min_samples_leaf are the ones to validate; n_estimators is raised until the curve flattens.

Real-world pitfalls

• Trusting the training accuracy: it will be almost 1.0 by construction. Look at OOB or validation.
• OOB with non-independent data: if you have duplicate rows or several rows from the same user, the example "out" of a tree has a twin inside, and the OOB estimate is optimistic. A split by groups is needed.
• Time series: bootstrap and random validation mix past and future (leakage). Use temporal splits.
• Imbalanced classes: use class_weight="balanced_subsample" and evaluate with recall/precision; accuracy lies (see below).
• Importance from impurity alone: as seen above, it favors useless high-cardinality columns.
• Uncalibrated probabilities: the vote fractions of the trees are not true probabilities; if you need them, calibrate them.

Four of the six pitfalls concern evaluation: training, OOB with non-independent data, time series, accuracy on imbalanced classes. The other two concern how to read importance and probabilities. None produces a runtime error: the forest runs and returns plausible numbers, and it is the evaluation method that has to be corrected. For duplicates and time series see the process of an ML project.

The hyperparameters that matter: max_features and min_samples_leaf

The forest has few critical hyperparameters, and two of them control the two ingredients of the lesson. max_features regulates the decorrelation of the trees: with 1 each split chooses among a single feature (almost random trees), with "sqrt" the default (about 4-5 out of 20), with 0.5 half, with None all of them (back to pure bagging). min_samples_leaf regulates the depth: with 1 the trees grow until purity, with larger values they are pruned upstream and each leaf averages more examples. We evaluate a grid of twelve combinations with three accuracy estimates: the OOB (free, computed during training), the cross-validation with 5 folds on the training set and the test. If the OOB is reliable, the first two must be close.

from sklearn.model_selection import cross_val_score, StratifiedKFold

cv = StratifiedKFold(5, shuffle=True, random_state=0)
print("max_features   min_samples_leaf   OOB     CV (5 folds)   test")
for max_features in (1, "sqrt", 0.5, None):
    for leaf in (1, 5, 20):
        forest = RandomForestClassifier(200, max_features=max_features, min_samples_leaf=leaf, oob_score=True, random_state=0, n_jobs=-1).fit(Xtr, ytr)
        cv_score = cross_val_score(RandomForestClassifier(200, max_features=max_features, min_samples_leaf=leaf, random_state=0, n_jobs=-1), Xtr, ytr, cv=cv).mean()
        print("%-12s   %15d    %.3f   %.3f          %.3f" % (max_features, leaf, forest.oob_score_, cv_score, forest.score(Xte, yte)))

max_features   min_samples_leaf   OOB     CV (5 folds)   test
1                            1    0.828   0.841          0.837
1                            5    0.824   0.829          0.828
1                           20    0.804   0.804          0.823
sqrt                         1    0.864   0.873          0.863
sqrt                         5    0.861   0.856          0.857
sqrt                        20    0.848   0.839          0.833
0.5                          1    0.883   0.879          0.870
0.5                          5    0.877   0.881          0.868
0.5                         20    0.858   0.849          0.842
None                         1    0.881   0.878          0.878
None                         5    0.884   0.876          0.877
None                        20    0.864   0.854          0.865

For each of the twelve combinations, the listing trains a forest of 200 trees with the OOB estimate enabled, computes its cross-validation (with a new forest for each fold) and measures the test.

Read three things. First, the OOB is a good estimate of generalization: its difference from cross-validation is at most 1.3 points (for example 0.864 against 0.873, 0.881 against 0.878) and it costs a fraction of the time, because it does not require retraining.

Second, max_features matters: with 1 the accuracy is 0.82-0.84, with the square root 0.83-0.86 and with 0.5 or with all the features 0.84-0.88. Here the default value is not the best: pure bagging (None) and 0.5 win, because in this problem there are 5 informative and 3 redundant features out of 20 and the random choice of only 4 or 5 features at each split often leaves out the informative ones. The advantage of "sqrt" shows when a few features dominate and many are useless or nearly duplicated: it depends on the dataset, and it is a hyperparameter to search.

Third, min_samples_leaf: raising it from 1 to 20 worsens things systematically (OOB from 0.864 to 0.848 with sqrt, from 0.881 to 0.864 with all the features): the forest already averages many trees and does not need pruned trees; excessive pruning increases the bias. These trends hold for this problem: the lesson is the method, not the values.

Extra Trees: even more randomness

The Extremely Randomized Trees (P. Geurts, D. Ernst and L. Wehenkel, 2006) push the idea further: at each node, besides choosing a random subset of features, they also choose the threshold at random (one per feature, between the observed minimum and maximum) and keep the best of the random thresholds, instead of searching for the optimal one. Each tree individually is worse, but the correlation between trees collapses (the ρ term of the variance formula) and the computation is much faster, because it does not sort the values. In scikit-learn they do not bootstrap by default (they use all the examples), as the ensemble guide indicates.

import time
from sklearn.ensemble import ExtraTreesClassifier

for name, cls in (("RandomForest", RandomForestClassifier), ("ExtraTrees", ExtraTreesClassifier)):
    t0 = time.time()
    model = cls(300, random_state=0, n_jobs=1).fit(Xtr, ytr)
    print("%-13s test accuracy %.3f   training time %.2f s (single process)" % (name, model.score(Xte, yte), time.time() - t0))

# Output of one run; the times change at every run
# RandomForest  test accuracy 0.865   training time 2.36 s (single process)
# ExtraTrees    test accuracy 0.873   training time 0.87 s (single process)

The listing trains the two families with 300 trees, everything else being equal, and measures accuracy and time with a single process (n_jobs=1, for a fair comparison). Extra Trees score 0.873 against 0.865 for the Random Forest, a difference that is within the noise (with 600 test examples it is five examples). On time, in three runs on a PC that at that moment was also being used for other things, the Random Forest took from 2.18 to 3.55 s and the Extra Trees from 0.87 to 1.92 s, that is between 0.4 and 0.8 times the time. The result is common because the optimal threshold is not searched but chosen among random thresholds; if you have a lot of data or little time, Extra Trees are an alternative to try and to compare with validation.

The probabilities of a forest: are they reliable?

The predict_proba method of a forest returns the mean of the class fractions in the leaves of the trees. In the case of Naive Bayes we saw probabilities that were too sure of themselves. What happens here? We build the same reliability table: groups of predictions by probability, mean declared confidence and actual fraction of positives.

forest = RandomForestClassifier(300, min_samples_leaf=3, random_state=0, n_jobs=-1).fit(Xtr, ytr)
p_positive = forest.predict_proba(Xte)[:, 1]

print("predicted probability     predictions   declared mean   actual fraction of positives")
edges = [0, 0.1, 0.3, 0.5, 0.7, 0.9, 1.0001]
for lo, hi in zip(edges[:-1], edges[1:]):
    in_bin = (p_positive >= lo) & (p_positive < hi)
    print("   [%.1f, %.1f)         %8d        %.3f              %.3f" % (lo, min(hi, 1.0), in_bin.sum(), p_positive[in_bin].mean(), yte[in_bin].mean()))

predicted probability     predictions   declared mean   actual fraction of positives
   [0.0, 0.1)               49        0.069              0.041
   [0.1, 0.3)              195        0.189              0.113
   [0.3, 0.5)               91        0.389              0.374
   [0.5, 0.7)               54        0.584              0.704
   [0.7, 0.9)               96        0.828              0.938
   [0.9, 1.0)              115        0.946              0.974

The listing groups the 600 probability predictions of the forest into six intervals and compares, in each, the mean declared probability with the actual fraction of positives. The profile is the opposite of Naive Bayes: the forest is cautious. It must be said that the caution is not uniform: in the lowest interval (0-0.1) it declares 0.069 and the reality is 0.041, that is the forest overestimates the positives; in the 0.7-0.9 interval it declares 0.828 but the positives are 93.8%; in the 0.9-1 one it declares 0.946 and the reality is 0.974. So in the upper part the probabilities are too low, in the lower part too high.

The average over many trees, which are never all in agreement even when the example is "easy", compresses the probabilities toward the center: it cannot give 0 or 1 if even a single tree dissents. It is a less serious defect than that of Naive Bayes (the ordering is right and the distortion is moderate), and it is corrected with calibration (CalibratedClassifierCV). The probabilities of a forest are therefore better than those of Naive Bayes, but should not be taken as true probabilities without a check.

Imbalanced classes: weights and threshold

With 5% positives, a forest that predicts "by majority" (threshold 0.5) tends to ignore the rare class. Two levers: the class weights (class_weight), which make errors on the positives count more in the impurity computation, and the threshold on the probabilities. We compare them on the same imbalanced problem as KNN and Naive Bayes (here 5,000 examples).

from sklearn.metrics import recall_score, precision_score

X_imb, y_imb = make_classification(n_samples=5000, n_features=10, n_informative=4, weights=[0.95, 0.05], class_sep=1.0, flip_y=0.02, random_state=1)
X_imb_tr, X_imb_te, y_imb_tr, y_imb_te = train_test_split(X_imb, y_imb, test_size=0.3, stratify=y_imb, random_state=0)

print("class_weight          threshold   recall   precision")
for class_weight in (None, "balanced_subsample"):
    f = RandomForestClassifier(300, class_weight=class_weight, min_samples_leaf=3, random_state=0, n_jobs=-1).fit(X_imb_tr, y_imb_tr)
    p = f.predict_proba(X_imb_te)[:, 1]
    for threshold in (0.5, 0.3, 0.2, 0.1):
        pred = p >= threshold
        print("%-20s  %.1f         %.3f    %.3f" % (class_weight, threshold, recall_score(y_imb_te, pred), precision_score(y_imb_te, pred, zero_division=0)))

class_weight          threshold   recall   precision
None                  0.5         0.149    0.867
None                  0.3         0.322    0.737
None                  0.2         0.414    0.429
None                  0.1         0.621    0.233
balanced_subsample    0.5         0.230    0.800
balanced_subsample    0.3         0.425    0.446
balanced_subsample    0.2         0.598    0.295
balanced_subsample    0.1         0.759    0.177

The listing trains two forests, without weights and with balanced_subsample (the class weights are recomputed for each bootstrap sample), and for each measures recall and precision at four thresholds. With threshold 0.5 and no weights the forest finds only 15% of the positives (precision 0.87); with the weights 23% (precision 0.80): a modest improvement. The most powerful lever is the threshold: lowering it to 0.2 brings the recall to 41% without weights (precision 0.43) and to 60% with the weights (precision 0.30), and at 0.1 to 62% and 76% (with precision of 0.23 and 0.18). The two tools are partly interchangeable, but they do not produce new information: they move the operating point along the same precision-recall curve. The gain from the weights, at the same threshold, is real (at 0.2: 60% against 41% recall with precision 0.30 against 0.43) but it is mostly a shift along the curve. As always, the right threshold depends on the costs, as the process of an ML project explains.

The cost in production: size and prediction time

A forest is a heavy model. Each tree stores its nodes (threshold, feature, children, values), and a hundred or three hundred deep trees easily become tens of megabytes. We measure the effect of depth on the size of the model and on the prediction time, comparing it with a logistic regression that predicts the same class.

import pickle
from sklearn.linear_model import LogisticRegression

def measure(model, name):
    size_kb = len(pickle.dumps(model)) / 1024
    t0 = time.time()
    for _ in range(20):
        model.predict(Xte)
    ms_per_1000 = (time.time() - t0) / 20 * 1000 * 1000 / len(Xte)
    nodes = sum(t.tree_.node_count for t in model.estimators_) if hasattr(model, "estimators_") else 0
    print("%-34s %9d %12.0f %20.2f" % (name, nodes, size_kb, ms_per_1000))

print("model                              total nodes   size (KB)   ms per 1000 predictions")
for leaf in (1, 5, 20):
    measure(RandomForestClassifier(300, min_samples_leaf=leaf, random_state=0, n_jobs=1).fit(Xtr, ytr), "RandomForest 300 trees, min. leaf %d" % leaf)
measure(LogisticRegression(max_iter=1000).fit(Xtr, ytr), "logistic regression")

# Output of one run; nodes and sizes are deterministic, the times are not
# model                              total nodes   size (KB)   ms per 1000 predictions
# RandomForest 300 trees, min. leaf 1     85832         6790               110.98
# RandomForest 300 trees, min. leaf 5     48946         3907                86.55
# RandomForest 300 trees, min. leaf 20     18502         1528               108.21
# logistic regression                        0            1                 0.22

The measure function serializes the model with pickle to estimate its size, times 20 predictions on the whole test set and reports the milliseconds per thousand predictions; for the forests it also counts the total nodes. The effect of min_samples_leaf on size is very strong and deterministic (same numbers at every run): with one-example leaves the forest has 85,832 nodes and takes 6.8 MB; raising the minimum to 20 it drops to 18,502 nodes and 1.5 MB, about a quarter, at the expense of some accuracy (measured above: from 0.864 to 0.848 with the OOB).

The prediction time is less stable. In three runs, on a PC that was also used for other things, the forests took from 73 to 158 ms per thousand predictions, with no consistent order among the three depths, because the time depends mainly on the number of trees to traverse (always 300). The logistic regression takes about 1 KB (twenty-one numbers: twenty coefficients and an intercept) and predicts in 0.22-0.98 ms per thousand examples: about eighty times faster or more in every run. This is the real cost of the forest in production: if latency, memory or loading time are constraints, min_samples_leaf and the number of trees are sizing levers, and you choose the smallest model whose accuracy is not worse by an acceptable margin.

Three things to try

Exercise 1. Repeat the grid with 20 and with 1,000 trees: does the OOB get closer to the cross-validation as the trees increase? Exercise 2. Add 30 columns of pure noise and redo the max_features grid: is the "sqrt" default now the best? Exercise 3. Apply CalibratedClassifierCV to the forest and redo the reliability table: how close do the two columns get?

When the forest is not the right tool

A random forest is an excellent starting point for tabular data, but it has structural limits worth knowing. First, it does not extrapolate: each tree predicts the mean of the leaf values, so outside the range of the training data the prediction stays flat at the last value seen, unlike a linear model that continues along the line; for quantities with a trend (prices over time, growing demand) this is a systematic error, not noise.

Second, feature importance computed with impurity decrease favors variables with many distinct values and is split among correlated variables: read it as an indication, not as a causal ranking, and check it with permutation on the validation set. Third, the model is large and slow to serve compared with a linear one: hundreds of deep trees must be kept in memory and traversed at each prediction, which matters when latency is a constraint. Fourth, the probabilities of a forest are averages of votes, so they tend to avoid the extremes and must be calibrated if they are used as true probabilities. None of these limits makes it inadvisable: they serve to decide when to replace it or pair it with a linear model or a boosting model.

To practise on the domain of predictive maintenance of irrigation pumps there is an exercise skeleton in Python with three parts to complete: the out-of-bag share and the 1/e limit, the variance of the mean of correlated estimators with the number of trees needed, and Condorcet's jury theorem. The comments in the file are in Italian. As long as the points to complete are empty, it prints which exercises are missing.

What you can do now

You can derive why the average of many trees is wrong less often and where it stops (the floor ρσ²), you can build a forest with bootstrap, feature subset, voting and out-of-bag estimate, and you can measure the correlation among the trees instead of taking it for granted. You can choose max_features and min_samples_leaf with a grid, recognize the limits of probabilities and of impurity importance, and size the model for production. In the next lesson, boosting, the trees are no longer independent: each one corrects the errors of the previous one.

Sources