Random Forest: Algorithm and Regression

Bagging trains independent trees on bootstrap samples. Random Forest adds one more source of randomness: at each split, only a random subset of features is considered. This decorrelates the trees and reduces the ensemble variance beyond what bootstrapping alone achieves.

Anchor: 6-sample house price dataset (2 features for hand trace), then California Housing for code.

import numpy as np

# Hand-trace anchor: sq_ft + bedrooms → price
X_2feat = np.array([
    [650, 2], [850, 2], [1100, 3],
    [1400, 3], [1600, 4], [1900, 4]
])
y_2feat = np.array([180, 220, 280, 340, 370, 430])

Bagging vs Random Forest: The One Difference

Bagging: bootstrap n samples → train a full decision tree on all p features.

Random Forest = Bagging + at each split, randomly select $\sqrt{p}$ features (classification) or $p/3$ features (regression) and search only those.

For $p=2$ features (sq_ft, bedrooms), $\sqrt{2}\approx 1$ → at each split, only 1 randomly chosen feature is tested. Trees can only grow using whichever feature the coin flip selected — making them different from each other even beyond the bootstrap sample differences.

Two Sources of Randomness

Both together make trees much less correlated than pure bagging. Less correlated trees → bigger variance reduction when averaging.

Building 3 Trees on the Anchor

With p=2 features, each split coin-flips between sq_ft and bedrooms (max_features=1 in this case).

Bootstrap 1 (indices: 0,0,2,3,4,5)

Samples: [650,2]×2, [1100,3], [1400,3], [1600,4], [1900,4]. OOB: indices 1, skipped.

Unique samples: 5. Feature at root: sq_ft (coin flip: heads).

Best sq_ft split: t=1250 (separates the two [650,2] and [1100,3] from [1400,3],[1600,4],[1900,4]).

Bootstrap 2 (indices: 0,1,2,2,4,4)

Samples: [650,2],[850,2],[1100,3]×2,[1600,4]×2. OOB: indices 3, 5.

Feature at root: bedrooms (coin flip: tails).

Best bedrooms split: t=2.5 (bedrooms ≤ 2 → indices 0,1; bedrooms > 2 → indices 2,2,4,4).

For x_new=[1250,3]: bedrooms=3 > 2 → predict 325.0.

Bootstrap 3 (indices: 1,2,3,4,5,5)

Samples: [850,2],[1100,3],[1400,3],[1600,4],[1900,4]×2. OOB: index 0.

Feature at root: sq_ft (coin flip: heads).

Best sq_ft split: t=1250 (similar to Tree 1 but with different sample weights due to [1900,4]×2).

For x_new=[1250,3]: sq_ft=1250 ≤ 1250 → predict 250.0.

Ensemble Prediction for x_new=[1250, 3]

Why Feature Subsampling Helps: The Dominant Feature Problem

Imagine features are [sq_ft, bedrooms, day_of_week, neighborhood_id].

Without feature subsampling (pure bagging): sq_ft is always the best split at the root of every tree → all trees start identically → high pairwise correlation ρ.

Recall the variance formula: $\text{Var}(\text{ensemble})=\rho \cdot {\sigma }^2+\frac{1-\rho }{T}\cdot {\sigma }^2$

As $T\to \infty$, $\text{Var}\to \rho {\sigma }^2$. The correlation ρ sets a floor on how much ensembling helps.

With feature subsampling ($\sqrt{4}=2$ features per split): sq_ft is excluded from ~50% of splits → different trees explore bedrooms and other features at the root → lower ρ → lower ensemble variance floor.

This is the entire justification for Random Forest's feature subsampling: not to find better splits, but to build less correlated trees.

sklearn Implementation

from sklearn.ensemble import RandomForestRegressor
from sklearn.tree import DecisionTreeRegressor
from sklearn.datasets import fetch_california_housing
from sklearn.model_selection import train_test_split
from sklearn.metrics import mean_squared_error, r2_score
import numpy as np

data = fetch_california_housing()
X, y = data.data, data.target

X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=42)

# Single tree baseline
dt = DecisionTreeRegressor(random_state=42)
dt.fit(X_train, y_train)

# Random Forest
rf = RandomForestRegressor(n_estimators=100, max_features='sqrt', random_state=42, n_jobs=-1)
rf.fit(X_train, y_train)

for name, model in [('Decision Tree', dt), ('Random Forest', rf)]:
    y_pred = model.predict(X_test)
    rmse = np.sqrt(mean_squared_error(y_test, y_pred))
    r2 = r2_score(y_test, y_pred)
    print(f"{name:15s}: RMSE={rmse:.4f}, R²={r2:.4f}")

100 trees reduce RMSE by 0.228 — a 31% improvement over the single tree — purely from ensemble averaging.

n_estimators: Diminishing Returns

import time

print(f"{'n_trees':>8} | {'RMSE':>8} | {'R²':>8} | {'time(s)':>8}")
for n in [1, 5, 10, 50, 100, 200, 500]:
    rf_n = RandomForestRegressor(n_estimators=n, max_features='sqrt', random_state=42, n_jobs=-1)
    t0 = time.time()
    rf_n.fit(X_train, y_train)
    elapsed = time.time() - t0
    rmse = np.sqrt(mean_squared_error(y_test, rf_n.predict(X_test)))
    r2 = r2_score(y_test, rf_n.predict(X_test))
    print(f"{n:>8} | {rmse:>8.4f} | {r2:>8.4f} | {elapsed:>8.2f}")

Adding trees from 1→10 is high-value (RMSE drops 0.19). Adding trees from 100→500 gains only 0.003 at 4.7× more training time. The default n=100 sits on the diminishing-returns knee.

max_features: Feature Subsampling Sensitivity

print(f"{'max_features':>14} | {'n_features':>10} | {'RMSE':>8}")
for mf in ['sqrt', 0.3, 0.5, 0.7, 1.0]:
    rf_mf = RandomForestRegressor(n_estimators=100, max_features=mf, random_state=42, n_jobs=-1)
    rf_mf.fit(X_train, y_train)
    rmse = np.sqrt(mean_squared_error(y_test, rf_mf.predict(X_test)))
    n_feat = int(mf * 8) if isinstance(mf, float) else int(8**0.5)
    print(f"{str(mf):>14} | {n_feat:>10} | {rmse:>8.4f}")

max_features=1.0 (pure bagging) gives the worst RMSE among the options — confirming that feature subsampling helps by decorrelating trees. sqrt (the default) is near-optimal.

OOB Score: Free Validation

rf_oob = RandomForestRegressor(n_estimators=100, oob_score=True,
                                max_features='sqrt', random_state=42)
rf_oob.fit(X_train, y_train)
print(f"OOB  R²: {rf_oob.oob_score_:.4f}")
print(f"Test R²: {r2_score(y_test, rf_oob.predict(X_test)):.4f}")

OOB R² (0.765) is within 0.005 of the test R² (0.770). The OOB estimate is a reliable substitute for a held-out validation set — no CV needed. The small underestimate is because each sample is evaluated by fewer trees (only those for which it was OOB), but the approximation converges quickly with n_estimators ≥ 50.

Test Your Understanding

1. For the hand-trace anchor (n=6 samples, p=2 features), each bootstrap draws 6 samples with replacement. The expected number of unique samples is $6\times (1-(5/6{)}^6)\approx 4.0$. Verify this: compute $(1-(5/6{)}^6)$ and multiply by 6. How many OOB samples does each tree typically have?

2. In the max_features sweep, max_features=1.0 (pure bagging) was worse than max_features='sqrt'. But in theory, considering more features at each split should find better individual splits. Why does giving each split access to all features make the ensemble worse, even though it makes individual trees better?

3. The variance formula says $\text{Var}\to \rho {\sigma }^2$ as $T\to \infty$. In the sweep, RMSE at n=500 (0.5001) was barely better than n=100 (0.5031). What does this tell you about the minimum achievable RMSE from this approach — regardless of how many trees you add?

4. OOB R²=0.765 underestimates test R²=0.770 slightly. Each OOB estimate uses only ~36.8% of trees. A sample that is OOB for k trees gets a prediction averaged over those k trees. At n_estimators=10, each sample is OOB for ~3–4 trees on average. At n_estimators=100, it's ~37 trees. How does this explain the OOB underestimate at small n, and why does OOB converge to the true generalization error as n_estimators increases?

5. Random Forest's feature importance is computed as the total decrease in node impurity (MSE for regression) attributable to each feature, averaged over all trees. For the California Housing dataset (8 features: MedInc, HouseAge, AveRooms, AveBedrms, Population, AveOccup, Latitude, Longitude), which feature would you expect to have the highest importance, and why? (MedInc = median income in the block group.)