Facteurs croisés et modèles à effets mixtes

Un billet pour vous décrire un joli problème statistique auquel les modèles à effets mixtes sont censés répondre, mais qui malheureusement, n’aident pas trop ; sauf si on les utilise vraiment très bien.

Il s’agit d’une étude portant sur l’évaluation d’une intelligence artificielle (https://pubs.rsna.org/doi/10.1148/radiol.2021203886) assistant le travail de lecture radiographique. L’intelligence artificielle détecte de potentielles fractures osseuses sur des radiographies conventionnelles et les surligne avec un carré blanc, aidant ainsi le lecteur à les identifier. L’objectif de l’étude était de montrer que les lecteurs (radiologues et urgentistes de divers niveaux d’expérience) gagnaient en performance diagnostique grace à l’intelligence artificielle. L’analyse principale avait pour objectif de montrer une supériorité sur la sensibilité et une non-infériorité sur la spécificité (au seuil -3%). Pour être complètement transparent, il existe un conflit d’intérêt puisque la plupart des auteurs appartiennent à l’entreprise qui commercialise l’intelligence artificielle. En tant que statisticien j’étais néanmoins indépendant et n’ai reçu aucune rémunération.

La méthodologie a consisté à sélectionner 600 examens, dont 300 avec fracture et 300 sans fracture. Le gold standard a été créé par la lecture de trois experts du domaine. Les 600 examens ont été subdivisés aléatoirement en deux jeux de données (A et B) équilibrés de 300 examens, dont 150 avec et 150 sans fracture. Ensuite, 12 lecteurs de différents niveaux d’expérience (indépendants des experts du gold standard), ont lu les deux jeux de données, dont un des deux avec aide de l’IA et l’autre sans l’aide de l’IA. Sur les 12 lecteurs, 6 ont lu le jeu A avec l’aide de l’IA et le jeu B sans aide et les 6 autres ont lu le jeu B avec l’aide de l’IA et le jeu A sans aide.

Une nouvelle étude, avec la version suivante de l’IA (en cours de reviewing) utilise une méthodologie similaire, mais cette fois-ci avec 24 lecteurs et 480 examens. Cette fois-ci les 480 examens sont tous lus deux fois, avec et sans IA, avec une période de wash-out évitant une mémorisation des images. C’est plutôt de cette étude dont nous parlerons car certains phénomènes sont plus visibles avec celle-ci. Il n’y avait aucune donnée manquante. Tout était parfaitement équilibré (480 examens pour chacun des 24 lecteurs et 24 lecteurs pour chacun des examens).

On souhaite généraliser les résultats à l’ensemble des examens possibles, en supposant que l’échantillon d’examens radiologiques est représentatif, ainsi qu’à l’ensemble des médecins susceptibles de lire ces examens. En considérant qu’on s’intéresse aux moyennes et différences de moyennes des sensibilités et spécificités, on s’aperçoit qu’il est extrêmement simple de faire une estimation ponctuelle de l’effet de l’aide par l’IA. Considérons la sensibilité, égale à la moyenne de la variable binaire « vrai positif » (qui vaut 0 ou 1) dans les 240 examens avec fracture. Toutes les manières intuitives de calculer l’effet conduisaient à la même estimation ponctuelle, du fait de l’équilibre de taille des groupes:

  1. Calculer l’effet de l’IA pour chacun des lecteurs, puis faire la moyenne des 24 effets obtenus
  2. Calculer l’effet de l’IA pour chacun des examens, puis faire la moyenne des 240 effets obtenus
  3. Calculer la sensibilité avec l’IA sur les 240×24 observations sans aide de l’IA, puis calculer la sensibilité sur les 240×24 observations avec aide de l’IA, puis faire la différence de ces deux sensibilités, en ignorant complètement les colonnes « numéro d’examen » et « numéro de lecteur »

Ces trois estimateurs conduisent à la même estimation, parce que tous les lecteurs lisent le même nombre d’examens et tous les examens sont lus le même nombre de fois. Par ailleurs, ils ont tous des propriétés de pondération uniforme souhaitables. On n’a pas vraiment raison de donner plus de poids à certains lecteurs qu’à d’autres (même si après réflexion, on pourrait imaginer donner plus de poids à ceux qui ont un plus grand volume d’activité, du fait d’un impact de santé publique plus important). On n’a pas non plus raison de donner un plus grand poids à certains examens qu’à d’autres. Le design équilibré est intéressant parce qu’il fournit une précision statistique optimale pour cette moyenne non pondérée, en négligeant les problèmes d’hétéroscédasticité. En réalité, il existe des modèles foireux (modèles linéaires généralisés à effets mixtes) qui vont exploiter cette hétéroscédasticité, mais vu qu’ils ne convergent pas, je ne m’étendrai pas trop dessus.

Donc, la question ne porte pas vraiment sur la statistique en elle même, mais seulement sur la manière d’estimer son incertitude.

Description des fluctuations : modèle N°1

Il existe alors deux sources de fluctuations d’échantillonnage : les examens (n=480) et les lecteurs (n=24). Ces deux sources sont croisées, avec les problèmes suivants à prendre en compte:

  1. Il existe une variance de la difficulté d’un examen modélisable comme la sensibilité moyenne pour cet examen sur un nombre virtuellement infini de lecteurs. Notation (1|examen) dans le package R ‘lme4’.
  2. Il existe une variance de la compétence de base (sans IA) d’un lecteur, modélisable comme la sensibilité moyenne de ce lecteur sur un nombre virtuellement infini d’examens. Notation (1|lecteur)
  3. Il existe une interaction statistique entre les deux, c’est-à-dire, qu’un examen facile pour un lecteur peut être difficile pour un autre et vice versa, parce que différents lecteurs ne tombent pas dans les mêmes « pièges » et qu’un examen est difficile parce qu’il est piège. Notation (1|paste(lecteur,examen))
  4. Il existe une variance de l’effet de l’IA selon le lecteur ; c’est-à-dire, que l’IA peut faire gagner beaucoup de sensibilité à certains lecteurs alors qu’elle n’apportera peu, voire rien à d’autres. Notamment, un lecteur dont la sensibilité est très élevée ne pourra presque rien gagner avec l’IA alors qu’un lecteur qui a une sensibilité mauvaise pourra progresser beaucoup, ou pas. Notation (aide|lecteur).
  5. Il existe une interaction entre l’effet moyen (sur une infinité de lecteurs) de l’IA et l’examen ; notamment, il est évident qu’une IA n’apportera rien du tout à un examen sur lequel elle ne détecte aucune fracture réelle alors qu’elle pourra apporter beaucoup sur un examen pour lequel elle détecte une fracture très difficile à voir. Notation (aide|examen).
  6. Il existe bien sûr, une variance résiduelle, c’est-à-dire, que le même examen, lu par le même lecteur, deux fois à un mois d’intervalle, aura un résultat différent. Il y a des effets contextuels dans ce résidu.
  7. Il existe une interaction entre l’effet de l’IA, le lecteur et l’examen ; c’est-à-dire que l’effet de l’IA (virtuellement moyenné sur des centaines de lectures à plusieurs mois d’intervalles pour neutraliser la variance résiduelle) pourra être particulièrement élevé pour un lecteur donné sur un examen donné ; notamment on peut imaginer que l’IA aura un effet maximal si, à la fois, le lecteur était fortement susceptible de rater la fracture sur l’examen, que l’IA surligne la fracture considérée et qu’en plus, le lecteur est fortement susceptible de faire confiance à l’IA ; en effet, beaucoup de lecteurs (notamment pour les fractures vertébrales) ne font pas confiance à l’IA quand bien même elle a raison sur la fracture qu’elle surligne. Notation (aide|paste(lecteur,examen))
  8. Il existe une covariance entre le N°2 et le N°4, a priori, fortement négative. C’est-à-dire que ceux qui sont les moins bons vont s’améliorer plus fortement en moyenne ; à l’opposé, un lecteur dont la sensibilité est de 100% ne pourra voir qu’une stabilité ou une dégradation avec l’aide de l’IA. Implicitement inclus dans la notation (aide|lecteur).
  9. Il existe une covariance entre le N°1 et le N°5 (même concept que pour la covariance entre le N°2 et le N°4 mais s’appliquant aux examens plutôt qu’aux lecteurs). Implicitement inclus dans la notation (examen|lecteur).
  10. Il existe une covariance entre le N°3 et le N°7 (même concept que la covariance entre le N°1 et le N°5 mais conditionnel au lecteur). Implicitement inclus dans la notation (aide|paste(lecteur,examen))

Toutes ces sources de fluctuations d’échantillonnage ne sont pas pure spéculation de ma part : j’ai même pu évaluer, de manière directe ou indirecte les huits premiers éléments que je décris, notamment l’existence de l’interaction IA×lecteur×examen (effet N°7).

Avant de parler des modèles à effets mixtes, je vais décrire des estimateurs d’intervalles de confiance beaucoup plus simples et génériques et, fortement ou faiblement biaisés.

Simplification : modélisation des deltas (modèle N°2)

Pour le sous-groupes des examens avec fracture, le jeu de données contient 240 examens × 24 lecteurs × 2 observations (avec et sans IA) soit 11520 observations. Pour chaque observation, on a une variable binaire « fracture détectée par le lecteur » 1 ou 0. Comme on s’intéresse à la différence de moyennes de cette variable avec IA versus sans IA, on peut commencer par calculer les différences appariées pour chacune des 240×24 combinaisons (examen, lecteur), conduisant ainsi à un échantillon de 240×24=5760 observations. Cette différence appariée (delta(IA)) aura trois valeurs possibles : -1, 0 ou 1. La moyenne de la variable delta(IA) sur ces 5760 observations est alors l’estimateur de l’effet moyen de l’IA. Ces différences appariées vont simplifier beaucoup les choses, puisqu’au lieu de prendre en compte à chaque fois trois effets (la sensibilité sans IA, le gain de sensibilité par l’IA et la covariance entre les deux) pour trois facteurs (le lecteur, l’examen et la combinaison lecteur×examen) on n’a plus qu’à prendre en compte un seul effet, toujours pour les trois facteurs considérés.

Il reste alors à prendre en compte:

  • A : Une variance inter-lecteur du delta(IA) (prenant implicitement en compte les effets N°2, N°4 et N°8 du modèle précédent)
  • B : Une variance inter-examen du delta(IA) (effets N°1, N°5 et N°9 du modèle précédent)
  • C : Une interaction lecteur-examen sur le delta(IA) (effets N°3, N°7 et N°10 du modèle précédent)
  • D : Un résidu

En réalité, seuls les effets N°4, N°5 et N°7 du modèle précédent ont réellement à être modélisés, puisque les autres n’ont aucun impact sur la statistique.

Méthode N°1 : Student

Pour commencer, on peut calculer, pour chacun des 24 lecteurs, le gain de sensibilité attribuable à l’intelligence artificielle, comme la différence de sensibilité avec vs sans IA. On peut ensuite utiliser la méthode de Student à un seul échantillon, pour en estimer la moyenne. C’est équivalent à un test de Student sur séries appariées, avec n=24 paires d’observations avec vs sans IA. Cette méthode néglige la covariance entre les gains des différents lecteurs attribuable au fait que les mêmes examens sont lus. Néanmoins, cette méthode est acceptable si on considère que l’échantillon de 480 examens est suffisamment grand pour que l’erreur d’estimation du gain par l’IA de chacun des lecteurs (erreur type sur le delta(IA)), est négligeable en comparaison à la variance inter-lecteur de ce gain. On observait alors un effet moyen de l’IA sur la sensibilité à +10,4% (IC 95%: 7,7% à 13,1%), avec donc, une largeur d’intervalle de confiance à 5,4%.

Ce test de Student est nettement moins mauvais qu’il n’y paraît naïvement. Déjà, il bénéficie des différences appariées, et donc a tous les bénéfices du modèle N°2. Il prend évidemment en compte le A du modèle simplifié (variance inter-lecteur des deltas), mais aussi le B (variance inter-examen des deltas) car cette variance va être incluse dans la variance observée par le Student. En effet, pour un lecteur donné, avec n=240 deltas, cette variance inter-lecteur conduira à une légère fluctuation d’échantillonnage de la moyenne de ces 240 deltas, égale à la variance inter-examen divisée par 240, car les examens sont indépendants. Cela va indirectement se ressentir dans la variance du Student avec n=24 lecteurs, car le Student ne fait pas de shrinkage : la variance empirique du Student va à la fois contenir la vraie variance inter-lecteur mais aussi l’erreur de mesure de l’effet lecteur attribuable aux variances inter-examens qui va s’y additionner. Le Student prendra aussi en compte la covariance entre les lecteurs, attribuable au fait qu’ils ont observé les mêmes examens. En effet, cette covariance positive va se traduire par une diminution de la variance empirique du Student. Le problème, c’est que la variance empirique du Student va TROP fortement diminuer à cause de cette covariance. Pour caricaturer, si les lecteurs étaient parfaitement stéréotypés (variance inter-lecteur nulle et résidu nul), alors avec n=5 examens et n=100 lecteurs, la variance empirique de Student serait nulle bien qu’en réalité il existe de véritables fluctuations d’échantillonnage sur les n=5 examens. Enfin, l’effet C du modèle simplifié va tendre à diminuer la covariance entre les lecteurs attribuables au fait qu’ils évaluent les mêmes examens, et donc, va atténuer le problème de surcorrection de la covariance dont je viens de parler.

Si on se réfère au premier modèle, non simplifié, le Student fait quand même totalement disparaître les effets N°1, N°2, N°3, N°8, N°9 et N°10) grâce à l’usage des différences appariées du modèle simplifié. Ça prend en compte les effets N°4 et N°5 (A et B du modèle simplifié) mais pas bien la covariance qu’elles engendrent, même si, l’effet N°7 tendrait à atténuer le problème.

Méthode N°2 : Double bootstrap

Globalement, quand je ne sais pas comment estimer une incertitude, j’ai le réflexe de penser au boostrap. Néanmoins, le schéma du boostrap doit suivre le schéma d’échantillonnage ! Ceux qui m’avaient commandé les analyses avaient boostrappé les lignes d’observation (n=480×24×2 = 23040 observations) en ignorant complètement que les lecteurs étaient soumis aux fluctuations d’échantillonnage et que les 48 observations d’un même examen étaient corrélées ! Cela conduit à un intervalle de confiance nettement trop étroit : (IC95% : 8,8% à 12,0%) d’une largeur de 3,3%. L’intervalle n’est pas aussi étroit qu’on pourrait craindre parce qu’en ignorant la covariance intra-lecteur et intra-examen, ça surestime aussi la variance sur certains aspects ; en bref, ça biaise fortement dans des directions opposées.

Ici, on peut proposer un boostrap respectant les règles d’échantillonnage : on tire au hasard avec remise 240 examens parmi 240 (on se limite aux examens avec fracture pour le calcul de la sensibilité), puis on les renumérote de N°1 à N°240 (puisque les « doublons » artificiellement créés par la méthode de rééchantillonnage sont censés représenter des examens différents si on tirait au hasard depuis la population). Ensuite, dans cet échantillon, on tire au hasard 24 lecteurs parmi 24 avec remise, et on les renumérote de 1 à 24. On peut alors calculer la statistique de gain moyen de sensibilité attribuable à l’IA.

On peut faire ce double boostrap sur le jeu des 240×24×2 observations ou sur le jeu des différences appariées de 240×24 observations, proposé dans la section « Simplification : modélisation des deltas », mais les deux sont strictement équivalents. En effet, les 2 observations d’un même couple (lecteur,examen) resteront toujours soudés dans un même bloc durant les deux rééchantillonnages.

Le double bootstrap (le terme correct serait plutôt bootstrap de 2 facteurs croisés, mais j’aime bien double bootstrap) est particulièrement simple et rapide à réaliser informatiquement sous R sur 240×24 différences appariées. Il suffit de faire une matrice de 24 lignes et 240 colonnes (ou 240 lignes et 24 colonnes, peu importe). On peut rééchantillonner simultanément les lignes et les colonnes avec la syntaxe du logiciel R mat[indices_lignes, indices_colonnes]. Ce rééchantillonnage simultané est équivalent à un rééchantillonnage successif des lignes puis des colonnes (ou des colonnes puis des lignes). Enfin, en appliquant la fonction mean() à la matrice entière, on calcule la moyenne de toutes les cellules de la matrice : c’est l’estimation ponctuelle de l’effet de l’IA.

Le rôle des lecteurs et des examens est parfaitement symétrique dans cette méthodologie de double bootstrap. Le bootstrap va prendre empiriquement en compte la structure de corrélation dans toute sa complexité. Cela semble simple intuitif et non biaisé au premier abord, mais en réalité, c’est un plus biaisé qu’il n’y paraît. Il y a une tendance à surestimer les fluctuations d’échantillonnage. Même si le rôle des examens et des lecteurs est symétrique, considérons la matrice de 24 lignes × 240 colonnes, contenant dans chaque cellule, un delta (gain de sensibilité sur cet examen pour ce lecteur par l’IA) et considérons des scenarii virtuels qui vont beaucoup nous aider à comprendre ce qu’il se passe.

Comme précisé précédemment, les 24 moyennes des 24 lignes sont corrélées parce que les 24 lecteurs observent les mêmes examens. Supposons dans un premier temps que cette corrélation soit nulle parce que tous les examens auraient la même difficulté et l’IA fournirait le même gain pour tous les examens : il y aurait juste un effet lecteur et un résidu mais pas d’effet examen. Même si c’est probablement faux, on peut imaginer qu’on puisse être plus ou moins proche de cette situation. Dans ce contexte d’absence d’effet examen, le simple bootstrap des lignes (ou l’intervalle de confiance de Student sur les n=24 lignes) ne sera pas biaisé du tout : il prendra à la fois en compte la variance inter-lecteur (varL) et la variance résiduelle déjà atténuée par n=240 mesures (varR/240) : ces deux variances (varL et varR/240) s’additionneront sous forme de la variance empirique des 24 moyennes des 24 lignes. Comme il n’y a aucune covariance entre ces lignes, l’estimateur de la variance n’est pas biaisé. Si on double boostrappe, on va surestimer la variance, parce que la variance résiduelle (divisée par 240) sera ajoutée une seconde fois lors du rééchantillonnage des colonnes ! La cause profonde de ce double comptage de la variance résiduelle est explicable par l’absence de shrinkage fait par le boostrap des lignes, c’est-à-dire qu’elle n’enlève pas la varR/240. Si on échantillonnait les 24 lignes à partir de la population plutôt que de faire un rééchantillonnage, on aurait la variance inter-lecteur varL, mais il ne s’y ajouterait pas varR/240 ; c’est alors que le rééchantillonnage des colonnes ajouterait de manière tout à fait juste ce varR/240 qui participe quand même aux fluctuations d’échantillonnage. De la même manière, comme il n’y a pas de shrinkage de la variance inter-examen, la varR/24 lui est artificiellement additionnée et cette varR/24 est comptée deux fois.

L’intérêt principal du double bootstrap réside néanmoins dans le scenarii où on aurait une forte variance de l’effet de l’IA selon l’examen, conduisant à une covariance non négligeable des 24 moyennes correspondant aux 24 lignes. Le double bootstrap prend en compte les effets les plus fins générant cette covariance, y compris l’interaction entre l’effet de l’IA, l’examen et le lecteur (effet N°7 du premier modèle et effet C du modèle simplifié).

Le double boostrap percentile conduisait à l’intervalle de confiance 10,4% (IC 95%: 6,8 à 14,1%) alors que le double bootstrap normal avec approximation à la loi de Student à 23 degrés de liberté conduisait à l’intervalle de confiance 10,4% (6,6% à 14,2%) avec une largeur d’intervalle à 7,57%, plus élevée que celle du Student naïf mais peut-être trop élevée.

Méthode N°3 : hybride Student/boostrap

La méthode de Student est biaisée parce qu’elle ne prend pas en compte la covariance entre les moyennes des lignes de la matrice de 24 lignes et 240 colonnes. Peut-on estimer et prendre en compte cette covariance ? Oui ! On peut estimer cette covariance empiriquement en rééchantillonnant les colonnes, sans toucher aux lignes. Une opération de boostrap correspond à un rééchantillonnage des 240 colonnes, suivi d’un calcul de 24 moyennes (de chacune 240 observations). Avec dix mille boostrap, on calcule la matrice de variance-covariance des 24 moyennes, c’est-à-dire, une matrice 24×24, qui donne la covariance de chaque paire de lecteurs. Il existe ensuite deux stratégies totalement équivalentes. Une première stratégie qui consiste à estimer la covariance entre deux lecteurs comme la moyenne de toutes les cellules de la matrice 24×24 après exclusion de la diagonale. C’est cette fameuse diagonale qu’il ne faut pas additionner à la variance du Student sur les 24 lignes pour ne pas doublement compter varR/240. Une fois qu’on a cet estimateur de la covariance, en appliquant les formules VAR(X+Y)=VAR(X)+VAR(Y)+2×COV(X,Y) aux 24 lecteurs, on s’aperçoit que l’erreur type de la moyenne sur les 24 lecteurs se calcule comme sqrt(n×V+C*(n^2-n))/n où n=24 est le nombre de lecteurs, V est la variance empirique inter-lecteur (qui est une estimation de varL+varR/240), et C est la covariance empirique qu’on a calculé par boostrap. On s’aperçoit aussi que le calcul est identique si on considère que les covariances de toutes les paires de lecteurs sont toutes différentes et qu’on n’en fait pas la moyenne. Cette méthode d’addition de la covariance est donc valable, quand bien même les lecteurs ne seraient pas indépendants les uns des autres. En pratique, on repose quand même sur l’hypothèse d’indépendance des lecteurs pour le calcul de varL. De manière intéressante, on obtient exactement la même erreur type si on inverse le rôle des lignes et des colonnes, c’est-à-dire si on fait un Student avec les 240 examinateurs et qu’on corrige le Student par la covariance empiriquement obtenue en rééchantillonnant les 24 lecteurs. Il existera éventuellement une différence mineure entre les deux approches : le nombre de degrés de liberté de la procédure de Student. Dans un cas, il paraît logique d’utiliser 23 ddl, en se disant que la covariance est assez faible et donc, devrait peu altérer la loi de Student, alors que dans l’autre cas, on est tenté d’utiliser 239 ddl. Si on veut un intervalle de confiance qui ne soit ni conservatif ni libéral, il faudrait un nombre de ddl intermédiaire qu’on devrait pouvoir calculer par des approximations de la somme de plusieurs lois de Student de variance différente et qui covarient. En l’absence de covariance, il y a l’approximation de Satterthwaite. En présence de covariance (modeste), on doit pouvoir généraliser cette approximation, mais c’est le cadet de nos soucis. Avec 23 ddl, la loi de Student est déjà proche de la loi normale (quantile 0,975 à 2,07 vs 1,96). On a des problèmes plus sévères que le calcul du nombre de degrés de libertés pour l’approximation à une statistique t de Student, tels qu’une imperfection de l’approximation normale avec n=24 lecteurs dont l’un d’entre eux s’avère d’ailleurs être un outlier.

Avec cette méthode, on trouve un intervalle de confiance à 10,4% (IC 95%: 6,9 à 13.9%), avec une largeur de 7,0%, soit un petit peu moins que le double bootstrap, parce qu’on n’a pas compté doublement les variances.

Méthode N°4 : double boostrap rétréci (shrinked)

On pouvait aussi se baser sur du double boostrap, mais enlever à la variance double-bootstrappée, les variances qui sont comptées doublement. Avec cette stratégie, on trouve une erreur type très proche de celle de la méthode hybride Student/boostrap soit un intervalle de confiance à 10,4% (IC 95% : 6,9 à 13,9%) ; pour lequel j’ai utilisé la loi de Student à 23 ddl, afin d’avoir le moins de différences possibles pour comparer aux autres méthodes. L’intervalle de confiance est de largeur équivalente avec ce double bootstrap rétréci (6,95%) et avec la méthode hybride (7,04%) ; la différence étant principalement explicable par la correction du biais de la variance (24/23) dans la méthode hybride ; la méthode hybride sans le facteur correcteur de la variance fournit une largeur d’intervalle de confiance à 6,96%. Les deux stratégies semblent donc très proches.

Méthode N°5 : modèle à effets mixtes sur le jeu de données simplifié

On va essayer de modéliser ce système de deux facteurs croisés, indépendamment échantillonnés, avec un modèle linéaire à effets mixtes (linear mixed effect model ou LMM). On va commencer à partir du modèle simplifié (les deltas), qui efface quand même la majeure partie de la complexité du problème. A priori, on peut modéliser les effets A, B et D du modèle simplifié, mais l’effet C (interaction examen×lecteur) est impossible à estimer car il n’y aura qu’une seule observation pour estimer chacun des effets, alors qu’on doit en même temps évaluer une variance résiduelle. Ce modèle aurait autant d’effets aléatoires que d’observations, en plus d’un effet fixe (intercept) ; c’est impossible à estimer. On doit alors faire un modèle ignorant l’interaction examen×lecteur. Le modèle suivant (syntaxe lme4) a été estimé :

deltaIA ~ (1|lecteur)+(1|examen)

Cela a permis d’estimer un intervalle de confiance (IC95% : 7,0% à 13,8%) d’une largeur de 6,790% selon la méthode du rapport de vraisemblance ; avec Wald, on a une largeur à 6,739%. Il semble que le modèle fournisse des résultats à peu près cohérents avec les méthodes de double bootsrap rétréci et hybride. En creusant il semblerait que la différence soit principalement explicable par une différence du calcul du nombre de degrés de liberté.

Méthode N°6 : modèle à effet mixte sur le jeu de données complet

Malheureusement, je crains que la plupart des usagers du modèle à effets mixtes ne passent pas par l’étape de simplification modélisatoire engendré par l’usage des deltas sur chaque couple (examen,lecteur).

On arrive à modéliser un maximum d’effets tout en conservant un modèle qui converge, on arrive à:

vrai_positif ~ IA+(IA|examen)+(IA|lecteur)+(1|couple_examen_lecteur)

On arrive alors à modéliser les effets N°1 à N°9 sauf l’effet N°7. On s’en sort donc aussi bien que pour le modèle simplifié. Les largeurs des intervalles de confiance sont presque les mêmes : 6,739% pour Wald et 6,788% pour le rapport de vraisemblance. On a juste des problèmes de performances computationnelles, qui restent encore tolérables (39 secondes pour calculer l’IC95% du rapport de vraisemblance).

Méthode N°7 : modèle à effets mixtes naïf

Tout le monde ne passe pas par la phase de réflexion de modélisation des 10 effets du modèle le plus complexe que j’ai décrit. Souvent, un modèle à effets mixtes est juste proposé avec un ou deux intercepts aléatoires. Qu’est-ce que ça donne ?

vrai_positif ~ IA+(1|examen)+(1|lecteur)

La largeur de l’intervalle de confiance est de 2,20%, ridiculement petit. Il faut dire qu’on fait l’hypothèse que les gains attribuables à l’IA sont strictement identiques pour tous les examens, tous les lecteurs et toutes les combinaisons examen×lecteur. Par rapport au modèle simplifié, il ne modélise, ni l’effet A, ni l’effet B, ni l’effet C. Ça modélise peut-être un peu la sensibilité sans IA, mais très mal l’effet de l’IA. C’est tellement mauvais que même le bootstrap des 23040 observations, ignorant totalement toute structure de corrélation est plus proche de la réalité (erreur type 3,3%).

Méthode N°8 : modèle identité-binomial à effets mixtes

Pour des raisons qui m’échappent, les modèles linéaires gaussiens sont mal vus sur les variables binaires : les modèles linéaires généralisés leurs sont souvent préférés, même si ça ajoute plus de biais que ça n’en enlève. Que donne le modèle linéaire généralisé à effets mixtes identité-binomial ? Ça ne converge pas, même avec le modèle très simplifié :

vrai_positif ~ IA+(1|examen)+(1|lecteur)

Ça ne converge toujours pas quand on enlève l’effet lecteur. Il faut enlever l’effet examen si on veut le faire converger. Le modèle le plus élaboré qui converge est:

vrai_positif ~ IA+(IA|lecteur)

Qui n’apporte alors rien par rapport au Student (methode N°1). Avec ça, on obtient un intervalle de confiance d’une largeur de 5,1% ; un peu plus étroit que celui du Student (méthode N°1).

Méthode N°9 : modèle logistique à effets mixtes

Pour le modèle logistique à effets mixtes, on arrivait à modéliser plus de termes :

vrai_positif ~ IA+(IA|examen)+(IA|lecteur)

C’est-à-dire qu’on est arrivé à mettre tous les effets les importants. Sauf que les coefficients ressortis par le modèle ne sont pas interprétables car conditionnels aux effets aléatoires inconnus. Ici, on observe un odds ratio conditionnel (dans le modèle logistique à effets mixtes) à 6,78 (IC95% : 4,41 à 10,44) alors que l’estimation ponctuelle de l’odds ratio marginal est à 1,65. En bref, c’est complètement bidon.

Si on veut comparer l’incertitude à celle d’autres méthodes, on va calculer le petit z d’un test de Wald de comparaison à l’effet nul (ici, odds ratio à 1). Le petit z était égal à 8,7 contre 6,1 pour la méthode hybride N°3.

Méthode N°10 : GEE

Les équations d’estimation généralisées (GEE) servent à estimer des effets marginaux dans des modèles linéaires généralisés en prenant en compte une structure de corrélation. Malheureusement, ils ne gèrent pas les facteurs croisés : il faut donc choisir si le cluster est défini par l’examen ou le lecteur. Le reviewer qui proposa les GEE suggérait l’usage de l’examen radiologique comme cluster, ignorant totalement les variances attribuables au lecteur. Ainsi, avec un modèle GEE identité-gaussien, on observait un intervalle de confiance d’une largeur 4,72%, avec l’estimateur sandwich de la variance ; la structure de corrélation (exchangeable, unstructured) n’a pas d’importance dans un design équilibré à partir du moment où on utilise l’estimateur sandwich, à ceci près que certains ont du mal à converger. En creusant un peu, on s’aperçoit qu’en cas de design équilibré et avec l’usage de la variance sandwich le modèle est équivalent à faire la moyenne de l’effet brut dans chaque cluster, puis faire un intervalle de confiance de Student grace à la variance empirique inter-cluster ; la seule différence porte sur les ddl ; la variance du GEE/sandwich n’est pas débiaisée par la multiplication par 24/23 et l’intervalle de confiance est calculé avec une loi normale plutôt qu’un Student à 23 ddl ; autrement, c’est exactement la même chose. À partir du moment où on utilise la variance sandwich, le GEE est invariant à la transformation : que l’on utilise un modèle identité-gaussien, logit-binomial, identité-binomial, on arrive au même intervalle de confiance, après usage de la méthode de delta pour tout exprimer sur l’échelle de la différence absolue de pourcentages. Attention, ce n’est pas une vérité générale sur les GEE, mais une vérité contextuelle au fait que (1) le design est équilibré et (2) les modèles fournissent exactement les mêmes prédictions parce qu’il n’y a a qu’un seul effet binaire ; ce serait différent si on intégrait plusieurs effets fixes, avec hypothèse d’absence d’interaction qui conduiraient à des contraintes différentes selon la fonction de lien utilisée. Donc, au final, un GEE ne fait pas mieux qu’un simple Student (méthode N°1).

Synthèse des résultats et discussion

Pour commencer, il faut oublier les GLMM. Ça converge rarement et dès qu’on met une fonction de lien différente de l’identité, ça fournit des résultats conditionnels ininterprétables et complètement gonflés.

Pour le problème des deux facteurs croisés équilibrés, les GEE n’apportent rien par rapport aux méthodes les plus naïves : le Student tout bête après avoir calculé la moyenne dans chaque cluster.

Les modèles linéaires à effets mixtes (LMM) fonctionnent probablement correctement s’ils sont bien spécifiés ; ils peuvent modéliser les principaux effets sans trop de problème de convergence (ils en avaient un peu plus sur l’étude sur 600 radiographies × 12 lecteurs × 1 lecture). Néanmoins, je ne suis vraiment pas à l’aise avec ces modèles, parce que je ne vois pas en détail leur mécanique s’exécuter, et donc je ne peux pas de manière complètement intuitive percevoir leurs conditions de validité ; la relative fiabilité de l’estimateur hybride (méthode N°2) me paraît bien plus facile à percevoir. Sans Monte Carlo de validation, je n’arrive pas à donner du crédit au LMM même dans des conditions aussi bonnes que celle de deux facteurs croisés équilibrés.

Le bootstrap est très intuitif et compréhensible, mais il faut respecter le schéma d’échantillonnage dans le rééchantillonnage et bien concevoir son fonctionnement pour percevoir les inflations de variance dès qu’il y a des schémas non triviaux ; il aide aussi beaucoup à percevoir les limites dues aux tailles d’échantillon : même si on a 23040 observations, on a que 24 lecteurs.

Les bidouilles comme l’estimateur hybride ne sont probablement pas parfaites, mais sont très compréhensibles et transparentes dans leur fonctionnement ; on se doute que ce qu’on fait n’est pas parfait, mais on sait que c’est largement suffisant et on comprend ce qu’on fait. Plutôt que de chercher à modéliser très finement les paramètres d’un modèle, on peut se contenter d’estimer et corriger des biais, en allant toujours du plus gros au plus fin. On peut aussi estimer et corriger un biais sans forcément le décomposer complètement dans un modèle.

Les effets aléatoires sont généralement nombreux si on veut une modélisation à peu près correcte. L’intercept aléatoire seul, peut être une solution catastrophique à un problème qui nécessite une ou plusieurs pentes aléatoires. L’usage d’un modèle à effets mixtes ne doit pas non plus vous empêcher de pré-traiter les données avant le calcul, avec la méthode des deltas, par exemple.

L’usage de moyennes et modèles linéaires estimés par les moindres carrés simplifie énormément la vie, surtout avec un design équilibré : tout est beaucoup plus cohérent parce que la plupart des manières d’attaquer le problème sont équivalents. En bonus, ça fournit des informations interprétables en absolu (différences absolues de pourcentages dont l’inverse est le Number Needed to Treat). C’est simple, compréhensible, efficace. Dès qu’on utilise des odds ratio, on s’aperçoit que c’est bien plus difficile de comprendre ce qui se passe en présence d’interactions ou d’ajustements. D’ailleurs, l’estimateur des moindres carrés offre de nombreuses propriétés extraordinaires, telles qu’un résidu moyen nul, même en cas d’écart aux conditions de validité (interactions, effets non linéaires), comme décrit dans un billet précédent.

Ai-je insisté sur le fait qu’il fallait comprendre ce qu’on faisait ? Je suis effrayé par les usages délirants des modèles à effets mixtes que je vois : résidus ultra-corrélés considérés comme indépendants, moyennes faites entre des objets de nature différente, confusion entre la position de variable réponse et de covariable, absence totale de questionnement sur la nature des effets qu’on mesure : on parle « d’effet » plutôt que de différences de moyennes pondérées ; on ne sait même plus quels sont les poids. On a aussi un gros souci avec la notion d’hypothèse. On fait des tonnes « d’hypothèses » manifestement fausses sans préciser, ni savoir, ce qu’il se passe lorsque ces hypothèses sont violées. Je suis bien plus à l’aise quand je définis précisément LA statistique qui m’intéresse avec une formule simple (p.e. différences de moyennes arithmétiques) et que je ne parle pas d’hypothèse mais de biais d’estimateur ; ces biais étant en partie évaluables et contrôlés comme ne dépassant pas un seuil déraisonnable. Même si ça peut être plus fortement biaisé qu’un modèle avec hypothèses, au moins, on sait à quoi s’attendre au pire des cas.

Laisser un commentaire

Votre adresse e-mail ne sera pas publiée. Les champs obligatoires sont indiqués avec *