Interprétation des modèles PDF

Interprétation des modèles Introduction Subsection 1 Introduction Denis Larocque © Techniques d’exploitation de données 545 / 954  Interprétation des modèles Introduction Introduction Il...

Interprétation des modèles Introduction Subsection 1 Introduction Denis Larocque © Techniques d’exploitation de données 545 / 954  Interprétation des modèles Introduction Introduction Il existe de nombreuses méthodes pour faire de l’apprentisage supervisé. Il y a les modèles de régression classiques (linéaire, logistique). Dans les chapitres précédents, nous avons vu les arbres et les combinaisons d’arbres (forêt aléatoire et boosting), ainsi que les réseaux de neurones (RN). Un RN, une forêt ou un boosting d’arbres sont des modèles difficilement interprétables. Nous avons vu qu’il est possible d’avoir une idée de l’importance des variables (VIMP), mais cela reste incomplet. On voudrait aussi en savoir plus sur la forme des liens entre les variables explicatives et la variable cible. Nous avons vu qu’un arbre seul peut parfois être utile pour l’interprétation, mais que cela dépend fortement de sa taille. Denis Larocque © Techniques d’exploitation de données 546 / 954  Interprétation des modèles Introduction Un modèle de régression classique est généralement interprétable. Mais même ici il y a des nuances. Si le modèle utilise des transformations complexes des variables explicatives et que nous voulons voir les effets des variables originales, avant les transformations, ce n’est pas toujours simple. De plus, si le modèle contient des interactions d’ordres élevés, alors l’interprétation peut aussi devenir plus complexe. Mais s’il est basé sur les variables qu’on veut interpréter et qu’il y a seulement des interactions d’ordre 2, alors c’est généralement facile de l’interpréter. Denis Larocque © Techniques d’exploitation de données 547 / 954  Interprétation des modèles Introduction Notez que ce n’est pas le fait d’être paramétrique qui rend un modèle facile à interpréter. Tout comme une régression linéaire, un RN est un modèle paramétrique mais qui est impossible à interpréter directement, à cause qu’il contient un grand nombre de paramètres et de transformations complexes des variables. Denis Larocque © Techniques d’exploitation de données 548 / 954  Interprétation des modèles Introduction Interpréter un modèle et comprendre pourquoi il fournit telle ou telle prédiction est important pour plusieurs raisons. En voici quelques-unes. Une entreprise veut généralement bien comprendre tout le processus décisionnel avant de déployer un nouveau système. Cela inclut comprendre les modèles qui fournissent des prédictions, s’ils font partie du processus, Avoir de bonnes prédictions est très bien, mais on veut souvent pouvoir comprendre quelles variables influencent notre variable cible et de quelles manières. Par exemple, dans un problème de désabonnement, c’est bien d’avoir une estimation de la probabilité de désabonnement, mais on voudra généralement agir en amont avant même que le client ne soit rendu à considérer se désabonner. Pour ce faire, il faut bien comprendre les facteurs qui influencent le désabonnement. Denis Larocque © Techniques d’exploitation de données 549 / 954  Interprétation des modèles Introduction Dans certaines industries, notamment en finance, les instances règlementaires imposent souvent des balises qui demandent plus de transparence, d’où la nécessité d’avoir des méthodes pour interpréter les modèles. Il y a également des questions éthiques, notammnent dû aux biais que peuvent avoir certains modèles, et investiguer ces aspects nécessite de pouvoir interpréter les modèles. Finalement, dans certains cas, les personnes affectées par les décisions impliquant des modèles doivent leur faire confiance et les accepter, ce qui fait encore entrer en jeu l’aspect interprétation. Denis Larocque © Techniques d’exploitation de données 550 / 954  Interprétation des modèles Introduction Un domaine de recherche très actif actuellement se nomme Interpretable AI ou Explainable AI (XAI) et cherche à développer des méthodes pour faciliter l’interprétation de modèles complexes. Nous allons introduire certains principes de base, et présenter quelques méthodes, parmi les plus connues. Mais comme il y a beaucoup de développements dans ce domaine, notez que le panorama risque de changer rapidement au cours des prochaines années. Denis Larocque © Techniques d’exploitation de données 551 / 954  Interprétation des modèles Introduction Il y a deux dimensions à considérer: 1) interprétation globale ou locale et 2) méthode spécialisées à un modèle ou agnostique Interprétation globale: On veut interpréter le modèle de manière générale. Quelles sont les variables importantes? Comment les variables affectent-elles le 𝑌 ?. Interprétation locale: On veut ici comprendre ce qui influence une prédiction spécifique pour une observation. Méthode spécialisée à un modèle: Méthode qui est développée pour un type de modèle (par exemple une forêt) et qui utilise ses caratéristiques pour guider l’interprétation. Méthode agnostique: Méthode qui peut s’appliquer aux prédictions de n’importe quel modèle. Typiquement, elle prend les prédictions du modèle et les variables explicatives en entrées et fournit l’interprétation en sortie, sans se soucier d’où viennent les prédictions. Denis Larocque © Techniques d’exploitation de données 552 / 954  Interprétation des modèles Introduction Un aspect important doit être compris dès le départ. Un modèle est une approximation de la réalité, notamment, de la vraie fonction qui relie les variables explicatives à la variable cible. Les méthodes que nous allons voir servent à interpréter le modèle. Si le modèle n’est pas adéquat, son interprétation ne sera d’aucune utilité ou, pire encore, pourra donner de fausses indications. Les exemples simples qui suivent illustrent cela. Denis Larocque © Techniques d’exploitation de données 553 / 954  Interprétation des modèles Introduction Dans ce premier exemple, il y une seule variable explicative 𝑥 et une variable cible 𝑦 continue. Le graphe de la page suivante montre les données. Denis Larocque © Techniques d’exploitation de données 554 / 954  Interprétation des modèles Introduction 1.0 0.8 0.6 y 0.4 0.2 0.0 0.0 0.2 0.4 0.6 0.8 1.0 x Denis Larocque © Techniques d’exploitation de données 555 / 954  Interprétation des modèles Introduction Nous allons utiliser un modèle qui est directement interprétable, une régression linéaire simple, c’est-à-dire: 𝐸[𝑦|𝑥] = 𝛽0 + 𝛽1 𝑥. Les estimations des paramètres sont 𝛽0̂ = −0, 016 et 𝛽1̂ = 1, 01. Le graphe de la page suivante montre le modèle ajusté. Son interprétation est simple, lorsque 𝑥 augmente de 1, 𝑦 augmente de 1,01, en moyenne. De plus le modèle est adéquat alors cette interprétation est utile et reflète bien la réalité. Denis Larocque © Techniques d’exploitation de données 556 / 954  Interprétation des modèles Introduction ## (Intercept) x ## -0.01589033 1.0 0.8 0.6 1.01145502 y 0.4 0.2 0.0 0.0 0.2 0.4 0.6 0.8 1.0 Denis Larocque © Techniques d’exploitation x de données 557 / 954  Interprétation des modèles Introduction Voici un deuxième exemple, toujours avec une seule variable explicative 𝑥 et une variable cible 𝑦 continue. Le graphe de la page suivante montre les données. Denis Larocque © Techniques d’exploitation de données 558 / 954  Interprétation des modèles Introduction 0.05 0.00 −0.05 −0.10 y −0.15 −0.20 −0.25 0.0 0.2 0.4 0.6 0.8 1.0 x Denis Larocque © Techniques d’exploitation de données 559 / 954  Interprétation des modèles Introduction Nous allons encore utiliser une régression linéaire simple, c’est-à-dire: 𝐸[𝑦|𝑥] = 𝛽0 + 𝛽1 𝑥. Évidemment ce modèle n’est pas adéquat car il y a un lien non-linéaire entre 𝑥 et 𝑦. Les estimations des paramètres sont 𝛽0̂ = −0, 076 et 𝛽1̂ = −0, 019. Le graphe de la page suivante montre le modèle ajusté, qui confirme que le modèle n’est pas adéquat. L’interprétation du modèle est toujours simple. Lorsque 𝑥 augmente de 1, 𝑦 diminue de 0,019 , en moyenne. Mais comme le modèle n’est pas adéquat, cette interprétation n’est pas utile. Pire encore, si on ne se rend pas compte que le modèle n’est pas adéquat, on va tirer des conclusions fausses de cette interprétation. De plus, le paramètre de la variable 𝑥 n’est pas significativement différent de 0, ce qui donne une indication faussée de son importance. Denis Larocque © Techniques d’exploitation de données 560 / 954  Interprétation des modèles Introduction ## (Intercept) x ## -0.07646241 -0.01942221 0.05 0.00 −0.05 −0.10 y −0.15 −0.20 −0.25 0.0 0.2 0.4 0.6 0.8 1.0 Denis Larocque © Techniques d’exploitation x de données 561 / 954  Interprétation des modèles Introduction Évidemment, ici, on pourrait essayer un modèle de régression avec un terme quadratique pour 𝑥, c’est-à-dire: 𝐸[𝑦|𝑥] = 𝛽0 + 𝛽1 𝑥 + 𝛽2 𝑥2. Le graphe de la page suivante montre le modèle ajusté, qui est à présent adéquat. Donc, on pourrait l’interpréter. Comme il est non-linéaire, le plus simple est d’utiliser des graphiques, comme celui de la page suivante, pour interpréter l’effet de 𝑥 sur le 𝑦. La plupart des méthodes que nous verrons utilisent justement des graphiques pour l’interprétation. Denis Larocque © Techniques d’exploitation de données 562 / 954  Interprétation des modèles Introduction 0.05 0.00 −0.05 −0.10 y −0.15 −0.20 −0.25 0.0 0.2 0.4 0.6 0.8 1.0 x Denis Larocque © Techniques d’exploitation de données 563 / 954  Interprétation des modèles Introduction Ces exemples simples servent de mise en garde. Les méthodes que nous allons voir servent à interpréter le modèle. *Et pas la réalité Moins le modèle est adéquat, moins son interprétation est utile. Pire encore, l’interpétation d’un modèle inadéquat pourrait donner de fausses indications quant à la forme des liens entre les variables explicatives et la variable cible et quant à l’importance des variables explicatives. Denis Larocque © Techniques d’exploitation de données 564 / 954  Interprétation des modèles Introduction Idéalement, nous aurions un modèle performant pour les prédictions et interprétable. Par exemple, supposons qu’une régression linéaire ou logistique soit aussi bon, en termes de prédictions, que des modèles de type boites noires (forêt, RN). De plus, supposons que ce modèle est basé sur les variables que nous voulons interpréter, et non pas sur des transformations trop obscures, et ne contient pas d’interactions trop complexes. Nous aurions alors le meilleur des deux mondes. Ainsi, lorsque c’est possible, il faut essayer le plus possible d’être dans cette situation. Denis Larocque © Techniques d’exploitation de données 565 / 954  Interprétation des modèles Introduction Certains chercheurs sont très critiques par rapport au courant Interpretable AI. Voir Rudin (2019) pour une discussion à ce sujet. Je crois qu’il faut effectivement tenter le plus possible d’être dans la situation idéale d’avoir un bon modèle prédictif qui est aussi interprétable. Mais si ce n’est pas possible, alors les méthodes d’interprétation sont utiles si elles sont bien utilisées et si on comprend bien leurs limites. Denis Larocque © Techniques d’exploitation de données 566 / 954  Interprétation des modèles Introduction Dans ce qui suit, nous allons utiliser des exemples simples, avec des données simulées, de telle sorte que nous connaissons la réalité. Cela nous permettra de voir comment interpréter les résultats des méthodes d’interprétations. Mais aussi de voir comment elles réagissent dans certaines situations ce qui nous permettra de fournir certaines mises en garde. Denis Larocque © Techniques d’exploitation de données 567 / 954  Interprétation des modèles Implémentations de méthodes d’interprétation Subsection 2 Implémentations de méthodes d’interprétation Denis Larocque © Techniques d’exploitation de données 568 / 954  Interprétation des modèles Implémentations de méthodes d’interprétation Implémentations de méthodes d’interprétation Il y a de nombreux packages en R. Une comparaison de 27 d’entre eux est faite dans Maksymiuk, Gosiewska, and Biecek (2020) et la liste continue de s’allonger. Deux, parmi les plus complets, sont DALEX et iml mais certaines méthodes sont seulement disponibles dans des packages plus spécialisés. Il y a aussi de nombreuses librairies en Python, telles ELI5, SHAP, Alibi. Denis Larocque © Techniques d’exploitation de données 569 / 954  Interprétation des modèles Mesures d’importance des variables basées sur les permutations Subsection 3 Mesures d’importance des variables basées sur les permutations Denis Larocque © Techniques d’exploitation de données 570 / 954  Interprétation des modèles Mesures d’importance des variables basées sur les permutations Mesures d’importance des variables basées sur les permutations Le but ici est d’attribuer un score d’importance à une variable explicative ou, plus généralement, un groupe de variables explicatives, qui reflète son rôle pour prédire la variable cible. Il y a certaines subtilités dans cet exercice, surtout dues au fait que les variables travaillent habituellement conjointement pour expliquer le 𝑌 et qu’elles sont aussi typiquement corrélées (ou plus généralement dépendantes). Comme d’habitude, nous avons 𝑝 variables explicatives 𝑥1 , … , 𝑥𝑝. Denis Larocque © Techniques d’exploitation de données 571 / 954  Interprétation des modèles Mesures d’importance des variables basées sur les permutations Une première idée naturelle qu’on pourrait avoir est de mesurer l’impact de retirer une variable (ou un groupe de variables). C’est-à-dire, on pourrait ajuster le modèle avec toutes les variables et calculer sa performance (une mesure d’erreur), possiblement sur un échantillon de validation. Ensuite, on pourrait enlever une variable, réajuster le modèle sans elle, et calculer à nouveau sa performance. La différence entre les deux erreurs représente la perte de performance due au fait de ne pas avoir inclus cette variable dans le modèle. On pourrait refaire cela pour chaque variable et ainsi attribuer des scores de performances basés sur ces pertes de performance. Plus l’erreur augmente en retirant la variable, plus elle est importante. Denis Larocque © Techniques d’exploitation de données 572 / 954  Interprétation des modèles Mesures d’importance des variables basées sur les permutations En régression linéaire, le test d’hypothèse pour une variable est essentiellement basé sur ce principe. En effet, la statistique de test est une fonction de 𝑆𝑆𝐸(réduit) − 𝑆𝑆𝐸(complet) où, 𝑆𝑆𝐸(complet) est la somme des carrés des erreurs du modèle avec toutes les variables et 𝑆𝑆𝐸(réduit) est la somme des carrés des erreurs du modèle en ayant retiré la variable à tester. Cette version pour faire le test utilise les données d’entrainement mais on pourrait utiliser des données de validation pour avoir un score d’importance. Denis Larocque © Techniques d’exploitation de données 573 / 954  Interprétation des modèles Mesures d’importance des variables basées sur les permutations On pourrait donc considérer faire cela pour n’importe quel modèle. Mais c’est très rarement fait en pratique. Une des raisons est que ça peut être très lourd computationnellement. En effet, le modèle doit être réajusté chaque fois qu’on retire une variable. Avec un grand 𝑝 et/ou des modèles lourds à ajuster (par exemple un RN), cela devient rapidement impraticable. Mais il y a une raison plus subtile liée à ce que l’on veut réellement mesurer. Nous y reviendrons plus loin avec un exemple. Denis Larocque © Techniques d’exploitation de données 574 / 954  Interprétation des modèles Mesures d’importance des variables basées sur les permutations C’est pourquoi, l’approche standard est plutôt basée sur l’idée de perturber les données, que nous avons vue plus tôt avec les forêts aléatoires (les VIMP). La méthode de perturbation la plus courante consiste à permuter les observations pour une variable. Ce type de méthodes permet de juger globalement de l’importance relative des variables dans un modèle donné. C’est Breiman (2001) qui a introduit ce concept dans l’article original sur les RF. Mais cette idée peut-être utilisée avec n’importe quel modèle. La page suivante présente un algorithme générique pour calculer un VIMP. Denis Larocque © Techniques d’exploitation de données 575 / 954  Interprétation des modèles Mesures d’importance des variables basées sur les permutations Algorithme générique pour VIMP basé sur les permutations 0 On a un modèle ajusté avec des données d’entrainement qui contient 𝑝 variables explicatives. On a un critère de performance 𝐿 (plus il est petit, mieux c’est). On calcule 𝐿 avec des données non-modifiées, appelons cela 𝐿𝑜. 1 Pour 𝑗 = 1, … , 𝑝, répéter les étapes 2 à 4. 2 Permuter les valeurs de la variable 𝑥𝑗 dans l’échantillon. 3 Calculer 𝐿 pour cet échantillon modifié, appelons cela 𝐿𝑗. 4 Calculer la différence 𝐷𝑗 = 𝐿𝑗 − 𝐿𝑜 ou le ratio 𝑅𝑗 = 𝐿𝑗 /𝐿𝑜 du critère calculé avec l’échantillon permuté et celui calculé avec l’échantillon original. 5 Retourner les valeurs 𝐷𝑗 ou 𝑅𝑗 ordonnées de la plus grande à la plus petite. Denis Larocque © Techniques d’exploitation de données 576 / 954  Interprétation des modèles Mesures d’importance des variables basées sur les permutations Voici quelques explications supplémentaires. Cet algorithme est très général. Il peut s’appliquer avec n’importe quel modèle (pas seulement une RF) et type de variable cible. Le modèle est ajusté une seule fois au début. Le critère de performance peut-être n’importe quoi de raisonnable. Pour un 𝑌 continu, on peut prendre le MSE, le MAE, le MAPE. Pour un 𝑌 binaire, on peut prendre le taux de mauvaise classification, le AUC, le score de Brier. À la fin plus la différence 𝐷𝑗 (ou le ratio 𝑅𝑗 ) est grand, plus la variable est importante. Afin de diminuer la variabilité (car les permutations sont aléatoires), on peut répéter le processus plusieurs fois et faire la moyenne des résultats. Denis Larocque © Techniques d’exploitation de données 577 / 954  Interprétation des modèles Mesures d’importance des variables basées sur les permutations On peut utiliser l’échantillon d’entrainement ou un échantillon de validation pour calculer les 𝐿. On pourrait aussi le faire par validation-croisée. C’est plus simple de le faire avec les données d’entrainement et ce sera généralement suffisant. Mais si c’est possible, il est préférable de le faire en calculant les 𝐿 hors-échantillon (out of sample). La méthode originale de Breiman avec les RF utilise les données OOB pour calculer 𝐿 et est dont faite hors-échantillon. Il faut bien vérifier ce que font les différentes implémentations des VIMP. Lorsqu’on utilise les données d’entraînement, on espère que les biais vont s’annuler. Denis Larocque © Techniques d’exploitation de données 578 / 954  Interprétation des modèles Mesures d’importance des variables basées sur les permutations La méthode originale de Breiman avec les RF que nous avons décrite plus tôt est légèrement différente. La différence de performance entre l’échantillon original et l’échantillon permuté est calculé séparément pour chaque arbre (en prenant la prédiction de l’arbre). Ensuite, on fait la moyenne sur tous les arbres pour avoir la mesure finale. Avec l’approche générique qu’on vient de présenter, la mesure de performance est calculée sur toute la RF d’un seul coup (en prenant la prédiction de la RF). Comme on utilise la prédiction de la RF, il est préférable de calculer les VIMP de cette façon (c’est-à-dire avec l’approche générique), mais ça ne devrait pas faire trop de différence la plupart du temps. Le package randomForestSRC offre les 2 options. Denis Larocque © Techniques d’exploitation de données 579 / 954  Interprétation des modèles Mesures d’importance des variables basées sur les permutations Lorsque les variables explicatives sont corrélées (ce qui est le cas presque toujours), il n’est pas évident de séparer les effets des variables et d’attribuer des importances individuelles. C’est pourquoi il peut être intéressant de calculer l’importance d’un groupe de variables. Cela est aussi pertinent si on utilise le hot one encoding pour les variables catégorielles. En effet, un VIMP individuel nous donne seulement l’importance d’une indicatrice à la fois, alors qu’on est plutôt intéressé à avoir l’importance de la variable globalement. Il y a peu d’implémentations de VIMP par groupes, randomForestSRC étant une rare exception qui permet de le faire dans le cadre des RF. Denis Larocque © Techniques d’exploitation de données 580 / 954  Interprétation des modèles Exemple de régression avec les données de jeux en ligne (suite) Subsection 4 Exemple de régression avec les données de jeux en ligne (suite) Denis Larocque © Techniques d’exploitation de données 581 / 954  Interprétation des modèles Exemple de régression avec les données de jeux en ligne (suite) Exemple de régression avec les données de jeux en ligne (suite) Nous poursuivons l’exemple de prédiction du montant dépensé avec le jeu mobile. Nous avions déjà obtenu les VIMP classiques avec ranger. Il s’agit de la version où la différence de performance entre l’échantillon original et l’échantillon permuté est calculé séparément pour chaque arbre, tel qu’expliqué plus haut. Revoici le graphe de ces VIMP. Denis Larocque © Techniques d’exploitation de données 582 / 954  Interprétation des modèles Exemple de régression avec les données de jeux en ligne (suite) totpurchases trendpurchase country numpurchases numsessions numelements totplaytime totscore acquis skill1 trendsession numdays numlives opsys skill2 difflevel 0 20 40 60 80 Importance Denis Larocque © Techniques d’exploitation de données 583 / 954  Interprétation des modèles Exemple de régression avec les données de jeux en ligne (suite) Nous allons à présent calculer les VIMP pour le boosting avec lightgbm avec le package iml. Nous prenons le MSE comme critère et la différence (au lieu du ratio) afin de comparer plus facilement avec les résultats de ranger. Nous faisons une seule répétitions mais en pratique, il serait préférable de faire quelques répétitions. Denis Larocque © Techniques d’exploitation de données 584 / 954  Interprétation des modèles Exemple de régression avec les données de jeux en ligne (suite) # R # Pour faire une graphe des VIMP avec iml pour le modèle lightgbm library(iml) # On doit d'abord définir la fonction qui sera utilisée pour calculer les # prédictions iml_predict = function(X.model, newdata) {predict(X.model, as.matrix(newdata))} # Rappels: # modelmobile135 est l'objet qui contient le modèle que nous avons ajusté plus tôt. # lgbmobiletrain1 est le data frame avec les variables explicatives. # mobiletrain$y est le vecteur avec la variable cible exp_iml = Predictor$new(modelmobile135, data = lgbmobiletrain1, y = mobiletrain$y, predict.function = iml_predict) # Calcul des VIMP. On prend une seule répétitions pour que ça aille plus # vite. La fonction de perte L est le MSE. set.seed(45343) imp = FeatureImp$new(exp_iml, loss = "mse",compare="difference",n.repetitions =1) # Graphe des VIMP plot(imp) Denis Larocque © Techniques d’exploitation de données 585 / 954  Interprétation des modèles Exemple de régression avec les données de jeux en ligne (suite) trendpurchase totpurchases country numpurchases numsessions acquis totplaytime trendsession skill2 opsys skill1 totscore numlives numelements numdays difflevel 0 25 50 75 100 Feature Importance (loss: mse) Denis Larocque © Techniques d’exploitation de données 586 / 954  Interprétation des modèles Exemple de régression avec les données de jeux en ligne (suite) Les résultats sont très similaires. Les deux variables trendpurchase et totpurchases se démarquent clairement des autres. Les trois suivantes sont aussi les mêmes qu’avec ranger. Il s’agit de country, numpurchases et numsessions. Plus le score d’importance diminue, plus l’ordre devient variable car il s’agit essentiellement de bruit. Les VIMP sont des mesures intéressantes, mais elles ne nous donnent pas d’idées quant à la nature des liens entre les variables explicatives et la variable cible. Les outils qui suivent nous permettront de faire cela. Mais avant voici un exemple afin de bien comprendre ce que mesurent les VIMP. Denis Larocque © Techniques d’exploitation de données 587 / 954  Mesures d’importance des variables basées sur les permutations Interprétation des modèles (suite) Subsection 5 Mesures d’importance des variables basées sur les permutations (suite) Denis Larocque © Techniques d’exploitation de données 588 / 954  Mesures d’importance des variables basées sur les permutations Interprétation des modèles (suite) Mesures d’importance des variables basées sur les permutations (suite) Nous avons mentionné plus tôt qu’une autre approche consiste à quantifier l’impact de retirer une variable. Nous allons faire un exemple pour illustrer que cette approche peut donner des résultats différents. Dans cet exemple simple, il y a 3 variables explicatives 𝑥1 , 𝑥2 , 𝑥3. Les 2 premières sont directement liées au 𝑦 mais pas 𝑥3. Mais elle l’est indirectement car elle est fortement corrélée avec 𝑥1. 𝑥1 et 𝑥2 sont indépendantes et uniformément distribuées dans l’intervalle [0, 1]. 𝑥3 est également normalisée pour être dans [0, 1]. Le vrai modèle qui génère le 𝑦 est 𝑦 = 200 ∗ 𝑥1 + 100 ∗ 𝑥2 + 𝜖, où 𝜖 est une erreur de loi normale de moyenne 0. Voir le code pour plus de détails. Denis Larocque © Techniques d’exploitation de données 589 / 954  Mesures d’importance des variables basées sur les permutations Interprétation des modèles (suite) Comme 𝑥1 et 𝑥2 sont sur la même échelle, on peut imaginer que 𝑥1 est plus “importante” que 𝑥2 car son paramètre est 2 fois plus grand (200 contre 100). Voici d’abord le code qui génère les données, la matrice de corrélation et les coefficients du modèle de régression linéaire avec les 3 variables. Denis Larocque © Techniques d’exploitation de données 590 / 954  Mesures d’importance des variables basées sur les permutations Interprétation des modèles (suite) # R set.seed(233) n=100 x1=runif(n) x2=runif(n) x3=x1/10+rnorm(n,sd=.01) x3=(x3-min(x3))/(max(x3)-min(x3)) y=200*x1 + 100*x2 + 100*rnorm(n,sd=.05) datall=data.frame(y,x1,x2,x3) cor(cbind(y,x1,x2,x3)) ## y x1 x2 x3 ## y 1.0000000 0.87486265 0.45013987 0.82830342 ## x1 0.8748626 1.00000000 -0.03094521 0.95516480 ## x2 0.4501399 -0.03094521 1.00000000 -0.05307288 ## x3 0.8283034 0.95516480 -0.05307288 1.00000000 modall=lm(y~.,data=datall) Denis Larocque © Techniques d’exploitation de données 591 / 954  Mesures d’importance des variables basées sur les permutations Interprétation des modèles (suite) summary(modall) ## ## Call: ## lm(formula = y ~., data = datall) ## ## Residuals: ## Min 1Q Median 3Q Max ## -12.6558 -3.5658 -0.1701 3.4840 12.7756 ## ## Coefficients: ## Estimate Std. Error t value Pr(>|t|) ## (Intercept) -3.132 1.597 -1.961 0.0528. ## x1 192.964 6.342 30.428 𝑐. On va donc travailler avec des modèles qui sont en mesure de fournir des probabilités estimées, et non pas seulement des prédictions binaires (0-1). En principe, c’est le cas pour n’importe quel modèle mais peut-être que certaines implémentations fournissent seulement le 𝑌 ̂ = 0 ou 1. Denis Larocque © Techniques d’exploitation de données 735 / 954  Classification avec classes déséquilibrées et avec coûts Trois types de méthodes Voici trois types de méthodes. 1 On ajuste le seuil 𝑐 à la fin afin d’optimiser un critère adéquat, de préférence minimiser le CME (ou maximiser le GME) avec une matrice de coût (gain) judicieusement choisie. 2 On utilise une méthode d’échantillonnage, par exemple sous-échantillonner les ”0” ou sur-échantillonner les ”1” au début, avant d’entrainer les modèles. 3 On modifie l’algorithme du modèle de prédiction à l’interne afin de tenir compte des classes déséquilibrées. Ces trois types de méthodes ne sont pas exclusifs. On peut les utiliser conjointement. Nous allons donner plus de détails dans ce qui suit. Denis Larocque © Techniques d’exploitation de données 736 / 954  Classification avec classes déséquilibrées et avec coûts Trois types de méthodes Ajustement du seuil Débutons avec l’ajustement du seuil pour une matrice de coût. Comme on l’a vu, on peut supposer que 𝑐00 = 𝑐11 = 0. Dans ce cas, on a vu que la règle optimale théorique consiste à classer une observation à 1 si 𝑐10 (1 − 𝑝(𝑥)) < 𝑐01 𝑝(𝑥) où 𝑝(𝑥) = 𝑃 (𝑌 = 1|𝑋 = 𝑥). Une simple manipulation algébrique montre que cela est équivalent à classer une observation à 1 si 𝑐10 𝑝(𝑥) >. 𝑐01 + 𝑐10 Le seuil optimal théorique pour transformer une probabilité en décision binaire est donc 𝑐 = 𝑐10 /(𝑐01 + 𝑐10 ). Denis Larocque © Techniques d’exploitation de données 737 / 954 Classification avec classes déséquilibrées et avec coûts Trois types de méthodes Si le coût des 2 types d’erreur est le même, c’est-à-dire 𝑐01 = 𝑐10 , alors le seuil optimal est 𝑐 = 0, 5. Si par exemple 𝑐10 = 1 et 𝑐01 = 9, c’est-à-dire, il est 9 fois plus coûteux de rater un 1 que de rater un 0, alors le seuil optimal est 𝑐 = 0, 1. On voit bien que le seuil est déplacé vers le bas et va faire en sorte qu’on va classer plus d’observations à 1, quitte à se tromper plus souvent, pour être certain de ne pas rater trop de 1. Denis Larocque © Techniques d’exploitation de données 738 / 954 Classification avec classes déséquilibrées et avec coûts Trois types de méthodes L’approche la plus simple consiste donc à utiliser directement le seuil optimal théorique pour notre matrice de coût pour convertir les probabilités estimées afin d’obtenir les prédictions 0-1. Par contre, il est possible que les probabilités soient biaisées. Une méthode qui est préférable dans ce cas est d’optimiser le seuil avec les données; Sheng and Ling (2006). C’est-à-dire, on estime le CME pour une grille de valeurs pour le seuil 𝑐 et on choisit la valeur qui minimise le CME. Mais il faut le faire correctement. Une manière simple serait d’estimer le CME avec un échantillon de validation. Denis Larocque © Techniques d’exploitation de données 739 / 954 Classification avec classes déséquilibrées et avec coûts Trois types de méthodes Mais avec des classes déséquilibrées, il arrive souvent qu’il n’y ait pas assez d’observations (de la classe rare) pour faire cela. Il est alors préférable de le faire par validation-croisée. L’idée est d’obtenir les probabilités estimées pour toutes les observations de l’échantillon d’entrainement dans un schéma de validation-croisée. Plus précisément, on divise les données en 𝐾 groupes. Pour 𝑘 = 1, … , 𝐾 , on enlève le groupe 𝑘 et on entraine le modèle avec le reste des données. Toutes les étapes de l’entrainement doivent être faites ici. Il est donc possible qu’il y ait une autre boucle de validation-croisée d’imbriquée. On obtient ensuite les probabilités estimées pour le groupe 𝑘. Denis Larocque © Techniques d’exploitation de données 740 / 954 Classification avec classes déséquilibrées et avec coûts Trois types de méthodes À la fin, on a une probabilité estimée pour toutes les observations de l’échantillon d’entrainement. On peut alors estimer le CME pour chaque seuil dans la grille et trouver le seuil optimal empirique. Denis Larocque © Techniques d’exploitation de données 741 / 954 Classification avec classes déséquilibrées et avec coûts Trois types de méthodes Notez que si on utilise une méthode d’ensemble où un nouvel échantillon est prélevé à chaque itération, comme une forêt aléatoire, on peut s’en tirer sans faire de validation-croisée. En effet, prenons l’exemple d’une forêt. Il est possible d’utiliser les probabilités estimées OOB (out of bag) pour l’échantillon d’entrainement pour optimiser le seuil, au lieu de le faire par validation-croisée. Il faut juste se rappeler qu’il est parfois souhaitable d’augmenter le nombre d’arbres dans la forêt car les prédictions OOB sont basées sur environ 1/3 des arbres de la forêt. Denis Larocque © Techniques d’exploitation de données 742 / 954  Classification avec classes déséquilibrées et avec coûts Trois types de méthodes Échantillonnage Le deuxième type de méthode est l’échantillonnage. L’idée est de tenter d’aider le modèle à apprendre en rééquilibrant les classe 0 et 1. Il existe un très grand nombre de variantes mais il y a essentiellement deux approches: 1) sous-échantillonner la classe majoritaire 0 ; 2) sur-échantillonner la classe minoritaire 1. Les deux peuvent être faites en même temps. Il faut bien comprendre qu’on ne crée pas de nouvelles informations en faisant cela. On espère simplement que le modèle aura plus de facilité à discriminer les 0 des 1 si on rééquilibre les classes. Denis Larocque © Techniques d’exploitation de données 743 / 954 Classification avec classes déséquilibrées et avec coûts Trois types de méthodes Lorsque le modèle est déjà basé sur une méthode de rééchantillonnage, comme par exemple une forêt aléatoire qui utilise des échantillons bootstrap ou du sous-échantillonage pour construire chaque arbre, alors cette méthode peut facilement s’y intégrer. En effet, on peut alors simplement ajuster les poids d’échantillonnage. Par exemple, avec une forêt, au lieu de préléver des échantillons bootstrap où chaque observation a la même probabilité d’être choisie à chaque tirage, on peut ajuster ces probabilités afin de rééquilibrer les échantillons. Il y a plusieurs façons de le faire. Denis Larocque © Techniques d’exploitation de données 744 / 954 Classification avec classes déséquilibrées et avec coûts Trois types de méthodes L’une d’elles se nomme balanced random forest; Chen et al. (2004). L’idée est très simple. L’échantillon pour construire un arbre de la forêt est obtenu de la manière suivante. On prélève un échantillon bootstrap de la classe minoritaire. On prélève un échantillon de même taille, avec remise, de la classe majoritaire. Ainsi, s’il y a 𝑛1 observations de la classe minoritaire, l’échantillon pour construire l’arbre contient 2𝑛1 observations. Il s’agit donc d’une variante de l’idée de sous-échantillonner la classe majoritaire, mais en le faisant par stratification. Nous utiliserons cette approche plus loin car elle est directement implémentée dans le package R randomForestSRC. Denis Larocque © Techniques d’exploitation de données 745 / 954 Classification avec classes déséquilibrées et avec coûts Trois types de méthodes Il y a des approches plus complexes qui vont “créer” des nouvelles observations au lieu de faire un échantillonnage. La méthode la plus connue est possiblement SMOTE (Synthetic Minority Oversampling TEchnique). Il s’agit d’une méthode de sur-échantillonnage de la classe minoritaire. Plusieurs variantes sont disponibles, mais voici l’idée de base. Denis Larocque © Techniques d’exploitation de données 746 / 954 Classification avec classes déséquilibrées et avec coûts Trois types de méthodes On veut créer des “nouvelles” observations de la classe minoritaire afin de rééquilibrer les données. Le nombre d’observations à créer peut être vu comme un hyper-paramètre. Pour créer un “nouveau” 1, voici comment faire. Pour simplifier supposons que les variables explicatives sont continues. On sélectionne un 1 au hasard. Supposons que c’est l’observation 𝑗 avec variables explicatives 𝑥𝑗. On trouve ses 𝑘 voisins les plus proches dans l’espace des variables explicatives (ce qui suppose que les variables ont été standardisées au préalable) parmi les observations avec 𝑌 = 1. 𝑘 est donc un autre hyper-paramètre. Denis Larocque © Techniques d’exploitation de données 747 / 954 Classification avec classes déséquilibrées et avec coûts Trois types de méthodes On choisit un de ces 𝑘 voisins au hasard. Supposons que c’est l’observation 𝑣 avec variables explicatives 𝑥𝑣. On tire un nombre 𝛼 au hasard entre 0 et 1. La nouvelle observation est une moyenne pondérée de ces 2 observations, 𝛼𝑥𝑗 + (1 − 𝛼)𝑥𝑣. On répète le processus jusqu’à ce qu’on ait créé le nombre d’observations voulus. On entraine le modèle avec les données originales plus les nouvelles observations. Denis Larocque © Techniques d’exploitation de données 748 / 954 Classification avec classes déséquilibrées et avec coûts Trois types de méthodes Il faut être prudent si on utilise de l’échantillonnage pour rééquilibrer les classes de quelque manière que ce soit. Si on considère que l’échantillon d’entrainement original est représentatif de la population, et aussi de la population des observations futures pour lesquelles nous voudrons des prédictions, alors les probabilités estimées produites par un modèle qui a utilisé un rééquilibrage des classes ne sont pas valides pour estimer les probabilités de cette polulation. Cela n’est pas un probème si on optimise le seuil comme on l’a vu plus tôt (par validation-croisée) pour convertir ces probabilités en prédictions binaires. Dans ce cas, seulement l’ordre des probabilités est important et il est toujours valide. Denis Larocque © Techniques d’exploitation de données 749 / 954 Classification avec classes déséquilibrées et avec coûts Trois types de méthodes Mais si on veut une estimation de la probabilité elle-même pour la population cible, il faut corriger les probabilités. Voici comment faire en général. Soit 𝑝𝑣 la proportion de 1 dans la population cible pour laquelle nous voulons obtenir des estimations de la probabilité. Soit 𝑝𝑢 la proportion de 1 dans l’échantillon utilisé pour entrainer le modèle. Il s’agit de la proportion après le rééquilibrage et non pas la proportion dans l’échantillon d’entrainement original. Soit 𝑝𝑢̂ (𝑥) l’estimation de 𝑃 (𝑌 = 1|𝑋 = 𝑥) provenant du modèle entrainé avec les données rééquilibrées. Soit 𝑝𝑣̂ (𝑥) l’estimation de 𝑃 (𝑌 = 1|𝑋 = 𝑥) pour la population cible que nous voulons calculer. Denis Larocque © Techniques d’exploitation de données 750 / 954 Classification avec classes déséquilibrées et avec coûts Trois types de méthodes Pour convertir 𝑝𝑢̂ (𝑥) en 𝑝𝑣̂ (𝑥), il suffit de faire la transformation suivante. (𝑝𝑣 /𝑝𝑢 )𝑝𝑢̂ (𝑥) 𝑝𝑣̂ (𝑥) =. ((1 − 𝑝𝑣 )/(1 − 𝑝𝑢 )(1 − 𝑝𝑢̂ (𝑥)) + (𝑝𝑣 /𝑝𝑢 )𝑝𝑢̂ (𝑥) Mais, encore une fois, cela n’est pas nécessaire si on veut seulement obtenir une prédiction 0-1. Dans ce cas, optimiser le seuil par validation croisée avec les 𝑝𝑢̂ (𝑥) est suffisant. La transformation est monotone et ne change pas l’ordre relatif lorsqu’on l’applique à un ensemble de probabilités. Denis Larocque © Techniques d’exploitation de données 751 / 954 Classification avec classes déséquilibrées et avec coûts Trois types de méthodes Voici un exemple concret. Supposons qu’il y a 10% de 1 dans la population d’intérêt où des estimations de la probabilités seront calculées. Supposons que le modèle a été entrainé en rééquilibrant, par sous ou sur-échantillonnage, sur un échantillon qui contient 50% de 1. Alors 𝑝𝑣 = 0, 1 et 𝑝𝑢 = 0, 5. Si le modèle produit une probabilité de 𝑝𝑢̂ (𝑥) = 0, 4, alors la conversion de cette probabilité pour la population cible est (0, 1/0, 5)0, 4 = 0, 069. (0, 9/0, 5)(0, 6) + (0, 1/0, 5)0, 4 Denis Larocque © Techniques d’exploitation de données 752 / 954  Classification avec classes déséquilibrées et avec coûts Trois types de méthodes Modification de l’algorithme à l’interne Les deux premiers types de méthodes agissent avant ou après l’entrainement du modèle. L’ajustement du seuil agit à la fin. L’échantillonnage agit avant. Le troisième type de méthode modifie le modèle à l’interne pour tenter d’améliorer sa performance dans le cas de classes déséquilibrées. Habituellement, cela est fait en changeant le critère de performance ou la fonction de perte. Notez que si on intègre l’échantillonnage à une méthode qui fait déjà du rééchantillonnage (par exemple une forêt), la distinction est moins claire. On peut voir cela comme une modification à l’interne et/ou comme de l’échantillonnage. Ce troisième type de méthode est très général et peut prendre différentes formes. Denis Larocque © Techniques d’exploitation de données 753 / 954 Classification avec classes déséquilibrées et avec coûts Trois types de méthodes Par exemple, avec une forêt, on pourrait utiliser un critère de split qui serait mieux adapté aux classes déséquilibrées que le critère habituel, qui est Gini. Une possibilité est d’utiliser le AUC comme critère de split. Nous verrons une illustration de cela plus loin car le package randomForestSRC implémente ce critère. S’il est possible d’utiliser un critère de split basé sur la perte moyenne pondérée, cela peut être une autre stratégie. En fait, il est possible de modifier les probabilités à priori pour arriver à cette fin. Avec un RN ou le boosting, utiliser une fonction de perte plus adaptée peut aussi être bénéfique. Encore ici, cela pourrait être une perte moyenne pondérée. Denis Larocque © Techniques d’exploitation de données 754 / 954 Classification avec classes déséquilibrées et avec coûts Quelques approches pratiques Subsection 4 Quelques approches pratiques Denis Larocque © Techniques d’exploitation de données 755 / 954  Classification avec classes déséquilibrées et avec coûts Quelques approches pratiques Quelques approches pratiques Une première approche concrète extrêmement simple consiste à entrainer les modèles comme d’habitude, en ne se préoccupant pas du fait d’avoir des classes déséquilibrées. On va se préoccuper des classes déséquilibrées seulement à la fin pour calculer les prédictions binaires. Pour ce faire, on fixe une matrice de coût pertinente et on utilise le seuil théorique optimal selon cette matrice pour convertir les probabilités estimées en prédictions binaires. Une approche tout juste un peu plus complexe consiste à optimiser le seuil par validation-croisée, au lieu d’utiliser le seuil théorique optimal. Même si elles sont très simples, ces approches peuvent souvent bien performer. Denis Larocque © Techniques d’exploitation de données 756 / 954 Classification avec classes déséquilibrées et avec coûts Quelques approches pratiques La raison est que plusieurs modèles seront en mesure de bien séparer les “0” des “1,” même avec des classes déséquilibrées. Les probabilités estimées pour la classe minoritaire “1” auront tendance à être plus élevées que celles de la classe majoritaire “0,” mais elles seront quand mêmes faibles. Par exemple, en moyenne le modèle pourrait fournir une probabilité estimée de 0,01 pour un vrai “0” et de 0,02 pour un vrai “1.” Dans ce cas, le modèle fait bien son travail. Mais il faut ensuite trouver le bon seuil pour convertir les probabilités estimées en 0-1. De plus, comme la détection des “1” est généralement plus importante, l’utilisation d’une matrice de coût (gain) est naturel. Denis Larocque © Techniques d’exploitation de données 757 / 954 Classification avec classes déséquilibrées et avec coûts Quelques approches pratiques Pour aider à comprendre cela, prenons l’exemple d’un arbre, mais l’idée peut être étendue à un ensemble d’arbres. Même s’ils sont relativement peu nombreux, les “1” vont quand même se retrouver dans des noeuds terminaux. Ils seront peut-être être “noyés” dans les nombreux “0” qui s’y trouvent mais ils seront quand même là. Donc, la proportion de “1” (qui est la probabilité estimée), dans un noeud terminal qui contient au moins un “1” sera peut-être très faible, mais elle sera quand même plus grande que dans un noeud qui n’en contient pas du tout. Donc en optimisant le seuil à la fin, selon une matrice de coût, on pourrait séparer ces petites probabilités des autres encore plus petites et avoir un modèle performant. Denis Larocque © Techniques d’exploitation de données 758 / 954 Classification avec classes déséquilibrées et avec coûts Quelques approches pratiques Ces approches sont très simples car on ne fait rien de différent pour l’entrainement des modèles. Une autre approche assez simple consiste à ajouter une étape d’échantillonnage au tout début; on va toujours optimiser le seuil à la fin. S’il s’agit d’un modèle qui n’utilise habituellement pas de rééchantillonnage, par exemple un RN, alors on lui donne seulement un ensemble d’entrainement différent et le reste se fait comme d’habitude. Comme il y a un élément aléatoire dans l’échantillonnage, si c’est possible, on peut répéter cela à quelques reprises afin de diminuer la variance. Les probabilités estimées pourraient être la moyenne des probabilités estimées sur les répétitions. Il s’agit donc en quelque sorte d’une méthode d’ensemble. Denis Larocque © Techniques d’exploitation de données 759 / 954 Classification avec classes déséquilibrées et avec coûts Quelques approches pratiques S’il s’agit d’un modèle qui utilise le rééchantillonnage à la base, comme une forêt, alors intégrer le rééquilibrage des classes par échantillonnage peut être fait naturellement. Nous avons déjà présenté une variante de cela plus tôt, balanced random forest. Dans ce cas particulier, comme chaque arbre est construit sur un échantillon différent, on atténue déjà la partie aléatoire de l’échantillonnage à même le modèle. Denis Larocque © Techniques d’exploitation de données 760 / 954 Classification avec classes déséquilibrées et avec coûts Quelques approches pratiques Finalement, une approche un peu moins accessible est d’utiliser un algorithme qui utilise un critère de performance ou une fonction de perte à l’interne qui est plus adaptée aux classes déséquilibrées. Cela dépend de la disponibilité des implémentations. Mais lorsqu’elles sont disponibles, ces approches ne sont pas vraiment très difficiles à utiliser. Denis Larocque © Techniques d’exploitation de données 761 / 954 Classification avec classes déséquilibrées et avec coûts Exemple modifié avec les données churn Subsection 5 Exemple modifié avec les données churn Denis Larocque © Techniques d’exploitation de données 762 / 954  Classification avec classes déséquilibrées et avec coûts Exemple modifié avec les données churn Exemple modifié avec les données churn Dans l’exemple de désabonnement que nous utilisons depuis le début du cours, nous avons environ 26% de “1,” c’est-à-dire de désabonnement, dans l’échantillon. Ici, nous allons utiliser les mêmes données mais en les modifiant afin que la proportion de “1” soit beaucoup plus faible, pour qu’on soit clairement dans une situation de classes déséquilibrées. On va garder tous les “0.” On traite donc cela comme un nouvel exemple. Le code qui suit prépare les données que nous allons utiliser, à partir des données existantes. Denis Larocque © Techniques d’exploitation de données 763 / 954  Classification avec classes déséquilibrées et avec coûts Exemple modifié avec les données churn # R set.seed(76424) # churntrain et churntest contiennent les données originales # de l'exemple de churn. churnrare=rbind(churntrain,churntest) # Indicatrice temporaire pour l'échantillon train churnrare$train=c(rep(1,nrow(churntrain)),rep(0,nrow(churntest))) # On garde seulement un petit nombre de 1 churnrare=churnrare[as.logical( rbinom(nrow(churnrare),1,.2*(churnrare$Churn.Value==1)+(churnrare$Churn.Value==0))),] Denis Larocque © Techniques d’exploitation de données 764 / 954  Classification avec classes déséquilibrées et avec coûts Exemple modifié avec les données churn # Convertir les variables "character" en "factor" pour randomForestSRC sapply(churnrare, class)[1:5] ## Referred.a.Friend Number.of.Referrals Tenure.in.Months Offer ## "character" "integer" "integer" "character" ## Phone.Service ## "character" churnrare = as.data.frame(unclass(churnrare),stringsAsFactors = TRUE) sapply(churnrare, class)[1:5] ## Referred.a.Friend Number.of.Referrals Tenure.in.Months Offer ## "factor" "integer" "integer" "factor" ## Phone.Service ## "factor" Denis Larocque © Techniques d’exploitation de données 765 / 954  Classification avec classes déséquilibrées et avec coûts Exemple modifié avec les données churn # Création des données d'entrainement et de test churnraretrain = churnrare[churnrare$train==1,] churnraretest = churnrare[churnrare$train==0,] # On enlève l'indicatrice de train churnraretrain$train=NULL churnraretest$train=NULL # Nombre d'observations, de 1 et proportion de 1. dim(churnraretrain) ## 2375 39 table(churnraretrain$Churn.Value) ## ## 0 1 ## 2217 158 mean(churnraretrain$Churn.Value) ## 0.06652632 dim(churnraretest) ## 3184 39 table(churnraretest$Churn.Value) ## ## 0 1 ## 2989 195 mean(churnraretest$Churn.Value) ## 0.06124372 Denis Larocque © Techniques d’exploitation de données 766 / 954 Classification avec classes déséquilibrées et avec coûts Exemple modifié avec les données churn Les échantillons d’entrainement et de tests comptent 2375 et 3184 observations respectivement. Il y a 158 “1” (désabonnement) dans l’échantillon d’entrainement ce qui correspond à 6,7% des observations. Il y a 195 “1” dans l’échantillon de test ce qui correspond à 6,1% des observations. Voici d’abord une fonction pour calculer le coût (gain) moyen. Denis Larocque © Techniques d’exploitation de données 767 / 954  Classification avec classes déséquilibrées et avec coûts Exemple modifié avec les données churn # R # Fonction pour calculer le coût ou gain moyen coutgain=function(pred01,y,gcmat=cbind(c(0,1),c(1,0))) { ## INPUT: # pred01 = vecteur des prédictions (0-1) # y = vecteur des vraies valeurs de Y (0-1) # gcmat = matrice de gain ou de coût (2X2) # (on veut maximiser le gain moyen ou minimiser le coût moyen) # (1,1) = coût (gain) si pred=0 et Y=0 # (1,2) = coût (gain) si pred=0 et Y=1 # (2,1) = coût (gain) si pred=1 et Y=0 # (2,2) = coût (gain) si pred=1 et Y=1 ## OUTPUT: # le coût (ou gain) moyen mean(gcmat[1,1]*(pred01==0)*(y==0)+gcmat[1,2]*(pred01==0)*(y==1)+ gcmat[2,1]*(pred01==1)*(y==0)+gcmat[2,2]*(pred01==1)*(y==1)) } Denis Larocque © Techniques d’exploitation de données 768 / 954 Classification avec classes déséquilibrées et avec coûts Exemple modifié avec les données churn La fonction suivante nous servira à optimiser le seuil afin de minimiser le coût moyen ou maximiser le gain moyen. Habituellement, les probabilités estimées qu’on lui fournit doivent être obtenues hors-échantillon, tel qu’expliqué plus haut. Cela peut être fait avec un échantillon de validation, par validation-croisée ou avec les estimations OOB. Denis Larocque © Techniques d’exploitation de données 769 / 954  Classification avec classes déséquilibrées et avec coûts Exemple modifié avec les données churn # R # fonction pour optimiser le seuil selon une matrice de coût (ou gain) seuiloptim=function(pred,y,gcmat=cbind(c(0,1),c(1,0)),cout=TRUE,cutp=seq(0,1,.005),graphe=FALSE) { ## INPUT: # pred = vecteur des probabilités estimées. Idéalement provenant d'un échantillon # de validation, d'une validation-croisée ou par OOB. # y = vecteur des vraies valeurs de Y (0-1) # gcmat = matrice de gain ou de coût (2X2) # (on veut maximiser le gain moyen ou minimiser le coût moyen) # (1,1) = coût (gain) si pred01=0 et Y=0 # (1,2) = coût (gain) si pred01=0 et Y=1 # (2,1) = coût (gain) si pred01=1 et Y=0 # (2,2) = coût (gain) si pred01=1 et Y=1 # cout = TRUE si gcmat est une matrice de coût et # FALSE si c'est une matrice de gain # cutp = vecteur des seuils à essayer # graphe = TRUE si on veut un graphe des résultats ## OUTPUT # liste avec 2 éléments: # 1) data frame avec les seuils et estimations des gains (coûts) moyens # 2) le seuil avec le gain moyen maximum (ou le coût moyen minimum) # et le gain moyen (coût moyen) associé. Denis Larocque © Techniques d’exploitation de données 770 / 954  Classification avec classes déséquilibrées et avec coûts Exemple modifié avec les données churn nc=length(cutp) res=rep(0,nc) for(i in 1:nc) { predi=as.numeric(pred>cutp[i]) res[i]=mean(gcmat[1,1]*(predi==0)*(y==0)+gcmat[1,2]*(predi==0)*(y==1)+ gcmat[2,1]*(predi==1)*(y==0)+gcmat[2,2]*(predi==1)*(y==1)) } if(cout==TRUE) { if(graphe==TRUE){plot(cutp,res,type="l",xlab="seuil",ylab="coût moyen espéré")} out=list(NULL,NULL) out[]=data.frame("seuil"=cutp,"cout"=res) out[]=out[][which.min(res),] } else { if(graphe==TRUE){plot(cutp,res,type="l",xlab="seuil",ylab="gain moyen espéré")} out=list(NULL,NULL) out[]=data.frame("seuil"=cutp,"gain"=res) out[]=out[][which.max(res),] } out } Denis Larocque © Techniques d’exploitation de données 771 / 954 Classification avec classes déséquilibrées et avec coûts Exemple modifié avec les données churn Comme nous avons vu plus tôt, il y a un très grand nombre d’approches que nous pourrions utiliser. Nous allons illustrer certaines d’entre elles avec randomForestSRC qui est le package le plus complet pour les forêts aléatoires en R. En utilisant toujours le même algorithme, ce sera plus facile de comparer l’impact des différentes approches en faisant varier un seul élément à la fois. Notre critère de performance sera le CME avec une matrice de coût qu’on va spécifier nous-même. Denis Larocque © Techniques d’exploitation de données 772 / 954  Classification avec classes déséquilibrées et avec coûts Exemple modifié avec les données churn Nous allons utiliser les 6 variantes suivantes: 1 Une RF de base avec le critère de split de Gini en utilisant le seuil optimal théorique pour convertir les probabilités. 2 Une RF de base avec le critère de split de Gini en utilisant le seuil optimisé empiriquement avec les prédictions OOB pour convertir les probabilités. 3 Une RF avec le critère de split AUC en utilisant le seuil optimal théorique pour convertir les probabilités. Denis Larocque © Techniques d’exploitation de données 773 / 954  Classification avec classes déséquilibrées et avec coûts Exemple modifié avec les données churn 4 Une RF avec le critère de split AUC en utilisant le seuil optimisé empiriquement avec les prédictions OOB pour convertir les probabilités. 5 Une balanced random forest (BRF) avec le critère de split de Gini en utilisant le seuil optimisé empiriquement avec les prédictions OOB pour convertir les probabilités. 6 Une balanced random forest (BRF) avec le critère de split AUC en utilisant le seuil optimisé empiriquement avec les prédictions OOB pour convertir les probabilités. Denis Larocque © Techniques d’exploitation de données 774 / 954 Classification avec classes déséquilibrées et avec coûts Exemple modifié avec les données churn La méthode 1 est la plus simple. On construit une RF comme d’habitude et on fait la conversion des probabilités estimées avec le seuil optimal théorique obtenu de la matrice de coût. La méthode 2 utilise la même forêt, mais on trouve ensuite le seuil optimal empiriquement. Pour ce faire, on va utiliser la fonction seuiloptim vue plus haut avec les probabilités OOB pour l’échantillon d’entrainement. Ainsi, on a besoin d’entrainer une seule forêt. La méthode 3 est un exemple d’adaptation à l’interne de l’algorithme de base. Au lieu d’utiliser le critère de split habituel (Gini), on va utiliser le critère AUC. Comme le AUC est essentiellement un critère qui évalue si on est en mesure de bien ordonner les observations, il est potentiellement plus efficace avec des données déséquilibrées. Denis Larocque © Techniques d’exploitation de données 775 / 954 Classification avec classes déséquilibrées et avec coûts Exemple modifié avec les données churn La méthode 4 utilise la même forêt que la méthode 3 mais on trouve ensuite le seuil optimal empiriquement, comme pour la méthode 2. La méthode 5 est un exemple d’utilisation de l’échantillonnage. Tel qu’expliqué plus tôt une BRF utilise un échantillonnage bootstrap stratifié afin d’avoir autant de 1 que de 0. Chaque arbre est construit avec un tel échantillon différent. Dans cette version, on utilise le critère de split de Gini. Comme les probabilités estimées ne correspondent plus à la population originales, on va optimiser empiriquement le seuil à utiliser pour la conversion avec les prédictions OOB. La méthode 6 est encore une BRF mais cette fois on utilise le critère de split AUC. Cette méthode constitue en quelque sorte le cocktail complet car elle utilise l’échantillonnage (bootstrap stratifié), une modification à l’interne (AUC), et un ajustement empirique du seuil à la fin. Denis Larocque © Techniques d’exploitation de données 776 / 954 Classification avec classes déséquilibrées et avec coûts Exemple modifié avec les données churn Notre critère de performance pour évaluer les méthodes, mais aussi pour optimiser le seuil empiriquement est le CME avec la matrice de coût suivante: vérité 0 1 prédiction 0 0 20 1 1 0 Ainsi, il est 20 fois plus coûteux de rater un client qui va se désabonner que de rater un client qui va rester. Mais afin d’illustrer les autres critères que nous avons présentés, nous allons également les calculer. La fonction manycrit présentée à la page suivante servira à cette fin. Denis Larocque © Techniques d’exploitation de données 777 / 954  Classification avec classes déséquilibrées et avec coûts Exemple modifié avec les données churn # R # Fonction pour calculer plusieurs critères de performance manycrit=function(pred,y,seuil,gcmat=cbind(c(0,1),c(1,0)),nom="modele1") { ## INPUT: # pred = vecteur des probabilités estimées. # y = vecteur des vraies valeurs de Y (0-1) # seuil = seuil pour convertir les probabilités # (si pred>seuil, alors pred01=1) # gcmat = matrice de gain ou de coût (2X2) # (on veut maximiser le gain moyen ou minimiser le coût moyen) # (1,1) = coût (gain) si pred=0 et Y=0 # (1,2) = coût (gain) si pred=0 et Y=1 # (2,1) = coût (gain) si pred=1 et Y=0 # (2,2) = coût (gain) si pred=1 et Y=1 # nom = nom du modèle ## OUTPUT: # data frame avec # coutmoyen = cout moyen (ou gain moyen) selon la matrice de cout (gain) # tmc = taux de mauvaise classification # auc = aire sous la courbe ROC # sens = sensitivité # spec = spécificité # pred = precision # bac = balanced accuracy # gmean = G-mean # F1 = F1 # gmeas = G-measure # seuil = seuil utilisé Denis Larocque © Techniques d’exploitation de données 778 / 954  Classification avec classes déséquilibrées et avec coûts Exemple modifié avec les données churn library(ROCR) if(is.factor(y)){y=as.numeric(y)-1} pred01=as.numeric(pred>seuil) procr=ROCR::prediction(pred,y) auc=performance(procr,"auc")@y.values[] tmc=mean((y-pred01)^2) sens=sum((pred01==1)*(y==1))/sum(y==1) spec=sum((pred01==0)*(y==0))/sum(y==0) prec=sum((pred01==1)*(y==1))/sum(pred01==1) bac=(sens+spec)/2 gmean=sqrt(sens*spec) F1=2*(sens*prec)/(sens+prec) gmeas=sqrt(sens*prec) coutmoyen=coutgain(pred01=as.numeric(pred>seuil),y=y,gcmat=gcmat) out=data.frame("coutmoyen"=coutmoyen,"tmc"=tmc,"auc"=auc,"sens"=sens, "spec"=spec,"prec"=prec,"bac"=bac,"gmean"=gmean,"F1"=F1,"gmeas"=gmeas, "seuil"=seuil) row.names(out)=nom out } Denis Larocque © Techniques d’exploitation de données 779 / 954 Classification avec classes déséquilibrées et avec coûts Exemple modifié avec les données churn Voici une dernière fonction, fitrare, afin d’automatiser l’entrainement des modèles car plusieurs parties se répètent. Elle nous servira à entrainer les 6 variantes présentées plus haut et à calculer les critères de performance avec la fonction manycrit. Denis Larocque © Techniques d’exploitation de données 780 / 954  Classification avec classes déséquilibrées et avec coûts Exemple modifié avec les données churn # R fitrare=function(modele,ntree,gcmat,nom,setseed) { set.seed(setseed) if(modele==1){mod=rfsrc(Churn.Value~.,data=churnraretrain,ntree=ntree)} if(modele==2){mod=rfsrc(Churn.Value~.,data=churnraretrain,ntree=ntree,splitrule="auc")} if(modele==3){mod=imbalanced(Churn.Value~.,data=churnraretrain,ntree=ntree,method="brf")} if(modele==4){mod=imbalanced(Churn.Value~.,data=churnraretrain,ntree=ntree,method="brf",splitrule="auc")} seuilth=gcmat[2,1]/(gcmat[1,2]+gcmat[2,1]) predmod=predict(mod,newdata=churnraretest)$predicted[,2] seuilo=seuiloptim(mod$predicted.oob[,2],churnraretrain$Churn.Value,gcmat=gcmat, cout=TRUE,cutp=seq(0,1,.001),graphe=TRUE) out=manycrit(predmod,churnraretest$Churn.Value,seuilo[]$seuil,gcmat,paste(nom,"seuil optimisé",sep=" ")) if(modele 0, ce critère pénalise les splits déséquilibrés. Denis Larocque © Techniques d’exploitation de données 870 / 954 Modélisation incrémentale et effets de traitements individuels Estimation du CATE Une autre idée proposée par Rzepakowski and Jaroszewicz (2012) est de maximiser le critère de split suivant 𝑛𝐺 2 𝑛𝐷 2 𝜏̂ + 𝜏̂ , 𝑛 𝐺 𝑛 𝐷 où 𝑛 = 𝑛𝐺 + 𝑛𝐷 est le nombre total d’observations dans les deux noeuds. Il s’agit donc d’une moyenne pondérée des carrés des effets de traitement estimés à gauche et à droite. L’idée ici est d’utiliser un critère de split qui va chercher les plus grands effets de traitement possibles. Denis Larocque © Techniques d’exploitation de données 871 / 954 Modélisation incrémentale et effets de traitements individuels Estimation du CATE Nous avons brièvement mentionné dans le chapitre sur les arbres qu’une forêt aléatoire peut être utilisée pour estimer à peu près n’importe quelle quantité d’intérêt. La méthode generalized random foretst, ou GRF, proposée dans Athey, Tibshirani, and Wager (2019) est une méthode très générale qui permet de faire cela. Elle est implémentée dans le package R grf. Elle peut entre autre servir à estimer le CATE. Denis Larocque © Techniques d’exploitation de données 872 / 954 Modélisation incrémentale et effets de traitements individuels Estimation du CATE Sans entrer dans les détails, dans ce contexte, cette approche est basée sur des arbres construits avec le critère de split suivant: √ 𝑛𝐺 𝑛𝐷 |𝜏𝐺 ̂ − 𝜏𝐷 ̂ |. Il s’agit donc d’un critère similaire à ceux que l’on vient de présenter. Denis Larocque © Techniques d’exploitation de données 873 / 954 Modélisation incrémentale et effets de traitements individuels Estimation du CATE Pour être précis, l’approche GRF utilise une approximation de ce critère afin d’améliorer la vitesse de calcul. De plus, par défaut, la méthode utilise le concept d’honnêteté (honesty), nécessaire pour que les résultats théoriques dérivés dans l’article soient valides. Essentiellement ce concept demande que la détermination des splits et l’estimation finale de la quantité d’intérêt dans les noeuds terminaux ne soient pas faits avec les mêmes observations. Une façon simple de faire cela est de diviser l’échantillon en deux dès le départ. La première partie sert à construire les arbres (trouver les splits) et la deuxième sert à estimer la quantité d’intérêt dans les noeuds terminaux. Cela est fait par défaut avec le package grf. Denis Larocque © Techniques d’exploitation de données 874 / 954 Modélisation incrémentale et effets de traitements individuels Estimation du CATE 𝑅-learner Une autre approche récente proposée par Nie and Wager (2021) est un peu similaire aux approches utilisant des transformations vues plus tôt mais elle utilise plutôt une fonction de perte specialisée. Seulement une description sommaire est donnée ici. Soit 𝑚(𝑋) = 𝐸[𝑌 |𝑋], la moyenne conditionnelle de 𝑌 sachant 𝑋. Notez que l’on ne conditionne pas sur le traitement ici. Il s’avère que 𝑌 − 𝑚(𝑋) = (𝑊 − 𝑝𝑊 (𝑋))𝜏 (𝑋) + 𝜖 où 𝜖 est un terme d’erreur de moyenne 0. Denis Larocque © Techniques d’exploitation de données 875 / 954 Modélisation incrémentale et effets de traitements individuels Estimation du CATE Si on connaissait les fonctions 𝑚 et 𝑝𝑊 , on pourrait estimer la fonction 𝜏 en trouvant la fonction qui minimise la fonction de perte suivante 𝑛 2 ∑ (𝑌𝑖 − 𝑚(𝑋𝑖 ) − (𝑊𝑖 − 𝑝𝑊 (𝑋𝑖 ))𝜏 (𝑋𝑖 )) + Λ(𝜏 ) 𝑖=1 où Λ(𝜏 ) est un terme de régularisation qui contrôle la complexité de la fonction 𝜏. Denis Larocque © Techniques d’exploitation de données 876 / 954 Modélisation incrémentale et effets de traitements individuels Estimation du CATE L’idée est alors de procéder en deux étapes. Étape 1: on estime les fonctions 𝑚 et 𝑝𝑊 par validation-croisée. Étape 2: on utilise ces estimations dans la fonction de perte de la page précédente afin de l’approximer et on l’optimise pour trouver 𝜏 ̂. Denis Larocque © Techniques d’exploitation de données 877 / 954 Modélisation incrémentale et effets de traitements individuels Implémentations Subsection 4 Implémentations Denis Larocque © Techniques d’exploitation de données 878 / 954 Modélisation incrémentale et effets de traitements individuels Implémentations Implémentations Les approches à un modèle (𝑆 -learner), à deux modèles (𝑇 -learner), l’approche 𝑋 -learner et l’approche par transformation peuvent être implémentées facilement car elles utilisent directement les modèles existants. Mais il existe tout de même des implémentations pour nous simplifier encore plus la vie. En R, il y a les packages causalToolbox et rlearner sur GitHub. En Python, il y a la librairie causalML. Denis Larocque © Techniques d’exploitation de données 879 / 954 Modélisation incrémentale et effets de traitements individuels Implémentations L’approche 𝑅-learner est dans le package R rlearner sur GitHub et aussi dans la librairie causalML en Python. Des arbres avec des critères de splits spécialisés et des forêts qui utilisent ces arbres sont disponibles dans la librairie causalML en Python. En R, il y a le package grf et aussi le package uplift mais ce dernier est archivé sur CRAN car il n’a pas été mis à jour depuis longtemps. La courbe et le coeffient Qini (voir la prochaine section) sont disponibles dans le package R tools4uplift et dans la librairie causalML en Python. Comme ce domaine est en pleine effervescence, gardez l’oeil ouvert pour de nouvelles méthodes, packages et librairies. Denis Larocque © Techniques d’exploitation de données 880 / 954 Modélisation incrémentale et effets de traitements individuels Évaluation de la performance d’un modèle Subsection 5 Évaluation de la performance d’un modèle Denis Larocque © Techniques d’exploitation de données 881 / 954 Modélisation incrémentale et effets de traitements individuels Évaluation de la performance d’un modèle Évaluation de la performance d’un modèle Tel que nous l’avons mentionné, un défi important en analyse incrémentale est l’évaluation de la performance d’un modèle. Dans le contexte d’apprentissage supervisé habituel, la réponse d’intérêt est le 𝑌 et nous avons des observations directes. Par conséquent, on peut calculer plusieurs mesures de performance pour un modèle, généralement avec un échantillon de validation ou bien par validation-croisée. Par exemple, pour un 𝑌 continu, on peut calculer l’erreur quadratique moyenne 1 𝑛 ∑(𝑦 − 𝑦𝑖̂ )2 , 𝑛 𝑖=1 𝑖 où les 𝑦𝑖 sont les vraies valeurs dans l’échantillon de validation et où les 𝑦𝑖̂ sont les prédictions pour ces valeurs obtenues du modèle qui a été ajusté sur l’échantillon d’entrainement. Denis Larocque © Techniques d’exploitation de données 882 / 954 Modélisation incrémentale et effets de traitements individuels Évaluation de la performance d’un modèle En analyse incrémentale on s’intéresse à un effet de traitement. On peut voir cela du point de vue de l’effet de traitement individuel 𝑌 (1) − 𝑌 (0), où, rappelez-vous, 𝑌 (1) est la réponse qu’on observerait si le sujet recevait le traitement et 𝑌 (0) est la réponse qu’on observerait si le sujet était dans le groupe contrôle. Il s’agit de réponses potentielles. 𝑌 (1) − 𝑌 (0) est une variable aléatoire. En pratique, on observe seulement l’une d’elles. On ne connait donc pas la “vérité” pour vérifier si notre prédiction de l’effet de traitement pour ce sujet est bonne ou non. Denis Larocque © Techniques d’exploitation de données 883 / 954 Modélisation incrémentale et effets de traitements individuels Évaluation de la performance d’un modèle On peut aussi voir cela du point de vue du CATE, 𝜏 (𝑥) = 𝐸[𝑌 (1) − 𝑌 (0)|𝑋 = 𝑥]. Cette fois-ci, il s’agit d’une moyenne théorique. Mais on a le même problème, on n’a pas de vrai CATE. On ne peut donc pas vérifier si notre estimation du CATE, 𝜏 ̂(𝑥), est bonne ou non. Pour être bien clair, si on connaissait les vrai 𝜏 pour nos observations, on pourrait calculer des mesures telles 1 𝑛 ∑(𝜏 − 𝜏𝑖̂ )2 , 𝑛 𝑖=1 𝑖 comme dans le problème d’apprentissage supervisé habituel. Denis Larocque © Techniques d’exploitation de données 884 / 954 Modélisation incrémentale et effets de traitements individuels Évaluation de la performance d’un modèle Notez que lorsqu’on fait des simulations et qu’on génère nous-même les données, on connait généralement les vrais 𝜏𝑖 et on peut alors utiliser les mesures de performances habituelles. Mais avec des vraies données, ce n’est pas possible. Ceci à plusieurs impacts. L’évaluation d’une mesure de performance avec un échantillon de validation (ou par validation-croisée) intervient à plusieurs endroits. On en a besoin pour optimiser les hyper-paramètres pour un modèle donné. On en a besoin pour comparer divers modèles afin de choisir le meilleur. On en a besoin pour avoir une évaluation juste de la performance du modèle retenu à la fin. Denis Larocque © Techniques d’exploitation de données 885 / 954 Modélisation incrémentale et effets de traitements individuels Évaluation de la performance d’un modèle Il existe des manières d’évaluer la performance mais elles sont toutes indirectes et aucune n’est pleinement satisfaisante. La plus connue, qui est pour le cas où 𝑌 est binaire, provient d’une idée similaire à celle de la courbe des gains cumulés que nous avons vu dans le contexte de prédiction d’un 𝑌 binaire. Il s’agit de la courbe Qini, aussi appelée (relative incremental uplift), et du coefficient Qini qui lui est associé. Comme pour n’importe quelle mesure de peformance, il est préférable de la calculer avec un échantillon de validation. Soient 𝜏1̂ , … , 𝜏𝑛̂ les estimations de 𝜏 (𝑥) provenant du modèle pour les observation de cet échantillon. Denis Larocque © Techniques d’exploitation de données 886 / 954 Modélisation incrémentale et effets de traitements individuels Évaluation de la performance d’un modèle On débute en ordonnant les observations en ordre décroissant selon les 𝜏𝑖̂. N’oubliez pas que plus 𝜏𝑖̂ est grand, plus le traitement risque d’être bénéfique pour le sujet, selon le modèle du moins. Pour 𝜙 ∈ [0, 1], soit 𝑆(𝜙) l’ensemble des observations de taille 𝑛𝜙 (arrondi à l’entier vers le haut) qui ont les valeurs de 𝜏𝑖̂ les plus élévées. Par exemple, pour 𝜙 = 1/10, 𝑆(1/10) contient le top 10% des observations avec les 𝜏 (𝑥) estimés les plus grands. Denis Larocque © Techniques d’exploitation de données 887 / 954 Modélisation incrémentale et effets de traitements individuels Évaluation de la performance d’un modèle (𝜙) Soit 𝑅𝑡 le nombre de réponses positives (𝑌 = 1) parmi les sujets du groupe traitement qui sont dans 𝑆(𝜙). (𝜙) Soit 𝑛𝑡 le nombre de sujets du groupe traitement qui sont dans 𝑆(𝜙). (𝜙) Soit 𝑅𝑐 le nombre de réponses positives (𝑌 = 1) parmi les sujets du groupe contrôle qui sont dans 𝑆(𝜙). (𝜙) Soit 𝑛𝑐 le nombre de sujets du groupe contrôle qui sont dans 𝑆(𝜙). Soit 𝑛𝑡 le nombre de sujets du groupe traitement dans tous (1) l’échantillon. Donc, selon la notation, 𝑛𝑡 = 𝑛𝑡. Soit (𝜙) 1 (𝜙) (𝜙) 𝑛 𝑔(𝜙) = (𝑅𝑡 − 𝑅𝑐 𝑡(𝜙) ). 𝑛𝑡 𝑛𝑐 Denis Larocque © Techniques d’exploitation de données 888 / 954 Modélisation incrémentale et effets de traitements individuels Évaluation de la performance d’un modèle La courbe Qini est le graphe de 𝑔(𝜙) (axe des 𝑦) lorsque 𝜙 (axe des 𝑥) varie de 0 à 1. Plus la courbe est élévée, meilleur est le modèle. (𝜙) (𝜙) (𝜙) (𝜙) La quantité 𝑅𝑡 − 𝑅𝑐 𝑛𝑡 /𝑛𝑐 représente, pour les sujets dans 𝑆(𝜙), le nombre de réponses positives parmi les traités moins le nombre de réponse positives parmi les contrôles, mais corrigé comme s’il y avait le même nombre de contrôles et de traités dans 𝑆(𝜙). Donc, plus c’est grand, plus le traitement est bénéfique pour les sujets dans 𝑆(𝜙). Mais comme les sujets ont été ordonnés selon leurs estimations de l’effet de traitement, cela veut dire que plus cette valeur est grande, meilleur est le modèle car il a bien identifié les sujets pour lesquels le traitement est bénéfique. Denis Larocque © Techniques d’exploitation de données 889 / 954 Modélisation incrémentale et effets de traitements individuels Évaluation de la performance d’un modèle Notez que (1) (1) 𝑅𝑡 𝑅𝑐 𝑔(1) = (1) − (1) 𝑛𝑡 𝑛𝑐 est la proportion de réponses positives parmi les traités moins la proportion de réponses positives parmi les contrôles, pour tout l’échantillon. Le coefficient de Qini est l’aire sous la courbe Qini. Plus il est grand, meilleur est le modèle. Denis Larocque © Techniques d’exploitation de données 890 / 954 Modélisation incrémentale et effets de traitements individuels Exemple Subsection 6 Exemple Denis Larocque © Techniques d’exploitation de données 891 / 954 Modélisation incrémentale et effets de traitements individuels Exemple Exemple Nous allons illustrer certaines des méthodes que nous venons de voir avec un exemple. Avec des vraies données, le vrai CATE est inconnu. Nous allons donc utiliser des données simulées où le CATE est connu. Cela nous permettra de voir comment les méthodes performent réellement. Mais nous allons aussi utiliser le coefficient de Qini. Ceci nous permettra de voir si ce coefficient est un guide raisonnable pour mesurer la performance des modèles. Denis Larocque © Techniques d’exploitation de données 892 / 954 Modélisation incrémentale et effets de traitements individuels Exemple Il s’agit d’un contexte fictif où nous avons un 𝑌 binaire et un traitement 𝑊 binaire. Vous pouvez penser qu’il s’agit du contexte présenté dans l’exemple au début du chapitre. C’est-à-dire, 𝑌 indique si le client a acheté le produit ou non et 𝑊 indique si on lui a offert un rabais ou non. Les données d’entrainement et de test se trouvent dans les fichiers uplifttrain.txt et uplifttest.txt. On suppose que 10% des clients de l’échantillon d’entrainement ont été assignés au hasard au groupe traitement. Il y a 8 variables explicatives 𝑋1 , … , 𝑋8 , toutes continues. Denis Larocque © Techniques d’exploitation de données 893 / 954 Modélisation incrémentale et effets de traitements individuels Exemple Comme dans les exemples précédents, on doit tout faire avec le fichier d’entrainement (uplifttrain). le fichier test (uplifttest) n’existe pas et sert seulement à évaluer la performance des modèles à la fin. Lisons d’abord les données. On voit qu’il y a 10000 observations dans chacun des fichiers. Denis Larocque © Techniques d’exploitation de données 894 / 954 Modélisation incrémentale et effets de traitements individuels Exemple # R # répertoire où se trouvent les fichiers pa="d:/11000347/Desktop/nouvelle_version_6600_A2022/data/" # Importations des données uplifttrain=read.table(paste(pa,"uplifttrain.txt",sep=""),header=TRUE) uplifttest=read.table(paste(pa,"uplifttest.txt",sep=""),header=TRUE) dim(uplifttrain) ## 10000 11 dim(uplifttest) ## 10000 11 Denis Larocque © Techniques d’exploitation de données 895 / 954 Modélisation incrémentale et effets de traitements individuels Exemple # Pour voir les premières lignes du fichier d'entrainement head(uplifttrain) ## cate x1 x2 x3 x4 x5 x6 ## 1 0.2181113 0.8758560 5.479076e-01 0.2161522 0.1598646 0.7159878 0.82847070 ## 2 0.6636012 1.0268702 2.850826e-05 0.4345003 0.8665823 0.3797207 0.09788352 ## 3 0.2416642 2.0834525 1.121645e-02 0.3308756 0.2061300 0.9188507 0.07389726 ## 4 0.2972226 0.8782224 7.138470e-02 0.3511046 0.4883497 0.3416479 0.31043139 ## 5 0.3025411 0.1614190 6.758491e-02 0.1639272 0.5189354 0.2606775 0.21473742 ## 6 0.4355378 -0.9119912 3.028674e-01 0.0514517 0.4939535 0.5020920 0.03398525 ## x7 x8 w y ## 1 0.13505591 -0.7916774 1 0 ## 2 0.30388948 0.4882207 1 1 ## 3 0.01877887 -0.8784694 1 1 ## 4 0.13826250 -0.3477786 1 0 ## 5 0.13830260 -1.2601037 1 0 ## 6 0.27145513 0.7300980 1 0 Denis Larocque © Techniques d’exploitation de données 896 / 954 Modélisation incrémentale et effets de traitements individuels Exemple Le code qui suit montre comment obtenir des résumés détaillés des données avec le package summarytools que nous avons déjà utilisé dans les exemples précédents. Tel que mentionné, on voit qu’il y a bien 10% (1000) sujets dans le groupe traitement dans l’échantillon d’entrainement. Il ne semble y avoir rien d’étrange dans les variables explicatives. Les huit sont bien continues. Certaines ont une distribution symétrique et d’autres une distribution asymétrique. Denis Larocque © Techniques d’exploitation de données 897 / 954 Modélisation incrémentale et effets de traitements individuels Exemple Nous obtenons également la distribution de 𝑌 pour le groupe traitement et le groupe contrôle dans l’échantillon d’entrainement. Il y a 12,6% (1130/9000) des clients qui ont acheté le produit dans le groupe contrôle et 50.5% (505/1000) dans le groupe traitement. Le traitement (offir le rabais) a clairement un effet global sur la probabilité d’achat. Mais on veut cibler plus précisément les clients pour lesquels le rabais aura le plus grand effet sur la probabilité d’achat. Denis Larocque © Techniques d’exploitation de données 898 / 954 Modélisation incrémentale et effets de traitements individuels Exemple # R # Il est possible d'avoir des résumés détaillés avec ce package (non exécuté) # library(summarytools) # Pour obtenir un résumé des données directement en format html (non exécuté) # view(dfSummary(uplifttrain),file="d:/11000347/Desktop/summaryuplifttrain.html") # view(dfSummary(uplifttest),file="d:/11000347/Desktop/summaryuplifttest.html") # Les résultats se trouvent dans les fichiers "summaryuplifttrain.html" et # "summaryuplifttest.html"" fournis avec le matériel du cours # Distribution de Y dans les groupes traitement et contrôle table(uplifttrain$y[uplifttrain$w==0]) ## ## 0 1 ## 7870 1130 table(uplifttrain$y[uplifttrain$w==1]) ## ## 0 1 ## 495 505 Denis Larocque © Techniques d’exploitation de données 899 / 954 Modélisation incrémentale et effets de traitements individuels Exemple Pour que ce soit très clair, la variable cate, ce qu’on cherche à estimer, ne serait pas disponible dans une vraie application. On ne peut l’utiliser aucunement pour entrainer les modèles. Nous allons donc la retirer des fichiers et la sauvegarder dans des vecteurs à part. On voit que l’effet de traitement cate varie de 0,028 à 0,723 avec une moyenne de 0,374 dans l’échantillon d’entrainement. De plus, en inspectant son histogramme, dans summaryuplifttrain.html, on voit qu’il varie beaucoup d’un sujet à l’autre. Tenter de trouver les clients avec les grandes valeurs de cate est donc potentiellement profitable. Mais, encore une fois, on ne connaitrait pas cette distribution dans une vraie application. Denis Larocque © Techniques d’exploitation de données 900 / 954 Modélisation incrémentale et effets de traitements individuels Exemple # R # Retirer "cate" des data frame et la mettre dans des fichiers à part uplifttraincate=uplifttrain[,c("cate")] uplifttestcate=uplifttest[,c("cate")] uplifttrain=uplifttrain[,c("x1","x2","x3","x4","x5","x6","x7","x8","w","y")] uplifttrain=uplifttrain[,c("x1","x2","x3","x4","x5","x6","x7","x8","w","y")] # Résumé de "cate" dans l'échantillon d'entrainement summary(uplifttraincate) ## Min. 1st Qu. Median Mean 3rd Qu. Max. ## 0.02774 0.22569 0.37341 0.37436 0.52513 0.72307 Denis Larocque © Techniques d’exploitation de données 901 / 954 Modélisation incrémentale et effets de traitements individuels Exemple Nous allons utiliser deux critères pour mesurer la performance des modèles. Les deux seront calculés sur les données test. Le premier est l’erreur absolue moyenne 1 𝑛 𝑀 𝐴𝐸𝜏 = ∑ |𝜏𝑖̂ − 𝜏𝑖 | 𝑛 𝑖=1 où les 𝜏𝑖̂ et les 𝜏𝑖 sont les estimations et les vraies valeurs du CATE. Le deuxième critère sera le coefficient Qini que nous avons vu plus tôt. La raison principale est de voir si les deux critères tendent à ranger les modèles dans le même ordre. Denis Larocque © Techniques d’exploitation de données 902 / 954 Modélisation incrémentale et effets de traitements individuels Exemple Dans une vraie application, 𝑀 𝐴𝐸𝜏 ne sera jamais disponible, même pour l’échantillon d’entrainement. Le coefficient de Qini quant à lui sera estimable avec l’échantillon d’entrainement. Voici d’abord une fonction pour calculer ces deux critères. Le coefficient de Qini sera obtenu avec le package tools4uplift. Denis Larocque © Techniques d’exploitation de données 903 / 954 Modélisation incrémentale et effets de traitements individuels Exemple # R # Fonction pour calculer la performance (MAE et Qini) perflift=function(w,y,pred,cate,nom) { # w = traitement (0-1) # y = vrai y (0-1) # pred = estimation du cate selon un modèle # cate = vraie valeur du cate (seulement disponible avec des # données simulées) # nom = nom du modèle library(tools4uplift) mae=mean(abs(cate-pred)) qini=NA if(sd(pred)>0) { perfm=PerformanceUplift(data = data.frame("w"=w,"y"=y,"pred"=pred), treat = "w", outcome = "y", prediction = "pred", equal.intervals = TRUE, nb.group = 50) qini=QiniArea(perfm) } out=data.frame(mae,qini) names(out)=c("MAE","QINI") row.names(out)=nom out } Denis Larocque © Techniques d’exploitation de données 904 / 954 Modélisation incrémentale et effets de traitements individuels Exemple On va se donner un benchmark simple. Il s’agit d’utiliser la proportion de clients ayant acheté le produit dans le groupe traitement moins la proportion de clients ayant acheté le produit dans le groupe contrôle, pour l’ensemble de l’échantillon d’ent