[Python] Why shap's explainer.model.predict() and model.predict don't match?

I have a machine learning model and I calculated SHAP on it using following code:

import shap

background = shap.kmeans(X_dev, k=100)
explainer = shap.TreeExplainer(model, feature_perturbation="interventional", model_output='probability', data=background.data)

shap_values = explainer.shap_values(X_val)


I observed that when I do explainer.model.predict() and model.predict() both predictions don't match.

model.predict_proba(X_val)[:10]
>>
array([[0.90563095, 0.09436908],
[0.675441 , 0.324559 ],
[0.7728198 , 0.22718018],
[0.00906086, 0.99093914],
[0.5687084 , 0.4312916 ],
[0.146478 , 0.853522 ],
[0.21917653, 0.78082347],
[0.871528 , 0.12847197],
[0.34144473, 0.65855527],
[0.25084436, 0.74915564]], dtype=float32)

explainer.model.predict(X_val)[:10]
>>>
array([0.09436912, 0.32455891, 0.22718007, 0.99093911, 0.43129163,
0.85352203, 0.7808235 , 0.12847196, 0.65855535, 0.74915555])


Question: Why is there a difference in predictions after a few decimal points?? Why aren't they exactly similar?

What I observed: I have also observed that this difference increases if you train on bigger dataset. Like if you use larger dataset to train a model, like with 1 millon rows and 100s of columns, this differece in prediction increases. The example which I showed above is created on Titanic dataset (which is very small) but if you try it on larger dataset, the difference will be massive.

What I tried:

• Tried both TPD and Interventional SHAP.
• Tried comparing raw scores.
• Tried using xgb.DMatrix to predict on the model.
• Tried making prediction using XGBClassifier API on pd.DataFrame data (not DMatrix).

Minimal reproducible code:

import pandas as pd
import numpy as np
from sklearn.model_selection import train_test_split
import xgboost as xgb
import shap

# Load Titanic dataset
data = pd.read_csv('https://raw.githubusercontent.com/datasciencedojo/datasets/master/titanic.csv')

# Preprocessing
data = data.drop(['Name', 'Ticket', 'Cabin'], axis=1)
data = pd.get_dummies(data, columns=['Sex', 'Embarked'], drop_first=True)
data = data.fillna(data.median())

X = data.drop('Survived', axis=1)
y = data['Survived']

#----------
# If you want to use large data to see big difference
#n_samples = 1000000
#n_features = 150
#X_large = pd.DataFrame(np.random.rand(n_samples, n_features), columns=[f'feature_{i}' for i in range(n_features)])
#y_large = pd.Series(np.random.randint(0, 2, size=n_samples), name='target')
#X_dev_large, X_val_large, y_dev_large, y_val_large = train_test_split(X_large, y_large, test_size=0.2, random_state=42)
#----------

X_dev, X_val, y_dev, y_val = train_test_split(X, y, test_size=0.2, random_state=42)
dtrain = xgb.DMatrix(X_dev, label=y_dev)
dval = xgb.DMatrix(X_val, label=y_val)

params = {
'objective': 'binary:logistic',
'eval_metric': 'logloss',
'max_depth': 4,
'eta': 0.1
}
model = xgb.train(params, dtrain, num_boost_round=100)
model.save_model('xgboost_model.json')

from xgboost import XGBClassifier

loaded_model = XGBClassifier()
loaded_model.load_model('xgboost_model.json')

background = shap.kmeans(X_dev, k=100)
explainer = shap.TreeExplainer(loaded_model, feature_perturbation="interventional", model_output='probability', data=background.data)
shap_values = explainer.shap_values(X_val)

print(loaded_model.predict_proba(X_val)[:10])
print(explainer.model.predict(X_val)[:10])


Other details:

• 
  Shap version: 0.47.2 (I also tried 0.37 something)
  
  
• 
  Python version: 3.9.21
  
  
• 
  XGBoost version: 2.1.1 (tried 1.7 as well)
  

Continue reading...