R: Co widzę na wykresach częściowej zależności GBM i RandomForest?

Właściwie myślałem, że zrozumiałem, co można pokazać z częściową fabułą zależności, ale używając bardzo prostego hipotetycznego przykładu, byłem dość zdziwiony. W poniższym fragmencie kodu wygenerować trzy zmienne niezależne ( , b , c ) i jedną zmienną zależną ( y ) z c pokazującym zbliżenie liniową zależność y , a i b są nieskorelowane, przy czym Y . Wykonuję analizę regresji za pomocą wzmocnionego drzewa regresji przy użyciu pakietu R :gbm

a <- runif(100, 1, 100)
b <- runif(100, 1, 100)
c <- 1:100 + rnorm(100, mean = 0, sd = 5)
y <- 1:100 + rnorm(100, mean = 0, sd = 5)
par(mfrow = c(2,2))
plot(y ~ a); plot(y ~ b); plot(y ~ c)
Data <- data.frame(matrix(c(y, a, b, c), ncol = 4))
names(Data) <- c("y", "a", "b", "c")
library(gbm)
gbm.gaus <- gbm(y ~ a + b + c, data = Data, distribution = "gaussian")
par(mfrow = c(2,2))
plot(gbm.gaus, i.var = 1)
plot(gbm.gaus, i.var = 2)
plot(gbm.gaus, i.var = 3)

Nic dziwnego, że dla zmiennych a i b wykresy częściowej zależności dają linie poziome wokół średniej a . Intryguje mnie fabuła dla zmiennej c . Otrzymuję linie poziome dla zakresów c <40 ic > 60, a oś y jest ograniczona do wartości bliskich średniej y . Od i b są całkowicie niezwiązane z y (a więc tam zmienne znaczenie w modelu 0), to oczekuje się, że cwykazywałby częściową zależność w całym zakresie zamiast tego sigmoidalnego kształtu dla bardzo ograniczonego zakresu jego wartości. Próbowałem znaleźć informacje w Friedmanie (2001) „Przybliżenie funkcji chciwości: maszyna zwiększająca gradient” oraz w Hastie i in. (2011) „Elementy uczenia statystycznego”, ale moje umiejętności matematyczne są zbyt niskie, aby zrozumieć wszystkie zawarte w nich równania i wzory. Zatem moje pytanie: co determinuje kształt wykresu zależności częściowej dla zmiennej c ? (Proszę wyjaśnić słowami zrozumiałymi dla nie-matematyka!)

DODANO 17 kwietnia 2014 r .:

Podczas oczekiwania na odpowiedź, to stosuje się takie same przykładowe dane dla analiz R-pakietu randomForest. Wykresy częściowej zależności randomForest bardziej przypominają to, czego oczekiwałem po wykresach gbm: częściowa zależność zmiennych objaśniających a i b zmienia się losowo i blisko około 50, podczas gdy zmienna objaśniająca c wykazuje częściową zależność w całym zakresie (i prawie w całym cały zakres y ). Jakie mogą być przyczyny tych różnych kształtów wykresów częściowej zależności w gbmi randomForest?

Oto zmodyfikowany kod, który porównuje wykresy:

a <- runif(100, 1, 100)
b <- runif(100, 1, 100)
c <- 1:100 + rnorm(100, mean = 0, sd = 5)
y <- 1:100 + rnorm(100, mean = 0, sd = 5)
par(mfrow = c(2,2))
plot(y ~ a); plot(y ~ b); plot(y ~ c)
Data <- data.frame(matrix(c(y, a, b, c), ncol = 4))
names(Data) <- c("y", "a", "b", "c")

library(gbm)
gbm.gaus <- gbm(y ~ a + b + c, data = Data, distribution = "gaussian")

library(randomForest)
rf.model <- randomForest(y ~ a + b + c, data = Data)

x11(height = 8, width = 5)
par(mfrow = c(3,2))
par(oma = c(1,1,4,1))
plot(gbm.gaus, i.var = 1)
partialPlot(rf.model, Data[,2:4], x.var = "a")
plot(gbm.gaus, i.var = 2)
partialPlot(rf.model, Data[,2:4], x.var = "b")
plot(gbm.gaus, i.var = 3)
partialPlot(rf.model, Data[,2:4], x.var = "c")
title(main = "Boosted regression tree", outer = TRUE, adj = 0.15)
title(main = "Random forest", outer = TRUE, adj = 0.85)

Poświęciłem trochę czasu na napisanie własnego „częściowego plotera funkcyjnego”, zanim zdałem sobie sprawę, że jest już dołączony do biblioteki R randomForest.

[EDYTUJ ... ale potem spędziłem rok na tworzeniu pakietu CRAN forestFloor , który moim zdaniem jest znacznie lepszy niż klasyczne wykresy zależności częściowej]

Wykres częściowej funkcji jest świetny w instancjach, ponieważ pokazany tutaj przykład symulacji, w którym zmienna objaśniająca nie wchodzi w interakcje z innymi zmiennymi. Jeśli każda zmienna objaśniająca przyczynia się dodatkowo do docelowego Y przez jakąś nieznaną funkcję, metoda ta świetnie pokazuje szacowaną funkcję ukrytą. Często widzę takie spłaszczenie w granicach funkcji cząstkowych.

Kilka powodów: randomForsest ma argument o nazwie „nodesize = 5”, co oznacza, że ​​żadne drzewo nie podzieli grupy 5 członków lub mniej. Dlatego każde drzewo nie może odróżnić się z większą precyzją. Warstwa zestawu do pakowania / ładowania początkowego wygładza się poprzez głosowanie wielu funkcji krokowych poszczególnych drzew - ale tylko w środku obszaru danych. Zbliżając się do granic danych reprezentujących przestrzeń, „amplituda” funkcji częściowej spadnie. Ustawienie wielkości węzła = 3 i / lub uzyskanie większej liczby obserwacji w porównaniu z hałasem może zmniejszyć ten efekt spłaszczenia granicy ... Gdy stosunek sygnału do szumu ogólnie spada w losowym lesie, skala prognoz ulega kondensacji. Zatem prognozy nie są absolutnie dokładne, ale jedynie liniowo skorelowane z celem. Możesz zobaczyć wartości aib jako przykłady bardzo niskiego stosunku sygnału do szumu, i dlatego te funkcje częściowe są bardzo płaskie. Jest to miła cecha losowego lasu, że już z zakresu prognoz zestawu treningowego możesz zgadnąć, jak dobrze model działa. OOB.predictions jest również świetny ..

spłaszczanie częściowej działki w regionach bez danych jest uzasadnione: Ponieważ losowy las i CART są modelowaniem opartym na danych, osobiście podoba mi się koncepcja, że ​​modele te nie ekstrapolują. Zatem przewidywanie c = 500 lub c = 1100 jest dokładnie takie samo jak c = 100 lub w większości przypadków również c = 98.

Oto przykład kodu ze zmniejszonym spłaszczaniem krawędzi:

Nie próbowałem pakietu GBM ...

oto przykładowy kod oparty na twoim przykładzie ...

#more observations are created...
a <- runif(5000, 1, 100)
b <- runif(5000, 1, 100)
c <- (1:5000)/50 + rnorm(100, mean = 0, sd = 0.1)
y <- (1:5000)/50 + rnorm(100, mean = 0, sd = 0.1)
par(mfrow = c(1,3))
plot(y ~ a); plot(y ~ b); plot(y ~ c)
Data <- data.frame(matrix(c(y, a, b, c), ncol = 4))
names(Data) <- c("y", "a", "b", "c")
library(randomForest)
#smaller nodesize "not as important" when there number of observartion is increased
#more tress can smooth flattening so boundery regions have best possible signal to             noise, data specific how many needed

plot.partial = function() {
partialPlot(rf.model, Data[,2:4], x.var = "a",xlim=c(1,100),ylim=c(1,100))
partialPlot(rf.model, Data[,2:4], x.var = "b",xlim=c(1,100),ylim=c(1,100))
partialPlot(rf.model, Data[,2:4], x.var = "c",xlim=c(1,100),ylim=c(1,100))
}

#worst case! : with 100 samples from Data and nodesize=30
rf.model <- randomForest(y ~ a + b + c, data = Data[sample(5000,100),],nodesize=30)
plot.partial()

#reasonble settings for least partial flattening by few observations: 100 samples and nodesize=3 and ntrees=2000
#more tress can smooth flattening so boundery regions have best possiblefidelity
rf.model <- randomForest(y ~ a + b + c, data = Data[sample(5000,100),],nodesize=5,ntress=2000)
plot.partial()

#more observations is great!
rf.model <- randomForest(y ~ a + b + c,
 data = Data[sample(5000,5000),],
 nodesize=5,ntress=2000)
plot.partial()

Jak wspomniano w komentarzach powyżej, model gbm byłby lepszy z pewnym dostrajaniem parametrów. Łatwym sposobem wykrycia problemów w modelu i potrzeby takich parametrów jest wygenerowanie pewnych wykresów diagnostycznych. Na przykład dla powyższego modelu GBM z domyślnymi parametrami (i używając pakietu plotmo do tworzenia wykresów) mamy

gbm.gaus <- gbm(y~., data = Data, dist = "gaussian")
library(plotmo)   # for the plotres function
plotres(gbm.gaus) # plot the error per ntrees and the residuals

co daje

Na wykresie po lewej stronie widzimy, że krzywa błędu nie osiągnęła najniższego poziomu. A na prawej wykresie resztki nie są tym, czego byśmy chcieli.

Jeśli przebudujemy model z większą liczbą drzew

gbm.gaus1 <- gbm(y~., data = Data, dist = "gaussian",
                 n.trees=5000, interact=3)
plotres(gbm.gaus1)

dostajemy

Widzimy krzywą błędu na dole z dużą liczbą drzew, a wykres resztek jest zdrowszy. Możemy również wykreślić wykresy częściowej zależności dla nowego modelu GBM i losowego modelu lasu

library(plotmo)
plotmo(gbm.gaus1, pmethod="partdep", all1=TRUE, all2=TRUE)
plotmo(rf.model,  pmethod="partdep", all1=TRUE, all2=TRUE)

co daje

Gbm i losowe wykresy modeli lasu są teraz podobne, zgodnie z oczekiwaniami.

Musisz zaktualizować interaction.depthparametr podczas budowania modelu wzmocnionego. Domyślnie przyjmuje wartość 1, co spowoduje, że wszystkie drzewa, które gbmbuduje algorytm , zostaną podzielone tylko raz. Oznaczałoby to, że każde drzewo dzieli się na zmienne ci w zależności od próbki obserwacji, z których korzysta, podzieli się gdzieś w okolicach 40 - 60.

Oto częściowe wykresy z interaction.depth = 3