Improves random forest to use in EMR.

parent 9b1b8b33
......@@ -119,11 +119,13 @@ predict(m0)
# Valores preditos (médias em cada região).
unique(sort(predict(m0)))
table(sort(predict(m0)))
#-----------------------------------------------------------------------
# Deixar a árvore crescer mais.
m1 <- rpart(prod ~ ., data = teca,
m1 <- rpart(prod ~ .,
data = teca,
control = list(minsplit = 5,
cp = 0.001))
......@@ -143,9 +145,23 @@ summary(m0)
# Visualização da árvore de regressão.
rpart.plot(m0)
# Visualização alternativa.
plot(m0)
text(m0)
# Valores preditos (médias em cada região).
unique(sort(predict(m0)))
# Importância das variáveis.
cbind(m0$variable.importance)
# Soma de quadrados residual do modelo nulo.
m0$frame$dev[1] # SQres do ~1.
(m0$frame$n[1] - 1) * var(ap$lpreco) # SQres do ~1.
# R².
1 - sum(residuals(m0)^2)/(m0$frame$dev[1])
# Criando um grid nas variáveis consideradas pela árvore.
grid <- with(ap,
expand.grid(larea = seq(min(larea, na.rm = TRUE),
......@@ -179,6 +195,27 @@ cloud(lpreco ~ larea + vagas, data = ap, col = yp)
# Os patamares.
wireframe(y ~ larea + vagas, data = grid, drape = TRUE)
# Deixar a árvore crescer mais.
m0 <- rpart(lpreco ~ .,
data = ap,
control = list(cp = 0.0025))
# Visualização da árvore de regressão.
rpart.plot(m0)
# Visualização alternativa.
plot(m0)
text(m0)
# Valores preditos (médias em cada região).
unique(sort(predict(m0)))
# Importância das variáveis.
cbind(m0$variable.importance)
# R².
1 - sum(residuals(m0)^2)/(m0$frame$dev[1])
```
## Árvores de regressão com bagging
......@@ -225,14 +262,18 @@ set.seed(102030)
B <- 200
j <- 1
frac <- numeric(B)
fits <- replicate(B, simplify = FALSE, {
i <- sample(s, size = n, replace = TRUE)
frac[j] <<- length(unique(i))/n
j <<- j + 1
ap_bs <- ap[i, ]
m_bs <- rpart(lpreco ~ ., data = ap_bs)
return(m_bs)
})
fits <- replicate(B,
simplify = FALSE,
expr = {
# Reamostra com reposição.
i <- sample(s, size = n, replace = TRUE)
frac[j] <<- length(unique(i))/n
j <<- j + 1
ap_bs <- ap[i, ]
# Ajuste da árvore aos dados de treino.
m_bs <- rpart(lpreco ~ ., data = ap_bs)
return(m_bs)
})
# A proporção de valores usados nas amostras bootstrap.
mean(frac)
......@@ -241,27 +282,36 @@ mean(frac)
pred$y <- sapply(fits, FUN = predict, newdata = pred)
str(pred)
# O predito médio.
# Predição para os primeiros casos.
head(pred)
# Estatísticas para o B valores preditos para alguns casos.
mean(as.vector(pred[1, "y"]))
var(as.vector(pred[1, "y"]))
# O predito médio ("a sabedoria das multidões").
pred$ym <- rowMeans(pred$y)
# Predito contra observado.
# x11()
xyplot(pred$ym ~ ap$lpreco, aspect = "iso") +
layer(panel.abline(a = 0, b = 1))
# Qual o predito para o imóvel mediano?
new <- lapply(pred, FUN = function(x) {
if (is.numeric(x)) {
median(x, na.rm = TRUE)
} else {
levels(x)[1]
}
})
new <- lapply(pred,
FUN = function(x) {
if (is.numeric(x)) {
median(x, na.rm = TRUE)
} else {
levels(x)[1]
}
})
new
# Predito por cada árvore.
y <- sapply(fits, FUN = predict, newdata = new)
# Distribuição dos valores preditos e valor médio.
# Distribuição dos B valores preditos e valor médio.
plot(density(y))
rug(y)
m <- mean(y)
......@@ -275,7 +325,7 @@ library(ipred)
# help(package = "ipred", help_type = "html")
# Fazendo bagging.
bg <- bagging(lpreco ~ ., data = ap, nbagg = 50, coob = TRUE)
bg <- bagging(lpreco ~ ., data = ap, nbagg = 200, coob = TRUE)
bg
# Predito contra observado.
......@@ -312,18 +362,26 @@ rpart.plot(m0)
# Replicar.
set.seed(302010)
B <- 1000
rf <- replicate(B, simplify = FALSE, {
v <- sample(xvars, size = nv, replace = FALSE)
i <- sample(1:n, size = n, replace = TRUE)
m0 <- rpart(prod ~ .,
data = teca[i, c(v, "prod")],
control = list(minsplit = 3,
cp = 0.01))
return(m0)
})
rf <- replicate(B,
simplify = FALSE,
expr = {
v <- sample(xvars, size = nv, replace = FALSE)
i <- sample(1:n, size = n, replace = TRUE)
m0 <- rpart(prod ~ .,
data = teca[i, c(v, "prod")],
control = list(minsplit = 3,
cp = 0.01))
return(m0)
})
# ATTENTION: na árvore de regressão, o sorteio das variáveis é feito
# após cada split e não uma única vez como o que está neste código.
# Obtenção dos preditos.
y_rf <- sapply(rf, FUN = predict, newdata = teca)
head(y_rf[, 1:6])
# Cálculo da média.
ym <- rowMeans(y_rf)
xyplot(ym ~ teca$prod,
......@@ -338,7 +396,37 @@ cor(ym, teca$prod)
library(randomForest)
rf <- randomForest(prod ~ .,
data = teca,
ntree = 3,
mtry = 2,
keep.inbag = TRUE,
keep.forest = TRUE)
rf
1 - sum(rf$oob.times == 0)/length(rf$oob.times)
str(rf)
# Número de vezes que cada variável foi usada em cada bag.
head(rf$inbag)
# Valores preditos.
sort(unique(predict(rf)))
# Inspecionando uma das árvores da floresta.
t <- 1
one_tree <- getTree(rf, k = t, labelVar = TRUE)
nrow(one_tree)
rf$forest$ndbigtree[t]
one_tree
# Outra forma de acessar as variáveis usadas.
rf$forest$bestvar[, t]
names(rf$forest$xlevels)[rf$forest$bestvar[, t]]
#-----------------------------------------------------------------------
# help(package = "randomForest", help_type = "html")
rf <- randomForest(prod ~ .,
data = teca,
ntree = B,
......@@ -357,6 +445,17 @@ xyplot(ym + yp ~ teca$prod,
type = c("p", "smooth")) +
layer(panel.abline(a = 0, b = 1))
# ATTENTION: a randomForest() criou árvores mais profundas que o código
# didático feito algumas linhas acima. Além do mais, as `m` variáveis
# são sorteadas a cada split e não uma vez apenas.
# o <- order(yp)
# xyplot(ym[o] + yp[o] ~ seq_along(ym), auto.key = TRUE)
```
Para construir o gráfico de uma árvore, leia esse post: <https://stats.stackexchange.com/questions/41443/how-to-actually-plot-a-sample-tree-from-randomforestgettree>.
```{r}
#-----------------------------------------------------------------------
# Dados de imóveis.
......@@ -382,10 +481,22 @@ xyplot(yp ~ ap2$lpreco,
auto.key = TRUE,
type = c("p", "smooth")) +
layer(panel.abline(a = 0, b = 1))
importance(rf, type = 2)
im <- importance(rf, type = 2)
100 * im/sum(im)
varImpPlot(rf)
```
## Árvore de regressão com boosting
Veja as ilustrações aqui para entender:
<https://medium.com/mlreview/gradient-boosting-from-scratch-1e317ae4587d>.
Visite também [gradiente descendente boosting]("./02-gradient-methods.html#gradiente-descendente-boosting").
```{r}
library(gbm)
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment