Random Forests (R) — ranger¶
Fast random forests with built-in OOB and permutation importance¶
ranger is the standard fast R implementation of Breiman's random forest — the package cross-check for the from-scratch rf.py. It gives the out-of-bag error for free, both impurity and permutation importance, and Extremely Randomized Trees via one argument. Same data: Taiwan credit-card default (predict default next month), out of sample on a 30% held-out test set. (ROC-AUC = probability the model ranks a random defaulter above a random non-defaulter; 0.5 = chance, 1 = perfect; threshold-free and imbalance-robust — defined in the CART notebook.)
options(repr.plot.width=8, repr.plot.height=4.6)
.libPaths(c("C:/Users/user/R/win-library/4.6", .libPaths()))
for(p in c("ranger")) if(!requireNamespace(p,quietly=TRUE)) install.packages(p, repos="https://cloud.r-project.org", quiet=TRUE)
suppressMessages(library(ranger))
d <- read.csv("credit_default.csv"); d$default <- factor(d$default)
set.seed(0); i <- sample(nrow(d), 0.7*nrow(d)); tr <- d[i,]; te <- d[-i,]
auc <- function(y,p){ r<-rank(p); n1<-sum(y==1); n0<-sum(y==0); (sum(r[y==1])-n1*(n1+1)/2)/(n1*n0) }
yte <- as.integer(te$default==1)
cat(sprintf("credit default: %d train / %d test; default rate %.1f%%\n", nrow(tr), nrow(te), 100*mean(d$default==1)))
Warning message: "package 'ranger' was built under R version 4.6.1"
credit default: 21000 train / 9000 test; default rate 22.1%
1. Fit a forest, read the OOB error¶
ranger(default ~ ., probability=TRUE) grows the forest; its out-of-bag predictions ($predictions) score each training row using only the trees that never saw it — a free validation estimate that we compare to the held-out test AUC.
set.seed(1)
rf <- ranger(default ~ ., data=tr, num.trees=500, probability=TRUE, importance="impurity")
oob_auc <- auc(as.integer(tr$default==1), rf$predictions[,2])
te_auc <- auc(yte, predict(rf, te)$predictions[,2])
cat(sprintf("OOB AUC %.4f vs test AUC %.4f -> OOB is an honest free estimate\n", oob_auc, te_auc))
cat(sprintf("(OOB prediction error reported by ranger: %.4f)\n", rf$prediction.error))
OOB AUC 0.7725 vs test AUC 0.7737 -> OOB is an honest free estimate
(OOB prediction error reported by ranger: 0.1361)
2. Tuning mtry (the de-correlation knob)¶
mtry is R's name for max_features — how many features each split may consider. Small mtry de-correlates the trees (good) but weakens each one; the test AUC peaks at an intermediate value. We sweep it using the OOB estimate.
p <- ncol(tr)-1; grid <- c(1,2,3,round(sqrt(p)),6,10,p); grid <- sort(unique(grid))
res <- sapply(grid, function(m){ set.seed(1); f<-ranger(default~.,data=tr,num.trees=300,probability=TRUE,mtry=m); auc(as.integer(tr$default==1), f$predictions[,2]) })
plot(grid, res, type="b", pch=19, col="#2f855a", xlab="mtry (features per split)", ylab="OOB AUC", main="mtry: intermediate is best (de-correlation vs strength)")
abline(v=round(sqrt(p)), lty=2, col="#c53030"); legend("bottomright", sprintf("sqrt(p)=%d (default)",round(sqrt(p))), lty=2, col="#c53030", bty="n")
cat(sprintf("best mtry = %d (OOB AUC %.4f); using all %d features is worse -- de-correlation matters.\n", grid[which.max(res)], max(res), p))
best mtry = 1 (OOB AUC 0.7780); using all 23 features is worse -- de-correlation matters.
3. Importance — impurity is biased, permutation is honest¶
As in Python, we inject two pure-noise features (one continuous, one binary) and compare ranger's impurity importance (biased toward the many-valued continuous noise) with its permutation importance (which correctly ignores it).
set.seed(0); trn <- tr; ten <- te
trn$NOISE_cont <- rnorm(nrow(trn)); trn$NOISE_bin <- rbinom(nrow(trn),1,.5)
fi <- ranger(default~.,data=trn,num.trees=400,importance="impurity")
fp <- ranger(default~.,data=trn,num.trees=400,importance="permutation")
gi <- sort(fi$variable.importance, decreasing=TRUE); gp <- sort(fp$variable.importance, decreasing=TRUE)
par(mfrow=c(1,2), mar=c(4,7,3,1))
b<-head(gi,10); barplot(rev(b),horiz=TRUE,las=1,col=ifelse(grepl("NOISE",names(rev(b))),"#c53030","#2b6cb0"),main="Impurity importance (inflates NOISE_cont)",xlab="")
b<-head(gp,10); barplot(rev(b),horiz=TRUE,las=1,col=ifelse(grepl("NOISE",names(rev(b))),"#c53030","#2f855a"),main="Permutation importance (noise ~ 0)",xlab="")
par(mfrow=c(1,1))
cat(sprintf("impurity importance ranks NOISE_cont #%d of %d; permutation importance puts it near zero -- prefer permutation.\n",
which(names(gi)=="NOISE_cont"), length(gi)))
impurity importance ranks NOISE_cont #2 of 25; permutation importance puts it near zero -- prefer permutation.
4. Extremely Randomized Trees, and regression¶
splitrule="extratrees" picks split thresholds at random — more de-correlation, faster fits (Geurts et al. 2006). And ranger regresses by swapping the response for a numeric one: California housing median value, out of sample.
set.seed(1); et <- ranger(default~.,data=tr,num.trees=500,probability=TRUE,splitrule="extratrees")
cat(sprintf("Random Forest test AUC %.4f | Extra-Trees test AUC %.4f\n", te_auc, auc(yte,predict(et,te)$predictions[,2])))
h <- read.csv("cali_housing.csv"); set.seed(0); j <- sample(nrow(h),0.7*nrow(h)); htr<-h[j,]; hte<-h[-j,]
set.seed(1); rr <- ranger(MedHouseVal~.,data=htr,num.trees=500)
cat(sprintf("California housing regression forest: out-of-sample RMSE %.4f ($100k)\n", sqrt(mean((predict(rr,hte)$predictions-hte$MedHouseVal)^2))))
Random Forest test AUC 0.7737 | Extra-Trees test AUC 0.7762
California housing regression forest: out-of-sample RMSE 0.4906 ($100k)
5. How the forest improves — graphically¶
As in the Python notebook, two one-feature views over continuous predictors show what the forest buys. Classification: $P(\text{default})$ vs the credit limit — empirical proportions (quantile bins), logistic regression, one rpart tree, and the ranger forest. Regression: California value vs median income — linear fit, one tree, forest. Averaging hundreds of bootstrapped trees turns the single tree's coarse steps into a smooth curve, while still bending where the straight-line logit/linear baseline cannot.
suppressMessages(library(rpart))
par(mfrow=c(1,2), mar=c(4,4,3,1))
# classification: P(default) vs credit limit
lim <- tr$LIMIT_BAL/1000; yb <- as.integer(tr$default==1)
qb <- unique(quantile(lim, seq(0,1,length=16))); b <- findInterval(lim, qb, all.inside=TRUE)
ctr <- tapply(lim,b,mean); emp <- tapply(yb,b,mean); g <- seq(min(lim), quantile(lim,0.99), length=300)
lo <- predict(glm(yb~lim, family=binomial), data.frame(lim=g), type="response")
t1 <- predict(rpart(default~LIMIT_BAL, tr, method="class", control=rpart.control(maxdepth=4)), data.frame(LIMIT_BAL=g*1000), type="prob")[,2]
f1 <- predict(ranger(default~LIMIT_BAL, tr, probability=TRUE, num.trees=400, min.node.size=200), data.frame(LIMIT_BAL=g*1000))$predictions[,2]
plot(ctr, emp, pch=19, xlab="credit limit (NT$ thousands)", ylab="P(default)", main="Classification: P(default) vs credit limit", ylim=c(0,max(emp)*1.1))
lines(g,lo,col="#dd6b20",lwd=2); lines(g,t1,col="#a0aec0",lwd=2); lines(g,f1,col="#2f855a",lwd=2)
legend("topright", c("empirical","logistic","single tree","random forest"), col=c("black","#dd6b20","#a0aec0","#2f855a"), pch=c(19,NA,NA,NA), lty=c(NA,1,1,1), lwd=2, bty="n", cex=0.8)
# regression: California value vs median income
h <- read.csv("cali_housing.csv"); set.seed(0); jj <- sample(nrow(h),0.7*nrow(h)); htr <- h[jj,]
xi <- htr$MedInc; gh <- seq(min(xi), quantile(xi,0.99), length=300)
ols <- predict(lm(MedHouseVal~MedInc, htr), data.frame(MedInc=gh))
trg <- predict(rpart(MedHouseVal~MedInc, htr, control=rpart.control(maxdepth=4)), data.frame(MedInc=gh))
frg <- predict(ranger(MedHouseVal~MedInc, htr, num.trees=400, min.node.size=40), data.frame(MedInc=gh))$predictions
plot(xi, htr$MedHouseVal, pch=".", col="grey", xlab="median income", ylab="median house value ($100k)", main="Regression: value vs income", xlim=c(min(xi),quantile(xi,0.99)))
lines(gh,ols,col="#dd6b20",lwd=2); lines(gh,trg,col="#a0aec0",lwd=2); lines(gh,frg,col="#2f855a",lwd=2)
legend("topleft", c("linear","single tree","random forest"), col=c("#dd6b20","#a0aec0","#2f855a"), lty=1, lwd=2, bty="n", cex=0.8)
par(mfrow=c(1,1))
cat("The forest curve is smooth (variance averaged over 400 bootstrapped trees) and bent (nonlinearity kept) --\n")
cat("smoother than the single tree, more flexible than the straight-line logit/linear baseline.\n")
The forest curve is smooth (variance averaged over 400 bootstrapped trees) and bent (nonlinearity kept) --
smoother than the single tree, more flexible than the straight-line logit/linear baseline.
Proportions vs predictions — the regression forest¶
Bin the held-out block groups by the forest's predicted value and plot the mean actual value in each bin, the same diagnostic the rest of the collection uses. Averaging hundreds of trees should remove the step-function coarseness a single tree shows — and should also, by construction, shrink predictions toward the mean, which this view makes visible where RMSE does not.
options(repr.plot.width=7, repr.plot.height=5); par(mar=c(4,4,3,1))
prd <- predict(rr, hte)$predictions; act <- hte$MedHouseVal
q <- unique(quantile(prd, seq(0,1,length=11))); b <- cut(prd, q, include.lowest=TRUE)
pm <- as.numeric(tapply(prd,b,mean)); am <- as.numeric(tapply(act,b,mean))
plot(pm, am, type="b", pch=19, col="#2f855a", lwd=2, xlab="random-forest predicted ($100k)",
ylab="mean actual value", main=sprintf("Proportions vs predictions (RMSE %.3f)", sqrt(mean((prd-act)^2))))
abline(0,1,lty=2); legend("topleft", c("perfect","binned means"), lty=c(2,1), pch=c(NA,19),
col=c("black","#2f855a"), lwd=2, bty="n", cex=0.8)
cat(sprintf("distinct predicted values: %d (a single tree emits far fewer)\n", length(unique(round(prd,6)))))
cat(sprintf("binned gaps (actual - predicted) run %+.3f to %+.3f; slope of actual on predicted %.3f\n",
min(am-pm), max(am-pm), coef(lm(act~prd))[2]))
cat("\nThe forest predicts a continuum rather than a handful of constants -- the coarseness of a single tree is gone.\n")
cat("The slope is the quantity to read, and the benchmark is 1.0: an optimal forecast has cov(actual, prediction)\n")
cat("equal to var(prediction), which makes this regression's slope exactly 1. Above 1 means the predictions are\n")
cat("COMPRESSED toward the mean -- they move less than the outcome does -- which is precisely what averaging many\n")
cat("trees produces, since each additional tree pulls the ensemble back toward the sample mean.\n")
cat("\nThat is the cost side of bagging's variance reduction, and it is invisible in RMSE: the forest buys stability\n")
cat("by under-reacting at the extremes. Compare the single tree in the CART notebook, which is coarse but close to\n")
cat("unbiased in slope -- the ensemble trades one defect for the other.\n")
distinct predicted values: 6189 (a single tree emits far fewer)
binned gaps (actual - predicted) run -0.103 to +0.182; slope of actual on predicted 1.092
The forest predicts a continuum rather than a handful of constants -- the coarseness of a single tree is gone.
The slope is the quantity to read, and the benchmark is 1.0: an optimal forecast has cov(actual, prediction)
equal to var(prediction), which makes this regression's slope exactly 1. Above 1 means the predictions are
COMPRESSED toward the mean -- they move less than the outcome does -- which is precisely what averaging many
trees produces, since each additional tree pulls the ensemble back toward the sample mean.
That is the cost side of bagging's variance reduction, and it is invisible in RMSE: the forest buys stability
by under-reacting at the extremes. Compare the single tree in the CART notebook, which is coarse but close to
unbiased in slope -- the ensemble trades one defect for the other.
6. Summary¶
ranger reproduces the from-scratch random forest: a large out-of-sample AUC on credit default (well above the single tree of the CART notebook), an out-of-bag estimate that tracks the test error for free, an intermediate-mtry optimum from de-correlation, the impurity-vs-permutation importance bias (and its fix), and Extra-Trees for extra de-correlation. This is the package mirror of rf_python.ipynb.
Random forests cut error by averaging independent low-bias trees to kill variance. The next notebook takes the opposite tack — boosting adds dependent small trees in sequence to cut bias (Friedman 2001), leading to XGBoost / LightGBM / CatBoost.