Trees, Forests and Boosting

Week 6 of 16 · Supervised learning · 7 days

Week 06 · Supervised learning

Trees, Forests and Boosting

Week 06 · Day 1 of 7

How a Decision Tree Decides

Impurity, greedy splitting, and the only model you can read end to end

Every model so far draws a single boundary through the whole space. A decision tree does something different: it asks a question, splits the data, and asks a different question of each half.

How a split gets chosen

The tree considers every feature and every possible cut point, and keeps the one that most reduces impurity.

Gini impurity: The probability that two rows drawn at random from a node have different labels. Zero means the node is pure, every row the same class. For two classes the maximum is 0.5, at a 50/50 mix. The tree picks the split that lowers the weighted average impurity of the children the most.

import numpy as np

def gini(labels):
    labels = np.asarray(labels)
    if len(labels) == 0:
        return 0.0
    p1 = labels.mean()
    return 1 - (p1 ** 2 + (1 - p1) ** 2)

print('all one class      %.4f' % gini([1, 1, 1, 1]))
print('three to one       %.4f' % gini([1, 1, 1, 0]))
print('even mix           %.4f' % gini([1, 1, 0, 0]))

all one class      0.0000
three to one       0.3750
even mix           0.5000

import numpy as np
import pandas as pd

df = pd.read_csv('customers.csv').drop_duplicates().copy()
df['contract'] = df['contract'].str.strip().str.title()
df['total_charges'] = pd.to_numeric(df['total_charges'], errors='coerce')

NUM = ['tenure_months', 'support_calls', 'monthly_charges']
CAT = ['contract', 'internet_service', 'payment_method', 'has_dependents']
TARGET = 'churned'
import numpy as np

def gini(labels):
    labels = np.asarray(labels)
    if len(labels) == 0:
        return 0.0
    p1 = labels.mean()
    return 1 - (p1 ** 2 + (1 - p1) ** 2)

y = df[TARGET].to_numpy()
x = df['tenure_months'].to_numpy()
parent = gini(y)

print('impurity before any split: %.4f' % parent)
print('\n%12s %10s %10s %12s' % ('split at', 'left n', 'right n', 'gain'))
best = None
for cut in [6, 12, 18, 24, 36, 48]:
    left, right = y[x <= cut], y[x > cut]
    weighted = (len(left) * gini(left) + len(right) * gini(right)) / len(y)
    gain = parent - weighted
    if best is None or gain > best[1]:
        best = (cut, gain)
    print('%12d %10d %10d %12.5f' % (cut, len(left), len(right), gain))

print('\nbest of these: tenure <= %d, gain %.5f' % best)

impurity before any split: 0.3924

    split at     left n    right n         gain
           6        345       2655      0.00723
          12        949       2051      0.02653
          18       1545       1455      0.03244
          24       1918       1082      0.03664
          36       2420        580      0.02537
          48       2686        314      0.01562

best of these: tenure <= 24, gain 0.03664

That is the entire algorithm, applied recursively. Split, then repeat on each side, until a stopping rule fires.

Let scikit-learn do it and read the result

import numpy as np
import pandas as pd

df = pd.read_csv('customers.csv').drop_duplicates().copy()
df['contract'] = df['contract'].str.strip().str.title()
df['total_charges'] = pd.to_numeric(df['total_charges'], errors='coerce')

NUM = ['tenure_months', 'support_calls', 'monthly_charges']
CAT = ['contract', 'internet_service', 'payment_method', 'has_dependents']
TARGET = 'churned'
from sklearn.compose import ColumnTransformer
from sklearn.pipeline import Pipeline
from sklearn.impute import SimpleImputer
from sklearn.preprocessing import StandardScaler, OneHotEncoder

def preprocessor(scale=True):
    num_steps = [('i', SimpleImputer(strategy='median'))]
    if scale:
        num_steps.append(('s', StandardScaler()))
    return ColumnTransformer([
        ('num', Pipeline(num_steps), NUM),
        ('cat', Pipeline([('i', SimpleImputer(strategy='most_frequent')),
                          ('o', OneHotEncoder(handle_unknown='ignore'))]), CAT),
    ])
from sklearn.model_selection import train_test_split

X, y = df[NUM + CAT], df[TARGET]
X_tr, X_te, y_tr, y_te = train_test_split(X, y, test_size=0.25,
                                          stratify=y, random_state=42)
from sklearn.tree import DecisionTreeClassifier, export_text

tree = Pipeline([('prep', preprocessor(scale=False)),
                 ('clf', DecisionTreeClassifier(max_depth=3,
                                                random_state=42))]).fit(X_tr, y_tr)
names = list(tree.named_steps['prep'].get_feature_names_out())
print(export_text(tree.named_steps['clf'], feature_names=names))

|--- cat__contract_Month-To-Month <= 0.50
|   |--- num__tenure_months <= 17.50
|   |   |--- cat__contract_One Year <= 0.50
|   |   |   |--- class: 0
|   |   |--- cat__contract_One Year >  0.50
|   |   |   |--- class: 0
|   |--- num__tenure_months >  17.50
|   |   |--- cat__contract_Two Year <= 0.50
|   |   |   |--- class: 0
|   |   |--- cat__contract_Two Year >  0.50
|   |   |   |--- class: 0
|--- cat__contract_Month-To-Month >  0.50
|   |--- num__monthly_charges <= 51.02
|   |   |--- num__tenure_months <= 9.50
|   |   |   |--- class: 0
|   |   |--- num__tenure_months >  9.50
|   |   |   |--- class: 0
|   |--- num__monthly_charges >  51.02
|   |   |--- num__tenure_months <= 32.50
|   |   |   |--- class: 1
|   |   |--- num__tenure_months >  32.50
|   |   |   |--- class: 0

This is why people like trees

You can read the whole model. Every prediction is a path from the root to a leaf, and you can print that path for any customer and hand it to someone who has never heard of machine learning. No other model in this course offers that.

Trees do not need scaling, and barely care about outliers

import numpy as np
import pandas as pd

df = pd.read_csv('customers.csv').drop_duplicates().copy()
df['contract'] = df['contract'].str.strip().str.title()
df['total_charges'] = pd.to_numeric(df['total_charges'], errors='coerce')

NUM = ['tenure_months', 'support_calls', 'monthly_charges']
CAT = ['contract', 'internet_service', 'payment_method', 'has_dependents']
TARGET = 'churned'
from sklearn.compose import ColumnTransformer
from sklearn.pipeline import Pipeline
from sklearn.impute import SimpleImputer
from sklearn.preprocessing import StandardScaler, OneHotEncoder

def preprocessor(scale=True):
    num_steps = [('i', SimpleImputer(strategy='median'))]
    if scale:
        num_steps.append(('s', StandardScaler()))
    return ColumnTransformer([
        ('num', Pipeline(num_steps), NUM),
        ('cat', Pipeline([('i', SimpleImputer(strategy='most_frequent')),
                          ('o', OneHotEncoder(handle_unknown='ignore'))]), CAT),
    ])
from sklearn.model_selection import train_test_split

X, y = df[NUM + CAT], df[TARGET]
X_tr, X_te, y_tr, y_te = train_test_split(X, y, test_size=0.25,
                                          stratify=y, random_state=42)
from sklearn.tree import DecisionTreeClassifier
from sklearn.metrics import roc_auc_score

for scale in [True, False]:
    m = Pipeline([('prep', preprocessor(scale=scale)),
                  ('clf', DecisionTreeClassifier(max_depth=5,
                                                 random_state=42))]).fit(X_tr, y_tr)
    print('scaled=%-6s AUC %.6f'
          % (scale, roc_auc_score(y_te, m.predict_proba(X_te)[:, 1])))

scaled=True   AUC 0.774524
scaled=False  AUC 0.774524

Identical to six decimal places. A split asks “is this value above the threshold”, and scaling preserves order, so it changes nothing. The same reasoning is why the 97-call outlier in week 2 moved a tree not at all.

Gini or entropy

import numpy as np
import pandas as pd

df = pd.read_csv('customers.csv').drop_duplicates().copy()
df['contract'] = df['contract'].str.strip().str.title()
df['total_charges'] = pd.to_numeric(df['total_charges'], errors='coerce')

NUM = ['tenure_months', 'support_calls', 'monthly_charges']
CAT = ['contract', 'internet_service', 'payment_method', 'has_dependents']
TARGET = 'churned'
from sklearn.compose import ColumnTransformer
from sklearn.pipeline import Pipeline
from sklearn.impute import SimpleImputer
from sklearn.preprocessing import StandardScaler, OneHotEncoder

def preprocessor(scale=True):
    num_steps = [('i', SimpleImputer(strategy='median'))]
    if scale:
        num_steps.append(('s', StandardScaler()))
    return ColumnTransformer([
        ('num', Pipeline(num_steps), NUM),
        ('cat', Pipeline([('i', SimpleImputer(strategy='most_frequent')),
                          ('o', OneHotEncoder(handle_unknown='ignore'))]), CAT),
    ])
from sklearn.model_selection import train_test_split

X, y = df[NUM + CAT], df[TARGET]
X_tr, X_te, y_tr, y_te = train_test_split(X, y, test_size=0.25,
                                          stratify=y, random_state=42)
from sklearn.tree import DecisionTreeClassifier
from sklearn.metrics import roc_auc_score

for criterion in ['gini', 'entropy', 'log_loss']:
    m = Pipeline([('prep', preprocessor(scale=False)),
                  ('clf', DecisionTreeClassifier(criterion=criterion,
                                                 max_depth=5,
                                                 random_state=42))]).fit(X_tr, y_tr)
    print('%-9s AUC %.4f'
          % (criterion, roc_auc_score(y_te, m.predict_proba(X_te)[:, 1])))

gini      AUC 0.7745
entropy   AUC 0.7684
log_loss  AUC 0.7684

Practically indistinguishable, which is the usual result. Do not spend time tuning the criterion; spend it on depth, which is tomorrow.

Day 1 takeaway

A tree greedily picks the split that most reduces impurity, then recurses. It is the only model here you can read end to end. It needs no scaling, ignores monotone transformations and shrugs off outliers, because it only compares values against thresholds.

Week 06 · Day 2 of 7

Depth, Pruning and Instability

Stopping a tree memorising, and the weakness that becomes tomorrow's strength

By Pawan Rama Mali1377 words

Left alone, a tree will keep splitting until every leaf is pure. That means it memorises the training set perfectly and generalises terribly.

Watch it happen

import numpy as np
import pandas as pd

df = pd.read_csv('customers.csv').drop_duplicates().copy()
df['contract'] = df['contract'].str.strip().str.title()
df['total_charges'] = pd.to_numeric(df['total_charges'], errors='coerce')

NUM = ['tenure_months', 'support_calls', 'monthly_charges']
CAT = ['contract', 'internet_service', 'payment_method', 'has_dependents']
TARGET = 'churned'
from sklearn.compose import ColumnTransformer
from sklearn.pipeline import Pipeline
from sklearn.impute import SimpleImputer
from sklearn.preprocessing import StandardScaler, OneHotEncoder

def preprocessor(scale=True):
    num_steps = [('i', SimpleImputer(strategy='median'))]
    if scale:
        num_steps.append(('s', StandardScaler()))
    return ColumnTransformer([
        ('num', Pipeline(num_steps), NUM),
        ('cat', Pipeline([('i', SimpleImputer(strategy='most_frequent')),
                          ('o', OneHotEncoder(handle_unknown='ignore'))]), CAT),
    ])
from sklearn.model_selection import train_test_split

X, y = df[NUM + CAT], df[TARGET]
X_tr, X_te, y_tr, y_te = train_test_split(X, y, test_size=0.25,
                                          stratify=y, random_state=42)
from sklearn.tree import DecisionTreeClassifier
from sklearn.metrics import roc_auc_score

print('%8s %10s %10s %8s %8s' % ('depth', 'train', 'test', 'gap', 'leaves'))
for depth in [1, 2, 3, 5, 8, 12, None]:
    m = Pipeline([('prep', preprocessor(scale=False)),
                  ('clf', DecisionTreeClassifier(max_depth=depth,
                                                 random_state=42))]).fit(X_tr, y_tr)
    tr = roc_auc_score(y_tr, m.predict_proba(X_tr)[:, 1])
    te = roc_auc_score(y_te, m.predict_proba(X_te)[:, 1])
    print('%8s %10.4f %10.4f %8.4f %8d'
          % (depth, tr, te, tr - te, m.named_steps['clf'].get_n_leaves()))

   depth      train       test      gap   leaves
       1     0.7057     0.7056   0.0001        2
       2     0.7660     0.7733  -0.0072        4
       3     0.7980     0.7774   0.0206        8
       5     0.8409     0.7745   0.0664       29
       8     0.9047     0.7309   0.1738      116
      12     0.9779     0.6763   0.3016      311
    None     1.0000     0.6291   0.3709      539

Unlimited depth gives a perfect training score and the worst test score of the lot. Hundreds of leaves, most holding a single customer. The tree has recorded the training set rather than learned from it.

The four ways to stop it

Parameter	Controls	Sensible starting point
`max_depth`	How many questions deep	3 to 8 for a single tree
`min_samples_leaf`	Smallest allowed leaf	1 to 5 percent of rows
`min_samples_split`	Smallest node worth splitting	Twice the leaf minimum
`max_leaf_nodes`	Total leaves, grown best-first	An alternative to depth

import numpy as np
import pandas as pd

df = pd.read_csv('customers.csv').drop_duplicates().copy()
df['contract'] = df['contract'].str.strip().str.title()
df['total_charges'] = pd.to_numeric(df['total_charges'], errors='coerce')

NUM = ['tenure_months', 'support_calls', 'monthly_charges']
CAT = ['contract', 'internet_service', 'payment_method', 'has_dependents']
TARGET = 'churned'
from sklearn.compose import ColumnTransformer
from sklearn.pipeline import Pipeline
from sklearn.impute import SimpleImputer
from sklearn.preprocessing import StandardScaler, OneHotEncoder

def preprocessor(scale=True):
    num_steps = [('i', SimpleImputer(strategy='median'))]
    if scale:
        num_steps.append(('s', StandardScaler()))
    return ColumnTransformer([
        ('num', Pipeline(num_steps), NUM),
        ('cat', Pipeline([('i', SimpleImputer(strategy='most_frequent')),
                          ('o', OneHotEncoder(handle_unknown='ignore'))]), CAT),
    ])
from sklearn.model_selection import train_test_split

X, y = df[NUM + CAT], df[TARGET]
X_tr, X_te, y_tr, y_te = train_test_split(X, y, test_size=0.25,
                                          stratify=y, random_state=42)
from sklearn.tree import DecisionTreeClassifier
from sklearn.metrics import roc_auc_score

print('%18s %10s %10s' % ('min_samples_leaf', 'test AUC', 'leaves'))
for leaf in [1, 5, 20, 50, 150]:
    m = Pipeline([('prep', preprocessor(scale=False)),
                  ('clf', DecisionTreeClassifier(min_samples_leaf=leaf,
                                                 random_state=42))]).fit(X_tr, y_tr)
    print('%18d %10.4f %10d'
          % (leaf, roc_auc_score(y_te, m.predict_proba(X_te)[:, 1]),
             m.named_steps['clf'].get_n_leaves()))

  min_samples_leaf   test AUC     leaves
                 1     0.6291        539
                 5     0.7260        232
                20     0.7630         75
                50     0.7915         31
               150     0.7974         12

min_samples_leaf is the better knob

Depth limits every branch equally, whether or not that branch had enough data to justify going deeper. A minimum leaf size adapts: dense regions of the data get to split further, sparse ones stop early. If you only tune one parameter on a tree, tune this one.

Cost-complexity pruning

Rather than guessing a limit in advance, grow the tree fully and then cut back the branches that do not pay for themselves.

import numpy as np
import pandas as pd

df = pd.read_csv('customers.csv').drop_duplicates().copy()
df['contract'] = df['contract'].str.strip().str.title()
df['total_charges'] = pd.to_numeric(df['total_charges'], errors='coerce')

NUM = ['tenure_months', 'support_calls', 'monthly_charges']
CAT = ['contract', 'internet_service', 'payment_method', 'has_dependents']
TARGET = 'churned'
from sklearn.compose import ColumnTransformer
from sklearn.pipeline import Pipeline
from sklearn.impute import SimpleImputer
from sklearn.preprocessing import StandardScaler, OneHotEncoder

def preprocessor(scale=True):
    num_steps = [('i', SimpleImputer(strategy='median'))]
    if scale:
        num_steps.append(('s', StandardScaler()))
    return ColumnTransformer([
        ('num', Pipeline(num_steps), NUM),
        ('cat', Pipeline([('i', SimpleImputer(strategy='most_frequent')),
                          ('o', OneHotEncoder(handle_unknown='ignore'))]), CAT),
    ])
from sklearn.model_selection import train_test_split

X, y = df[NUM + CAT], df[TARGET]
X_tr, X_te, y_tr, y_te = train_test_split(X, y, test_size=0.25,
                                          stratify=y, random_state=42)
from sklearn.tree import DecisionTreeClassifier
from sklearn.metrics import roc_auc_score
import numpy as np

prep = preprocessor(scale=False).fit(X_tr)
Xt_tr, Xt_te = prep.transform(X_tr), prep.transform(X_te)

path = DecisionTreeClassifier(random_state=42).cost_complexity_pruning_path(
    Xt_tr, y_tr)
alphas = path.ccp_alphas[:max(1, len(path.ccp_alphas) // 8)][:8]

print('%12s %10s %10s %8s' % ('ccp_alpha', 'train', 'test', 'leaves'))
for a in alphas:
    t = DecisionTreeClassifier(random_state=42, ccp_alpha=a).fit(Xt_tr, y_tr)
    print('%12.6f %10.4f %10.4f %8d'
          % (a, roc_auc_score(y_tr, t.predict_proba(Xt_tr)[:, 1]),
             roc_auc_score(y_te, t.predict_proba(Xt_te)[:, 1]),
             t.get_n_leaves()))

   ccp_alpha      train       test   leaves
000000     1.0000     0.6291      539
000111     1.0000     0.6308      537
000221     0.9999     0.6071      533
000254     0.9999     0.6093      530
000267     0.9999     0.6141      521
000267     0.9999     0.6141      521
000274     0.9998     0.6143      518
000279     0.9998     0.6199      515

As alpha rises the tree shrinks and the training score falls, while the test score improves up to a point. Pruning finds a smaller and better tree than growing an unrestricted one ever could.

A single tree is unstable

The measure that matters is not whether the printed tree looks different, but whether it predicts differently. Fit trees on overlapping samples and count how often two of them disagree about the same customer.

import numpy as np
import pandas as pd

df = pd.read_csv('customers.csv').drop_duplicates().copy()
df['contract'] = df['contract'].str.strip().str.title()
df['total_charges'] = pd.to_numeric(df['total_charges'], errors='coerce')

NUM = ['tenure_months', 'support_calls', 'monthly_charges']
CAT = ['contract', 'internet_service', 'payment_method', 'has_dependents']
TARGET = 'churned'
from sklearn.compose import ColumnTransformer
from sklearn.pipeline import Pipeline
from sklearn.impute import SimpleImputer
from sklearn.preprocessing import StandardScaler, OneHotEncoder

def preprocessor(scale=True):
    num_steps = [('i', SimpleImputer(strategy='median'))]
    if scale:
        num_steps.append(('s', StandardScaler()))
    return ColumnTransformer([
        ('num', Pipeline(num_steps), NUM),
        ('cat', Pipeline([('i', SimpleImputer(strategy='most_frequent')),
                          ('o', OneHotEncoder(handle_unknown='ignore'))]), CAT),
    ])
from sklearn.model_selection import train_test_split

X, y = df[NUM + CAT], df[TARGET]
X_tr, X_te, y_tr, y_te = train_test_split(X, y, test_size=0.25,
                                          stratify=y, random_state=42)
from sklearn.tree import DecisionTreeClassifier
from sklearn.model_selection import train_test_split
import numpy as np
from itertools import combinations

def disagreement(depth, n_models=6):
    probas = []
    for seed in range(n_models):
        Xs, _, ys, _ = train_test_split(X_tr, y_tr, train_size=0.7,
                                        random_state=seed)
        m = Pipeline([('prep', preprocessor(scale=False)),
                      ('clf', DecisionTreeClassifier(max_depth=depth,
                                                     random_state=42))]).fit(Xs, ys)
        probas.append(m.predict_proba(X_te)[:, 1])
    label = np.mean([((a > 0.5) != (b > 0.5)).mean()
                     for a, b in combinations(probas, 2)])
    prob = np.mean([np.abs(a - b).mean() for a, b in combinations(probas, 2)])
    return label, prob

print('%10s %16s %16s' % ('max_depth', 'labels differ', 'mean |p1 - p2|'))
for depth in [2, 4, 8, 12, None]:
    label, prob = disagreement(depth)
    print('%10s %15.1f%% %16.4f' % (depth, 100 * label, prob))

 max_depth    labels differ   mean |p1 - p2|
         2            19.9%           0.0342
         4            15.2%           0.0899
         8            18.8%           0.1741
        12            23.5%           0.2371
      None            23.6%           0.2366

Two trees on 70 percent of the same data, disagreeing on a quarter of customers

Read the probability column first: it rises steadily with depth, from 0.03 to 0.24. The greedy algorithm commits at each step with no way to reconsider, so once two candidate splits are close, a handful of rows decides which wins and every branch below inherits that coin toss. The deeper the tree, the more of those coin tosses accumulate.

The label column is not monotone, and the exception is worth understanding. A depth-2 tree has four leaves, so its probabilities are four numbers, extremely stable. But if one of those numbers sits near 0.5, a small shift flips an entire leaf's worth of customers at once. Stable probabilities, unstable labels. It is a compact reminder that thresholding throws information away.

Either way, this variance is what tomorrow exploits. If a deep tree disagrees with itself depending on which rows it saw, build hundreds on different samples and average them: the disagreement cancels, the agreement survives.

Day 2 takeaway

An unrestricted tree memorises its training data. Control it with min_samples_leaf in preference to max_depth, because it adapts to how much data each branch has. Cost-complexity pruning finds a better tree than guessing limits. And a single tree is unstable, which is a weakness that becomes a strength tomorrow.

Week 06 · Day 3 of 7

Bagging and Random Forests

Averaging away instability, and the validation set you get for free

By Pawan Rama Mali971 words

Yesterday ended on instability. Today that becomes the mechanism: average many unstable models and the noise cancels while the signal survives.

Bagging

Bootstrap aggregating: Draw many random samples of the training data with replacement, fit an unrestricted tree on each, and average their predictions. Each tree overfits its own sample differently, and averaging cancels the individual errors.

import numpy as np

rng = np.random.default_rng(0)
n = 1000
rows = np.arange(n)
sample = rng.choice(rows, size=n, replace=True)

unique = len(np.unique(sample))
print('bootstrap sample of %d drawn from %d rows' % (n, n))
print('distinct rows included: %d (%.1f%%)' % (unique, 100 * unique / n))
print('left out entirely:      %d (%.1f%%)' % (n - unique, 100 * (n - unique) / n))
print('\ntheory says 1 - 1/e = %.3f are included' % (1 - np.exp(-1)))

bootstrap sample of 1000 drawn from 1000 rows
distinct rows included: 622 (62.2%)
left out entirely:      378 (37.8%)

theory says 1 - 1/e = 0.632 are included

About 63 percent of rows appear in any given bootstrap sample. The 37 percent left out are that tree's out-of-bag rows, and they give you a free validation set.

Random forests add a second source of randomness

Bagging alone leaves the trees correlated: if one feature is much stronger than the rest, every tree splits on it first. A random forest fixes that by offering each split only a random subset of features.

import numpy as np
import pandas as pd

df = pd.read_csv('customers.csv').drop_duplicates().copy()
df['contract'] = df['contract'].str.strip().str.title()
df['total_charges'] = pd.to_numeric(df['total_charges'], errors='coerce')

NUM = ['tenure_months', 'support_calls', 'monthly_charges']
CAT = ['contract', 'internet_service', 'payment_method', 'has_dependents']
TARGET = 'churned'
from sklearn.compose import ColumnTransformer
from sklearn.pipeline import Pipeline
from sklearn.impute import SimpleImputer
from sklearn.preprocessing import StandardScaler, OneHotEncoder

def preprocessor(scale=True):
    num_steps = [('i', SimpleImputer(strategy='median'))]
    if scale:
        num_steps.append(('s', StandardScaler()))
    return ColumnTransformer([
        ('num', Pipeline(num_steps), NUM),
        ('cat', Pipeline([('i', SimpleImputer(strategy='most_frequent')),
                          ('o', OneHotEncoder(handle_unknown='ignore'))]), CAT),
    ])
from sklearn.model_selection import train_test_split

X, y = df[NUM + CAT], df[TARGET]
X_tr, X_te, y_tr, y_te = train_test_split(X, y, test_size=0.25,
                                          stratify=y, random_state=42)
from sklearn.ensemble import BaggingClassifier, RandomForestClassifier
from sklearn.tree import DecisionTreeClassifier
from sklearn.metrics import roc_auc_score

models = {
    'single tree ': DecisionTreeClassifier(random_state=42),
    'bagging x200': BaggingClassifier(DecisionTreeClassifier(random_state=42),
                                      n_estimators=200, random_state=42),
    'forest x200 ': RandomForestClassifier(n_estimators=200, random_state=42),
}
for name, clf in models.items():
    m = Pipeline([('prep', preprocessor(scale=False)), ('clf', clf)]).fit(X_tr, y_tr)
    print('%s AUC %.4f'
          % (name, roc_auc_score(y_te, m.predict_proba(X_te)[:, 1])))

single tree  AUC 0.6291
bagging x200 AUC 0.7699
forest x200  AUC 0.7788

The out-of-bag score

import numpy as np
import pandas as pd

df = pd.read_csv('customers.csv').drop_duplicates().copy()
df['contract'] = df['contract'].str.strip().str.title()
df['total_charges'] = pd.to_numeric(df['total_charges'], errors='coerce')

NUM = ['tenure_months', 'support_calls', 'monthly_charges']
CAT = ['contract', 'internet_service', 'payment_method', 'has_dependents']
TARGET = 'churned'
from sklearn.compose import ColumnTransformer
from sklearn.pipeline import Pipeline
from sklearn.impute import SimpleImputer
from sklearn.preprocessing import StandardScaler, OneHotEncoder

def preprocessor(scale=True):
    num_steps = [('i', SimpleImputer(strategy='median'))]
    if scale:
        num_steps.append(('s', StandardScaler()))
    return ColumnTransformer([
        ('num', Pipeline(num_steps), NUM),
        ('cat', Pipeline([('i', SimpleImputer(strategy='most_frequent')),
                          ('o', OneHotEncoder(handle_unknown='ignore'))]), CAT),
    ])
from sklearn.model_selection import train_test_split

X, y = df[NUM + CAT], df[TARGET]
X_tr, X_te, y_tr, y_te = train_test_split(X, y, test_size=0.25,
                                          stratify=y, random_state=42)
from sklearn.ensemble import RandomForestClassifier
from sklearn.metrics import roc_auc_score, accuracy_score

m = Pipeline([('prep', preprocessor(scale=False)),
              ('clf', RandomForestClassifier(n_estimators=300, oob_score=True,
                                             random_state=42))]).fit(X_tr, y_tr)
forest = m.named_steps['clf']
print('out-of-bag accuracy %.4f' % forest.oob_score_)
print('test accuracy       %.4f' % accuracy_score(y_te, m.predict(X_te)))

out-of-bag accuracy 0.7560
test accuracy       0.7587

A validation score for free

Each tree is scored on the rows it never saw, and those scores are pooled. You get something close to cross-validation without fitting anything extra. Useful when a full cross-validation would be expensive, though it is slightly pessimistic because each row is judged by only about a third of the forest.

What to tune

import numpy as np
import pandas as pd

df = pd.read_csv('customers.csv').drop_duplicates().copy()
df['contract'] = df['contract'].str.strip().str.title()
df['total_charges'] = pd.to_numeric(df['total_charges'], errors='coerce')

NUM = ['tenure_months', 'support_calls', 'monthly_charges']
CAT = ['contract', 'internet_service', 'payment_method', 'has_dependents']
TARGET = 'churned'
from sklearn.compose import ColumnTransformer
from sklearn.pipeline import Pipeline
from sklearn.impute import SimpleImputer
from sklearn.preprocessing import StandardScaler, OneHotEncoder

def preprocessor(scale=True):
    num_steps = [('i', SimpleImputer(strategy='median'))]
    if scale:
        num_steps.append(('s', StandardScaler()))
    return ColumnTransformer([
        ('num', Pipeline(num_steps), NUM),
        ('cat', Pipeline([('i', SimpleImputer(strategy='most_frequent')),
                          ('o', OneHotEncoder(handle_unknown='ignore'))]), CAT),
    ])
from sklearn.model_selection import train_test_split

X, y = df[NUM + CAT], df[TARGET]
X_tr, X_te, y_tr, y_te = train_test_split(X, y, test_size=0.25,
                                          stratify=y, random_state=42)
from sklearn.ensemble import RandomForestClassifier
from sklearn.metrics import roc_auc_score

print('%12s %10s' % ('n_estimators', 'test AUC'))
for n in [5, 25, 100, 300, 800]:
    m = Pipeline([('prep', preprocessor(scale=False)),
                  ('clf', RandomForestClassifier(n_estimators=n,
                                                 random_state=42))]).fit(X_tr, y_tr)
    print('%12d %10.4f'
          % (n, roc_auc_score(y_te, m.predict_proba(X_te)[:, 1])))

print('\n%18s %10s' % ('min_samples_leaf', 'test AUC'))
for leaf in [1, 5, 20, 50]:
    m = Pipeline([('prep', preprocessor(scale=False)),
                  ('clf', RandomForestClassifier(n_estimators=300,
                                                 min_samples_leaf=leaf,
                                                 random_state=42))]).fit(X_tr, y_tr)
    print('%18d %10.4f'
          % (leaf, roc_auc_score(y_te, m.predict_proba(X_te)[:, 1])))

n_estimators   test AUC
           5     0.7205
          25     0.7671
         100     0.7728
         300     0.7776
         800     0.7779

  min_samples_leaf   test AUC
                 1     0.7776
                 5     0.7989
                20     0.8098
                50     0.8158

More trees never overfit; deeper trees still can

Adding trees to a forest cannot make it worse, the score plateaus and stays there, so n_estimators is a budget decision, not a tuning decision. The parameters that genuinely matter are min_samples_leaf and max_features. This is the opposite of boosting, where more rounds absolutely can overfit.

Day 3 takeaway

A bootstrap sample contains about 63 percent of the rows, leaving the rest as a free validation set. Bagging averages away the instability of individual trees; a random forest goes further by decorrelating them through random feature subsets. More trees never hurt, so tune leaf size instead.

Week 06 · Day 4 of 7

Boosting

Fitting each model to the last one's mistakes, implemented then industrialised

By Pawan Rama Mali1034 words

A forest builds hundreds of independent trees and averages them. Boosting builds them in sequence, each one working on what the previous ones got wrong.

The idea, implemented

import numpy as np
from sklearn.tree import DecisionTreeRegressor

rng = np.random.default_rng(0)
x = rng.uniform(0, 6, 300).reshape(-1, 1)
y = np.sin(1.5 * x.ravel()) + 0.3 * x.ravel() + rng.normal(0, 0.2, 300)

prediction = np.full(len(y), y.mean())    # start with the mean
lr = 0.3

print('%8s %12s' % ('round', 'train MSE'))
print('%8s %12.5f' % ('start', ((y - prediction) ** 2).mean()))
for r in range(1, 41):
    residual = y - prediction              # what we still get wrong
    stump = DecisionTreeRegressor(max_depth=2).fit(x, residual)
    prediction = prediction + lr * stump.predict(x)
    if r in (1, 5, 10, 20, 40):
        print('%8d %12.5f' % (r, ((y - prediction) ** 2).mean()))

   round    train MSE
   start      0.86697
       1      0.49189
       5      0.10329
      10      0.04882
      20      0.03073
      40      0.02490

Gradient boosting: Fit a weak model to the errors of everything built so far, add a fraction of it to the running prediction, and repeat. Each round is a step down the gradient of the loss, which is where the name comes from, and why week 3's gradient descent was worth implementing.

Thirty lines, and it is the core of the algorithm that wins most tabular competitions.

The library versions

import numpy as np
import pandas as pd

df = pd.read_csv('customers.csv').drop_duplicates().copy()
df['contract'] = df['contract'].str.strip().str.title()
df['total_charges'] = pd.to_numeric(df['total_charges'], errors='coerce')

NUM = ['tenure_months', 'support_calls', 'monthly_charges']
CAT = ['contract', 'internet_service', 'payment_method', 'has_dependents']
TARGET = 'churned'
from sklearn.compose import ColumnTransformer
from sklearn.pipeline import Pipeline
from sklearn.impute import SimpleImputer
from sklearn.preprocessing import StandardScaler, OneHotEncoder

def preprocessor(scale=True):
    num_steps = [('i', SimpleImputer(strategy='median'))]
    if scale:
        num_steps.append(('s', StandardScaler()))
    return ColumnTransformer([
        ('num', Pipeline(num_steps), NUM),
        ('cat', Pipeline([('i', SimpleImputer(strategy='most_frequent')),
                          ('o', OneHotEncoder(handle_unknown='ignore'))]), CAT),
    ])
from sklearn.model_selection import train_test_split

X, y = df[NUM + CAT], df[TARGET]
X_tr, X_te, y_tr, y_te = train_test_split(X, y, test_size=0.25,
                                          stratify=y, random_state=42)
from sklearn.ensemble import (AdaBoostClassifier, GradientBoostingClassifier,
                              HistGradientBoostingClassifier,
                              RandomForestClassifier)
from sklearn.metrics import roc_auc_score
import time

models = {
    'random forest': RandomForestClassifier(n_estimators=300, random_state=42),
    'adaboost     ': AdaBoostClassifier(n_estimators=200, random_state=42),
    'grad boosting': GradientBoostingClassifier(random_state=42),
    'hist grad    ': HistGradientBoostingClassifier(random_state=42),
}
for name, clf in models.items():
    t = time.perf_counter()
    m = Pipeline([('prep', preprocessor(scale=False)), ('clf', clf)]).fit(X_tr, y_tr)
    el = time.perf_counter() - t
    print('%s AUC %.4f   %.2fs'
          % (name, roc_auc_score(y_te, m.predict_proba(X_te)[:, 1]), el))

random forest AUC 0.7776   2.00s
adaboost      AUC 0.8098   1.19s
grad boosting AUC 0.8034   0.74s
hist grad     AUC 0.7797   5.20s

Learning rate against number of rounds

import numpy as np
import pandas as pd

df = pd.read_csv('customers.csv').drop_duplicates().copy()
df['contract'] = df['contract'].str.strip().str.title()
df['total_charges'] = pd.to_numeric(df['total_charges'], errors='coerce')

NUM = ['tenure_months', 'support_calls', 'monthly_charges']
CAT = ['contract', 'internet_service', 'payment_method', 'has_dependents']
TARGET = 'churned'
from sklearn.compose import ColumnTransformer
from sklearn.pipeline import Pipeline
from sklearn.impute import SimpleImputer
from sklearn.preprocessing import StandardScaler, OneHotEncoder

def preprocessor(scale=True):
    num_steps = [('i', SimpleImputer(strategy='median'))]
    if scale:
        num_steps.append(('s', StandardScaler()))
    return ColumnTransformer([
        ('num', Pipeline(num_steps), NUM),
        ('cat', Pipeline([('i', SimpleImputer(strategy='most_frequent')),
                          ('o', OneHotEncoder(handle_unknown='ignore'))]), CAT),
    ])
from sklearn.model_selection import train_test_split

X, y = df[NUM + CAT], df[TARGET]
X_tr, X_te, y_tr, y_te = train_test_split(X, y, test_size=0.25,
                                          stratify=y, random_state=42)
from sklearn.ensemble import GradientBoostingClassifier
from sklearn.metrics import roc_auc_score

print('%8s %8s %10s %10s' % ('lr', 'rounds', 'train', 'test'))
for lr, n in [(1.0, 200), (0.3, 200), (0.1, 200), (0.05, 200), (0.05, 800)]:
    m = Pipeline([('prep', preprocessor(scale=False)),
                  ('clf', GradientBoostingClassifier(learning_rate=lr,
                                                     n_estimators=n,
                                                     random_state=42))]).fit(X_tr, y_tr)
    print('%8s %8d %10.4f %10.4f'
          % (lr, n, roc_auc_score(y_tr, m.predict_proba(X_tr)[:, 1]),
             roc_auc_score(y_te, m.predict_proba(X_te)[:, 1])))

      lr   rounds      train       test
0      200     0.9998     0.7410
3      200     0.9797     0.7774
1      200     0.9157     0.7968
05      200     0.8857     0.8082
05      800     0.9547     0.7843

Boosting can and will overfit

At a learning rate of 1.0 the training AUC runs away toward 1.0 while the test score falls. Unlike a forest, more rounds is not free. The two parameters trade against each other: halve the learning rate and you need roughly twice the rounds, which costs time but usually generalises better. A small rate with early stopping is the standard recipe.

Early stopping

import numpy as np
import pandas as pd

df = pd.read_csv('customers.csv').drop_duplicates().copy()
df['contract'] = df['contract'].str.strip().str.title()
df['total_charges'] = pd.to_numeric(df['total_charges'], errors='coerce')

NUM = ['tenure_months', 'support_calls', 'monthly_charges']
CAT = ['contract', 'internet_service', 'payment_method', 'has_dependents']
TARGET = 'churned'
from sklearn.compose import ColumnTransformer
from sklearn.pipeline import Pipeline
from sklearn.impute import SimpleImputer
from sklearn.preprocessing import StandardScaler, OneHotEncoder

def preprocessor(scale=True):
    num_steps = [('i', SimpleImputer(strategy='median'))]
    if scale:
        num_steps.append(('s', StandardScaler()))
    return ColumnTransformer([
        ('num', Pipeline(num_steps), NUM),
        ('cat', Pipeline([('i', SimpleImputer(strategy='most_frequent')),
                          ('o', OneHotEncoder(handle_unknown='ignore'))]), CAT),
    ])
from sklearn.model_selection import train_test_split

X, y = df[NUM + CAT], df[TARGET]
X_tr, X_te, y_tr, y_te = train_test_split(X, y, test_size=0.25,
                                          stratify=y, random_state=42)
from sklearn.ensemble import HistGradientBoostingClassifier
from sklearn.metrics import roc_auc_score

m = Pipeline([('prep', preprocessor(scale=False)),
              ('clf', HistGradientBoostingClassifier(
                  max_iter=1000, learning_rate=0.05,
                  early_stopping=True, validation_fraction=0.15,
                  n_iter_no_change=20, random_state=42))]).fit(X_tr, y_tr)
clf = m.named_steps['clf']
print('allowed up to 1000 rounds, stopped at %d' % clf.n_iter_)
print('test AUC %.4f' % roc_auc_score(y_te, m.predict_proba(X_te)[:, 1]))

allowed up to 1000 rounds, stopped at 80
test AUC 0.7924

It holds back 15 percent of the training data, watches the score on it, and stops when twenty rounds pass without improvement. You set a generous ceiling and let the model decide, which is more reliable than tuning the round count by hand.

	Random forest	Gradient boosting
Trees built	Independently, in parallel	Sequentially, each on the last errors
Individual trees	Deep, overfitting	Shallow, underfitting
More trees	Never hurts	Eventually overfits
Key parameters	`min_samples_leaf`, `max_features`	`learning_rate`, rounds, depth
Tuning effort	Low	Higher, and it pays
Typical winner	When you have no time to tune	When you do

Day 4 takeaway

Boosting fits each new model to the errors of the ensemble so far, which is gradient descent in function space. It generally beats a forest but must be tuned: learning rate and round count trade against each other, and unlike a forest it can overfit. Use a small rate with early stopping.

Week 06 · Day 5 of 7

XGBoost and LightGBM

Regularised, histogram-based boosting that handles gaps natively

By Pawan Rama Mali864 words

XGBoost and LightGBM are gradient boosting done properly, faster, regularised, and with native handling of the things that make real data awkward.

What they add

Regularisation on the trees themselves, so leaf values cannot grow unchecked.
Histogram-based splitting: bucket each feature into 256 bins instead of testing every value, which is dramatically faster.
Native missing-value handling: each split learns which way a missing value should go.
Serious engineering: multi-core, cache-aware, and able to work out of memory.

import numpy as np
import pandas as pd

df = pd.read_csv('customers.csv').drop_duplicates().copy()
df['contract'] = df['contract'].str.strip().str.title()
df['total_charges'] = pd.to_numeric(df['total_charges'], errors='coerce')

NUM = ['tenure_months', 'support_calls', 'monthly_charges']
CAT = ['contract', 'internet_service', 'payment_method', 'has_dependents']
TARGET = 'churned'
from sklearn.compose import ColumnTransformer
from sklearn.pipeline import Pipeline
from sklearn.impute import SimpleImputer
from sklearn.preprocessing import StandardScaler, OneHotEncoder

def preprocessor(scale=True):
    num_steps = [('i', SimpleImputer(strategy='median'))]
    if scale:
        num_steps.append(('s', StandardScaler()))
    return ColumnTransformer([
        ('num', Pipeline(num_steps), NUM),
        ('cat', Pipeline([('i', SimpleImputer(strategy='most_frequent')),
                          ('o', OneHotEncoder(handle_unknown='ignore'))]), CAT),
    ])
from sklearn.model_selection import train_test_split

X, y = df[NUM + CAT], df[TARGET]
X_tr, X_te, y_tr, y_te = train_test_split(X, y, test_size=0.25,
                                          stratify=y, random_state=42)
import xgboost as xgb
import lightgbm as lgb
from sklearn.ensemble import HistGradientBoostingClassifier
from sklearn.metrics import roc_auc_score
import time

models = {
    'sklearn hist': HistGradientBoostingClassifier(random_state=42),
    'xgboost     ': xgb.XGBClassifier(n_estimators=300, learning_rate=0.1,
                                      max_depth=4, random_state=42,
                                      eval_metric='logloss'),
    'lightgbm    ': lgb.LGBMClassifier(n_estimators=300, learning_rate=0.1,
                                       max_depth=4, random_state=42,
                                       verbose=-1),
}
for name, clf in models.items():
    t = time.perf_counter()
    m = Pipeline([('prep', preprocessor(scale=False)), ('clf', clf)]).fit(X_tr, y_tr)
    el = time.perf_counter() - t
    print('%s AUC %.4f   %.2fs'
          % (name, roc_auc_score(y_te, m.predict_proba(X_te)[:, 1]), el))

sklearn hist AUC 0.7797   4.45s
xgboost      AUC 0.7765   3.21s
lightgbm     AUC 0.7755   1.40s

Missing values need no imputation

import numpy as np
import pandas as pd

df = pd.read_csv('customers.csv').drop_duplicates().copy()
df['contract'] = df['contract'].str.strip().str.title()
df['total_charges'] = pd.to_numeric(df['total_charges'], errors='coerce')

NUM = ['tenure_months', 'support_calls', 'monthly_charges']
CAT = ['contract', 'internet_service', 'payment_method', 'has_dependents']
TARGET = 'churned'
from sklearn.model_selection import train_test_split

X, y = df[NUM + CAT], df[TARGET]
X_tr, X_te, y_tr, y_te = train_test_split(X, y, test_size=0.25,
                                          stratify=y, random_state=42)
import xgboost as xgb
from sklearn.metrics import roc_auc_score

# No imputer at all -- the raw numeric columns, gaps included.
Xn_tr = X_tr[NUM]
Xn_te = X_te[NUM]
print('missing values passed straight in: %d' % Xn_tr.isna().sum().sum())

m = xgb.XGBClassifier(n_estimators=300, max_depth=4, learning_rate=0.1,
                      random_state=42, eval_metric='logloss').fit(Xn_tr, y_tr)
print('AUC %.4f' % roc_auc_score(y_te, m.predict_proba(Xn_te)[:, 1]))

missing values passed straight in: 133
AUC 0.7237

Learned, not guessed

At each split the algorithm tries sending missing values left and sending them right, and keeps whichever reduces the loss more. That is strictly better than imputing a median, because it lets missingness carry its own meaning, the point week 2 made about informative gaps, handled automatically.

The parameters that matter

Parameter	Does	Try
`n_estimators`	Boosting rounds	300 to 2000 with early stopping
`learning_rate`	Step size per round	0.01 to 0.1
`max_depth`	Tree depth	3 to 8
`subsample`	Row fraction per tree	0.6 to 1.0
`colsample_bytree`	Feature fraction per tree	0.6 to 1.0
`reg_lambda`	L2 on leaf weights	1 to 10
`min_child_weight`	Minimum weight in a leaf	1 to 20

import numpy as np
import pandas as pd

df = pd.read_csv('customers.csv').drop_duplicates().copy()
df['contract'] = df['contract'].str.strip().str.title()
df['total_charges'] = pd.to_numeric(df['total_charges'], errors='coerce')

NUM = ['tenure_months', 'support_calls', 'monthly_charges']
CAT = ['contract', 'internet_service', 'payment_method', 'has_dependents']
TARGET = 'churned'
from sklearn.compose import ColumnTransformer
from sklearn.pipeline import Pipeline
from sklearn.impute import SimpleImputer
from sklearn.preprocessing import StandardScaler, OneHotEncoder

def preprocessor(scale=True):
    num_steps = [('i', SimpleImputer(strategy='median'))]
    if scale:
        num_steps.append(('s', StandardScaler()))
    return ColumnTransformer([
        ('num', Pipeline(num_steps), NUM),
        ('cat', Pipeline([('i', SimpleImputer(strategy='most_frequent')),
                          ('o', OneHotEncoder(handle_unknown='ignore'))]), CAT),
    ])
from sklearn.model_selection import train_test_split

X, y = df[NUM + CAT], df[TARGET]
X_tr, X_te, y_tr, y_te = train_test_split(X, y, test_size=0.25,
                                          stratify=y, random_state=42)
import xgboost as xgb
from sklearn.metrics import roc_auc_score

configs = [
    ('defaults      ', {}),
    ('shallow, slow ', dict(max_depth=3, learning_rate=0.03, n_estimators=600)),
    ('subsampled    ', dict(max_depth=4, learning_rate=0.05, n_estimators=500,
                            subsample=0.8, colsample_bytree=0.8)),
    ('regularised   ', dict(max_depth=5, learning_rate=0.05, n_estimators=500,
                            reg_lambda=10, min_child_weight=10)),
]
for name, kw in configs:
    clf = xgb.XGBClassifier(random_state=42, eval_metric='logloss', **kw)
    m = Pipeline([('prep', preprocessor(scale=False)), ('clf', clf)]).fit(X_tr, y_tr)
    tr = roc_auc_score(y_tr, m.predict_proba(X_tr)[:, 1])
    te = roc_auc_score(y_te, m.predict_proba(X_te)[:, 1])
    print('%s train %.4f  test %.4f  gap %.4f' % (name, tr, te, tr - te))

defaults       train 0.9940  test 0.7574  gap 0.2366
shallow, slow  train 0.8913  test 0.8020  gap 0.0893
subsampled     train 0.9507  test 0.7804  gap 0.1703
regularised    train 0.9071  test 0.7905  gap 0.1166

Read the gap column. The default configuration overfits hardest; the shallow, slow and regularised versions trade training score for a smaller gap. The best test score of the four, 0.8020, comes from the slowest and shallowest configuration, and tomorrow shows that even that does not settle the question.

Day 5 takeaway

XGBoost and LightGBM add tree-level regularisation, histogram splitting and learned missing-value directions. Feed them raw numerics with gaps intact. Tune the learning rate and depth first, then subsampling and the L2 penalty, and watch the train-test gap rather than the test score alone.

Week 06 · Day 6 of 7

Feature Importance, and Why the Default Lies

Impurity importance ranking pure noise first, and the honest alternative

By Pawan Rama Mali910 words

Every tree model offers a feature_importances_ attribute. It is the most-quoted and least-trustworthy number in applied machine learning.

What the default importance measures

import numpy as np
import pandas as pd

df = pd.read_csv('customers.csv').drop_duplicates().copy()
df['contract'] = df['contract'].str.strip().str.title()
df['total_charges'] = pd.to_numeric(df['total_charges'], errors='coerce')

NUM = ['tenure_months', 'support_calls', 'monthly_charges']
CAT = ['contract', 'internet_service', 'payment_method', 'has_dependents']
TARGET = 'churned'
from sklearn.compose import ColumnTransformer
from sklearn.pipeline import Pipeline
from sklearn.impute import SimpleImputer
from sklearn.preprocessing import StandardScaler, OneHotEncoder

def preprocessor(scale=True):
    num_steps = [('i', SimpleImputer(strategy='median'))]
    if scale:
        num_steps.append(('s', StandardScaler()))
    return ColumnTransformer([
        ('num', Pipeline(num_steps), NUM),
        ('cat', Pipeline([('i', SimpleImputer(strategy='most_frequent')),
                          ('o', OneHotEncoder(handle_unknown='ignore'))]), CAT),
    ])
from sklearn.model_selection import train_test_split

X, y = df[NUM + CAT], df[TARGET]
X_tr, X_te, y_tr, y_te = train_test_split(X, y, test_size=0.25,
                                          stratify=y, random_state=42)
from sklearn.ensemble import RandomForestClassifier
import pandas as pd

m = Pipeline([('prep', preprocessor(scale=False)),
              ('clf', RandomForestClassifier(n_estimators=300,
                                             random_state=42))]).fit(X_tr, y_tr)
names = m.named_steps['prep'].get_feature_names_out()
imp = pd.Series(m.named_steps['clf'].feature_importances_, index=names)
print(imp.sort_values(ascending=False).head(8).round(4).to_string())

num__monthly_charges                    0.3216
num__tenure_months                      0.3167
num__support_calls                      0.1009
cat__contract_Month-To-Month            0.0701
cat__contract_Two Year                  0.0359
cat__payment_method_Electronic check    0.0213
cat__contract_One Year                  0.0185
cat__internet_service_Fibre optic       0.0180

Impurity importance: The total reduction in impurity contributed by each feature, summed over every split in every tree. It is free to compute, which is why it is the default, and it is biased in a way that is easy to demonstrate.

The bias, demonstrated

import numpy as np
import pandas as pd
from sklearn.ensemble import RandomForestClassifier

rng = np.random.default_rng(0)
n = 2000
y = rng.integers(0, 2, n)

X = pd.DataFrame({
    'real_signal': y + rng.normal(0, 1.2, n),   # genuinely predictive
    'binary_noise': rng.integers(0, 2, n),      # pure noise, 2 values
    'id_noise': rng.permutation(n),             # pure noise, 2000 values
})

m = RandomForestClassifier(n_estimators=300, random_state=42).fit(X, y)
print(pd.Series(m.feature_importances_, index=X.columns).round(4).to_string())

real_signal     0.5739
binary_noise    0.0035
id_noise        0.4225

It rewards high-cardinality columns for being noisy

id_noise is a random permutation. It contains no information whatsoever, and impurity importance ranks it comparably to the column that genuinely predicts the target, far above the binary noise column. A feature with thousands of distinct values offers thousands of candidate split points, so by chance some of them separate the training rows well. The importance measures opportunity to overfit, not usefulness.

Permutation importance is the honest one

Permutation importance: Shuffle one column, re-score the model, and see how much the score drops. If nothing happens, the model was not relying on that column. Measured on held-out data, it answers the question you actually asked.

import numpy as np
import pandas as pd
from sklearn.ensemble import RandomForestClassifier
from sklearn.inspection import permutation_importance
from sklearn.model_selection import train_test_split

rng = np.random.default_rng(0)
n = 2000
y = rng.integers(0, 2, n)
X = pd.DataFrame({
    'real_signal': y + rng.normal(0, 1.2, n),
    'binary_noise': rng.integers(0, 2, n),
    'id_noise': rng.permutation(n),
})
X_tr, X_te, y_tr, y_te = train_test_split(X, y, test_size=0.3, random_state=0)
m = RandomForestClassifier(n_estimators=300, random_state=42).fit(X_tr, y_tr)

r = permutation_importance(m, X_te, y_te, n_repeats=20, random_state=0,
                           scoring='roc_auc')
print(pd.DataFrame({'importance': r.importances_mean,
                    'std': r.importances_std},
                   index=X.columns).round(4).to_string())

              importance     std
real_signal       0.1876  0.0209
binary_noise     -0.0009  0.0128
id_noise          0.0107  0.0107

Now both noise columns sit within one standard deviation of zero, where they belong, while the real signal is an order of magnitude above them. Note that id_noise has not vanished entirely, the model did latch onto it slightly, and permutation importance correctly reports that as small rather than pretending it is absent. Same model, same data, a different and better question.

On the churn model

import numpy as np
import pandas as pd

df = pd.read_csv('customers.csv').drop_duplicates().copy()
df['contract'] = df['contract'].str.strip().str.title()
df['total_charges'] = pd.to_numeric(df['total_charges'], errors='coerce')

NUM = ['tenure_months', 'support_calls', 'monthly_charges']
CAT = ['contract', 'internet_service', 'payment_method', 'has_dependents']
TARGET = 'churned'
from sklearn.compose import ColumnTransformer
from sklearn.pipeline import Pipeline
from sklearn.impute import SimpleImputer
from sklearn.preprocessing import StandardScaler, OneHotEncoder

def preprocessor(scale=True):
    num_steps = [('i', SimpleImputer(strategy='median'))]
    if scale:
        num_steps.append(('s', StandardScaler()))
    return ColumnTransformer([
        ('num', Pipeline(num_steps), NUM),
        ('cat', Pipeline([('i', SimpleImputer(strategy='most_frequent')),
                          ('o', OneHotEncoder(handle_unknown='ignore'))]), CAT),
    ])
from sklearn.model_selection import train_test_split

X, y = df[NUM + CAT], df[TARGET]
X_tr, X_te, y_tr, y_te = train_test_split(X, y, test_size=0.25,
                                          stratify=y, random_state=42)
from sklearn.ensemble import RandomForestClassifier
from sklearn.inspection import permutation_importance
import pandas as pd

m = Pipeline([('prep', preprocessor(scale=False)),
              ('clf', RandomForestClassifier(n_estimators=300,
                                             random_state=42))]).fit(X_tr, y_tr)
r = permutation_importance(m, X_te, y_te, n_repeats=15, random_state=0,
                           scoring='roc_auc')
print(pd.Series(r.importances_mean, index=X_te.columns)
      .sort_values(ascending=False).round(4).to_string())

contract            0.1478
tenure_months       0.0667
monthly_charges     0.0171
support_calls       0.0127
payment_method      0.0111
internet_service    0.0102
has_dependents      0.0082

Permuting whole original columns, before encoding, which is what you actually want to know: is contract useful, not is contract_Two Year useful.

Correlated features share the blame

Shuffle tenure_months while total_charges still carries most of the same information and the score barely moves, so both look unimportant. Permutation importance answers “how much does this model rely on this column”, never “how much does this column matter in the world”. Group correlated columns and permute them together if you need the second answer.

Day 6 takeaway

The default feature_importances_ rewards high-cardinality columns and can rank pure noise above real signal. Use permutation importance on held-out data instead, permute original columns rather than encoded ones, and remember that correlated features will divide the credit between them.

Week 06 · Day 7 of 7

Why the Simple Model Won

The bake-off, and matching an algorithm to the shape of the truth

By Pawan Rama Mali1000 words

Week 1 set you an assignment: swap logistic regression for a random forest and explain why the forest scored worse. Here is the answer, and it matters more than any algorithm this week.

The bake-off

import numpy as np
import pandas as pd

df = pd.read_csv('customers.csv').drop_duplicates().copy()
df['contract'] = df['contract'].str.strip().str.title()
df['total_charges'] = pd.to_numeric(df['total_charges'], errors='coerce')

NUM = ['tenure_months', 'support_calls', 'monthly_charges']
CAT = ['contract', 'internet_service', 'payment_method', 'has_dependents']
TARGET = 'churned'
from sklearn.compose import ColumnTransformer
from sklearn.pipeline import Pipeline
from sklearn.impute import SimpleImputer
from sklearn.preprocessing import StandardScaler, OneHotEncoder

def preprocessor(scale=True):
    num_steps = [('i', SimpleImputer(strategy='median'))]
    if scale:
        num_steps.append(('s', StandardScaler()))
    return ColumnTransformer([
        ('num', Pipeline(num_steps), NUM),
        ('cat', Pipeline([('i', SimpleImputer(strategy='most_frequent')),
                          ('o', OneHotEncoder(handle_unknown='ignore'))]), CAT),
    ])
from sklearn.model_selection import train_test_split

X, y = df[NUM + CAT], df[TARGET]
X_tr, X_te, y_tr, y_te = train_test_split(X, y, test_size=0.25,
                                          stratify=y, random_state=42)
from sklearn.dummy import DummyClassifier
from sklearn.linear_model import LogisticRegression
from sklearn.tree import DecisionTreeClassifier
from sklearn.ensemble import RandomForestClassifier, HistGradientBoostingClassifier
import xgboost as xgb
from sklearn.model_selection import cross_val_score, StratifiedKFold

cv = StratifiedKFold(n_splits=5, shuffle=True, random_state=42)
candidates = {
    'baseline     ': DummyClassifier(strategy='most_frequent'),
    'logistic     ': LogisticRegression(max_iter=1000, random_state=42),
    'tree leaf=50 ': DecisionTreeClassifier(min_samples_leaf=50, random_state=42),
    'forest       ': RandomForestClassifier(n_estimators=300, min_samples_leaf=5,
                                            random_state=42),
    'hist boosting': HistGradientBoostingClassifier(random_state=42),
    'xgboost      ': xgb.XGBClassifier(n_estimators=400, max_depth=3,
                                       learning_rate=0.05, subsample=0.8,
                                       random_state=42, eval_metric='logloss'),
}
for name, clf in candidates.items():
    s = cross_val_score(Pipeline([('prep', preprocessor()), ('clf', clf)]),
                        X_tr, y_tr, cv=cv, scoring='roc_auc')
    print('%s %.4f  +/- %.4f' % (name, s.mean(), s.std()))

baseline      0.5000  +/- 0.0000
logistic      0.8208  +/- 0.0137
tree leaf=50  0.7922  +/- 0.0078
forest        0.8203  +/- 0.0101
hist boosting 0.7973  +/- 0.0088
xgboost       0.8133  +/- 0.0083

A week of ensembles, and nothing beats the straight line

Logistic regression scores 0.8208. The best ensemble, a random forest with min_samples_leaf=5: scores 0.8203, which is a tie: the difference is a twentieth of the fold-to-fold spread. XGBoost is a point behind, and boosting with default settings is two points behind. Nothing here is a bug, and nothing you did is wrong. The honest summary is that a week of increasingly sophisticated machinery has managed to draw with the simplest model in the course.

Why

Look at how the data was generated, back in week 1:

# From make_dataset.py, week 1 day 3
#
#   logit = (-2.8
#            + 1.55 * (contract == 'Month-to-month')
#            - 0.85 * (contract == 'Two year')
#            - 0.055 * tenure
#            + 0.021 * monthly
#            + 0.30 * support
#            + 0.42 * (payment == 'Electronic check')
#            - 0.25 * (dependents == 'Yes'))
#   churned = sigmoid(logit) > uniform()

print('The truth is a weighted sum, passed through a sigmoid.')
print('That is exactly, precisely, the functional form of')
print('logistic regression -- and nothing else in this course.')

The truth is a weighted sum, passed through a sigmoid.
That is exactly, precisely, the functional form of
logistic regression -- and nothing else in this course.

Logistic regression is not approximating the truth here. It is the truth, with the right coefficients to be found. A tree can only approximate a smooth linear relationship with a staircase of splits, and no amount of boosting turns a staircase into a straight line.

Where trees win instead

import numpy as np
from sklearn.linear_model import LogisticRegression
from sklearn.ensemble import HistGradientBoostingClassifier
from sklearn.preprocessing import StandardScaler
from sklearn.pipeline import make_pipeline
from sklearn.model_selection import cross_val_score

rng = np.random.default_rng(0)
n = 3000
a = rng.uniform(-3, 3, n)
b = rng.uniform(-3, 3, n)

# An interaction with a threshold -- no linear model can express this.
y = (((a > 0) & (b > 0)) | ((a < -1.5) & (b < -1.5))).astype(int)
y = np.where(rng.uniform(size=n) < 0.05, 1 - y, y)      # 5% label noise
X = np.column_stack([a, b])

for name, clf in [('logistic', LogisticRegression(max_iter=1000)),
                  ('boosting', HistGradientBoostingClassifier(random_state=42))]:
    s = cross_val_score(make_pipeline(StandardScaler(), clf), X, y,
                        cv=5, scoring='roc_auc')
    print('%-9s AUC %.4f' % (name, s.mean()))

logistic  AUC 0.7184
boosting  AUC 0.9352

Same two features, a rule built from thresholds and interactions rather than a weighted sum, and the result reverses completely. Neither algorithm is better. They match different shapes of truth.

Prefer a linear model when	Prefer trees when
Effects are additive and smooth	Thresholds and interactions dominate
You need interpretable coefficients	You need raw predictive power
Data is wide and rows are few	Rows are plentiful
You must extrapolate beyond the training range	All prediction is inside the training range
Probabilities must be well calibrated	Ranking is what matters

Trees cannot extrapolate at all

A tree's prediction for any input beyond the training range is whatever the nearest leaf says, forever. Give a forest trained on tenures of 1 to 72 a customer with 200 months and it predicts as though they had 72. A linear model would continue the trend. Which behaviour is right depends entirely on your problem, but you must know which one you are getting.

import numpy as np
from sklearn.linear_model import LinearRegression
from sklearn.ensemble import RandomForestRegressor

x = np.linspace(0, 10, 200).reshape(-1, 1)
y = 3 * x.ravel() + 5

lin = LinearRegression().fit(x, y)
forest = RandomForestRegressor(n_estimators=100, random_state=0).fit(x, y)

print('%10s %14s %14s %12s' % ('x', 'truth', 'linear', 'forest'))
for v in [5.0, 10.0, 15.0, 50.0]:
    print('%10.1f %14.2f %14.2f %12.2f'
          % (v, 3 * v + 5, lin.predict([[v]])[0], forest.predict([[v]])[0]))

         x          truth         linear       forest
       5.0          20.00          20.00        19.97
      10.0          35.00          35.00        34.92
      15.0          50.00          50.00        34.92
      50.0         155.00         155.00        34.92

Your assignment

Add an interaction the linear model cannot see, for example a flag for month-to-month AND support_calls > 2: to the churn features, and rerun the bake-off. Does the gap between logistic regression and boosting narrow? Then engineer that same interaction as an explicit column for the linear model and see who wins. That is week 8.

Day 7 takeaway

The best algorithm is the one whose shape matches the truth in your data. Logistic regression wins here because the data was generated by a logistic function; boosting wins when thresholds and interactions dominate. Trees cannot extrapolate. Always run the simple model first. It is sometimes the answer, and it is always the baseline.