Métaanalyse avec R
Sommaire
Objectifs
La méta-analyse est une méthode statistique qui fût inventée en 1904 par Karl Pearson pour tenter de surmonter le problème de puissance statistique réduite dans les études d'échantillons de petites tailles. Aujourd'hui, elle est très utilisée dans le domaine de la santé, où les chercheurs sont confrontés à la multiplicité d'études et d'essais thérapeutiques portant sur un même sujet, mais dont les protocoles mis en place et plus généralement les conditions des études diffèrent. En effet, la méta-analyse consiste en une synthèse quantitative de ces différents résultats pour obtenir une valeur "corrigée" qui intègre plus d'informations qu'une valeur unique prise dans la littérature (avec plus ou moins de justifications) ou bien qu'une simple moyenne des valeurs de la littérature.
Dans un premier temps, il s'agit donc de collecter les valeurs de résultats des essais déjà réalisés, ainsi que les conditions dans lesquelles ils ont été obtenus. Pour celà, il est donc nécessaire de définir au préalable des critères d'inclusion ou non des études qui soient précisément établis. A noter que cette recherche doit être la plus exhaustive possible afin que la méta-analyse soit une synthèse pertinente de la littérature existante sur le sujet défini.
Package METAFOR
Il existe trois packages dans R permettant de réaliser des méta-analyses : META, RMETA et METAFOR. Ce dernier est cependant le plus perfectionné puisqu’il offre la possibilité de réaliser des méta-régressions, d’utiliser un modèle à effets fixes ou bien aléatoires, ou encore car il intègre le calcul de différentes mesures d’hétérogénéité.
Pour la suite, on prend l’exemple d’un traitement binaire (succès : Y=1 ou échec : Y=0). appliqué à une population test (X=1) mais pas à une population de contrôle (X=0).
A noter que X peut aussi être l'exposition à une maladie.
Jeu de données BCG
Pour illustrer, on utilisera le jeu de données « BCG » du package metafor.qui contient les données de 13 essais cliniques du vaccin Bacillus Calmette–Guérin (BCG). Ainsi, le traitement (variable X) consiste en l'administration ou non du vaccin, et la variable observée (Y) est l'infection ou non par la tuberculose. Les données contiennent également trois variables explicatives possibles :
- l'année de publication
- la latitude de la région où s'est déroulé l'essai
- la méthode d'allocation du traitement choisie pour l'essai.
#metafor package library(metafor) ### load BCG vaccine data data(dat.bcg)
Les résultats d’un essai peuvent être représentés sous la forme suivante :
vacciné | non vacciné | Total individus | |
---|---|---|---|
Tuberculeux | a = tpos | c = cpos | n1 |
Non tuberculeux | b = tneg | d = cneg | n2 |
Calcul de la taille d'effet ("effect size")
Par exemple : Odd ratio, Risque Relatif, Ecart de moyenne,..
Cette grandeur est nécessaire pour comparer les différentes études, et le but de la méta-analyse est précisément d'en obtenir une valeur corrigée.
Il est de convention d’assumer la normalité des effect size (vrai à partir du moment où l’échantillon est de taille suffisante, selon le Théorème Central Limite). Pour normaliser les effect size il est habituel d’utiliser une transformation logarithmique.(Tutorial: Performing Meta-Analysis With R - Stephanie A. Kovalchik, National Cancer Institute)
A partir des observations (fonction escalc)
Dans le package metafor, la fonction escalc permet de calculer directement les logESs et les variances correspondantes (SElogES2) à partir du nombre d'individus observés (tableau). (Pour la suite on prendra l'exemple des Odd Ratios)
### calculate log odds ratios and corresponding sampling variances dat <- escalc(measure="OR", ai=tpos, bi=tneg, ci=cpos, di=cneg, data=dat.bcg) cbind(dat$yi,dat$vi)
A partir des données brutes
Si les données recueuillies sur les études sont déjà des Odd ratio, il est possible de retrouver ces grandeurs à l'aide de l'intervalle de confiance de l'OR (à 95%)
<math>\begin{align} IC(OR)=[a,b] \\
\Rightarrow IC(lnOR)=[lna,lnb] \\ \\ Comme lnOR \longrightarrow \mathcal{N}(\mu,SElnOR^2) \\ IC(lnOR)=[\mu \pm 1,96 \times SElnOR] \\ \\ Resoudre : \left\{ \begin{array}{l}
\mu - 1,96 \times SElnOR \\ \mu + 1,96 \times SElnOR
\end{array} \right. \\
\Rightarrow \left\{ \begin{array}{l}
\mu = \frac{ln(a \times b)}{2} \\ SElnOR = \frac{ln(\frac{a}{b})}{2\times 1,96}
\end{array} \right. \\ \end{align} </math>
Méta-analyse
Plusieurs fonctions sont disponibles dans le package metafor, voir la documentation du package pour plus de détails. Deux fonctions intéressantes :
- rma.uni
- rma.mv
rma.uni (Random effect Meta-Analysis)
Permet de réaliser des méta-analyses avec un modèle linéaire :
- à effets fixes ou mixed/random
- avec ou sans modérateurs
Par défaut la méthode utilisée est RE REML (Random Effect avec REstricted Maximum Likehood).
Pour la suite, on renomme alpha les lnOR, et beta la variance correspondante (soit SElnOR2)
#Alpha=LnOR #Beta=VarLnOR dat$alpha<-dat$yi dat$beta<-dat$vi
exemple 1 : Première régression, sans inclure de modérateurs dans le modèle :
### random-effects model using rma.uni() resrma<-rma.uni(yi=alpha, vi=beta, measure="OR",data=dat) resrma$b #lnOR=-0.745 exp(-0.745) #OR corrigé = 0.4747
avec :
- yi : vecteur de longueur k (le nombre d’études) des lnOR observés (ou plus généralement les outcomes normalisés)
- vi (ou sei) : vecteur de longueur k contenant les variances (ou standard errors) correspondantes.
Il est important de noter que dans toutes les fonctions utilisées, les arguments requis (yi et vi) sont les effect size normalisés (que l'on a renommé précédemment alpha et beta). Plusieurs arguments peuvent être ajoutés :
- mods : permet d’inclure un ou plusieurs modérateurs dans le modèle.
- weights : par défaut le poids donné à chaque étude dans l’analyse est égal à l’inverse de sa variance, mais il peut être spécifié.
- method : permet de spécifier un modèle à effets fixes par exemple.
exemple 2 : fixed effect modèle
### fixed-effects model using rma.uni() resfe<-rma.uni(yi=yi, vi=vi, measure="OR",data=dat, method="FE") resfe resfe$b #lnOR=-0.4361 exp(-0.4361) #OR corrigé = 0.6465531
On peut vouloir inclure des covariates dans le modèle, devant la forte hétérogénéité observée observée entre les études (I2=99.92%)
exemple 3 : metaregression avec une seule variable :
### Meta-regression incluant la variable latitude res.lat<-rma.uni(yi=yi, vi=vi, measure="OR",data=dat, mods =ablat) res.lat #Comparaison des résultats : exp(c(resrma$b, res.lat$b[1])) #OR corrigé : 0.4587 vs. 1.3513 c(resrma$I2, res.lat$I2) #I^2 (resrma$I2-res.lat$I2)/resrma$I2 #Baisse de l'hétérogénéité de 42.6%
exemple 4 : metaregression incluant les trois variables explicatives :
### Meta-regression incluant les variables latitude, année et allocation du traitement res.meta<-rma.uni(yi=yi, vi=vi, measure="OR",data=dat, mods = ~ ablat + factor(alloc) + year ) res.meta
rma.mv (Random effect Meta-Analysis via Multivariate/Multilevel Linear Models
Présente l’avantage d’offrir la possibilité de préciser la structure de la corrélation entre les publications (par exemple issues d’une même étude), via les arguments struct et random.
exemple : données BCG Dans l’exemple, on considère que chaque essai est composée de deux études différentes : une sur la population vaccinée et l’autre sur la population de contrôle. Dans chaque cas, on observe le nombre d’individus atteints ou non par la tuberculose (les variables sont respectivement out1 et out2).
### bivariate random-effects model using rma.mv() ### change data into long format dat.long <- to.long(measure="OR", ai=tpos, bi=tneg, ci=cpos, di=cneg, data=dat.bcg) #chaque étude divisée en deux : vaccinés (group=1) ou non (group=2) #out1=nb d'individus atteints #out2=nb non atteints ### calculate log odd-ratios and corresponding sampling variances dat.long <- escalc(measure="PLO", xi=out1, mi=out2, data=dat.long) #yi=ln(out1/out2) ### set "2" (controlled/non-vaccinated) to reference level of group variable dat.long$group <- relevel(dat.long$group, ref=2) ### bivariate random-effects model using rma.mv() res <- rma.mv(yi, vi, mods = ~ group, random = ~ group | study, struct="UN", method="REML", data=dat.long) res
avec :
- xi : nombre d’observations de l’évènement observé (ici la tuberculose)
- mi : nombre d’observations de l’évènement complémentaire
ici : yi = ln(out1/out2) method : de nouveau : fixed or random-effect model, différentes méthodes d’estimation pour tau2. A noter que la méthode du Restricted Maximum Likehood semble être la plus performante.
- random : permet de spécifier la structure du "random effect" du modèle :
~ inner | outer (un seul possible) : les outcomes avec des niveaux différents de la variable outer sont considérés indépendants tandis que les outcomes avec un niveau identique de la variable outer sont corrélés suivant un "random effect" correspondant au niveaux de la variable inner.
~ 1 | id. Random effect sur la variable "id" : les outcomes avec des niveaux différents de la variable id sont considérés indépendants, mais pas si le niveau est le même (par exemple si les essais proviennent d'une même étude).
- struct : permet de spécifier la structure de la variance indiquée dans l'argument random ("UN" pour unstructured variance- covariance matrix).