36 - Random Forest: perché cento alberi sbagliano meno di uno

Chiedi a una persona di stimare il peso di un bue a una fiera: sbaglierà, in un senso o nell'altro. Chiedilo a cento persone indipendenti e fai la media: gli errori in eccesso e in difetto tendono a compensarsi, e la media è molto più vicina alla verità di quasi tutte le singole stime. È la «saggezza della folla». Ma funziona solo a una condizione: gli errori devono essere poco correlati. Se cento persone copiano dallo stesso vicino, la media sbaglia come il vicino.

La lezione sugli alberi decisionali finiva con una cattiva notizia: un singolo albero è instabile, cambia molto se cambia di poco il campione di dati. La Random Forest trasforma quel difetto in un vantaggio: addestra molti alberi diversi e ne fa la media. In questa lezione deriviamo perché funziona (con la formula della varianza della media), la costruiamo a mano e la confrontiamo con la libreria.

Due parole da chiarire: ensemble, bias e varianza

Un ensemble (insieme) è un modello costruito combinando molti modelli più semplici, detti base learner. Per capire perché combinarli aiuta, serve scomporre l'errore. Immagina di giocare a freccette: il bias è quanto, in media, i tuoi tiri sono spostati rispetto al centro del bersaglio (errore sistematico: la mira è storta); la varianza è quanto i tiri sono sparpagliati tra loro (errore casuale: la mano trema). Un albero profondo è come un tiratore con la mira dritta ma la mano tremante: in media colpisce il centro, ma ogni singolo tiro finisce dove capita. Fare la media di molti tiratori cancella il tremore.

Perché la media di molti modelli sbaglia meno

L'errore di un modello si scompone in bias (errore sistematico, il modello è troppo semplice) e varianza (quanto la previsione cambia al cambiare del campione di addestramento). Un albero profondo ha bias basso e varianza alta: è l'identikit di un modello che farebbe bene a essere mediato.

Prendi B previsioni X₁...X_B, ciascuna con varianza σ² e con correlazione a coppie ρ. La varianza della loro media è:

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

Il primo passaggio conta le B varianze e le B(B−1) covarianze (ognuna vale ρσ²). Da questa formula si leggono due leve:

• Il termine (1−ρ)σ²/B scompare aumentando B: più alberi, meno varianza. Questo è il motivo per cui aggiungere alberi non produce overfitting.
• Il termine ρσ² non dipende da B: se gli alberi sono molto correlati, mediare serve a poco. L'unico modo per abbassarlo è ridurre la correlazione ρ tra gli alberi.

Le due leve corrispondono a due modi diversi di sbagliare. Con pochi alberi sbagli per la dispersione del primo termine, e la curi con la quantità. Con alberi troppo simili sbagli per il secondo termine, che nessuna quantità cura. Ne segue che «più alberi» è una ricetta incompleta: la parte nuova di una Random Forest rispetto a una semplice media è tutto ciò che rende gli alberi diversi, e sono due ingredienti di casualità.

Ingrediente 1: bagging (bootstrap aggregating)

Ogni albero viene addestrato su un campione bootstrap: n esempi estratti con reinserimento dal training set di n esempi. Alcuni esempi compaiono più volte, altri mai. Quanti restano fuori? La probabilità che un dato esempio non venga estratto in una singola estrazione è 1 − 1/n; in n estrazioni indipendenti è (1 − 1/n)ⁿ, che per n grande tende a 1/e ≈ 0,368. Quindi ogni albero non vede circa il 37% degli esempi. Verifichiamolo:

import numpy as np

rng = np.random.default_rng(0)
for n in (10, 100, 1000, 100000):
    fuori = [len(np.setdiff1d(np.arange(n), rng.integers(0, n, size=n))) / n for _ in range(20)]
    print("n=%6d  quota esclusa dal bootstrap: %.3f   (teoria (1-1/n)^n = %.3f)" % (n, np.mean(fuori), (1 - 1 / n) ** n))
print("limite 1/e = %.3f" % np.exp(-1))

n=    10  quota esclusa dal bootstrap: 0.400   (teoria (1-1/n)^n = 0.349)
n=   100  quota esclusa dal bootstrap: 0.367   (teoria (1-1/n)^n = 0.366)
n=  1000  quota esclusa dal bootstrap: 0.371   (teoria (1-1/n)^n = 0.368)
n=100000  quota esclusa dal bootstrap: 0.368   (teoria (1-1/n)^n = 0.368)
limite 1/e = 0.368

Il listato per ogni dimensione n estrae un bootstrap (rng.integers(0, n, size=n): n indici con reinserimento), conta con np.setdiff1d quanti indici tra 0 e n−1 non sono stati estratti, e ripete l'esperimento 20 volte per fare la media. Il confronto con (1−1/n)ⁿ mostra la teoria: già con n = 100 la frazione esclusa è 0,367, vicinissima al limite 1/e = 0,368.

Questi esempi esclusi sono un regalo: per ogni albero formano un piccolo set di validazione già pronto, e si chiamano out-of-bag (OOB). Per valutare la foresta si fa votare ogni esempio solo dagli alberi che non lo hanno visto in addestramento: una stima di generalizzazione senza set di validazione separato.

Il flusso del bagging: il training set genera più campioni bootstrap, un albero per campione, e la previsione finale è la media (o il voto) di tutti gli alberi.

Ingrediente 2: sottoinsieme casuale di feature a ogni divisione

Il bagging da solo non basta. Se c'è una feature molto forte, tutti gli alberi la sceglieranno nelle prime divisioni e finiranno per assomigliarsi (ρ alto). La Random Forest aggiunge una regola: a ogni nodo, la migliore divisione viene cercata solo tra m feature estratte a caso (di solito m = √d per la classificazione). Così alberi diversi sono costretti a usare feature diverse, e la correlazione scende. Ogni albero individualmente diventa un po' peggiore, ma la media migliora: è il compromesso della formula di sopra.

Implementazione a mano

Un albero lo abbiamo già, in scikit-learn, con il parametro max_features che realizza proprio l'estrazione casuale a ogni nodo. Costruiamo la foresta con un ciclo, il voto a maggioranza e la stima OOB:

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

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

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

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

    def predict_proba(self, X):
        return np.mean([t.predict_proba(X) for t in self.alberi], axis=0)   # media delle probabilità

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

    def oob_score(self, X, y):
        somma = np.zeros((len(X), 2)); voti = np.zeros(len(X))
        for t, f in zip(self.alberi, self.fuori):
            somma[f] += t.predict_proba(X[f]); voti[f] += 1
        ok = voti > 0                                              # esempi fuori da almeno un albero
        return (somma[ok].argmax(axis=1) == y[ok]).mean()

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

print("albero singolo       test %.3f" % singolo.score(Xte, yte))
print("foresta a mano       test %.3f   OOB %.3f" % ((mia.predict(Xte) == yte).mean(), mia.oob_score(Xtr, ytr)))
print("foresta sklearn      test %.3f   OOB %.3f" % (lib.score(Xte, yte), lib.oob_score_))

albero singolo       test 0.778
foresta a mano       test 0.875   OOB 0.866
foresta sklearn      test 0.863   OOB 0.864

La classe ForestaDaZero ha tre ruoli. fit: per ognuno dei n_alberi estrae un bootstrap (idx), segna in mask_fuori quali esempi non sono stati estratti (gli out-of-bag), e addestra un DecisionTreeClassifier con max_features="sqrt" sui soli esempi del bootstrap; random_state diverso per ogni albero. predict_proba: media delle probabilità dei singoli alberi. oob_score: per ogni esempio fa votare solo gli alberi per cui era out-of-bag (somma[f] += ... usa la maschera f per sommare le probabilità solo sugli esempi esclusi di quell'albero) e confronta con la vera etichetta.

Il singolo albero fa 0.778 sul test; la foresta a mano 0.875 e quella di scikit-learn 0.863 (le piccole differenze tra le due dipendono dai semi casuali). Nota che la stima OOB (calcolata solo sul training) è vicina all'accuratezza sul test: funziona come validazione gratuita. La foresta di scikit-learn fa la media delle probabilità dei singoli alberi, come la nostra, invece del voto a maggioranza dell'articolo originale di Breiman: lo dice la guida agli ensemble.

Quanti alberi servono?

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

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

Il listato addestra foreste con un numero crescente di alberi e stampa training e test. Con B = 1 la foresta è un singolo albero su un bootstrap (peggiore di un albero addestrato sul training intero, perché vede solo il 63% degli esempi distinti).

Il test passa da 0.743 con un solo albero a 0.857 con 25, 0.863 con 100 e 0.865 con 400: l'accuratezza cresce e poi si stabilizza: oltre un certo numero di alberi il guadagno è trascurabile e si paga solo in tempo e memoria. Aggiungere alberi non fa overfitting, come anticipava la formula; ma l'accuratezza di training vicina a 1.0 non significa nulla da sola, perché ogni albero profondo memorizza il proprio bootstrap.

Misurare la correlazione tra gli alberi

La teoria dice che la feature subsampling serve a decorrelare. Misuriamolo: confrontiamo foreste con max_features=None (il bagging puro, tutte le feature a ogni nodo) e con sqrt, calcolando la correlazione media tra le previsioni di probabilità di coppie di alberi sul test set:

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

for nome, mf in (("bagging puro (tutte le feature)", None), ("random forest (sqrt)", "sqrt")):
    f = RandomForestClassifier(n_estimators=100, max_features=mf, random_state=0).fit(Xtr, ytr)
    print("%-32s correlazione media tra alberi %.3f   test %.3f" % (nome, correlazione_media(f, Xte), f.score(Xte, yte)))

bagging puro (tutte le feature)  correlazione media tra alberi 0.549   test 0.878
random forest (sqrt)             correlazione media tra alberi 0.440   test 0.863

La funzione correlazione_media: estimators_ è la lista degli alberi della foresta. Per ciascuno si ricava la probabilità della classe 1 sul test (matrice alberi × campioni); np.corrcoef restituisce la matrice di correlazione tra alberi, e con np.triu_indices_from(C, k=1) si prende il triangolo superiore senza la diagonale (ogni coppia una sola volta) per farne la media.

La correlazione media tra gli alberi scende da 0.549 a 0.440 quando le feature sono sottocampionate: è il termine ρ della formula che si abbassa, come previsto. Ma guarda l'accuratezza: in questo dataset il bagging puro fa 0.878 contro 0.863 della foresta. Con 600 campioni di test la differenza è di 9 esempi, quindi non basta a dichiarare un vincitore. Il messaggio onesto è questo: decorrelare è un mezzo, non una garanzia. Il vantaggio di max_features dipende dal dataset (in genere cresce quando poche feature dominano), e va scelto con la validazione come ogni iperparametro.

Importanza delle feature: due metodi, un'insidia

La foresta offre feature_importances_, calcolata come riduzione media di impurità portata da ogni feature. È veloce ma ha un difetto noto: favorisce le feature con molti valori distinti, anche se inutili, perché offrono più soglie da provare. L'alternativa è la permutation importance: mescola a caso i valori di una colonna nei dati di validazione e misura quanto cala l'accuratezza; se la colonna non serve, non cala. Aggiungiamo una colonna di rumore puro con molti valori distinti:

from sklearn.inspection import permutation_importance

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

imp_gini = f.feature_importances_
perm = permutation_importance(f, Xr_te, yte, n_repeats=10, random_state=0).importances_mean
print("colonna di rumore (ultima)  importanza da impurità: %.4f   (media delle altre 20: %.4f)" % (imp_gini[-1], imp_gini[:-1].mean()))
print("colonna di rumore (ultima)  permutation importance:  %.4f   (media delle 5 più alte: %.4f)" % (perm[-1], np.sort(perm[:-1])[-5:].mean()))

colonna di rumore (ultima)  importanza da impurità: 0.0204   (media delle altre 20: 0.0490)
colonna di rumore (ultima)  permutation importance:  -0.0017   (media delle 5 più alte: 0.0669)

Il listato aggiunge una colonna di rumore gaussiano continuo (valori tutti distinti, perciò con molte soglie possibili) e allena la foresta. feature_importances_ è l'importanza per impurità; permutation_importance mescola una colonna alla volta per 10 volte sul test e misura il calo medio di accuratezza.

L'importanza basata sull'impurità attribuisce alla colonna di rumore 0.0204, cioè circa il 2% dell'importanza totale (che somma a 1): un valore non nullo per una colonna che per costruzione non dice niente. La permutation importance sul test la colloca a −0.0017, cioè zero entro il rumore, mentre le cinque feature più importanti in media superano 0.066. Regola pratica: usa quella da permutazione per capire cosa conta davvero; ricorda che, con feature fortemente correlate, anche quella si spartisce l'importanza tra le colonne «gemelle». La guida di scikit-learn conferma entrambi i punti: l'importanza da impurità per gli alberi è fortemente sbilanciata verso le feature con molti valori distinti, e si calcola sul training, mentre la permutazione può essere calcolata su dati non visti.

Trade-off

• Pro: ottimo punto di partenza su dati tabellari, pochi iperparametri critici, robusto al rumore, non richiede scaling, stima OOB gratuita, parallelizzabile (n_jobs=-1), stima l'importanza delle feature.
• Contro: perde l'interpretabilità del singolo albero; modello pesante in memoria e più lento in previsione; non estrapola fuori dall'intervallo del training; spesso battuto dal boosting (prossima lezione) quando serve l'ultimo punto di accuratezza.
• Iperparametri che contano: max_features (decorrelazione), min_samples_leaf (alberi meno profondi, meno rumore), n_estimators (di più è meglio finché non si stabilizza).

Il filo comune ai pro e ai contro è che la foresta è una media di molti alberi profondi: la media la rende robusta e quasi priva di messa a punto, ma toglie la leggibilità del singolo albero e porta con sé il peso di centinaia di alberi. Dei tre iperparametri, max_features e min_samples_leaf sono quelli da validare; n_estimators si alza finché la curva non si appiattisce.

Trappole reali

• Fidarsi dell'accuratezza di training: sarà quasi 1.0 per costruzione. Guarda OOB o validazione.
• OOB con dati non indipendenti: se hai righe duplicate o più righe dello stesso utente, l'esempio «fuori» da un albero ha un gemello dentro, e la stima OOB è ottimista. Serve una divisione per gruppi.
• Serie temporali: bootstrap e validazione casuale mescolano passato e futuro (leakage). Usa divisioni temporali.
• Classi sbilanciate: usa class_weight="balanced_subsample" e valuta con recall/precision; l'accuratezza mente (vedi sotto).
• Importanza dalla sola impurità: come visto sopra, favorisce colonne inutili ad alta cardinalità.
• Probabilità non calibrate: le frazioni di voto degli alberi non sono probabilità vere; se ti servono, calibrale.

Quattro delle sei trappole riguardano la valutazione: training, OOB con dati non indipendenti, serie temporali, accuratezza su classi sbilanciate. Le altre due riguardano come leggere importanza e probabilità. Nessuna produce un errore di esecuzione: la foresta gira e restituisce numeri plausibili, ed è il metodo di valutazione a dover essere corretto. Per i duplicati e le serie temporali vedi il processo di un progetto di ML.

Gli iperparametri che contano: max_features e min_samples_leaf

La foresta ha pochi iperparametri critici, e due di essi controllano i due ingredienti della lezione. max_features regola la decorrelazione degli alberi: con 1 ogni divisione sceglie tra una sola feature (alberi quasi casuali), con "sqrt" il default (circa 4-5 su 20), con 0,5 metà, con None tutte (si ricade nel bagging puro). min_samples_leaf regola la profondità: con 1 gli alberi crescono fino alla purezza, con valori più grandi vengono potati a monte e ciascuna foglia media più esempi. Valutiamo una griglia di dodici combinazioni con tre stime di accuratezza: l'OOB (gratis, calcolato durante l'addestramento), la validazione incrociata a 5 gruppi sul training e il test. Se l'OOB è affidabile, le prime due devono essere vicine.

from sklearn.model_selection import cross_val_score, StratifiedKFold

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

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

Il listato, per ciascuna delle dodici combinazioni, addestra una foresta da 200 alberi con la stima OOB attivata, ne calcola la validazione incrociata (con una foresta nuova per gruppo) e misura il test.

Leggi tre cose. Primo, l'OOB è una buona stima di generalizzazione: la sua differenza dalla validazione incrociata è al più di 1,3 punti (per esempio 0,864 contro 0,873, 0,881 contro 0,878) e costa una frazione del tempo, perché non richiede di riaddestrare.

Secondo, max_features conta: con 1 l'accuratezza è 0,82-0,84, con la radice quadrata 0,83-0,86 e con 0,5 o con tutte le feature 0,84-0,88. Qui il valore predefinito non è il migliore: il bagging puro (None) e 0,5 vincono, perché in questo problema ci sono 5 feature informative e 3 ridondanti su 20 e la scelta casuale di sole 4 o 5 feature a ogni divisione spesso lascia fuori le informative. Il vantaggio di "sqrt" si vede quando poche feature dominano e molte sono inutili o quasi duplicate: dipende dal dataset, ed è un iperparametro da cercare.

Terzo, min_samples_leaf: aumentarlo da 1 a 20 peggiora in modo sistematico (OOB da 0,864 a 0,848 con sqrt, da 0,881 a 0,864 con tutte le feature): la foresta media già molti alberi e non ha bisogno di alberi potati; una potatura eccessiva aumenta il bias. Queste tendenze valgono per questo problema: la lezione è il metodo, non i valori.

Extra Trees: ancora più casualità

Le Extremely Randomized Trees (P. Geurts, D. Ernst e L. Wehenkel, 2006) spingono l'idea oltre: a ogni nodo, oltre a scegliere un sottoinsieme casuale di feature, scelgono anche la soglia a caso (una per feature, tra il minimo e il massimo osservati) e tengono la migliore delle soglie casuali, invece di cercare quella ottima. Ogni albero individualmente è peggiore, ma la correlazione tra alberi crolla (il termine ρ della formula della varianza) e il calcolo è molto più rapido, perché non ordina i valori. In scikit-learn non fanno bootstrap per default (usano tutti gli esempi), come indica la guida agli ensemble.

import time
from sklearn.ensemble import ExtraTreesClassifier

for name, cls in (("RandomForest", RandomForestClassifier), ("ExtraTrees", ExtraTreesClassifier)):
    t0 = time.time()
    model = cls(300, random_state=0, n_jobs=1).fit(Xtr, ytr)
    print("%-13s accuratezza sul test %.3f   tempo di addestramento %.2f s (un solo processo)" % (name, model.score(Xte, yte), time.time() - t0))

# Uscita di una esecuzione; i tempi cambiano a ogni esecuzione
# RandomForest  accuratezza sul test 0.865   tempo di addestramento 2.52 s (un solo processo)
# ExtraTrees    accuratezza sul test 0.873   tempo di addestramento 1.02 s (un solo processo)

Il listato addestra le due famiglie con 300 alberi, a parità di tutto il resto, e misura accuratezza e tempo con un solo processo (n_jobs=1, per un confronto equo). Le Extra Trees fanno 0,873 contro 0,865 della Random Forest, una differenza che sta nel rumore (con 600 esempi di test sono cinque esempi). Sul tempo, in tre esecuzioni su un PC che in quel momento era usato anche per altro, la Random Forest ha impiegato da 1,55 a 2,52 s e le Extra Trees da 0,65 a 1,28 s, cioè tra 0,4 e 0,6 volte il tempo. Il risultato è frequente perché non si cerca la soglia ottima ma si sceglie tra soglie casuali; se hai molti dati o poco tempo, le Extra Trees sono un'alternativa da provare e da confrontare con la validazione.

Le probabilità di una foresta: affidabili?

Il metodo predict_proba di una foresta restituisce la media delle frazioni di classi nelle foglie degli alberi. Nel caso di Naive Bayes abbiamo visto probabilità troppo sicure di sé. Che cosa succede qui? Costruiamo la stessa tabella di affidabilità: gruppi di previsioni per probabilità, confidenza media dichiarata e frazione reale di positivi.

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

print("probabilità prevista     previsioni   media dichiarata   frazione reale di positivi")
edges = [0, 0.1, 0.3, 0.5, 0.7, 0.9, 1.0001]
for lo, hi in zip(edges[:-1], edges[1:]):
    in_bin = (p_positive >= lo) & (p_positive < hi)
    print("   [%.1f, %.1f)         %8d        %.3f              %.3f" % (lo, min(hi, 1.0), in_bin.sum(), p_positive[in_bin].mean(), yte[in_bin].mean()))

probabilità prevista     previsioni   media dichiarata   frazione reale di positivi
   [0.0, 0.1)               49        0.069              0.041
   [0.1, 0.3)              195        0.189              0.113
   [0.3, 0.5)               91        0.389              0.374
   [0.5, 0.7)               54        0.584              0.704
   [0.7, 0.9)               96        0.828              0.938
   [0.9, 1.0)              115        0.946              0.974

Il listato raggruppa le 600 previsioni di probabilità della foresta in sei intervalli e confronta, in ciascuno, la probabilità media dichiarata con la frazione reale di positivi. Il profilo è l'opposto di quello di Naive Bayes: la foresta è prudente. Va detto che la prudenza non è uniforme: nell'intervallo più basso (0-0,1) dichiara 0,069 e la realtà è 0,041, cioè la foresta sovrastima i positivi; nell'intervallo 0,7-0,9 dichiara 0,828 ma i positivi sono il 93,8%; in quello 0,9-1 dichiara 0,946 e la realtà è 0,974. Quindi nella parte alta le probabilità sono troppo basse, in quella bassa troppo alte.

La media su molti alberi, che non sono mai tutti d'accordo anche quando l'esempio è «facile», comprime le probabilità verso il centro: non può dare 0 o 1 se anche un solo albero dissente. È un difetto meno grave di quello di Naive Bayes (l'ordinamento è giusto e la distorsione è moderata), e si corregge con la calibrazione (CalibratedClassifierCV). Le probabilità di una foresta sono quindi migliori di quelle di Naive Bayes, ma non vanno prese come probabilità vere senza un controllo.

Classi sbilanciate: pesi e soglia

Con il 5% di positivi, una foresta predetta «a maggioranza» (soglia 0,5) tende a ignorare la classe rara. Due leve: i pesi di classe (class_weight), che fanno contare di più gli errori sui positivi nel calcolo dell'impurità, e la soglia sulle probabilità. Le confrontiamo sullo stesso problema sbilanciato di KNN e di Naive Bayes (qui 5 000 esempi).

from sklearn.metrics import recall_score, precision_score

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

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

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

Il listato addestra due foreste, senza pesi e con balanced_subsample (i pesi di classe sono ricalcolati per ciascun campione bootstrap), e per ciascuna misura recall e precision a quattro soglie. Con soglia 0,5 e senza pesi la foresta trova solo il 15% dei positivi (precision 0,87); con i pesi il 23% (precision 0,80): un miglioramento modesto. La leva più potente è la soglia: abbassarla a 0,2 porta la recall al 41% senza pesi (precision 0,43) e al 60% con i pesi (precision 0,30), e a 0,1 al 62% e al 76% (con precision di 0,23 e 0,18). I due strumenti sono in parte sostituibili, ma non producono informazione nuova: spostano il punto di lavoro lungo la stessa curva precision-recall. Il guadagno dei pesi, a parità di soglia, è reale (a 0,2: 60% contro 41% di recall con precision 0,30 contro 0,43) ma è soprattutto uno spostamento lungo la curva. Come sempre, la soglia giusta dipende dai costi, come spiega il processo di un progetto di ML.

Il costo in produzione: dimensioni e tempo di previsione

Una foresta è un modello pesante. Ogni albero memorizza i suoi nodi (soglia, feature, figli, valori), e cento o trecento alberi profondi diventano facilmente decine di megabyte. Misuriamo l'effetto della profondità sulle dimensioni del modello e sul tempo di previsione, confrontandolo con una regressione logistica che prevede la stessa classe.

import pickle
from sklearn.linear_model import LogisticRegression

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

print("modello                            nodi totali   dimensione (KB)   ms per 1000 previsioni")
for leaf in (1, 5, 20):
    measure(RandomForestClassifier(300, min_samples_leaf=leaf, random_state=0, n_jobs=1).fit(Xtr, ytr), "RandomForest 300 alberi, foglia min. %d" % leaf)
measure(LogisticRegression(max_iter=1000).fit(Xtr, ytr), "regressione logistica")

# Uscita di una esecuzione; nodi e dimensioni sono deterministici, i tempi no
# modello                            nodi totali   dimensione (KB)   ms per 1000 previsioni
# RandomForest 300 alberi, foglia min. 1     85832         6790                71.79
# RandomForest 300 alberi, foglia min. 5     48946         3907                93.78
# RandomForest 300 alberi, foglia min. 20     18502         1528                52.86
# regressione logistica                      0            1                 0.46

La funzione measure serializza il modello con pickle per stimarne la dimensione, cronometra 20 previsioni sull'intero test e riporta i millisecondi per mille previsioni; per le foreste conta anche i nodi totali. L'effetto di min_samples_leaf sulla dimensione è fortissimo ed è deterministico (stessi numeri in ogni esecuzione): con foglie di un solo esempio la foresta ha 85 832 nodi e occupa 6,8 MB; alzando il minimo a 20 scende a 18 502 nodi e 1,5 MB, circa un quarto, a scapito di un po' di accuratezza (misurata sopra: da 0,864 a 0,848 con l'OOB).

Il tempo di previsione è meno stabile. In tre esecuzioni, su un PC usato anche per altro, le foreste hanno impiegato da 51 a 147 ms per mille previsioni, senza un ordine coerente fra le tre profondità, perché il tempo dipende soprattutto dal numero di alberi da attraversare (sempre 300). La regressione logistica occupa circa 1 KB (ventuno numeri: venti coefficienti e un'intercetta) e prevede in 0,15-0,46 ms per mille esempi: da ottanta a oltre cento volte più veloce, secondo il carico della macchina al momento della misura. È il vero costo della foresta in produzione: se latenza, memoria o tempo di caricamento sono vincoli, min_samples_leaf e il numero di alberi sono leve di dimensionamento, e si sceglie il modello più piccolo la cui accuratezza non è peggiore di un margine accettabile.

Tre prove da fare

Esercizio 1. Ripeti la griglia con 20 e con 1 000 alberi: l'OOB si avvicina alla validazione incrociata al crescere degli alberi? Esercizio 2. Aggiungi 30 colonne di rumore puro e rifai la griglia di max_features: ora il default "sqrt" è il migliore? Esercizio 3. Applica CalibratedClassifierCV alla foresta e rifai la tabella di affidabilità: quanto si avvicinano le due colonne?

Quando la foresta non è lo strumento giusto

Una foresta casuale è un ottimo punto di partenza per dati tabulari, ma ha limiti strutturali da conoscere. Primo, non estrapola: ogni albero predice la media dei valori delle foglie, quindi fuori dall'intervallo dei dati di addestramento la previsione resta piatta sull'ultimo valore visto, a differenza di un modello lineare che prosegue lungo la retta; per grandezze con trend (prezzi nel tempo, domanda crescente) è un errore sistematico, non rumore.

Secondo, l'importanza delle feature calcolata con la diminuzione dell'impurità favorisce le variabili con molti valori distinti e si divide fra variabili correlate: leggila come un'indicazione, non come una graduatoria causale, e controllala con la permutazione sul set di validazione. Terzo, il modello è grande e lento da servire rispetto a uno lineare: centinaia di alberi profondi vanno tenuti in memoria e attraversati a ogni previsione, cosa che pesa quando la latenza è un vincolo. Quarto, le probabilità di una foresta sono medie di voti, quindi tendono a evitare gli estremi e vanno calibrate se vengono usate come probabilità vere. Nessuno di questi limiti la rende sconsigliabile: servono a decidere quando sostituirla o affiancarla con un modello lineare o con un boosting.

Per allenarti sul dominio della manutenzione predittiva delle pompe di irrigazione c'è uno scheletro di esercizi in Python con tre parti da completare: la quota fuori borsa e il limite 1/e, la varianza della media di stimatori correlati con il numero di alberi necessario, e il teorema della giuria di Condorcet. Finché i punti da completare sono vuoti, stampa quali esercizi mancano.

Che cosa sai fare adesso

Sai derivare perché la media di molti alberi sbaglia meno e dove si ferma (il pavimento ρσ²), sai costruire una foresta con bootstrap, sottoinsieme di feature, voto e stima out-of-bag, e sai misurare la correlazione fra gli alberi invece di darla per scontata. Sai scegliere max_features e min_samples_leaf con una griglia, riconoscere i limiti delle probabilità e dell'importanza da impurità, e dimensionare il modello per la produzione. Nella lezione successiva, boosting, gli alberi non sono più indipendenti: ciascuno corregge gli errori del precedente.

Fonti