class: center, middle, title-slide .title[ # Gradiente Boosting ] .subtitle[ ## Em regressão ] .author[ ### Jaime Utria ] .institute[ ### Departamento de Estatística - UFF ] --- <!-- macro comandos matemáticos Latex --> <script type="text/x-mathjax-config"> MathJax.Hub.Config({ TeX: { Macros: { vet: ["{\\mathbf #1}",1], prodint: ["\\langle #1, #2 \\rangle",2], posto: ["{\\mathrm{posto} (#1)}",1], tr: ["{\\mathrm{tr} (#1)}",1] } } }); </script> ## Gradiente Boosting: Ideias principais. * Como antes queremos construir um preditor forte baseado na combinação de uma sequência de preditores fracos. * Aqui consideramos como preditores fracos árvores de decisão com tipicamente 8 a 32 nós terminais (não tocos, como no caso de AdaBoost). * O preditor fraco inicial será uma árvore com um único nó raiz, os demais preditores são árvores com o número de nós folhas predeterminado. * Cada preditor subsequente tenta corrigir o erro de seu predecessor minimizando o gradiente da função de perda (pseudo-resíduo). * Em problemas de regressão, usa-se geralmente a função de perda quadrática, ao passo que em problemas de classificação usa-se a perda logaritmica. --- ## Gradiente Boosting em Regressão ``` r Boston_toy ``` ``` ## medv crim nox lstat ## 142 14.4 1.62864 0.624 34.41 ## 51 19.7 0.08873 0.439 13.45 ## 208 22.5 0.25199 0.489 18.06 ## 218 28.7 0.07013 0.550 9.69 ## 220 23.0 0.11425 0.550 10.50 ## 152 19.6 1.49632 0.871 13.28 ## 314 21.6 0.26938 0.544 7.90 ## 93 22.9 0.04203 0.464 8.16 ## 75 24.1 0.07896 0.437 6.78 ## 352 24.1 0.07950 0.411 5.49 ``` --- Vamos construir um preditor gradiente boosting com 4 árvores em que cada árvore têm 3 nós folhas (*toy example*). **Passo 1**: Construímos um preditor inicial como sendo uma árvore com um nó. ``` r # predição 1 pred1 <- mean(Boston_toy$medv) # resíduos 1 r1 <- Boston_toy$medv - pred1 ``` .center[ <!-- --> ] --- **Passo 2**: Sucessivamente construímos uma árvore de decisão para prever os resíduos ``` r dados2 <- data.frame(Boston_toy[,-1],r1) dados2 ``` ``` ## crim nox lstat r1 ## 142 1.62864 0.624 34.41 -7.66 ## 51 0.08873 0.439 13.45 -2.36 ## 208 0.25199 0.489 18.06 0.44 ## 218 0.07013 0.550 9.69 6.64 ## 220 0.11425 0.550 10.50 0.94 ## 152 1.49632 0.871 13.28 -2.46 ## 314 0.26938 0.544 7.90 -0.46 ## 93 0.04203 0.464 8.16 0.84 ## 75 0.07896 0.437 6.78 2.04 ## 352 0.07950 0.411 5.49 2.04 ``` --- ``` r rpart2 <- rpart(r1~.,data = dados2,control = rpart.control(minsplit = 5)) ``` .center[ <!-- --> ] --- ``` r predictions2 <- predict(rpart2,dados2) r2 <- dados2$r1 - predictions2 dados3 <- cbind.data.frame(dados2[,-4],r2) dados3 ``` ``` ## crim nox lstat r2 ## 142 1.62864 0.624 34.41 -2.60 ## 51 0.08873 0.439 13.45 -2.00 ## 208 0.25199 0.489 18.06 0.80 ## 218 0.07013 0.550 9.69 3.75 ## 220 0.11425 0.550 10.50 1.30 ## 152 1.49632 0.871 13.28 2.60 ## 314 0.26938 0.544 7.90 -0.10 ## 93 0.04203 0.464 8.16 -2.05 ## 75 0.07896 0.437 6.78 -0.85 ## 352 0.07950 0.411 5.49 -0.85 ``` --- .center[ ``` r rpart3 <- rpart(r2~.,data = dados3,control = rpart.control(minsplit = 5)) rpart.plot(rpart3, type = 4, extra = 101) ``` <!-- --> ] --- ``` r predictions3 <- predict(rpart3,dados3) r3 <- dados3$r2 - predictions3 dados4 <- cbind.data.frame(dados3[,-4],r3) residuos <- data.frame(r1,r2,r3) residuos ``` ``` ## r1 r2 r3 ## 142 -7.66 -2.60 -1.7000 ## 51 -2.36 -2.00 -0.5625 ## 208 0.44 0.80 1.7000 ## 218 6.64 3.75 1.8625 ## 220 0.94 1.30 -0.5875 ## 152 -2.46 2.60 0.7125 ## 314 -0.46 -0.10 -1.9875 ## 93 0.84 -2.05 -0.6125 ## 75 2.04 -0.85 0.5875 ## 352 2.04 -0.85 0.5875 ``` --- ``` r rpart4 <- rpart(r3~.,data = dados4,control = rpart.control(minsplit = 5)) ``` .center[ ``` r rpart.plot(rpart4, type = 4, extra = 101) ``` <!-- --> ] --- ``` r predictions4 <- predict(rpart4,dados4) r4 <- dados4$r3 - predictions4 residuos <- data.frame(r1,r2,r3,r4) residuos ``` ``` ## r1 r2 r3 r4 ## 142 -7.66 -2.60 -1.7000 -0.7083333 ## 51 -2.36 -2.00 -0.5625 -0.5625000 ## 208 0.44 0.80 1.7000 0.7083333 ## 218 6.64 3.75 1.8625 0.8708333 ## 220 0.94 1.30 -0.5875 -1.5791667 ## 152 -2.46 2.60 0.7125 1.7041667 ## 314 -0.46 -0.10 -1.9875 -0.9958333 ## 93 0.84 -2.05 -0.6125 -0.6125000 ## 75 2.04 -0.85 0.5875 0.5875000 ## 352 2.04 -0.85 0.5875 0.5875000 ``` --- **Fazendo predições com Gradiente Boosting** <img src="data:image/png;base64,#Aula16_Gradiente_Boosting_files/figure-html/unnamed-chunk-14-1.svg" style="display: block; margin: auto;" /> Vamos prever, por exemplo, o `medv` para a primeira amostra, isto é, para [crim = 1,63; nox = 0,62, lstat = 34,41]. Logo, a previsão é dada por: 22,06 -5,1 - 0,9 - 0,99 = 15,07. --- ## Algoritmo 1. Iniciar `\(\hat{f}_0(\vet x) = \text{argmin}_\gamma\sum_{i=1}^n L(y_i,\gamma)\)`. 2. Para `\(b = 1, \ldots, B\)`: (a) Calcular os pseudo-resíduos `\(r_{ib} = - \left[\frac{\partial L(y_i,\hat{f}(\mathbf{x}_i))}{\partial f(\mathbf{x}_i)}\right]_{\hat f = \hat f_{b-1}}\)`. (b) Ajustar uma árvore de regressão aos pseudo-resíduos `\(r_{ib}\)`, `\(i=1,\ldots,n\)` com regiões terminais `\(R_{jb}\)`, `\(j=1,\ldots,J_b\)`. (c) Calcular `\(\gamma_{jb} =\text{argmin}_\gamma \sum_{\vet x_i \in R_{jb}} L(y_i,\hat f_{b-1}(\mathbf x_i)+\gamma)\)` (d) Atualizar `\(\hat f_b(\mathbf{x}) = \hat f_{b-1}(\mathbf{x}) + \eta \sum_{j=1}^{J_b} \gamma_{jb}I(\mathbf{x} \in R_{jb})\)` 3. Preditor final: `\(\hat f_{\text{GB}}(\mathbf{x}) = \hat f_B(\mathbf{x})\)`. --- - A função de perda usada para problemas de regressão é a perda qudrática definida como segue: `$$L(y,\hat y) = \frac{1}{2}(y - \hat y)^2,$$` já seu gradiente (negativo) é dado por: `$$- \frac{\partial L(y_i,\hat{f}(\mathbf{x}_i))}{\partial f(\mathbf{x}_i)} = y_i - \hat{f}(\mathbf{x}_i)$$` que coincide com o resíduo (erro) obtido ao prever `\(y\)` por `\(\hat f\)`. - Daí, vemos que a ideia do gradiente boosting consiste em que cada preditor busca minimizar o erro de previsão do seu predecessor. --- ## Boston ``` r data(Boston, package = "MASS") library(caret) library(gbm) # para gradiente boosting set.seed(12345) index_train <- createDataPartition(y = Boston$medv, p = 0.7, list = F) training <- Boston[index_train,] testing <- Boston[-index_train,] gbm_boston <- gbm(medv ~., data = training, distribution = "gaussian", # para regressão n.trees = 5000, interaction.depth = 1, shrinkage = 0.1, bag.fraction = 1, cv.folds = 10) ``` --- ``` r gbm_boston ``` ``` ## gbm(formula = medv ~ ., distribution = "gaussian", data = training, ## n.trees = 5000, interaction.depth = 1, shrinkage = 0.1, bag.fraction = 1, ## cv.folds = 10) ## A gradient boosted model with gaussian loss function. ## 5000 iterations were performed. ## The best cross-validation iteration was 418. ## There were 13 predictors of which 12 had non-zero influence. ``` --- **Importância das variáveis** ``` r summary(gbm_boston, plotit = FALSE) ``` ``` ## var rel.inf ## lstat lstat 42.3490491 ## rm rm 36.2975225 ## crim crim 5.5666396 ## nox nox 4.9455821 ## dis dis 4.3254268 ## ptratio ptratio 2.2172029 ## black black 1.4956989 ## tax tax 1.1463872 ## age age 0.7403533 ## indus indus 0.6017425 ## rad rad 0.1313249 ## chas chas 0.1231708 ## zn zn 0.0598993 ``` --- ``` r summary(gbm_boston) ``` <img src="data:image/png;base64,#Aula16_Gradiente_Boosting_files/figure-html/unnamed-chunk-18-1.png" style="display: block; margin: auto;" /> ``` ## var rel.inf ## lstat lstat 42.3490491 ## rm rm 36.2975225 ## crim crim 5.5666396 ## nox nox 4.9455821 ## dis dis 4.3254268 ## ptratio ptratio 2.2172029 ## black black 1.4956989 ## tax tax 1.1463872 ## age age 0.7403533 ## indus indus 0.6017425 ## rad rad 0.1313249 ## chas chas 0.1231708 ## zn zn 0.0598993 ``` --- ## Avaliando o desempenho do ajuste ``` r predictions <- predict(gbm_boston , newdata = testing, n.trees = 5000) postResample(predictions,testing$medv) ``` ``` ## RMSE Rsquared MAE ## 2.795457 0.901135 2.171551 ``` --- ## Número ideal de árvores ``` r gbm.perf(gbm_boston) ``` <img src="data:image/png;base64,#Aula16_Gradiente_Boosting_files/figure-html/unnamed-chunk-20-1.svg" style="display: block; margin: auto;" /> ``` ## [1] 418 ``` --- ## Avaliando o desempenho do ajuste (com o número de árvores escolhido) ``` r prediction_cv<-predict(gbm_boston,testing, n.trees=418) postResample(prediction_cv,testing$medv) ``` ``` ## RMSE Rsquared MAE ## 2.7660788 0.9019897 2.1701701 ```