Calculating performance scores by subtracting the Shapley feature effects from predictions

Question

One of my friends has done the following: He trained an XGBoost model let's call it Model 1 and then calculated the feature effects using Shapley for the different features, let's say Feature_1, Feature_2, Feature_3. He then calculated the performance scores using the following combinations:

• Prediction 1 = Effect Feature_1 + Effect Feature_2 + Effect Feature_3
• Prediction 2 = Effect Feature_1 + Effect Feature_2
• Prediction 3 = Effect Feature_1 + Effect Feature_3

He found out that Prediction 3 is the one that it is better than all in terms of performance. For me it's wrong what he did since XGBoost is a correlation model and by changing the training data you change the correlation among the features so if you train a model using only the Feature_1 and Feature_3 then the predictions will not be equal to Prediction 3.

I was planning to prepare a dummy example for him by training a model Model 2 with only the Feature_1 and Feature_3 and show that the predictions of Model 2 are not equal to Prediction 3 but I wanted to ask the following:

• Do you have any resource that can explain it further and in a more scientifically sound way?
• And lastly but most importantly, is there a case where what he did actually holds true?

Thank you!

Answer 1

Overall, agree that using Shapley values to choose a subset of features is not the same as estimating something like model 2 with just those features.

For example, assume I've got a model like this for model 1 (using R's datasets::mtcars)

model_1 <- lm(mpg ~ am + vs + cyl, data = mtcars)

I can get Shapley values (using method in the {shapr} package) by

set.seed(11042022)

lm_shapr <- shapr::shapr(mtcars[c("vs", "am", "cyl")], model_1)

lm_explain <- shapr::explain(mtcars |> subset(select = c(vs, am, cyl)), 
                             lm_shapr, approach = "empirical", prediction_zero = mean(mtcars$mpg))

shap_preds <- lm_explain$dt

Which produces a result like:

> shap_preds |> head()
       none        am        vs        cyl
1: 20.09063  2.330967 -1.556853  0.8572712
2: 20.09063  2.330967 -1.556853  0.8572712
3: 20.09063  2.322544  2.031903  2.8468138
4: 20.09063 -1.560581  2.144306 -0.6842069
5: 20.09062 -1.563329 -1.612503 -2.4945257
6: 20.09063 -1.560581  2.144306 -0.6842069

It so happens that the highest prediction is combining am and cyl like so:

shap_preds |> transform(comb = am + cyl) |> subset(select = comb) |> cor(mtcars$mpg)
         [,1]
cor 0.8298505

which ostensibly performs better than all three variables:

> shap_preds |> transform(comb = cyl + vs + cyl) |> subset(select = comb) |> 
    cor(mtcars$mpg)
         [,1]
cor 0.8085383

A similar situation as to that described above, but with a linear model.

One complication is that Shapley values are designed to approximate the effect of a variable's impact across combinations of 'different coalitions of players' or different features being included or excluded in terms of their impact.

That means each variable's effect is an approximation and doesn't necessarily sum up to the predicted value. The actual correlation between all 3 variables, from the full model, with the outcome is:

> predict(model_1) |> cor(mtcars$mpg)
[1] 0.8729131

And, in fact, when the model is refit with the two predictors considered best from the last model, as you had suggested above, the refit model (i.e., model 2) changes in a way that does not look like model 1's Shapley values and effectively recovers the same correlation with the outcome.

> lm(mpg ~ am + cyl, data = mtcars) |> predict() |> cor(mtcars$mpg)
[1] 0.8712138

Do you have any resource that can explain it further and in a more scientifically sound way?

The Shapley value literature is fairly complex but the idea is, again, that Shapley values try to approximate what the effect the model has when a feature is omitted compared to when it is included which, for most implementations, is conditional on the model as estimated. As I noted in the comment above, these values are intended for explanation purposes to understand model predictions for complex functions like xgboost and not for selection of features.

And lastly but most importantly, is there a case where what he did actually holds true?

I believe it would in conditions where the features are completely uncorrelated with one another. I suspect that in any other situation that some form of the above would occur.