Interprétabilité et explicabilité (formalisé) des modèles prédictifs

Dans Confessiones, Saint Augustin écrivait

quid est ergo tempus? si nemo ex me quaerat, scio; si quaerenti explicare velim, nescio

que l’on traduit

Qu’est-ce-que le temps ? Si personne ne me le demande, je le sais. Si je veux l’expliquer à qui me le demande, je ne le sais plus.

Pour aller un peu plus loin (car souvent, si on nous demande d’expliquer, on a quelques idées), dans Une étude en rouge de Sir Arthur Conan Doyle, paru en 1887, on a l’échange suivant, entre Sherlock Holmes et le docteur Watson,

– Je me demande ce que cherche ce type là-bas, demandai-je,
désignant un grand individu habillé simplement qui suivait
l’autre côté de la rue, en examinant anxieusement les numéros.
Il tenait à la main une grande enveloppe bleue et, de toute
évidence, portait un message.
– Vous parlez de ce sergent d’infanterie de marine ? dit Sherlock Holmes.

puis, comme il s’avère que la personne est effectivement sergent dans la marine (tout comme un autre personnage de l’histoire, un certain Arthur Charpentier), le docteur Holmes lui demande une explication, il veut savoir comment il est arrivé à cette conclusion

« Comment diable avez-vous pu deviner cela ? demandai-je.
– Deviner quoi ? fit-il sans aménité.
– Eh bien, qu’il était un sergent de marine en retraite ?
– Je n’ai pas de temps à perdre en bagatelles ! répondit-il
avec brusquerie avant d’ajouter dans un sourire : excusez ma rudesse ! Vous avez rompu le fil de mes pensées. Mais c’est peut-être aussi bien. Ainsi donc vous ne voyiez pas que cet homme était un sergent de marine ?
Non, certainement pas !
– Décidément, l’explication de ma méthode me coûte plus que son application ! Si l’on vous demandait de prouver que deux et deux font quatre, vous seriez peut-être embarrassé ; et cependant, vous êtes sûr qu’il en est ainsi. Malgré la largeur de la rue, j’avais pu voir une grosse ancre bleue tatouée sur le dos de la main du gaillard. Cela sentait la mer. Il avait la démarche militaire et les favoris réglementaires ; c’était, à n’en pas douter, un marin. Il avait un certain air de commandement et d’importance. Rappelez-vous son port de tête et le balancement de sa canne ! En outre, son visage annonçait un homme d’âge moyen, sérieux, respectable. Tous ces détails m’ont amené à penser qu’il était sergent.
– C’est merveilleux ! m’écriai-je.
– Peuh ! L’enfance de l’art ! dit Holmes, mais d’un air qui me
parut trahir sa satisfaction devant ma surprise et mon admiration manifestes.

(pour être honnête, c’est Liu Cixin qui en parle dans Le problème à trois corps). Pour l’anecdote, c’est la première histoire du couple Holmes-Watson, qui introduit la méthode de travail de Sherlock Holmes. Pour ceux qui sont familier avec les nouvelles, cette approche narrative sera largement reprise par la suite: Sherlock Holmes énonce un fait, le docteur Watson est étonné et demande une explication, et Sherlock Holmes explique, point par point, comment il est arrivé à cette conclusion. C’est un peu cette approche qu’on tente de mettre en place quand on va construire un modèle prédictif : sur la base des données du Titanic, si on prédit que telle personne va mourir, et que telle autre va survivre, on veut comprendre pourquoi le modèle arrive à cette conclusion.

Au delà des réflexions générales, allant du “les modèles linéaires ne sont pas aussi simples à interpréter et expliquer qu’il y paraît” à “on peut expliquer sans pouvoir prévoir, et prévoir sans pouvoir expliquer“, je voulais revenir un peu sur les notions mathématiques classiques utilisées quand on parle d’explicabilité de modèles prédictifs. Sur la base des données Titanic,

library(DALEX)
data("titanic")
titanic = DALEX::titanic
idx = which((is.na(titanic$age))+
(is.na(titanic$sibsp))+
(is.na(titanic$parch))==0)
titanicNA = titanic[idx,]

on va retenir 6 modèles (pour prédire le décès ou pas du passager)

  • un modèle logistique GLM et une version GAM où l’âge est lissée
  • un arbre de classification CART et une forêt aléatoire RF
  • un modèle de boosting GBM et un support vecteur machine SVM
library(splines)
titanic_glm = glm(survived == "yes" ~ gender + age + class +
  sibsp  + parch + embarked, titanicNA, family="binomial")
titanic_gam = glm(survived == "yes" ~ gender + bs(age) + class +
  sibsp  + parch + embarked, titanicNA, family="binomial")
library("rpart")
titanic_cart = rpart(survived == "yes" ~ class + gender + age +
  sibsp + parch  + embarked, data = titanicNA)
library("gbm")
set.seed(1234)
titanic_gbm = gbm(survived == "yes" ~ class + gender + age +
  sibsp + parch  + embarked, data = titanicNA, n.trees = 15000, 
  distribution = "bernoulli")
library("randomForest")
set.seed(1234)
titanic_rf = randomForest(survived ~ class + gender + age +
  sibsp + parch + embarked, data = titanicNA)
library("e1071")
titanic_svm = svm(survived == "yes" ~ class + gender + age +
  sibsp + parch + embarked, data = titanicNA, 
  type = "C-classification", probability = TRUE)

et pour tenter d’expliquer, on verra des méthodes locales, avec deux passagers fictifs, Kate et Leonardo (comme j’avais déjà pu le faire par le passé).

newbase = data.frame(
  class = factor(c("1st","3rd"), 
        levels = c("1st", "2nd", "3rd",
        "deck crew", "engineering crew",
        "restaurant staff", "victualling crew")), 
  gender = factor(c("female","male"), levels = c("female", "male")),
  age = c(17,20),
  sibsp = c(1,0),
  parch = c(2,0),
  embarked = factor(c("Southampton","Southampton"), 
        levels = c("Belfast","Cherbourg","Queenstown","Southampton")))
rownames(newbase) = c("Winslet, Miss. Kate","DiCaprio, Mr. Leonardo")

.On va utiliser le package DALEX, écrit par Przemyslaw Biecek, et largement détaillé dans le livre co-écrit avec Tomasz Burzykowski, Explanatory model analysis,

titanic_cart_exp = DALEX::explain(model = titanic_cart, 
          data = titanicNA[, -9],
          y = titanicNA$survived == "yes", 
          label = "cart", 
          type = "classification")
titanic_glm_exp = DALEX::explain(model = titanic_glm, 
          data = titanicNA[, -9],
          y = titanicNA$survived == "yes", 
          label = "glm", 
          type = "classification")
titanic_gam_= DALEX::explain(model = titanic_gam, 
          data = titanicNA[, -9],
          y = titanicNA$survived == "yes", 
          label = "gam", 
          type = "classification")
titanic_rf_exp = DALEX::explain(model = titanic_rf, 
          data = titanicNA[, -9],
          y = titanicNA$survived == "yes",
          label = "rf") 
titanic_gbm_exp = DALEX::explain(model = titanic_gbm,
          data = titanicNA[, -9],
          y = titanicNA$survived == "yes",
          label = "gbm")
titanic_svm_exp = DALEX::explain(model = titanic_svm, 
          data = titanicNA[, -9],
          y = titanicNA$survived == "yes", 
          label = "svm")

Mais tout d’abord, on va revenir sur la distinction qui existe entre ceteris paribus et mutatis mutandis. Ceteris paribus (ou plutôt ceteris paribus sic stantibus) est la locution latine qui se traduit
par “toutes choses étant égales par ailleurs”. 
Mutatis mutandis se traduit par “ce qui devait être changé ayant été changé” ou “une fois effectuées les modifications nécessaires”. Cette distinction sera importante par la suite. Formellement, supposons que l’on considère un modèle \text{surpoids}=m(\text{poids},\text{taille}). Dans le premier cas, pour comprendre l’impact du poids sur le fait d’être, ou pas, en surpoids, autorise à regarderm(\text{poids}=x+dx,\text{taille}=y)-m(\text{poids}=x,\text{taille}=y)comme sur la figure ci-dessous

Dans le second cas, on souhaite tenir compte du fait qu’un individu de poids différent serait aussi, probablement, de taille différente, et donc on devrait regarderm(\text{poids}=x+dx,\text{taille}=y+\epsilon)-m(\text{poids}=x,\text{taille}=y)comme sur la figure ci-dessous

\epsilon tiendrait compte de la dépendance entre les deux variables.

Pour illustrer un peu plus (dans un contexte probabiliste), supposons un modèle Gaussien, autrement dit on a deux variables (X_1,X_2) (comme dans l’exemple précédant), telles que \begin{pmatrix}X_1\\X_2\end{pmatrix}\sim\mathcal{N}\left(\begin{pmatrix}\mu_1\\\mu_2\end{pmatrix},\begin{pmatrix}\sigma_1^2 & r\sigma_1\sigma_2 \\r\sigma_1\sigma_2 & \sigma_2^2\end{pmatrix}\right). Dans ce cas, notons (classiquement) (X_1^\perp,X_2^\perp) une version indépendante de (X_1,X_2). L’espérance conditionnelle vaut ici \displaystyle{\mathbb{E}_{X_1}[X_2|x_1^*]=\mu_2+\frac{r\sigma_2}{\sigma_1}(x_1^*-\mu_1)}et on notera\displaystyle{\mathbb{E}_{X_1^\perp}[X_2|x_1^*]=\mathbb{E}[X_2]=\mu_2}la version où les variables sont supposées indépendantes. Cette notation avec un indice peut surprendre, mais je vais utiliser \mathbb{E}_{X_1} pour dire qu’on prend l’espérance non pas sous \mathbb{P}, mais sous la loi \mathbb{P}_{X_2|X_1}, parfois notée \mathbb{E}^{\mathbb{P}_{X_2|X_1}}. Cette notation pourra choquer les puristes, mais elle me permettra de simplifier beaucoup de choses par la suite. Si on regarde la version statistique (ou empirique), dans le second cas (avec la version indépendante) \mathbb{E}_{X_1^\perp}[h(X_1,X_2)|x_1^*] \approx \frac{1}{n}\sum_{i=1}^n h(x_{1}^*,x_{i,2})alors que la vraie espérance conditionnelle signifie qu’on se place localement au voisinaeg de x_1^*et donc \mathbb{E}_{X_1}[h(X_1,X_2)|x_1^*] \approx \frac{1}{\|\mathcal{V}_\epsilon(x_1^*)\|}\sum_{i\in\mathcal{V}_\epsilon(x_1^*)} h(x_{1}^*,x_{i,2})\mathcal{V}_\epsilon(x_1^*) est le voisinage dont on parlait à l’instant, autrement dit\mathcal{V}_\epsilon(x_1^*)=\big\lbrace i:|x_{i,1}-x_1^*|\leq \epsilon\big\rbrace

Avant de me lancer dans le formalisme, on parlera d’approche locale quand on veut expliquer la prévision donnée par un modèle m en un point bien spécifique (qu’on notera \boldsymbol{x}^*) alors que l’approche globale consiste à se demander plutôt quelles sont les variables importantes dans le modèle m (et éventuellement, comment le modèle évolue marginalement en fonction de x_j^*).

  • Approche globale: variable importance

L’importance des variables est un premier outils pour comprendre quelles sont les variables importantes dans un modèle prédictif mFisher et al. (2019) ont proposé la définition suivante : étant donnée une fonction de perte \ell,VI_{j}=\mathbb{E}\big[\ell(Y,m(\boldsymbol{X}_{-j},X_j))\big]-\mathbb{E}\big[\ell(Y,m(\boldsymbol{X}_{-j},X_j^\perp))\big]et la version empirique est\widehat{VI}_{j}=\frac{1}{n}\sum_{i=1}^n \ell(y_i,m(\boldsymbol{x}_{i,-j},x_{i,j}))-\ell(y_i,m(\boldsymbol{x}_{i,-j},\tilde{x}_{i,j}))pour une permutation \tilde{\boldsymbol{x}}_j de {\boldsymbol{x}}_j.

vip_glm_999 = model_parts(explainer = titanic_glm_exp,
       loss_function = 1-AUC,
       B = 999,
       type = "difference"))
plot(vip_glm_999)

  • Approche locale: ICE ou ceteris paribus

Goldstein et al. (2015) ont introduit la notion de ICE (individual conditional expectation), correspondant simplement à la fonctionnelle ceteris paribusz\mapsto m_{\boldsymbol{x}^*,j}(z) = m(\boldsymbol{x}^*_{-j},z)=m(x_1^*,\cdots,x_{j-1}^*,z,x_{j+1}^*,\cdots,x_p^*) en un point \boldsymbol{x}^*\in\mathcal{X}  (avec un petit abus de notation sur les indices, puisqu’on notera abusivement (\boldsymbol{x}^*_{-j},x_j^*) quelle que soit la position j).

On peut alors regarder \delta m_{\boldsymbol{x}^*,j}(z)=m_{\boldsymbol{x}^*,j}(z) - m_{\boldsymbol{x}^*,j}(x_j^*), ou plus intéressant, La déviation moyenne absolue de la jème variable, en \boldsymbol{x}^*, est dm_{j}(\boldsymbol{x}^*) \displaystyle{dm_{j}(\boldsymbol{x}^*) =\mathbb{E}\big[|\delta m_{\boldsymbol{x}^*,j}(X_j)|\big] =\mathbb{E}\big[|m(\boldsymbol{x}_{-j},X_j) - m(\boldsymbol{x}_{-j},x_j^*)|\big]}La déviation moyenne absolue empirique de la j-ème variable, en \boldsymbol{x}^*, est\displaystyle{\widehat{dm}{j}(\boldsymbol{x}^*)=\frac{1}{n}\sum_{i=1}^n| m(\boldsymbol{x}_{-j},x_{i,j})-m(\boldsymbol{x}_{-j},x_j^*)|}.

Par exemple, juste pour Kate, et juste pour la régression logistique, le code ressemble à

cp_titanic_glm = predict_profile(explainer = titanic_glm_exp, new_observation = newbase[1,])
plot(cp_titanic_glm, variables = "age")

Ici, pour nos 6 modèles (on remplace simplement le modèle à expliquer dans le code précédent), on regarde la courbe d’évolution de l’âge, pour Kate (que se serait-il passé si Kate n’avait pas un âge x^*_j=17 mais x – toutes autres choses étant égales, par ailleurs, ceteris paribus).

ou pour la classe (que se serait-il passé si Kate n’avait pas voyagé en classe x^*_j=1 (première classe) mais dans la classe x – toutes autres choses étant égales, par ailleurs, ceteris paribus).

plot(cp_titanic_glm, variables = "class")

  • Approche locale : décomposition additive (break-point)

Pour un modèle linéaire\widehat{m}(\boldsymbol{x}^*) = \widehat{\beta}_0 + \widehat{\boldsymbol{\beta}}^\top\boldsymbol{x}^*=\widehat{\beta}_0+\sum_{j=1}^p \widehat{\beta}_j x_j^*=\overline{y} + \sum_{j=1}^p \underbrace{\widehat{\beta}_j\big(x_j^*-\overline{x}_j\big)}_{=v_j(\boldsymbol{x}^*)}
v_j(\boldsymbol{x}^*) est la contribution de la variable j dans la prédiction pour \boldsymbol{x}^*.

Plus généralement, Robnik-Šikonja and Kononenk (1997, 2003 et 2008) définissent la contribution de la j-ème variable, en \boldsymbol{x}^*, sous la formev_j(\boldsymbol{x}^*)=m(x_1^*,\cdots,,x_{j-1}^*,x^*_j,x_{j+1}^*,\cdots,x^*_p) - \mathbb{E}_{\boldsymbol{X}_{-j}^\perp}[m(x^*_1,\cdots,x^*_{j-1},X_j,x^*_{j+1},\cdots,x^*_p)]de telle sorte quem(\boldsymbol{x}^*) = \mathbb{E}_{}\big[m(\boldsymbol{X})\big] + \sum_{j=1}^pv_j(\boldsymbol{x}^*)donc, pour un modèle linéaire (avec une notation pompeuse mais qui permettra de mieux comprendre la généralisation, je pense) \displaystyle{v_j(\boldsymbol{x}^*) = \beta_j\big( x_j^* - \mathbb{E}_{\boldsymbol{X}_{-j}^\perp} [X_j]\big)} et \widehat{v}_j(\boldsymbol{x}^*) =\widehat{\beta}_j \big( x^*_j-\overline{x}_j\big).

Mais plus généralement \displaystyle{v_j(\boldsymbol{x}^*)=m(\boldsymbol{x}^*) - \mathbb{E}_{\boldsymbol{X}_{-j}^\perp}[m(\boldsymbol{x}^*_{-j},X_j))]} où on peut aussi écrire m(\boldsymbol{x}^*) sous la forme \mathbb{E}_{\boldsymbol{X}}[m(\boldsymbol{x}^*)], i.e.
v_j(\boldsymbol{x}^*)=\begin{cases}\mathbb{E}_{\boldsymbol{X}}\big[m(\boldsymbol{X})\big\vert x^*_1,\cdots,x^*_p\big] -\mathbb{E}_{\boldsymbol{X}_{-j}^\perp}\big[m(\boldsymbol{X})\big\vert x^*_1,\cdots,x^*_{j-1},x^*_{j+1},\cdots,x^*_p\big]\\\mathbb{E}_{\boldsymbol{X}}\big[m(\boldsymbol{X})\big\vert \boldsymbol{x}^*\big] - \mathbb{E}_{\boldsymbol{X}_{-j}^\perp}\big[m(\boldsymbol{X})\big\vert\boldsymbol{x}^*_{-j}\big]\end{cases}

La contribution de la j-ème variable, en \boldsymbol{x}^*, est\gamma_j^{bd}(\boldsymbol{x}^*)=v_j(\boldsymbol{x}^*)=\mathbb{E}_{\boldsymbol{X}}\big[m(\boldsymbol{X})\big\vert \boldsymbol{x}^*\big] -\mathbb{E}_{\boldsymbol{X}_{-j}^\perp}\big[m(\boldsymbol{X})\big\vert\boldsymbol{x}^*_{-j}\big]On retrouve ici clairement un approche ceteris paribus, qu’on peut illustrer sur l’exemple ci-dessous (je vais me contenter de faire des dessins). Considérons un ensemble \mathcal{X} une peu particulier puisque 0\leq x_2\leq x_1\leq 1.

On veut comprendre la prévision donnée par la régression par moindres carrés, pour un des points

On commence par se débarasser de la valeur moyenne de y

Ensuite, on va sommer, toutes choses étant égales par ailleurs, avec \displaystyle{<br /> \mathbb{E}_{X_1^\perp}\big[m(x_1^*,X_2)\big] \approx \frac{1}{n}\sum_{i=1}^n m(x_1^*,x_{2,i})}en sommant à x_1^*fixé (pourquoi pas)

mais aussi \displaystyle{<br /> \mathbb{E}_{X_2^\perp}\big[m(X_1,x_2^*)\big] \approx \frac{1}{n}\sum_{i=1}^n m(x_{1,i},x_2^*)}[\latex]en sommant à [latex]x_2^*fixé, qui est un peu plus génant sur cet exemple parce qu'on génére des pseudo-observations (x_{1,i},x_2^*) dans la zone rouge, impossible à atteindre.

Dans les deux derniers cas, la contribution sera la moyenne conditionnelle, à laquelle on retranche la moyenne globale (c'est la hauteur qu'on peut visualiser en jaune). Et on peut montrer que la moyenne, à laquelle on ajoute les deux contributions, donne la valeur prédite.

Sur les données du Titanic, pour Kate et juste la régression logistique

 bd_glm_kate = DALEX::predict_parts(explainer = titanic_glm_exp, 
new_observation = newbase[1,],
type = "break_down",
order = c("gender","class", "age", 
"parch", "sibsp", "embarked"))
plot(bd_glm_kate)

pour Kate.

Pour aller un peu plus loin, notons qu'on peut définir la contribution de la j-ème variable conditionnelle à un groupe de variables S\subset\{1,\cdots,p\}\backslash\{j\}, en \boldsymbol{x}^*, est\Delta_{j|S}(\boldsymbol{x}^*)=\begin{cases}\mathbb{E}_{\boldsymbol{X}^\perp_S,X^\perp_j}\big[m(\boldsymbol{X})\big\vert\boldsymbol{x} ^*_S, x^*_j\big] - \mathbb{E}_{\boldsymbol{X}_S^\perp}\big[m(\boldsymbol{X})\big\vert\boldsymbol{x}^*_S\big]\\\mathbb{E}_{\boldsymbol{X}^\perp_{S \cup\{j\}}}\big[ m(\boldsymbol{X})\big\vert \boldsymbol{x}^*_{S\cup\{j\}}\big] -\mathbb{E}_{\boldsymbol{X}_S^\perp}\big[m(\boldsymbol{X})\big\vert\boldsymbol{x}^*_S\big]\end{cases}de telle sorte que v_j(\boldsymbol{x}^*)=\Delta_{j|\{1,2,\cdots,p\}\backslash\{j\}}=\Delta_{j|-j}.

On peut aussi définir \Delta_{K|S}(\boldsymbol{x}^*), ou \Delta_{i,j|S}(\boldsymbol{x}^*) (ce qui permettrait d'analyser davantage les possibles interractions, mais on va passer rapidement).

  • Approche locale : la décomposition de Shapley

Pour rappel, dans un contexte (très) général, \forall S\subseteq\{1,\ldots,p\}, on dispose d'une fonction \text{val}(S), et on cherche des contributions \phi_j(\text{val}), vérifiant quelques critères

  • efficiency: \displaystyle{\sum_{j=1}^p \phi_j(\text{val}) = \text{val}(\{1,\ldots,p\})}
  • symmetry: si \text{val}\left(S\cup\{j\}\right)=\text{val}\left(S\cup\{k\}\right), \forall S\subseteq\{1,\ldots,p\}\setminus\{j,k\}, alors \phi_j=\phi_k
  • dummy: si \text{val}\left(S\cup\{j\}\right)=\text{val}\left(S\right), \forall S\subseteq\{1,\ldots,p\}, alors \phi_j=0
  • additivity: si \text{val}^{(1)} et \text{val}^{(2)} ont les decompositions respectives \phi^{(1)} and \phi^{(2)}, alors \text{val}^{(1)}+\text{val}^{(2)} a pour decomposition \phi^{(1)}+\phi^{(2)}

Shapley (1953) a montré que la seule contribution qui vérifie ces conditions est\phi_j(\text{val})=\sum_{S\subseteq\{1,\ldots,p\}\setminus\{j\}}\frac{|S|!\left(p-|S|-1\right)!}{p!}\left(\text{val}\left(S\cup\{j\}\right)-\text{val}(S)\right)ou encore\phi_{j}(\text{val})=\frac{1}{p} \sum_{S \subseteq\{1,\ldots,p\}\setminus\{j\}}\binom{p-1}{|S|}^{-1}(\text{val}(S\cup\{j\})-\text{val}(S))

On va utiliser ici \text{val}(S)=\mathbb{E}_{\boldsymbol{X}_S^\perp}\big[m(\boldsymbol{X})\big\vert \boldsymbol{x}^*_S\big]

La contribution de la j-ème variable, en \boldsymbol{x}^*, est\gamma_j^{shap}(\boldsymbol{x}^*)=\frac{1}{p} \sum_{S \subseteq\{1,\ldots,p\}\setminus\{j\}}\binom{p-1}{|S|}^{-1} \Delta_{j|S}(\boldsymbol{x}^*)

  • local accuracy: \displaystyle{\sum_{j=1}^p \gamma_j^{shap}(\boldsymbol{x}^*)=m(\boldsymbol{x}^*)-\mathbb{E}\big[m(\boldsymbol{X})\big]}
  • symmétrie: si j et k sont interchangeables, \gamma_j^{shap}(\boldsymbol{x}^*)=\gamma_k^{shap}(\boldsymbol{x}^*)
    \item dummy: si X_j ne contribue jamais, \gamma_j^{shap}(\boldsymbol{x}^*)=0
  • additivity: si m=m_1+m_2, \gamma_j^{shap}(\boldsymbol{x}^*;m)=\gamma_j^{shap}(\boldsymbol{x}^*;m_1)+\gamma_j^{shap}(\boldsymbol{x}^*;m_2)

Si p=2, \gamma_1^{shap}(\boldsymbol{x}^*)=\Delta_{1|2}(\boldsymbol{x}^*)=\gamma_1^{bd}(\boldsymbol{x}^*)

Si p\gg2, les calculs peuvent vite devenir lourds. Štrumbelj and Kononenko (2014) ont proposé une méthode par simulations. À partir de \boldsymbol{x}^* et d'un individu \boldsymbol{x}_i, on construit\tilde x_{j}=\begin{cases}x^*_j\text{ avec probabilité }1/2\\x_{i,j}\text{ avec probabilité }1/2\\\end{cases}et\begin{cases}\boldsymbol{x}^{*+}_i = (\tilde x_{1},\cdots,x_{j}^*,\cdots,\tilde x_{p})\\ \boldsymbol{x}^{*-}_i = (\tilde x_{1},\cdots,x_{i,j},\cdots,\tilde x_{p})\end{cases}et on note que \gamma_j^{shap}(\boldsymbol{x}^*)\approx m(\boldsymbol{x}^{*+}_i )-m(\boldsymbol{x}^{*-}_i) , et donc\widehat{\gamma}_j^{shap}(\boldsymbol{x}^*) = \frac{1}{s}\sum_{i\in\{1,\cdots,n\}} m(\boldsymbol{x}^{*+}_i )-m(\boldsymbol{x}^{*-}_i)
(on tire à chaque étape un individu i dans la base d'apprentissage, s fois).

shap_glm = DALEX::predict_parts(explainer = titanic_glm_exp, 
new_observation = newbase[1,],
type = "shap",
order = c("gender","class", "age", 
"parch", "sibsp", "embarked"))
plot(shap_gLm, show_boxplots=FALSE)

ou, avec des boxplots de confiance

plot(shap_gLm, show_boxplots=TRUE)

Au lieu de prendre \boldsymbol{x}^*, on peut prendre tous les individus de la base

ou regarder plus préciséments les nuages de points \{(x_{i,j},\gamma_j^{shap}(\boldsymbol{x}_i))\} (construits sur la base d'apprentissage) qui sont dit "Shapley Dependence Plots".

ou encore

Au lieu d'une vision locale (en \boldsymbol{x}^*) on peut avoir une vision globale.

Štrumbelj and Kononenko (2014), puis Lundberg and Lee (2017) ont proposé d'utiliser cette décomposition pour expliquer la contribution de chaque variable. La Shapley Feature Importance est\gamma_j^{shap} = \frac{1}{n}\sum_{i=1}^n \gamma_j^{shap}(\boldsymbol{x}_i)

  • Approche locale : LIME (Local Interpretable Model-Agnostic Explanations)

Étant donné un modèle m sur \mathcal{X}, Ribeiro, Singh and Guestrin (2016) propose de résoudre \underset{m_e\in\mathcal{E}}{\text{argmin}} \big\lbrace\ell_{\boldsymbol{x}^*}(m,m_e)\big\rbrace+\mathcal{P}(m_e)

  • \mathcal{E} est un sous-ensemble de fonctions \mathcal{X}\to\mathbb{R} "explicables'' ou \tilde{\mathcal{X}}\to\mathbb{R}, où \tilde{\mathcal{X}} est un espace plus simple (space for interpretable representation)
  • \ell_{\boldsymbol{x}^*} est une distance entre deux modèles au voisinage de \boldsymbol{x}^*
  • \mathcal{P} est une fonction de pénalité, croissante en la complexité du modèle

(je vais renvoyer ici à mes slides utilisés pour parler d'analyse d'images, il y a quelques années).Les images sont ici des individus \boldsymbol{x} en dimension 30000 (on a des images de 100\times100 pixels, en couleurs, ce qui triple la dimension. On voit clairement en haut à droite que la première étape est de simplifier en découpant (intelligemment) en 6\times.

Classiquement, on peut utiliser un modèle linéaire (avec LASSO) ou des arbres. On parle d'explicabilité par substitution (surrogate models) locale

library(localModel)
library(DALEXtra)
localModel_1_glm = predict_surrogate(explainer = titanic_glm_exp, 
         new_observation = newbase[1,],
         size = 1000,
         type = "localModel")
plot_interpretable_feature(localModel_1_glm,"age")

 

  • Approche locale : local-diagnostic plots

A partir de \boldsymbol{x}^*, on cherche les plus proches voisins. On compare ici

  •  la distribution globale des résidus
  • la distribution des résidus des voisins de \boldsymbol{x}^*

pour la mesure de similarité de Gower d_G(\boldsymbol{x}_i,\boldsymbol{x}^*)=\frac{1}{p}\sum_{j=1}^p d_j({x}_{i,j},{x}_j^*),d_j({x}_{i,j},{x}_j^*)=\begin{cases}\displaystyle{\frac{|{x}_{i,j}-{x}_j^*|}{\max\{x_j\}-\min\{x_j\}}}, \text{ si }j\text{ continue}\\\boldsymbol{1}({x}_{i,j}\neq{x}_j^*), \text{ si }j\text{ qualitative}\end{cases}

res_glm_1 = predict_diagnostics(explainer = titanic_glm_exp, 
    new_observation = newbase[1,],
    neighbors = 25)
plot(res_glm_1)

avec les résidus sur l'ensemble de la base en bas en vert, et pour les plus proches voisins de \boldsymbol{x}^* en bleu, en haut.

Mais on peut aussi regarder la courbe ceteris paribus de \boldsymbol{x}^* et celle de ses plus proches voisins, par exemple pour Kate (juste sur l'âge)

res_glm_1_age = predict_diagnostics(explainer = titanic_glm_exp, 
    new_observation = newbase[1,],
    neighbors = 25,
    variables = "age")
plot(res_glm_1_age)

  • Approche globale : Partial dependence plot

Le Partial Dependence Plot de la variable j est la fonction \mathcal{X}_j\to\mathbb{R}p_j(x_j^*)=\mathbb{E}_{X_j^\perp}\big[ m(\boldsymbol{X})|x_j^*\big]et sa version empirique est\widehat{p}_j(x_j^*)=\frac{1}{n} \sum_{i=1}^n {m}(x_j^*,\boldsymbol{x}_{i,-j})=\frac{1}{n}\sum_{i=1}^n \underbrace{{m}_{\boldsymbol{x}_{i}}(x_j^*)}_{\text{ceteris paribus}}

pdp_glm = model_profile(explainer = explain_glm, 
variables = "age")
plot(pdp_glm)

ou avec quelques courbes ceteris paribus d'observations dans la base

Et on peut regarder la moyenne non pas globale, mais sur deux sous-groupes, par exemple les hommes, et les femmes

pdp_glm = model_profile(explainer = explain_glm, 
variables = "age", groups = "gender")
plot(pdp_glm)

ou regroupés en deux classes de manière non supervisée (les sous-groupes sont ici constitués par modèles, on n'a donc aucune cohérence entre les dessins)

pdp_glm = model_profile(explainer = explain_glm, 
variables = "age", k = 3)
plot(pdp_glm)

  • Approche globale : Local dependence plot

Notons que \displaystyle{p_j(x_j^*)=\mathbb{E}_{X_j^\perp}\big[m(\boldsymbol{X})|x_j^*\big]} ne tenait pas compte du faire que (\boldsymbol{X}_{-j}|X_j)\overset{\mathcal{L}}{\neq}\boldsymbol{X}_{-j}. Au lieu de faire du ceteris paribus comme avec les PDP, on a utiliser ici une approche mutandis mutatis

Apley and Zhu (2020) propose{\ell}_j(x_j^*)=\mathbb{E}_{X_j}\big[m(\boldsymbol{X})|x_j^*\big]pour l'approche théorique, et la version empirique est alors\widehat{{\ell}}_j(x_j^*)=\frac{1}{\text{card}(V(x_j^*))}\sum_{i\in V(x_j^*)} {m}(x_j^*,\boldsymbol{x}_{i,-j})V(x_j^*)=\big\lbrace i:d(x_{i,j},x_j^*)\leq \epsilon\big\rbrace est le voisinage, ou\widehat{{\ell}}_j(x_j^*)=\frac{1}{\sum_i\omega_i(x_j^*)}\sum_{i=1}^n\omega_i(x_j^*) {m}(x_j^*,\boldsymbol{x}_{i,-j})\omega_i(x_j^*)=K_h(x_j^*-x_{i,j}) pour une version lissée par noyau.

loc_glm = model_profile(explainer = explain_glm, variables = "age", type="conditional")

  • Approche globale : Effets locaux accumulés

On peut définira_j(x_j^*)=\int_{-{\infty}}^{x_j^*}\mathbb{E}_{X_j}\left[\frac{\partial m(x_j,\boldsymbol{X}_{-j})}{\partial{x}_j}\Big\vert x_j\right]dx_j\displaystyle{\frac{\partial m(x_j,\boldsymbol{X}_{-j})}{\partial{x}_j}} décrit le changement local du modèle m du à x_j (ceteris paribus),m(x_j+dx_j,\boldsymbol{x}_{-j})-m(x_j,\boldsymbol{x}_{-j}) \approx \displaystyle{\frac{\partial m(x_j,\boldsymbol{X}_{-j})}{\partial{x}_j}}dx_jet on regarde la valeur moyenne locale.

Pour la version empirique\widehat{a}_j(x_j^*)= \alpha+ \sum_{u=1}^{k_j^*}\frac{1}{n_u}\sum_{u:x_{i,j}\in(a_{u-1},a_u]} \left[m(a_k,\boldsymbol{x}_{i,-j})-m(a_{k-1},\boldsymbol{x}_{i,-j})\right](où \alpha est une constante de normalisation car \mathbb{E}[a_j(X_j)]=0)

Apley and Zhu (2020) propose une discrétisation de X_j (partition \{(a_{k-1},a_k]\}), ou une version lissée par noyau.

acc_glm = model_profile(explainer = explain_glm, variables = "age", type="accumulated")


Leave a Reply

Your email address will not be published.

This site uses Akismet to reduce spam. Learn how your comment data is processed.