Ensembler: bagging, random forest og boosting

Et beslutningstre er ustabilt: litt andre treningsdata gir et annet tre. Det er en svakhet når du har ett tre, men en mulighet når du kan lage mange. Hadde vi hatt mange uavhengige treningssett fra samme fordeling, kunne vi trent ett tre på hvert og tatt snittet. Det har vi ikke. Vi har ett datasett.

Løsningen er et bootstrap-utvalg: trekk $n$ rader fra de $n$ radene i treningssettet, med tilbakelegging. Hver trekning velger hvilken som helst rad med sannsynlighet $1/n$, uansett hva som ble trukket før. Noen rader kommer derfor med flere ganger, og andre kommer ikke med i det hele tatt. Hvert tre ser en litt annen versjon av dataene.

Hvor mange rader blir utelatt? En bestemt rad blir ikke valgt i én trekning med sannsynlighet $1-1/n$, og ikke i noen av de $n$ uavhengige trekningene med sannsynlighet

Hvert bootstrap-utvalg inneholder altså bare rundt $63,2\%$ av de unike radene. De rundt $36,8\%$ som ikke ble trukket, er out-of-bag for det treet. Det er en gave: for hver treningsrad finnes det trær som aldri har sett den. Lar vi bare de trærne stemme over raden, får vi en ærlig prediksjon. Gjør vi det for alle rader, har vi et OOB-estimat av testfeilen uten å sette av et eget valideringssett. I scikit-learn slår du det på med oob_score=True.

Bagging (bootstrap aggregating) trener én modell per bootstrap-utvalg og slår prediksjonene sammen. Ved klassifikasjon bruker ensemblet ƒflertallsstemme, ved regresjon gjennomsnittet. Hvert enkelt tre er litt dårligere enn et tre som er trent på alle dataene, siden det ser færre unike rader. Likevel blir ensemblet bedre. Hvorfor?

Tenk på feilen til hvert tre som et tilfeldig utslag rundt det riktige svaret. Hvis utslagene var helt uavhengige, ville snittet av $M$ trær hatt variansen ${\sigma }^2/M$: hundre trær ville gitt en hundredel av variansen. Men trærne er trent på nesten de samme dataene, så feilene deres henger sammen. Med korrelasjon $\rho$ mellom feilene til to trær gir ƒvarians for et snitt av modeller

Det andre leddet krymper når $M$ vokser. Det første gjør det ikke. Uansett hvor mange trær du legger til, blir variansen aldri lavere enn $\rho {\sigma }^2$, så korrelasjonen setter et gulv. Variasjon mellom modellene er derfor hele poenget: bootstrap gjør trærne forskjellige nok til at snittet hjelper.

To konsekvenser er verdt å merke seg. Bagging reduserer varians, ikke bias: gjør alle trærne den samme systematiske feilen, overlever den avstemningen. Derfor bagger vi dype, ustabile trær med lav bias, ikke grunne trær. Og flere trær gir ikke overtilpasning slik dypere trær gjør. De stabiliserer bare snittet, mot en kostnad i tid og minne.

import numpy as np
import matplotlib.pyplot as plt

rng = np.random.default_rng(10)
rho = 0.5  # @param 0..1 0.05
sigma2 = 4.0                                   # varians i feilen til én modell
M = np.arange(1, 101)                          # antall modeller i ensemblet
felles = rng.normal(size=(5000, 1))            # feil alle modellene deler
egen = rng.normal(size=(5000, 100))            # feil hver modell har alene
feil = np.sqrt(sigma2) * (np.sqrt(rho) * felles + np.sqrt(1 - rho) * egen)
simulert = (np.cumsum(feil, axis=1) / M).var(axis=0)   # varians til snittet av de M første
formel = rho * sigma2 + (1 - rho) * sigma2 / M
print(f"rho = {rho:.2f}: varians {simulert[0]:.2f} for én modell, {simulert[24]:.2f} for 25, {simulert[99]:.2f} for 100")
print(f"gulvet rho·σ² = {rho * sigma2:.2f} forsvinner ikke uansett hvor mange modeller du legger til")

plt.figure()
plt.plot(M, simulert, label="simulert")
plt.plot(M, formel, label="formel ρσ² + (1 − ρ)σ²/M")
plt.ylim(0, sigma2 * 1.05)
plt.xlabel("antall modeller M")
plt.ylabel("varians til snittet")
plt.legend()
plt.tight_layout()
plt.show()

Bootstrap alene gjør ikke trærne særlig ulike hvis én feature er mye sterkere enn resten. Da velger nesten hvert tre den samme featuren i roten, og de øverste splittene, som betyr mest, blir like. Trærne gjør samme slags feil, $\rho$ blir høy, og gulvet $\rho {\sigma }^2$ blir høyt.

En random forest legger til én tilfeldighet til: ved hver split får treet bare velge blant et tilfeldig utvalg på $k$ av de $p$ featurene. Et nytt utvalg trekkes i hver node. Den sterke featuren er da ikke tilgjengelig i alle splitter, og andre, svakere mønstre får komme fram. Trærne blir mindre like, korrelasjonen faller, og snittet blir bedre, selv om hvert enkelt tre blir litt svakere.

Et vanlig valg er $k\approx \sqrt{p}$ for klassifikasjon, i scikit-learn max_features="sqrt". Med $k=p$ er random forest det samme som bagging av trær. Med svært liten $k$ blir trærne veldig ulike, men hvert av dem så svakt at snittet lider. $k$ er derfor en hyperparameter som balanserer styrken til hvert tre mot korrelasjonen mellom dem.

Hold de to tilfeldighetene fra hverandre: bootstrap velger rader for hvert tre, og featureutvalget velger kolonner for hver split. Trærne trenes uavhengig av hverandre og kan bygges parallelt.

Bagging lager modeller side om side. AdaBoost (adaptive boosting) lager dem etter hverandre. Hver modell er en svak lærer, typisk en beslutningsstubbe med én split, som bare må være litt bedre enn tilfeldig gjetting. Kunsten er hvordan de neste lærerne får vite hva de forrige bommet på.

AdaBoost gir hvert treningspunkt en vekt, og algoritmen går slik:

1. 1.Start med like vekter, $w_i=1/n$.
2. 2.Tren en svak lærer $h_m$ som minimerer den vektede feilen ${\epsilon }_m$, summen av vektene til punktene den bommer på.
3. 3.Gi læreren stemmevekten ƒAdaBoost-stemmevekt, ${\alpha }_m=\frac{1}{2}\ln \frac{1-{\epsilon }_m}{{\epsilon }_m}$. Lav feil gir stor ${\alpha }_m$, og ${\epsilon }_m=0,5$ gir ${\alpha }_m=0$, altså ingen stemme.
4. 4.Oppdater vektene med ƒAdaBoost-vektoppdatering: feilklassifiserte punkter ganges med $e^{{\alpha }_m}$ og riktige med $e^{-{\alpha }_m}$, og så normaliseres alt slik at summen blir 1.
5. 5.Gjenta fra punkt 2 med de nye vektene.

Den endelige prediksjonen er en vektet avstemning. Med klassene kodet som $-1$ og $+1$ er $\hat{y}=\mathrm{sign}({\sum }_m{\alpha }_m{h}_m(x))$, så gode lærere teller mer enn svake.

Vektoppdateringen har en elegant egenskap: etter normaliseringen har de feilklassifiserte punktene til sammen nøyaktig vekt $0,5$. Den forrige læreren er dermed ikke bedre enn gjetting på de nye vektene, og neste lærer må finne noe nytt. Baksiden er at punkter som stadig feilklassifiseres, får enorm vekt. Feilmerkede punkter og uteliggere kan derfor dra AdaBoost langt av sted.

Gradient boosting bygger også modellen trinnvis, men i stedet for å vekte om punktene lar den hvert nytt tre lære feilen som gjenstår. Tenk deg at du skal treffe et mål med en rekke små puff: etter hvert puff måler du hvor langt unna du er, og neste puff sikter på resten.

Med kvadratisk tap er oppskriften konkret:

1. 1.Start med en konstant, $F_0(x)=\bar{y}$.
2. 2.Regn ut residualene $r_i=y_i-F_{m-1}(x_i)$, det som gjenstår for hvert punkt.
3. 3.Tilpass et lite regresjonstre $h_m$ til residualene, ikke til $y$.
4. 4.Oppdater med ƒgradient boosting, $F_m(x)=F_{m-1}(x)+\eta h_m(x)$.

Hvorfor heter det gradient? For tapet $L=\frac{1}{2}(y-F{)}^2$ er $-\partial L/\partial F=y-F$, akkurat residualen. Hvert tre tar altså et gradient descent-steg, men i funksjonsrommet: i stedet for å flytte en vekt legger vi til en funksjon som peker i retningen der tapet synker raskest. For andre tap, som absoluttfeil eller log-tap ved klassifikasjon, er minus gradienten noe annet enn den vanlige residualen, og da kalles målet for hvert tre en pseudo-residual.

Læringsraten $\eta$, også kalt shrinkage, gjør at hvert tre bare får rette opp en brøkdel av feilen. Med $\eta =0,1$ trengs mange flere trær, men resultatet generaliserer som regel bedre, fordi ingen enkelt tre får diktere modellen. Antall runder $M$ og $\eta$ virker sammen: liten $\eta$ krever stor $M$. For mange runder overtilpasser, siden residualene til slutt bare er støy, så antall runder velges med validering eller tidlig stopping.

Begge familiene kombinerer mange trær, men de angriper motsatte problemer:

• •Hva som reduseres. Bagging tar dype trær med lav bias og høy varians og reduserer variansen ved å ta snittet. Boosting tar svake, grunne trær med høy bias og reduserer biasen ved å legge korreksjoner oppå hverandre.
• •Rekkefølge. Baggingtrærne er uavhengige og kan trenes parallelt. Boosting er sekvensielt: modell $m$ trenger feilene til modellene før den.
• •Flere modeller. I bagging er flere trær nesten alltid ufarlig, det stabiliserer bare snittet. I boosting kan for mange runder overtilpasse, fordi hver runde tilpasser seg treningsfeilen tettere.
• •Følsomhet. Random forest fungerer ofte godt med standardinnstillinger. Gradient boosting gir ofte den beste ytelsen på tabelldata, men krever at læringsrate, antall runder og tredybde tunes sammen. AdaBoost er følsom for feilmerkede punkter, siden de får stadig større vekt.
• •Hyperparametre. Random forest styres av antall trær, max_features og dybden på trærne. Boosting styres av antall runder, læringsraten og dybden på hvert lille tre.

Felles for alle er at ingen ensemble redder dårlige features, feilmerkede data eller lekkasje fra testsettet.

Legg merke til at skogen har 100 % treningsaccuracy uten å være dårlig på test.

import numpy as np
from sklearn.datasets import make_classification
from sklearn.ensemble import GradientBoostingClassifier, RandomForestClassifier
from sklearn.model_selection import train_test_split

X, y = make_classification(n_samples=600, n_features=8, n_informative=4,
                           flip_y=0.08, random_state=3)
X_tr, X_te, y_tr, y_te = train_test_split(X, y, test_size=0.3, random_state=0)
bagging = 0  # @param 0..1
laeringsrate = 0.1  # @param 0.05..1 0.05

rf = RandomForestClassifier(n_estimators=100, max_features=None if bagging else "sqrt",
                            oob_score=True, random_state=0).fit(X_tr, y_tr)
print(f"{'bagging' if bagging else 'random forest'}: trening {rf.score(X_tr, y_tr):.2f}, "
      f"OOB {rf.oob_score_:.2f}, test {rf.score(X_te, y_te):.2f}")

gb = GradientBoostingClassifier(n_estimators=200, learning_rate=laeringsrate,
                                max_depth=2, random_state=0).fit(X_tr, y_tr)
tren = [np.mean(p == y_tr) for p in gb.staged_predict(X_tr)]   # accuracy etter hver runde
test = [np.mean(p == y_te) for p in gb.staged_predict(X_te)]
for m in [1, 10, 50, 200]:
    print(f"boosting, {m:>3} runder: trening {tren[m - 1]:.2f}, test {test[m - 1]:.2f}")