Implementation sanity check

[martinfleis]

Hi @TaylorOshan,

as promised, here's some stuff related to current state of gwlearn. At the moment, we have some classification models we wanted to play with but the idea is very much extensible. I need to refactor some bits to make it more modular but that is a detail. What I'd like to get now is some sanity check to make sure that what we're doing here is sensible. See the overview below.

An overview of the implementation of GW classification in gwlearn

• The BaseClassifier creates a libpysal Graph encoding the kernel and bandwidth - this is quite fast but it is to be noted that for large datasets and larger bandwidths, this can eat quite some memory.
  • here
• The fitting process is, in principle, a groupby over focal IDs in the graph adjacency. The fitting runs in parallel across groups (focal geometries).
  • [here]https://github.com/uscuni/gwlearn/blob/5af948905293e852d2cda4668d97904205ba5806/gwlearn/base.py#L328
• The models are currently restricted to binary classification case only. with multi-class, there's a likelihood that each model will see different classes due to their geographic distribution. All the performance metrics are then incomparable and prediction complicated (though we could inject 0 for unseen classes). In theory, we can eventually enable this but it requires deeper understanding of the whole process from the user perspective.
• Even in binary cases, we can get locations with only 0 or 1 in the whole bandwidth - we currently skip those and do not report anything in such locations. That yields outputs like here. To be precise, we skip all focals where the proportion of minority is below a set threshold. ([here]https://github.com/uscuni/gwlearn/blob/5af948905293e852d2cda4668d97904205ba5806/gwlearn/base.py#L422)
• In case of imbalanced data, we have an option to undersample majority class.
  • here
• When fitting, we use the leave-one-out approach of not including the values from focal in the training, only the neighbors. We could easily include focal if we want by assigning self weights in the Graph but that is not an option and the moment.

Performance is measured in multiple ways, which also depends on the model. Right now we have logistic regression, random forest and gradient boosted trees.

1. Focal prediction - thanks to the leave-one-out approach, we can use the model to predict class on focal as its values were not used in the training (beside spatial dependency leakage). We use this to measure score, precision, recall... (see here)
2. Out-of-bag score - For RF, we can get OOB predictions which are fetched from all individual models, concatenated and used as a single array of predictions against relevant true values. (here). Given this is pulled from individual models, we can report global metrics as well as local metrics. (here)
3. Feature importances - Each local model (if RF or GB) has an array of feature importances
4. For logistic regression, we are pulling all the arrays from individual models with modelled predictions and measure performance based on those in a similar way we do for OOB data. (here and then here). So both global (concat data) and local.
5. LR also returns all local coefficients and intercepts.

Prediction is implemented both as predict_proba and predict, though there might be some changes to how it works. Right now, we get a location, build a kernel to identify all the models within the bandwidth and their weights. Then we do the prediction using all of these models and report a weighted average. (here). Geogranos et al also combine this with a prediction from a global model, which we don't do at the moment.

@TaylorOshan do you think this makes sense or are we doing something we should not here? :D it is perfectly possible.

Also, I have not found anything on using classification models within the GW context, are you aware of anything? Thanks!

[TaylorOshan]

Hey @martinfleis , I have a few thoughts that might be useful (or perhaps not!).

The only classification work I am aware of is the current GWR implementation using a Binomial probability distribution to conduct logistic regression for binary classification, but the focus is on regression and not prediction. AFAIK, there is no MGR binomial regression, but it is under development by others and may be imminent.

This mostly seems pretty reasonable. Three things worth double checking:

1. You might need to be careful when creating the weights. At some point, there was a difference between the fault PySAL behavior and what we had set up in the MGWR module. The former I believe would truncate very small values in order to force the W to become more sparse and consume a smaller memory footprint. Though not an unreasonable strategy, it diverges from the canonical GWR behavior found not only in the MGWR module, but also other software and many of the publications. Even if quite far away, an outlier to could still become influential in this manner, and in any event, could make it tricky, too validate implementations against each other.

2. Typically the weights do include self neighbors as there is no real reason not to. In many cases, W is used to differentiate "neighbors" from "self", but since in the GW parody W is more so used for "data borrowing", excluding self neighbors means throwing away the most representative information at hand.

3. From what I can tell, it looks like you are passing X, y, and the focal weights on to a scikit class and that takes care of interacting those elements to create the local samples – is that right? I'm asking because there's a slight difference between applying W to X and y directly within the estimator, such least squares, and manually applying the weights to each X and y and then conducting model fitting as if the weighted X and y were simply exog/endog. The former allows you to use X and y as is but diagonalize the focal W before use within the matrix algebra, while latter requires that you take the square root of the focal W array before multiplying it by X and y individually before passing into local regression. I assume that since this is required for a local OLS regression to be considered correct, it's probably required for most local model fitting, but I could be wrong. And this may be happening already, but it wasn't immediately clear to me.

Hope this is somewhat helpful.

[martinfleis]

You might need to be careful when creating the weights. At some point, there was a difference between the fault PySAL behavior and what we had set up in the MGWR module. The former I believe would truncate very small values in order to force the W to become more sparse and consume a smaller memory footprint.

I am pretty certain that Graph is not doing any of it and all the values are preserved. It might have been or still is a case for W but I have not checked that.

Typically the weights do include self neighbors as there is no real reason not to. In many cases, W is used to differentiate "neighbors" from "self", but since in the GW parody W is more so used for "data borrowing", excluding self neighbors means throwing away the most representative information at hand.

I think I'll include a keyword to control this. The reasoning was coming from the application of RF within the geographically weighted context which means that if focal is included, the focal-based performance is useless as the model has seen the value and might be overfitted. That is not the case for logistic regression, so I suppose having this controllable is a sensible approach.

From what I can tell, it looks like you are passing X, y, and the focal weights on to a scikit class and that takes care of interacting those elements to create the local samples – is that right?

Yes, that is correct.

The former allows you to use X and y as is but diagonalize the focal W before use within the matrix algebra, while latter requires that you take the square root of the focal W array before multiplying it by X and y individually before passing into local regression. I assume that since this is required for a local OLS regression to be considered correct, it's probably required for most local model fitting, but I could be wrong. And this may be happening already, but it wasn't immediately clear to me.

Nah, I don't think any of this is currently happening. I would like to rely on sklearn as much as possible to avoid re-implementation of the regressors and also to allow the base classes to consume anything with sklearn compatible API. So if I am reading this correctly, we need to do the latter, what needs to happen is that we take square root of W before passing it to the sklearn estimator in here? https://github.com/uscuni/gwlearn/blob/5af948905293e852d2cda4668d97904205ba5806/gwlearn/base.py#L468-L472

So

 local_model.fit( 
     X=X, 
     y=y, 
-    sample_weight=data["_weight"], 
+    sample_weight=np.sqrt(data["_weight"]), 
 )

Or am I reading it completely wrong?

[martinfleis]

You might need to be careful when creating the weights. At some point, there was a difference between the fault PySAL behavior and what we had set up in the MGWR module. The former I believe would truncate very small values in order to force the W to become more sparse and consume a smaller memory footprint.

I am pretty certain that Graph is not doing any of it and all the values are preserved. It might have been or still is a case for W but I have not checked that.

Reading the chapter 2.2.3.1 in our MGWR book, I think I now understand what you meant by this and I suppose it is actually still the case contrary to what I wrote.

With the Gaussian kernel, computed in libpysal, we do actually clip the tails beyond bandwidth because those are not included in the Graph / W in the first place. Same thing will happen with exponential kernel. but fixing that would break the logic of implementation I have here. Completely. (see pysal/libpysal#47)