Ensemble models¶

gwlearn implements geographically weighted versions of ensemble models for both classification and regression tasks. For classification, there’s a caveat: given it is unlikely that all categories are present in all local models, fitting a non-binary would lead to inconsistent local models. Hence gwlearn currently supports only binary classification.

This notebook outlines how to fit the models and extract their attributes. For details on prediction, see the Prediction guide.

import geopandas as gpd
from geodatasets import get_path
from sklearn import metrics

from gwlearn.ensemble import (
    GWGradientBoostingClassifier,
    GWGradientBoostingRegressor,
    GWRandomForestClassifier,
    GWRandomForestRegressor,
)

Get sample data

gdf = gpd.read_file(get_path("geoda.south")).to_crs(5070)

gdf["point"] = gdf.representative_point()
gdf = gdf.set_geometry("point")

y = gdf["FH90"] > gdf["FH90"].median()  # binary target for classification
X = gdf.iloc[:, 9:15]

Random Forest Classifier¶

The implementation of geographically-weighted random forest classifier follows the logic of linear models, where each local neighborhood defined by a set bandwidth is used to fit a single local model.

gwrf = GWRandomForestClassifier(
    bandwidth=150,
    fixed=False,
)
gwrf.fit(
    X,
    y,
    geometry=gdf.geometry,
)

GWRandomForestClassifier(bandwidth=150)

Focal score¶

The performance of these models can be measured in three ways. The first one is using the focal prediction. Unlike in linear models, where the focal observation is part of the local model training, in this case it is excluded to ensure it can be used to evaluate the model. Otherwise, it would be a part of the data the local model has seen and would report unrealistically high performance. If you want to include the focal observation in the training nevertheless, just set include_focal=True.

Focal accuracy can be measured as follows. Note that some local models are not fitted due to imbalance rules, and report NA that needs to be filtered out.

na_mask = gwrf.pred_.notna()

metrics.accuracy_score(y[na_mask], gwrf.pred_[na_mask])

0.7595258255715496

Pooled out-of-bag score¶

Another option is to pull the out-of-bag predictions from individual local models and pool them together. This uses more data to evaluate the model but given the local models are tuned for their focal location, some predictions on locations far from the focal point may be artifically worse than the actual model prediction would be.

metrics.accuracy_score(gwrf.oob_y_pooled_, gwrf.oob_pred_pooled_)

0.7442901495907424

Local score¶

The final option is to take the local out-of-bag predicitons and measure performance per each local model. To do that, you can use a method local_metric().

local_accuracy = gwrf.local_metric(metrics.accuracy_score)
local_accuracy

array([0.66666667, 0.68666667, 0.72      , ..., 0.62      , 0.72      ,
       0.68666667], shape=(1412,))

gdf.set_geometry("geometry").plot(
    local_accuracy, legend=True, missing_kwds=dict(color="lightgray")
).set_axis_off()

Feature importance¶

Feature importances are reported for each local model.

gwrf.feature_importances_

gdf.set_geometry("geometry").plot(
    gwrf.feature_importances_["HC60"], legend=True
).set_axis_off()

You can compare all of this to values extracted from a global model, fitted alongside.

gwrf.global_model.feature_importances_

array([0.14155788, 0.15740365, 0.22588841, 0.25075606, 0.10676568,
       0.11762833])

gwrf.feature_importances_.mean()

HR60    0.154666
HR70    0.159607
HR80    0.202272
HR90    0.221558
HC60    0.139777
HC70    0.122120
dtype: float64

Random Forest Regressor¶

The GWRandomForestRegressor provides the same functionality as the classifier version but for continuous target variables. It uses local Random Forest models to capture spatial heterogeneity in regression relationships.

y_reg = gdf["FH90"]

gwrf_reg = GWRandomForestRegressor(
    bandwidth=150,
    fixed=False,
    random_state=0,
)
gwrf_reg.fit(
    X,
    y_reg,
    geometry=gdf.representative_point(),
)

GWRandomForestRegressor(bandwidth=150, random_state=0)

The regressor provides local R² values showing how well the model fits in each neighborhood:

gwrf_reg.local_r2_.head()

0    0.104145
1    0.035537
2   -0.039846
3   -0.048805
4    0.637768
dtype: float64

gdf.set_geometry("geometry").plot(gwrf_reg.local_r2_, legend=True).set_axis_off()

Out-of-bag performance¶

Like the classifier, the regressor supports out-of-bag predictions that can be used to evaluate model performance:

from sklearn.metrics import r2_score

r2_score(gwrf_reg.oob_y_pooled_, gwrf_reg.oob_pred_pooled_)

0.5014741368618677

Feature importance¶

Feature importances show the contribution of each feature in local models:

gwrf_reg.feature_importances_.head()

gdf.set_geometry("geometry").plot(
    gwrf_reg.feature_importances_["HC60"], legend=True
).set_axis_off()

Gradient Boosting Classifier¶

If you prefer to use gradient boosting, there is a minimal implementation of geographically weighted gradient boosting classifier, following the same model described for the random forest above.

gwgb = GWGradientBoostingClassifier(
    bandwidth=150,
    fixed=False,
)
gwgb.fit(
    X,
    y,
    geometry=gdf.geometry,
)

GWGradientBoostingClassifier(bandwidth=150)

Given the nature of the model, the outputs are a bit more limited. You can still extract focal predictions.

nan_mask = gwgb.pred_.notna()

metrics.accuracy_score(y[nan_mask], gwgb.pred_[nan_mask])

0.7502116850127011

And local feature importances.

gwgb.feature_importances_

Leave out samples¶

However, the pooled data are not available. In this case, you can use the leave_out keyword to leave out a fraction (when float) or a set number (when int) of random observations from each local model. For these, the local model does prediction and returns as left_out_proba_ and left_out_y_ arrays.

gwgb_leave = GWGradientBoostingClassifier(
    bandwidth=150,
    fixed=False,
    leave_out=0.2,
)
gwgb_leave.fit(
    X,
    y,
    geometry=gdf.geometry,
)

GWGradientBoostingClassifier(bandwidth=150)

metrics.accuracy_score(
    gwgb_leave.left_out_y_, gwgb_leave.left_out_proba_.argmax(axis=1)
)

0.7224103866779565

Gradient Boosting Regressor¶

In the GWGradientBoostingRegressor, gradient boosting is applied for regression, taking into consideration the spatial heterogeneity. Gradient boosting is a type of boosting that builds a model by adding trees on a sequentially iterative basis; every new tree attempts to predict the errors made by the preceding ensemble. This makes the gradient boosting regressor appropriately suited for modeling complex non-linear relationships.

Unlike Random Forest, out-of-bag predictions are not implemented for Gradient Boosting in gwlearn. If you want to evaluate your model without a separate test dataset, you can use the leave_out option (as shown in the classifier section) or consider using GWRandomForestRegressor, which provides built-in OOB scoring.

gwgb_reg = GWGradientBoostingRegressor(
    bandwidth=150,
    fixed=False,
    random_state=0,
)
gwgb_reg.fit(
    X,
    y_reg,
    geometry=gdf.representative_point(),
)

GWGradientBoostingRegressor(bandwidth=150, random_state=0)

The regressor provides local R² values showing how well the model fits in each neighborhood:

gwgb_reg.local_r2_.head()

0   -0.063802
1   -0.160997
2   -0.276301
3   -0.262151
4    0.596555
dtype: float64

gdf.set_geometry('geometry').plot(gwgb_reg.local_r2_, legend=True).set_axis_off()

Feature importance¶

Gradient boosting computes feature importances based on how much each feature reduces the loss function across all trees in the ensemble. These importances tend to be more focused on the most predictive features compared to Random Forest:

gwgb_reg.feature_importances_.head()

gdf.set_geometry('geometry').plot(gwgb_reg.feature_importances_["HC60"], legend=True).set_axis_off()

Comparison with global model¶

You can compare the local feature importances to values from the global model:

gwgb_reg.global_model.feature_importances_

array([0.06612286, 0.10260245, 0.20458955, 0.45924727, 0.10816795,
       0.05926993])

gwgb_reg.feature_importances_.mean()

HR60    0.127108
HR70    0.121543
HR80    0.222230
HR90    0.302311
HC60    0.123322
HC70    0.103485
dtype: float64