Quantile Regression Forests vs. Random Forests#
This example compares the predictions generated by a quantile regression forest (QRF) and a standard random forest regressor (RF) on a synthetic dataset with right-skewed, heavy-tailed noise. In a right-skewed distribution, the mean is to the right of the median, and a few large values pull it further. This example demonstrates how the median (quantile = 0.5) predicted by a quantile regressor (QRF) can be a more reliable estimator than the mean predicted by a standard random forest when dealing with such distributions: the median minimizes the absolute error and is less sensitive to the large values, so its mean absolute error is lower. The QRF can predict any other quantile too, selected with the slider.
import altair as alt
import numpy as np
import pandas as pd
from sklearn.ensemble import RandomForestRegressor
from sklearn.model_selection import train_test_split
from sklearn.utils.validation import check_random_state
from quantile_forest import RandomForestQuantileRegressor
random_state = np.random.RandomState(0)
n_samples = 5000
quantiles = np.linspace(0, 1, num=101, endpoint=True).round(2).tolist()
def make_skewed_dataset(n_samples, sigma=1, random_state=None):
"""Make a dataset with right-skewed, heavy-tailed (lognormal) noise."""
random_state = check_random_state(random_state)
X = random_state.uniform(0, 10, size=(n_samples, 2))
y = 2 * np.sin(X[:, 0]) + random_state.lognormal(sigma=sigma, size=n_samples)
return X, y
# Create a toy dataset with right-skewed noise.
X, y = make_skewed_dataset(n_samples, sigma=1, random_state=0)
regr_rf = RandomForestRegressor(min_samples_leaf=5, random_state=random_state)
regr_qrf = RandomForestQuantileRegressor(min_samples_leaf=5, random_state=random_state)
X_train, X_test, y_train, y_test = train_test_split(X, y, random_state=random_state)
regr_rf.fit(X_train, y_train)
regr_qrf.fit(X_train, y_train)
y_pred_rf = regr_rf.predict(X_test) # standard RF predictions (mean)
y_pred_qrf = regr_qrf.predict(X_test, quantiles=quantiles) # QRF predictions (quantiles)
legend = {
"Actual": "#c0c0c0",
"RF (Mean)": "#f2a619",
"QRF (Quantile)": "#006aff",
}
df = pd.DataFrame(
{
"actual": y_test,
"rf": y_pred_rf,
**{f"qrf_{q_i:.3g}": y_i.ravel() for q_i, y_i in zip(quantiles, y_pred_qrf.T)},
}
)
# Mean absolute error of the RF mean and of each QRF quantile.
df_mae = pd.DataFrame(
{
"quantile": quantiles,
"qrf": np.abs(y_pred_qrf - y_test[:, np.newaxis]).mean(axis=0),
"rf": np.abs(y_pred_rf - y_test).mean(),
}
)
def plot_predictions_and_errors(df, df_mae, legend, max_value=8):
"""Plot histogram of predictions by model, and the mean absolute errors."""
# Slider for varying the quantile value used for generating the QRF predictions.
slider = alt.binding_range(
name="Predicted Quantile: ",
min=0,
max=1,
step=0.5 if len(quantiles) == 1 else 1 / (len(quantiles) - 1),
)
quantile_val = alt.param(name="quantile", value=0.5, bind=slider)
click = alt.selection_point(bind="legend", fields=["label"], on="click")
histogram = (
alt.Chart(df)
.add_params(click)
.transform_calculate(qrf_col="'qrf_' + quantile")
.transform_calculate(qrf="datum[datum.qrf_col]")
.transform_calculate(calculate="round(datum.actual * 4) / 4", as_="Actual")
.transform_calculate(calculate="round(datum.rf * 4) / 4", as_="RF (Mean)")
.transform_calculate(calculate="round(datum.qrf * 4) / 4", as_="QRF (Quantile)")
.transform_fold(list(legend.keys()), as_=["label", "value"])
.transform_filter(f"datum.value <= {max_value}") # omit the long tail of large values
.mark_bar()
.encode(
x=alt.X(
"value:N",
axis=alt.Axis(
labelAngle=0,
labelExpr="datum.value % 1 == 0 ? datum.value : null",
),
title=f"Actual and Predicted Target Values (Up to {max_value})",
),
y=alt.Y("count():Q", axis=alt.Axis(format=",d", title="Counts")),
color=alt.condition(
click,
alt.Color("label:N", sort=list(legend.keys()), title=None),
alt.value("lightgray"),
),
opacity=alt.condition(click, alt.value(1), alt.value(0.5)),
xOffset=alt.XOffset("label:N"),
tooltip=[
alt.Tooltip("label:N", title="Label"),
alt.Tooltip("value:Q", title="Value (Binned)"),
alt.Tooltip("count():Q", format=",d", title="Counts"),
],
)
.properties(title="Distribution of Actual and Predicted Values", height=300, width=650)
)
tooltip = [
alt.Tooltip("quantile:Q", format=".2f", title="Quantile"),
alt.Tooltip("qrf:Q", format=",.3f", title="QRF Mean Absolute Error"),
alt.Tooltip("rf:Q", format=",.3f", title="RF Mean Absolute Error"),
]
# Show errors up to twice that of the RF; those of extreme quantiles are larger.
y_domain = [np.floor(df_mae["qrf"].min() * 10) / 10, np.ceil(df_mae["rf"].max() * 20) / 10]
line_qrf = (
alt.Chart(df_mae)
.mark_line(clip=True, color=legend["QRF (Quantile)"], size=3)
.encode(
x=alt.X("quantile:Q", title="Predicted Quantile"),
y=alt.Y("qrf:Q", scale=alt.Scale(domain=y_domain), title="Mean Absolute Error"),
tooltip=tooltip,
)
)
rule_rf = (
alt.Chart(df_mae)
.mark_rule(color=legend["RF (Mean)"], size=3, strokeDash=[4, 4])
.encode(y=alt.Y("mean(rf):Q"), tooltip=[tooltip[2]])
)
point_qrf = (
alt.Chart(df_mae)
.transform_filter("abs(datum.quantile - quantile) < 1e-6")
.mark_circle(clip=True, color=legend["QRF (Quantile)"], opacity=1, size=100)
.encode(x=alt.X("quantile:Q"), y=alt.Y("qrf:Q"), tooltip=tooltip)
)
errors = (line_qrf + rule_rf + point_qrf).properties(
title="Mean Absolute Error of Predictions", height=200, width=650
)
chart = (
alt.vconcat(histogram, errors)
.add_params(quantile_val)
.configure_range(category=alt.RangeScheme(list(legend.values())))
)
return chart
chart = plot_predictions_and_errors(df, df_mae, legend)
chart