Random Forest — 2.2 GB Model Crashes 1 GB Pod

Random Forest is one of the most widely deployed machine learning algorithms in production systems today. From detecting credit card fraud at banks to predicting patient readmission in hospitals, it quietly powers decisions that affect millions of people. If you've ever wondered why a seasoned ML engineer reaches for Random Forest before trying something fancier, this article will show you exactly why.

I've personally deployed Random Forest models in three different production systems — a fraud scoring pipeline processing 40,000 transactions per second at a payments company, a customer churn predictor at a subscription SaaS platform, and a demand forecasting model for a retail chain. In every case, Random Forest was the first model we built, and in two of those cases, it remained the production model even after we tried gradient boosting and neural networks. Not because it scored highest on the leaderboard, but because it was the most reliable, the easiest to debug at 2 AM when something broke, and the fastest to retrain when upstream data schemas changed.

The algorithm was introduced by Leo Breiman in 2001 — a statistician who had already co-invented CART (Classification and Regression Trees) two decades earlier. Breiman's insight was deceptively simple: a single decision tree is brittle because it commits to one specific way of splitting the data. But if you grow hundreds of trees, each on a slightly different view of the data, and average their predictions, the individual errors cancel out. The crowd wisdom of diverse trees beats any single expert tree.

The core problem Random Forest solves is overfitting. A single decision tree is like a very eager student who memorises the exam paper instead of learning the subject — it performs brilliantly on training data and falls apart on anything new. Random Forest fixes this by deliberately injecting two kinds of randomness: random subsets of training rows (bagging) and random subsets of features at each split. Those two tricks force each tree to be different, and different trees make different errors. When you average their predictions, the errors cancel out and the signal survives.

There's a reason Random Forest has survived 25 years of ML hype cycles while flashier algorithms have come and gone. It parallelises trivially across CPU cores — training 500 trees on 8 cores is nearly 8x faster than training sequentially, unlike gradient boosting where each tree depends on the previous one. It handles unscaled features, mixed data types, and moderate missingness without complaint. And it gives you built-in feature importance out of the box — critical when a compliance officer asks 'why did your model flag this transaction?'

By the end of this article you'll be able to build, tune, and interpret a Random Forest model in Python using scikit-learn. You'll understand what hyperparameters actually matter (and which ones are mostly noise), how to extract feature importances for stakeholder reports, and exactly when Random Forest is the right tool versus when you should reach for something else.

How Random Forest Actually Builds Its Trees (Bagging + Feature Randomness)

Random Forest is an ensemble method built on two independent randomisation strategies, and understanding both is the difference between using it like a black box and using it with confidence.

The first strategy is Bootstrap Aggregating, universally called bagging. For each tree, scikit-learn samples the training dataset with replacement — meaning the same row can appear multiple times in one tree's training set, while roughly 37% of rows never appear at all. That 37% figure isn't arbitrary: for a dataset of n rows, the probability that any single row is never selected in one bootstrap sample is (1 - 1/n)^n, which converges to 1/e ≈ 0.368 as n grows. So each tree trains on about 63% of the data and has never seen the other 37%. Those excluded rows are called the Out-of-Bag (OOB) samples, and they act as a free validation set for that tree. This is important: you get an honest error estimate without touching your test set.

The second strategy is feature randomness. At every single split in every single tree, the algorithm only considers a random subset of features — typically the square root of the total number of features for classification (so if you have 100 features, each split only sees about 10 candidates). This seems counterintuitive; why hide information from the tree? Because without this step, every tree would pick the same dominant feature at the top split, the trees would be correlated, and correlated trees don't cancel each other's errors — they amplify them. I learned this the hard way on a pricing model where two features (customer tenure and contract value) dominated every split. Without feature subsampling, all 500 trees looked nearly identical, and the ensemble performed barely better than a single tree. Turning on max_features='sqrt' immediately dropped our MAE by 15%.

Here's how a single split actually works under the hood. The tree considers each candidate feature and evaluates every possible threshold value in the sorted feature values. For each (feature, threshold) pair, it splits the data into two groups and computes a purity metric — either Gini impurity or entropy (information gain). The Gini impurity of a node is:

Gini = 1 - Σ(p_i)²

where p_i is the proportion of samples belonging to class i. A perfectly pure node (all one class) has Gini = 0. A maximally mixed node (equal proportions of all classes) has Gini approaching 1. The algorithm picks the (feature, threshold) split that produces the largest weighted decrease in Gini across the two child nodes.

Entropy uses a different formula:

Entropy = -Σ(p_i × log2(p_i))

In practice, Gini and entropy produce almost identical trees — the split points differ by tiny amounts on most real datasets. Don't waste tuning time on this choice; I've never seen it change model performance by more than 0.1% on any production dataset.

Combine these two strategies and you get an ensemble where each tree is both trained on different data and forced to explore different feature combinations. Their individual mistakes become uncorrelated noise, and the majority vote or average surfaces the true signal.

import numpy as np
from sklearn.ensemble import RandomForestClassifier
from sklearn.datasets import load_breast_cancer
from sklearn.model_selection import train_test_split
from sklearn.metrics import classification_report

# io.thecodeforge.ml.ensemble.RandomForestDemo
# Load a real medical dataset — predicting malignant vs benign tumours
cancer_data = load_breast_cancer()
feature_matrix = cancer_data.data      # 30 numeric features per tumour sample
target_labels  = cancer_data.target    # 0 = malignant, 1 = benign

# Hold out 20% of data — the model never sees this during training
X_train, X_test, y_train, y_test = train_test_split(
    feature_matrix, target_labels,
    test_size=0.20,
    random_state=42,          # fix seed so results are reproducible
    stratify=target_labels    # keep class proportions equal in both splits
)

# Build the forest
# n_estimators: number of trees — more is generally better up to a point
# max_features: 'sqrt' means each split considers sqrt(30) ≈ 5 features randomly
# oob_score: use the free out-of-bag rows to estimate generalisation error
# n_jobs=-1: use all CPU cores — forests are embarrassingly parallel
forest_model = RandomForestClassifier(
    n_estimators=200,
    max_features='sqrt',
    oob_score=True,
    random_state=42,
    n_jobs=-1
)

forest_model.fit(X_train, y_train)

# OOB score is computed from rows that were NOT used to train each tree
# It's a reliable estimate of generalisation without touching X_test
print(f"Out-of-Bag accuracy estimate : {forest_model.oob_score_:.4f}")

# Now evaluate on the truly held-out test set
y_predictions = forest_model.predict(X_test)
print("\n--- Test Set Performance ---")
print(classification_report(y_test, y_predictions,
                            target_names=cancer_data.target_names))

# Sanity check: verify the 37% OOB math
# Each tree trains on ~63% of data; the rest are OOB samples
n_total = len(X_train)
avg_oob_per_tree = np.mean([np.sum(estimator.indices_) for estimator in forest_model.estimators_])
print(f"\nAverage rows per tree     : {avg_oob_per_tree:.0f} / {n_total} ({avg_oob_per_tree/n_total:.1%})")
print(f"Average OOB rows per tree : {n_total - avg_oob_per_tree:.0f} / {n_total} ({1 - avg_oob_per_tree/n_total:.1%})")

Feature Importance — Turning the Model Into a Story Your Stakeholders Understand

One reason Random Forest dominates in industry — despite gradient boosting often scoring slightly higher on leaderboards — is interpretability. Every trained forest can tell you exactly which features drove its decisions. That matters enormously when you need to explain a fraud-detection model to a compliance team or a churn model to a product manager.

Scikit-learn computes Mean Decrease in Impurity (MDI) importance: for each feature, it sums how much Gini impurity dropped across all splits and all trees where that feature was used, then normalises the result to sum to 1.0. A high score means that feature consistently produced clean splits across the forest.

One important caveat: MDI can inflate the importance of high-cardinality numerical features. If you're working with features that have wildly different numbers of unique values — like 'age' (continuous) versus 'country' (5 categories) — consider using Permutation Importance instead. It measures how much the model's accuracy drops when a feature's values are randomly shuffled, which is a more honest reflection of real-world impact.

Always plot both and compare. If they broadly agree, you can trust the ranking. If they disagree significantly, dig into why — it's often a signal that two features are correlated and the model is leaning on whichever one it found first.

Beyond MDI and permutation importance, there's a third approach gaining traction in production systems: SHAP (SHapley Additive exPlanations) values. SHAP assigns each feature a contribution score for every individual prediction, based on game theory. The TreeExplainer variant is optimised for tree ensembles and runs fast even on large forests. I've started using SHAP in every model I ship to production because it answers the question permutation importance can't: 'why did the model make this specific decision for this specific customer?' That's what your customer support team actually needs when a user calls to dispute a fraud flag.

One practical workflow I use in production: run MDI first (it's instant), validate the top-10 with permutation importance (takes a few minutes), and generate SHAP summary plots for the final stakeholder presentation. If all three broadly agree, you have a robust feature ranking. If they diverge, you've found something interesting worth investigating — usually correlated features or a data leakage path you missed.

import numpy as np
import matplotlib.pyplot as plt
from sklearn.ensemble import RandomForestClassifier
from sklearn.inspection import permutation_importance
from sklearn.datasets import load_breast_cancer
from sklearn.model_selection import train_test_split

# io.thecodeforge.ml.interpretability.FeatureImportanceAnalysis

cancer_data = load_breast_cancer()
X_train, X_test, y_train, y_test = train_test_split(
    cancer_data.data, cancer_data.target,
    test_size=0.20, random_state=42, stratify=cancer_data.target
)

forest_model = RandomForestClassifier(
    n_estimators=200, max_features='sqrt',
    random_state=42, n_jobs=-1
)
forest_model.fit(X_train, y_train)

feature_names = cancer_data.feature_names

# --- Method 1: Mean Decrease in Impurity (MDI) ---
# Comes free with every fitted RandomForest — fast but can bias toward
# high-cardinality features
mdi_importances = forest_model.feature_importances_
mdi_sorted_idx  = np.argsort(mdi_importances)[::-1]

print("Top 5 Features by MDI Importance:")
for rank, idx in enumerate(mdi_sorted_idx[:5], start=1):
    print(f"  {rank}. {feature_names[idx]:<35} {mdi_importances[idx]:.4f}")

# --- Method 2: Permutation Importance ---
# Slower but more honest — shuffles one feature at a time and measures accuracy drop
# n_repeats=15 means shuffle each feature 15 times and average the result
perm_result = permutation_importance(
    forest_model, X_test, y_test,
    n_repeats=15,
    random_state=42,
    n_jobs=-1
)
perm_sorted_idx = np.argsort(perm_result.importances_mean)[::-1]

print("\nTop 5 Features by Permutation Importance:")
for rank, idx in enumerate(perm_sorted_idx[:5], start=1):
    print(f"  {rank}. {feature_names[idx]:<35} {perm_result.importances_mean[idx]:.4f}")

# --- Method 3: SHAP values (if shap is installed) ---
# Gives per-prediction explanations — the gold standard for stakeholder reports
try:
    import shap
    explainer = shap.TreeExplainer(forest_model)
    shap_values = explainer.shap_values(X_test[:200])  # sample for speed

    # For binary classification, shap_values is a list [class_0, class_1]
    # Use class 1 (benign) values
    shap_importances = np.abs(shap_values[1]).mean(axis=0)
    shap_sorted_idx = np.argsort(shap_importances)[::-1]

    print("\nTop 5 Features by SHAP Importance:")
    for rank, idx in enumerate(shap_sorted_idx[:5], start=1):
        print(f"  {rank}. {feature_names[idx]:<35} {shap_importances[idx]:.4f}")
except ImportError:
    print("\n[SHAP not installed — run: pip install shap]")
    shap_importances = None
    shap_sorted_idx = None

# --- Visual comparison side by side ---
n_methods = 3 if shap_importances is not None else 2
fig, axes = plt.subplots(1, n_methods, figsize=(6 * n_methods, 6))

top_n = 10

# MDI bar chart
axes[0].barh(
    range(top_n),
    mdi_importances[mdi_sorted_idx[:top_n]][::-1],
    color='steelblue'
)
axes[0].set_yticks(range(top_n))
axes[0].set_yticklabels([feature_names[i] for i in mdi_sorted_idx[:top_n]][::-1])
axes[0].set_title('MDI Feature Importance')
axes[0].set_xlabel('Mean Decrease in Impurity')

# Permutation bar chart
axes[1].barh(
    range(top_n),
    perm_result.importances_mean[perm_sorted_idx[:top_n]][::-1],
    color='darkorange'
)
axes[1].set_yticks(range(top_n))
axes[1].set_yticklabels([feature_names[i] for i in perm_sorted_idx[:top_n]][::-1])
axes[1].set_title('Permutation Feature Importance')
axes[1].set_xlabel('Mean Accuracy Drop on Test Set')

# SHAP bar chart (if available)
if shap_importances is not None:
    axes[2].barh(
        range(top_n),
        shap_importances[shap_sorted_idx[:top_n]][::-1],
        color='seagreen'
    )
    axes[2].set_yticks(range(top_n))
    axes[2].set_yticklabels([feature_names[i] for i in shap_sorted_idx[:top_n]][::-1])
    axes[2].set_title('SHAP Feature Importance')
    axes[2].set_xlabel('Mean |SHAP value|')

plt.tight_layout()
plt.savefig('feature_importance_comparison.png', dpi=150)
print("\nPlot saved to feature_importance_comparison.png")

# --- Correlation diagnostic: find features that disagree between MDI and permutation ---
print("\n--- Disagreement Diagnostic ---")
print("Features where MDI rank differs from permutation rank by >5 positions:")
for feat_idx in range(len(feature_names)):
    mdi_rank = np.where(mdi_sorted_idx == feat_idx)[0][0]
    perm_rank = np.where(perm_sorted_idx == feat_idx)[0][0]
    if abs(mdi_rank - perm_rank) > 5:
        print(f"  {feature_names[feat_idx]:<35} MDI rank: {mdi_rank+1}, Perm rank: {perm_rank+1}")

Hyperparameter Tuning — The 20% of Knobs That Do 80% of the Work

Random Forest has a reputation for working well out of the box, and that reputation is earned. But 'good enough out of the box' is not the same as 'optimised for your problem'. Knowing which hyperparameters actually move the needle — and which are mostly cosmetic — saves you hours of pointless grid search.

The three parameters that genuinely matter are: n_estimators (more trees reduces variance but with diminishing returns past ~300), max_depth (limiting tree depth is the single most powerful guard against overfitting on small datasets), and min_samples_leaf (requiring each leaf to contain at least N samples smooths the decision boundary and helps with noisy labels).

Parameters that matter less than people think: max_features almost always works well at 'sqrt' for classification and 'log2' is rarely better. criterion ('gini' vs 'entropy') barely changes outcomes on most real datasets — the shapes of the two functions are nearly identical for balanced distributions.

Use RandomizedSearchCV rather than GridSearchCV. With Random Forest you're exploring a large continuous space; random sampling finds good regions faster than exhaustive grid search, and you can control the compute budget directly with n_iter.

One pattern I've used in every production system: tune in two passes. First pass uses RandomizedSearchCV with n_iter=40 and wide ranges to find the neighbourhood. Second pass narrows the ranges around the best values and runs a finer search. This two-pass approach consistently finds better hyperparameters than a single large search with the same total compute budget.

Another production trick: use warm_start=True to incrementally grow your forest. Train with 100 trees, check the OOB score, add another 100, check again. When the OOB score stops improving, stop. This avoids the guesswork of picking n_estimators upfront and gives you the learning curve for free.

For imbalanced datasets — which is most production classification problems — go beyond class_weight='balanced'. Consider using class_weight as a dictionary where you manually set higher weights for the minority class based on business cost. If a missed fraud costs $500 and a false alarm costs $5 in customer support time, your class weights should roughly reflect that 100:1 ratio, not the inverse frequency ratio that 'balanced' computes.

import numpy as np
from sklearn.ensemble import RandomForestClassifier
from sklearn.model_selection import RandomizedSearchCV, StratifiedKFold, train_test_split
from sklearn.datasets import load_breast_cancer
from sklearn.metrics import classification_report
from scipy.stats import randint

# io.thecodeforge.ml.tuning.RandomForestTuning

cancer_data = load_breast_cancer()
X_train, X_test, y_train, y_test = train_test_split(
    cancer_data.data, cancer_data.target,
    test_size=0.20, random_state=42, stratify=cancer_data.target
)

# Define the search space — use distributions for continuous params
hyperparam_space = {
    'n_estimators': randint(100, 600),
    'max_depth': [None, 5, 10, 20, 30],
    'min_samples_leaf': randint(1, 20),
    'min_samples_split': randint(2, 20),
    'max_features': ['sqrt', 'log2', 0.3]
}

base_forest = RandomForestClassifier(random_state=42)
cv_strategy = StratifiedKFold(n_splits=5, shuffle=True, random_state=42)

random_search = RandomizedSearchCV(
    estimator=base_forest,
    param_distributions=hyperparam_space,
    n_iter=40,
    scoring='f1_weighted',
    cv=cv_strategy,
    verbose=1,
    random_state=42,
    n_jobs=-1
)

random_search.fit(X_train, y_train)

print("\nBest hyperparameters found:")
for param_name, param_value in random_search.best_params_.items():
    print(f"  {param_name:<22}: {param_value}")

print(f"\nBest cross-validated F1 (weighted): {random_search.best_score_:.4f}")

# Evaluate on test set
best_forest = random_search.best_estimator_
y_predictions = best_forest.predict(X_test)

print("\n--- Tuned Model — Test Set Performance ---")
print(classification_report(y_test, y_predictions, target_names=cancer_data.target_names))

# --- Warm Start: Incremental Tree Growth ---
print("\n--- Warm Start: Incremental Tree Growth ---")
warm_forest = RandomForestClassifier(
    n_estimators=0,
    warm_start=True,
    max_features='sqrt',
    oob_score=True,
    random_state=42,
    n_jobs=-1
)

oob_scores = []
for n_trees in range(50, 501, 50):
    warm_forest.n_estimators = n_trees
    warm_forest.fit(X_train, y_train)
    oob_scores.append((n_trees, warm_forest.oob_score_))
    print(f"  n_estimators={n_trees:>4d}  OOB score={warm_forest.oob_score_:.4f}")

for i in range(1, len(oob_scores)):
    improvement = oob_scores[i][1] - oob_scores[i-1][1]
    if improvement < 0.001:
        print(f"\n  OOB plateaued at {oob_scores[i][0]} trees (improvement: {improvement:.5f})")
        break

When to Use Random Forest — and When to Walk Away

Random Forest is not the right tool for every problem, and knowing when to walk away is just as important as knowing how to use it.

Use Random Forest when: your dataset has a mix of numerical and categorical features, you have moderate-to-high dimensional data (dozens to hundreds of features), you need a reliable baseline quickly with minimal preprocessing (no feature scaling required), you need built-in feature importance for stakeholder communication, or your dataset is moderately sized — say 1,000 to 1,000,000 rows.

Consider alternatives when: your data has millions of rows and inference latency matters (gradient boosting with LightGBM will be faster and often more accurate), you're working with sequential or spatial data where structure matters (tree ensembles ignore the ordering of features), interpretability must be ironclad for regulatory reasons (a single shallow decision tree or logistic regression is easier to audit), or your problem involves image or text data (neural networks handle raw pixels and tokens far better).

One underrated strength: Random Forest is almost impossible to catastrophically misconfigure. You can hand it unscaled features, a few NaN values (with some workarounds), and class imbalance, and it still produces a reasonable model. That robustness is why it's the go-to algorithm for early-stage data exploration and rapid prototyping in industry.

Let me give you a real example. At a fintech company I worked at, we needed a model to predict which loan applicants would default. The dataset had 200,000 rows, 45 features (mix of numeric financials and categorical demographics), and about 6% default rate. We built a Random Forest in 20 minutes — no feature scaling, no encoding gymnastics, just OrdinalEncoder on the categoricals and go. AUC: 0.87. The gradient boosting model we built the following week scored 0.89. Better, yes — but it took 3 days of tuning, and the marginal 0.02 AUC improvement didn't change the business outcome. The Random Forest went to production because it was good enough, fast to retrain weekly, and the compliance team could actually understand the feature importances.

That said, Random Forest has real limitations. It can't extrapolate — if your test data contains values outside the training range, the tree just predicts the nearest leaf value. I've seen this bite teams doing time-series forecasting where future values trend upward; the Random Forest flatlines at the maximum training value. It's also memory-hungry: each tree stores its full structure, and 500 deep trees can easily consume several gigabytes. For real-time inference at sub-millisecond latency (like ad bidding), you might need to reduce tree count or switch to a lighter model.

Random Forest also doesn't handle feature interactions as explicitly as gradient boosting. Boosting builds each new tree specifically to correct the residual errors of the ensemble so far, which lets it discover complex interactions naturally. Random Forest's trees are independent — they find interactions only if the random feature subsampling happens to include the right combination, which is less systematic.

# Random Forest isn't just for classification — regression is equally powerful.
# Here we predict house prices, a classic regression task with mixed feature types.

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

# io.thecodeforge.ml.regression.RandomForestRegression

# California housing: predict median house value from 8 features
housing_data   = fetch_california_housing()
feature_matrix = housing_data.data    # e.g. median income, avg rooms, latitude
house_prices   = housing_data.target  # median house value in units of $100,000

X_train, X_test, y_train, y_test = train_test_split(
    feature_matrix, house_prices,
    test_size=0.20, random_state=42
)

# For regression: RandomForestRegressor averages leaf node values instead of voting
# max_features='sqrt' still works well; some practitioners use 1/3 of features
forest_regressor = RandomForestRegressor(
    n_estimators=300,
    max_depth=None,
    min_samples_leaf=3,   # slightly higher than default helps smooth regression curves
    max_features='sqrt',
    oob_score=True,
    random_state=42,
    n_jobs=-1
)

forest_regressor.fit(X_train, y_train)

y_predicted_prices = forest_regressor.predict(X_test)

# R² score: 1.0 is perfect, 0.0 means the model is no better than predicting the mean
r2  = r2_score(y_test, y_predicted_prices)
# MAE: average absolute error in the target unit ($100,000 in this case)
mae = mean_absolute_error(y_test, y_predicted_prices)

print(f"R² Score (test set)     : {r2:.4f}")
print(f"Mean Absolute Error     : ${mae * 100_000:,.0f} per house")
print(f"OOB R² estimate         : {forest_regressor.oob_score_:.4f}")

# Show a few sample predictions vs actuals
print("\nSample Predictions vs Actual (first 6 test houses):")
print(f"  {'Actual':>12}  {'Predicted':>12}  {'Error':>10}")
for actual, predicted in zip(y_test[:6], y_predicted_prices[:6]):
    error = abs(actual - predicted) * 100_000
    print(f"  ${actual*100_000:>10,.0f}  ${predicted*100_000:>10,.0f}  ${error:>8,.0f}")

# --- Demonstrate the extrapolation problem ---
# Create synthetic data where y = 2x + noise, train on x in [0, 10],
# then predict on x in [15, 20] — RF can't extrapolate beyond training range
print("\n--- Extrapolation Problem Demo ---")
np.random.seed(42)
X_synth = np.random.uniform(0, 10, size=(500, 1))
y_synth = 2 * X_synth.ravel() + np.random.normal(0, 0.5, size=500)

rf_synth = RandomForestRegressor(n_estimators=100, random_state=42)
rf_synth.fit(X_synth, y_synth)

X_future = np.array([[12], [15], [18], [20]])  # outside training range
predictions = rf_synth.predict(X_future)
actual = 2 * X_future.ravel()

print(f"  {'x':>5}  {'True (2x)':>10}  {'RF Predicted':>12}  {'Error':>8}")
for x, true_val, pred in zip(X_future.ravel(), actual, predictions):
    print(f"  {x:>5.0f}  {true_val:>10.1f}  {pred:>12.1f}  {abs(true_val - pred):>8.1f}")
print("  → RF predictions plateau near the max training y value, not the linear trend")

Production Deployment and Monitoring — Avoiding OOM, Latency, and Drift

Shipping a Random Forest model to production is more than calling .predict(). You need to handle memory constraints, latency SLAs, and data drift — or the model will fail silently.

Memory is the biggest surprise. A forest with 500 trees and no depth limit can easily exceed 1 GB. Before deploying, always serialize with joblib and check size: sys.getsizeof(joblib.dump(model, '/dev/null')). If it's over 500 MB, reduce n_estimators or cap max_depth. In Kubernetes, set the memory limit to at least 2x the model size to account for peak allocation during prediction.

Inference latency grows linearly with n_estimators. For sub-10ms SLAs, consider ONNX export via sklearn-onnx. The ONNX runtime can be 2-5x faster than scikit-learn's predict() because of graph optimisations. Another trick: set n_jobs=1 in the model before exporting to avoid thread contention in web servers.

Data drift is the silent killer. After retraining on new data, compare OOB scores: a drop of more than 0.05 is a red flag. Set up a monitoring job that logs OOB score, feature distribution statistics (mean, std per feature), and a sample of SHAP values. When drift is detected, trigger an alert and automatically roll back to the previous model version.

One more gotcha: when you retrain with a different random_state, the bootstrapped samples change, and the model's predictions will shift slightly even on the same data. This is normal, but it can confuse stakeholders. Fix the random_state across all production retrain jobs to ensure reproducibility. If you need to change it (e.g., for hyperparameter search), version the random_state in the model metadata.

import joblib
import numpy as np
from sklearn.ensemble import RandomForestRegressor

# io.thecodeforge.ml.production.RandomForestProduction

# ---- Step 1: Serialize and check size ----
model = RandomForestRegressor(n_estimators=200, max_depth=20, random_state=42)
model.fit(X_train, y_train)

import sys
model_bytes = joblib.dump(model, '/tmp/model.joblib')
model_size_mb = sys.getsizeof(open('/tmp/model.joblib', 'rb').read()) / 1024**2
print(f"Model size: {model_size_mb:.1f} MB")

# ---- Step 2: Monitor OOB score after retraining ----
def retrain_and_alert(old_oob_score, X_new, y_new, threshold=0.05):
    new_model = RandomForestRegressor(n_estimators=200, max_depth=20, random_state=42, oob_score=True)
    new_model.fit(X_new, y_new)
    new_oob_score = new_model.oob_score_
    if old_oob_score - new_oob_score > threshold:
        print(f"ALERT: OOB score dropped from {old_oob_score:.4f} to {new_oob_score:.4f}. Possible drift.")
        # rollback logic here
    return new_model

# ---- Step 3: Track feature distributions ----
def feature_stats(X, feature_names):
    means = X.mean(axis=0)
    stds = X.std(axis=0)
    for name, mean, std in zip(feature_names, means, stds):
        print(f"  {name:<30} mean={mean:.2f}, std={std:.2f}")

# ---- Step 4: SHAP monitoring sample ----
try:
    import shap
    explainer = shap.TreeExplainer(model)
    shap_sample = shap_values = explainer.shap_values(X_test[:100])
    # Check if top features changed
    top_features_prev = ['worst concave points', 'worst radius']
    top_features_now = list(np.argsort(np.abs(shap_values).mean(axis=0))[::-1][:5])
    print(f"Previous top features: {top_features_prev}")
    print(f"Current top features: {feature_names[top_features_now]}")
except ImportError:
    print("SHAP not installed — install for production monitoring.")