Les analyses par sous-groupe
Auteur : Nicolas Stefaniak
Niveau : Licence - master - doctorat - chercheurs confirmés
Prérequis
- Avoir installé R, RStudio et le package ‘effectsize’ (Ben-Shachar et al., 2020). Pour rappel, on peut installer ce package et le charger de la manière suivante :
install.packages("metafor")
library("metafor")Introduction
Au chapitre précédent, nous avons vu, avec les analyses d’influences, que certaines études pouvaient causer de l’hétérogénéité. Cependant, nous étions uniquement guidé par les données sans se questionner sur les raisons qui permettent d’expliquer pourquoi ces études pouvaient augmenter l’hétérogénéité de la mesure des tailles d’effet.
Dans ce chapitre, notre attention se portera sur les raisons des variations et sur la manière dont les éta-analyses peuvent être utilisées pour comparer l’effet moyen de différents sous-groupes d’études.
Borenstein (2009) fournit deux illustrations qui rendent limpides les situations où les analyses en sous-groupes apparaissent évidentes.
Considérez que vous souhaitiez mesurer l’effet d’un traitement visant à diminuer le risque de décès chez des patients atteints d’arythmie cardiaque. On peut imaginer que l’effet ne sera pas le même pour les patients qui souffrent de troubles cardiaques de manière chronique par rapport à ceux qui souffrent d’arythmie de manière aigue. Il semble assez naturel et évident de vouloir examiner si l’effet diffère entre les deux.
Imaginez à présent que, parmi les études retenues, la moitié ont été randomisées selon les standards recommandés, ce qui n’est pas le cas pour l’autre moitié. Á nouveau, il est raisonnable de penser que la taille d’effet calculée ne sera pas la même pour les essais contrôlés randomisés par rapport aux études qui ne le sont pas.
On peut aussi réfléchir en termes plus larges : vous voulez montrer que les psychothérapies sont efficaces pour traiter la dépression, vous pouvez néanmoins faire l’hypothèse que toutes les psychothérapies ne se valent pas.
Pour rappel, dans les méta-analyses, les tailles d’effet que nous allons utiliser sont le g de Hedges (Hedges & Olkin, 1985) pour différences entre moyennes, la transformation z de Fisher pour les corrélations (Alexander et al., 1989; Fisher, 1921) et une transformation logarithmique pour les rapports de cotes (Borenstein et al., 2009).
Vouloir tester l’effet différentiel d’une taille d’effet en fonction des sous-groupes consiste à réaliser ce qu’on appelle une analyse de modération (Harrer et al., 2021). Une différence majeure par rapport au chapitre précédent est que les analyses par sous-groupe doivent avoir été décidées a priori. Nous ne sommes plus guides par les données. Évidemment, il existe une multitude raisons qui pourraient expliquer des différences de tailles d’effet. Il est donc nécessaire de se restreindre aux raisons qui sont le plus vraisemblables dans le contexte de la question de recherche que nous nous posons.
L’objectif de ces analyses en sous-groupe est donc double : il s’agit non seulement d’expliquer en quoi certaines caractéristiques d’une étude vont influencer l’estimation de la taille d’effet dans le cadre de notre méta-analyse mais aussi d’ouvrir des pistes de réflexions et de recherche pour les autres chercheurs du domaine.
Pour les praticiens, les analyses en sous-groupes peuvent aussi orienter la prise de décision des praticiens en fonction des caractéristiques favorables ou défavorables des usagers auxquels il va être confronter.
Dans ce chapitre, les modèles fixes seront abordés en premier, ensuite les effets aléatoires avec une estimation combinée ou non du \(\tau^2\)
Pour chaque type de modèles, il existe différentes méthodes pour comparer les tailles d’effets des sous-groupes : le test Z, une analyse de variance basée sur le test Q et le test Q d’hétérogénéité (Borenstein et al., 2009).
Description de l’exemple de travail
Pour illustrer les analyses en sous-groupe, nous nous appuierons sur l’exemple du chapitre 19 de Borenstein (2009). Dans cet exemple, 10 études cherchent à identifier l’impact du tutorat sur les compétences en mathématiques chez des élèves de \(3^{ème}\). Les cinq premières études ont utilisé une méthode A de tutorat et les 5 suivantes, une méthode B. On fait l’hypothèse que les deux méthodes sont efficaces mais que l’une est plus efficace que l’autre.
Les données sont présentées dans le Tableau 1 et le Forest plot de cette analyse est présenté en Figure 1. Elles sont également disponibles dans le fichier ‘Exemple.sous.groupe.xlsx’.
{{% attachments style="grey" pattern=".*\.xlsx$" %}}study | type | g | var |
|---|---|---|---|
Thornhill | A | 0.110 | 0.010 |
Kendall | A | 0.224 | 0.030 |
Vandamm | A | 0.338 | 0.020 |
Leonard | A | 0.451 | 0.015 |
Professor | A | 0.480 | 0.010 |
Jeffries | B | 0.440 | 0.015 |
Fremont | B | 0.492 | 0.020 |
Doyle | B | 0.651 | 0.015 |
Stella | B | 0.710 | 0.025 |
Thorwald | B | 0.740 | 0.012 |
Figure 1. Graphique en forêt des deux types de méthodes de tutorat pour le développement des compétences en mathématique
Les analyses en sous-groupes pour les modèles à effets fixes.
Si nous devions comparer la moyenne de deux groupes d’individus, les options qui s’offriraient à nous seraient de faire un t de Student pour échantillons indépendants ou une analyse de variance.
La méthode suivra une logique assez similaire dans le cas des méta-analyses. Pour pouvoir atteindre cet objectif, il est nécessaire de passer par deux étapes :
- calculer la taille d’effet moyenne de chacun des groupes, ainsi que sa variance ;
- comparer les moyennes de tailles d’effet entre les groupes.
Il est tout à fait possible d’obtenir ces informations directement
dans R avec la fonction rma du package “metafor” (Viechtbauer,
2010) en précisant l’argument subset. Ainsi, on
veut la moyenne et la variance pour le type = A et pour le type = B.
resA <- rma(g, var, data=dat, method="EE", subset=type=="A")
resA##
## Equal-Effects Model (k = 5)
##
## I^2 (total heterogeneity / total variability): 52.56%
## H^2 (total variability / sampling variability): 2.11
##
## Test for Heterogeneity:
## Q(df = 4) = 8.4316, p-val = 0.0770
##
## Model Results:
##
## estimate se zval pval ci.lb ci.ub
## 0.3241 0.0535 6.0633 <.0001 0.2193 0.4289 ***
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Exercice 1
Utilisez la fonction rma pour réaliser l’analyse sur le
tutorat de Type B et sur l’ensemble des données.
Cliquez pour la solution
resB <- rma(g, var, data=dat, method="EE", subset=type=="B")
resB##
## Equal-Effects Model (k = 5)
##
## I^2 (total heterogeneity / total variability): 11.95%
## H^2 (total variability / sampling variability): 1.14
##
## Test for Heterogeneity:
## Q(df = 4) = 4.5429, p-val = 0.3375
##
## Model Results:
##
## estimate se zval pval ci.lb ci.ub
## 0.6111 0.0571 10.7013 <.0001 0.4992 0.7230 ***
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1resT <- rma(g, var, data=dat, method="EE")
resT##
## Equal-Effects Model (k = 10)
##
## I^2 (total heterogeneity / total variability): 65.96%
## H^2 (total variability / sampling variability): 2.94
##
## Test for Heterogeneity:
## Q(df = 9) = 26.4371, p-val = 0.0017
##
## Model Results:
##
## estimate se zval pval ci.lb ci.ub
## 0.4581 0.0390 11.7396 <.0001 0.3816 0.5346 ***
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Á des fins pédagogiques (et à titre de rappel du chapitre sur les effet fixes), nous allons faire le développement du calcul.
On sait que les pondérations correspondent à l’inverse des variances. Nous allon commencer par calculer l’inverse des variances et ajouter la colonne à notre jeu de données.
pondérations<-1/dat$var
dat$pondérations<-pondérationsstudy | type | g | var | pondérations |
|---|---|---|---|---|
Thornhill | A | 0.11 | 0.01 | 100.00 |
Kendall | A | 0.22 | 0.03 | 33.33 |
Vandamm | A | 0.34 | 0.02 | 50.00 |
Leonard | A | 0.45 | 0.01 | 66.67 |
Professor | A | 0.48 | 0.01 | 100.00 |
Jeffries | B | 0.44 | 0.01 | 66.67 |
Fremont | B | 0.49 | 0.02 | 50.00 |
Doyle | B | 0.65 | 0.01 | 66.67 |
Stella | B | 0.71 | 0.02 | 40.00 |
Thorwald | B | 0.74 | 0.01 | 83.33 |
Nous pouvons à présent calculer la moyenne pondérée en multipliant les taille d’effet par les pondérations.
dat$g_x_pond<-dat$g * dat$pondérationsstudy | type | g | var | pondérations | g_x_pond |
|---|---|---|---|---|---|
Thornhill | A | 0.11 | 0.01 | 100.00 | 11.00 |
Kendall | A | 0.22 | 0.03 | 33.33 | 7.47 |
Vandamm | A | 0.34 | 0.02 | 50.00 | 16.90 |
Leonard | A | 0.45 | 0.01 | 66.67 | 30.07 |
Professor | A | 0.48 | 0.01 | 100.00 | 48.00 |
Jeffries | B | 0.44 | 0.01 | 66.67 | 29.33 |
Fremont | B | 0.49 | 0.02 | 50.00 | 24.60 |
Doyle | B | 0.65 | 0.01 | 66.67 | 43.40 |
Stella | B | 0.71 | 0.02 | 40.00 | 28.40 |
Thorwald | B | 0.74 | 0.01 | 83.33 | 61.67 |
Nous pouvons à présent calculer les moyennes pondérées pour chaque groupe en additionnant les tailles d’effet pondérées et en les divisant par la somme des pondérations. En l’occurrence, nous souhaitons le faire par sous-groupe.
Ainsi, la taille d’effet mioyenne pondérée pour le tutorat de type A vaut :
M_a<-sum(dat$g_x_pond[which(dat$type=="A")])/ sum(dat$pondérations[which(dat$type=="A")])
M_a## [1] 0.3240952et sa variance est l’inverse de la somme des pondérations :
var_a<-1/ sum(dat$pondérations[which(dat$type=="A")])
var_a ## [1] 0.002857143On peut à présent calculer l’erreur-type et l’intervalle de confiance autour de la taille d’effet moyenne pondérée.
se_a<-sqrt(var_a )
se_a## [1] 0.05345225M_a+1.96*se_a## [1] 0.4288616M_a-1.96*se_a## [1] 0.2193288On obtient la statistique z en divisant la taille d’effet moyenne pondérée par son erreur-type :
Nous avons ici les résultats principaux de la méta-analyse pour le tutorat de type A.
z_a<-M_a/se_a
z_a## [1] 6.063267En rappel du chapitre consacré à l’hétérogénéité, on peut également obtenir les résultats de la statistique Q
Pour cela, nous devons dans un premier temps calculer le produit des pondérations et des tailles d’effets élevées au carré.
dat$g2_x_pond<-dat$g^2 * dat$pondérationsstudy | type | g | var | pondérations | g_x_pond | g2_x_pond |
|---|---|---|---|---|---|---|
Thornhill | A | 0.11 | 0.01 | 100.00 | 11.00 | 1.21 |
Kendall | A | 0.22 | 0.03 | 33.33 | 7.47 | 1.67 |
Vandamm | A | 0.34 | 0.02 | 50.00 | 16.90 | 5.71 |
Leonard | A | 0.45 | 0.01 | 66.67 | 30.07 | 13.56 |
Professor | A | 0.48 | 0.01 | 100.00 | 48.00 | 23.04 |
Jeffries | B | 0.44 | 0.01 | 66.67 | 29.33 | 12.91 |
Fremont | B | 0.49 | 0.02 | 50.00 | 24.60 | 12.10 |
Doyle | B | 0.65 | 0.01 | 66.67 | 43.40 | 28.25 |
Stella | B | 0.71 | 0.02 | 40.00 | 28.40 | 20.16 |
Thorwald | B | 0.74 | 0.01 | 83.33 | 61.67 | 45.63 |
On peut à présent calculer la statistique Q.
Q_a<-sum(dat$g2_x_pond[which(dat$type=="A")])-(sum(dat$g_x_pond[which(dat$type=="A")])^2 / sum(dat$pondérations[which(dat$type=="A")]))
Q_a## [1] 8.431597Cette statistique suit une distribution \(\chi^2\) ayant comme degrés de liberté le nombre d’études moins 1.
pchisq(Q_a, df=4, lower.tail = F)## [1] 0.07698796Enfin, on peut calculer la statistique \(I^2\) :
I2_a<-100*(Q_a-4)/Q_a
I2_a## [1] 52.5594Nous disposons à présent de tous les résultats de modèle pour le tutorat de type A. Nous pouvons à présent réaliser la même opération pour le tutorat de type B et pour l’ensemble des études.
Exercice 2
Si vous le souhaitez, vous pouvez vous exercez en réalisant les mêmes calculs pour le tutorat de type B et pour l’ensemble du modèle.
Cliquez pour la solution
M_b<-sum(dat$g_x_pond[which(dat$type=="B")])/ sum(dat$pondérations[which(dat$type=="B")])
M_b## [1] 0.611087var_b<-1/ sum(dat$pondérations[which(dat$type=="B")])
var_b ## [1] 0.00326087se_b<-sqrt(var_b )
se_b## [1] 0.05710402M_b+1.96*se_b## [1] 0.7230108M_b-1.96*se_b## [1] 0.4991631z_b<-M_b/se_b
z_b## [1] 10.70129Q_b<-sum(dat$g2_x_pond[which(dat$type=="B")])-(sum(dat$g_x_pond[which(dat$type=="B")])^2 / sum(dat$pondérations[which(dat$type=="B")]))
Q_b## [1] 4.542904Cette statistique suit une distribution \(\chi^2\) ayant comme degrés de liberté le nombre d’études moins 1.
pchisq(Q_b, df=4, lower.tail = F)## [1] 0.3374904Enfin, on peut calculer la statistique \(I^2\) :
I2_b<-100*(Q_b-4)/Q_b
I2_b## [1] 11.9506Et pour le modèle complet :
M<-sum(dat$g_x_pond)/ sum(dat$pondérations)
M## [1] 0.4581218var<-1/ sum(dat$pondérations)
var ## [1] 0.001522843se<-sqrt(var )
se## [1] 0.03902362M+1.96*se## [1] 0.5346081M-1.96*se## [1] 0.3816355z<-M/se
z## [1] 11.7396Q<-sum(dat$g2_x_pond)-(sum(dat$g_x_pond)^2 / sum(dat$pondérations))
Q## [1] 26.43708Cette statistique suit une distribution \(\chi^2\) ayant comme degrés de liberté le nombre d’études moins 1.
pchisq(Q, df=9, lower.tail = F)## [1] 0.00173212Enfin, on peut calculer la statistique \(I^2\) :
I2<-100*(Q-9)/Q
I2## [1] 65.95691Ces résultats sont présentés dans le tableau 5
statistiques | Type.A | Type.B | Global |
|---|---|---|---|
moyenne | 0.3241 | 0.6111 | 0.4581 |
variance | 0.0029 | 0.0033 | 0.0015 |
Erreur-type | 0.0535 | 0.0571 | 0.0390 |
z | 6.0633 | 10.7013 | 11.7396 |
Q | 8.4316 | 4.5429 | 26.4371 |
I2 | 52.5594 | 11.9506 | 65.9569 |
Nous pouvons également faire une représentation graphique de ces résultats avec un graphique en forêt.
Quand on fait les graphique en forêt par groupe, il apparaît que a
première ligne du jeu de données correpond à la dernière ligne du
graphique et inversément. Bref, il y a eu inversion. Il est dès lors
indispensable de le préciser avec l’argument rows de la
fonction forest où on présente les études dans l’ordre
inverse 16à12 et 7 à 3.
Il sera également nécessaire de prévoir des lignes supplémentaires pour indiquer les informations relatives à chacun des sous-groupes.
Concrètement, la méta-analyse contient 10 études. Il faut 10 lignes.
Il faut ajouter à cela une ligne pour présenter les résultats de chacun des sous-groupes (c’est-à-dire 2 lignes). Afin de faire ressortir les les conditions, nous allons les écrires en plus grand. Par ailleurs, pour rendre le graphique clair, il faut créer un espace entre chaque condition. Nous allons donc ajouter 2 fois 3 lignes.
Cela nous fait un total 10 + 2 + 6 = 18. Il faut rajouter la ligne du titre et les titres. C’est pourquoi nous allons indiquer que les limite de y s’étendent entre -1.2 à 21.
Pour faire apparaître les résultats pour chacun des groupes, nous
allons utiliser la fonction addpoly. Les autres éléments du
graphique ont déjà été décrits dans le chapitre consacré aux
effets fixes et aléatoires
Lorsque nous faisons le graphique en forêt (Figure 5), on constate que la taille d’effet pour la méthode A de tutorat avoir un effet plus faible que la méthode B. En effet, la taille d’effet est deux fois plus grande dans ce dernier cas.
resT <- rma(g, var, data=dat, method="EE", slab=study)
resA <- rma(g, var, data=dat, method="EE", subset=type=="A")
resB <- rma(g, var, data=dat, method="EE", subset=type=="B")
par(mar=c(4,4,2,2))
forest(res, xlim=c(-1,2), cex=1, refline=0.5,
header=c("Type / Etude", "SMD [IC à 95%]"),
xlab="Différence de moyennes standardisée", mlab="Combiné",
at=seq(-0.2, 1.2, by=0.2), efac=c(0,1,1.5),
rows=c(16:12, 7:3),
ylim=c(-1.2,21))
text(-1, 17.5, "Type A", pos=4, font=2)
text(-1, 8.5, "Type B", pos=4, font=2)
addpoly(resA, row=10.5, mlab="Combiné")
addpoly(resB, row= 1.5, mlab="Combiné")
Figure 5. Graphique en forêt pour les analyses en
sous-groupe.
Comparer les tailles d’effets par groupe
Plusieurs méthodes existent pour faire une comparaison des tailles d’effet moyenne pondérées. La première consiste à comparer les deux méthodes avec un test Z. La seconde s’appuie sur l’idée d’analyse de variance avec le test Q et la troisième méthode consiste à utiliser le test Q pour tester l’hétérogénéité.
Le test Z
Cette méthode ne s’applique que lorsque nous avons que deux groupes. On obtient le z en utilisant l’équation 1 (Borenstein et al., 2009).
\[Z = \frac{M_A-M_B}{\sqrt{V_{M_A}+V_{M_B} }} \text{ (1)}\] où \(M_A\) est la moyenne pondérée du premier groupe, \(M_B\) est la moyenne pondérée du groupe B tandis que \(V_{M_A}\) et \(V_{M_B}\) renvoie à leur variance.
L’hypothèse nulle est que les deux moyennes pondérées des tailles d’effet sont équivalentes \[H_0 : \theta_A = \theta_B \text{ (2)}\] Tandis que l’hypothèe alternative est qu’il existe une différence entre les deux. \[H_1 : \theta_A \neq \theta_B \text{ (2)}\]
Concrètement les réultats dans notre cas valent :
\[Z = \frac{0.6111-0.3241}{\sqrt{0.0029+0.0033 }}=3.6449 \]
On peut à présent calculer la probabilité associée à cette statistique en s’appuyant sur une distribution normale. Il faut remarquer qu’on multiplie la probabilité obtenue par 2 pour que la différence soit examinée de manière bilatérale.
2*pnorm( 3.6449, lower.tail=F)## [1] 0.0002674958Ainsi, nous pouvons rejeter l’hypothèse selon laquelle les deux tailles d’effet sont égales.
Utiliser le test Q en s’appuyant sur les principes de décomposition de la variance
Lorsqu’il y a plus que deux groupes, la méthode que nous venons de décrire ne fonctionne plus. Il est alors nécesaire de partionner la variance totale entre une variance entre les sous-groupe et une variance à l’intérieur de chaque sous-groupe.
Concrètement, il va falloir calculer :
- \(Q_A\) : la somme des carrés pondérée de toutes les études utilisant la méthode A par rapport à la moyenne de A.
- \(Q_B\) : la somme des carrés pondérée de toutes les études utilisant la méthode B par rapport à la moyenne de B.
- \(Q_{intra}\) : la somme de \(Q_A\) et \(Q_B\).
- \(Q_{bet}\) : la somme des carrés pondérées des moyennes des sous-groupes à la grande moyenne.
- \(Q\) : la somme des carrés pondérées de tous les effets à la grande moyenne.
Nous avons déjà calculé \(Q_A\) et \(Q_B\). Nous pouvons calculer facilement \(Q_{intra}\) en faisant la somme.
Q_within<-Q_a+Q_b
Q_within## [1] 12.9745Sachant que nous avons également déjà calculé le Q, et que les sommes des carrés s’additionnent, on peut obtenir aisément le \(Q_{bet}\) en soustrayant le \(Q_{intra}\) du Q total.
Q_bet<-Q-Q_within
Q_bet## [1] 13.46258Pour nous rapprocher le plus possible d’un table d’ANOVA, nous allons présenter les résultats sous la forme d’un tableau (voir Tableau 6)
Effet | Q | ddl | p |
|---|---|---|---|
A | 8.4316 | 4.0000 | 0.0770 |
B | 4.5429 | 4.0000 | 0.3375 |
intra | 12.9745 | 8.0000 | 0.1127 |
inter | 13.4626 | 1.0000 | 0.0002 |
total | 26.4371 | 9.0000 | 0.0017 |
Ce tableau nous indique quels sont les effets pour lesquels le Q est significatif est ceux pour lesquels ce n’est pas le cas. Il apparaît qu’il y a un effet global nous amenant à rejeter l’hypothèse nulle selon laquelle la distribution des tailles d’effet est homogène (la dernière ligne du tableau) mais on constate que cet effet est essentiellement dû à la variance entre les condition pour la \(Q_{inter}\) nous amène à rejeter l’hypothèse nulle selon laquelle les tailles d’effet moyennes pondérées entre les groupes est homogène.
Il est à noter que 13.46 est le carré de la statistique Z (3.6449) que nous venons de calculer, tout comme c’est le cas entre le t et le F dans une anova où il n’y a que deux groupes.
Le test Q pour tester l’hétérogénéité
La dernière méthode consiste à faire une méta-analyse sur les moyennes pondérées de chacun des sous-groupes et de calculer la statistique Q pour identifier s’il y a hétérogénéité entre les sous-groupes.
Ainsi, cela nous amène à faire une méta-analyse sur deux données ayant comme valeur la moyenne et la variance que nous présentée dans le tableau 5. Nous ne développerons pas à nouveau les calculs qui ont été présentés plus tôt dans le chapitre.
df<-data.frame(si = c("Type A", "Type B"),
yi = c(results[1,2],results[1,3] ),
vi =c(results[2,2],results[2,3] ))
df## si yi vi
## 1 Type A 0.3240952 0.002857143
## 2 Type B 0.6110870 0.003260870df<-data.frame(si = c("Type A", "Type B"),
yi = c(results[1,2],results[1,3] ),
vi =c(results[2,2],results[2,3] ))
rma(yi, vi, method="EE", slab=si, data = df)##
## Equal-Effects Model (k = 2)
##
## I^2 (total heterogeneity / total variability): 92.57%
## H^2 (total variability / sampling variability): 13.46
##
## Test for Heterogeneity:
## Q(df = 1) = 13.4626, p-val = 0.0002
##
## Model Results:
##
## estimate se zval pval ci.lb ci.ub
## 0.4581 0.0390 11.7396 <.0001 0.3816 0.5346 ***
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Ainsi, nous pouvons rejeter l’hypothèse d’homogénéité des tailles d’effets entre les deux sous-groupes.
Il faut remarquer que pour des effets fixes, cette méthode renvoie strictement les mêmes résultats que ceux de la méthode précédente.
Tester des différences entre sous-groupes pour des modèles fixes avec R.
Évidemment, il n’est pas nécessaire de calculer manuellement ces valeurs. On les obtient assez aisément dans R.
Nous obtenons ces résultats avec la fonction rma.mv du
package ‘metafor’ (Viechtbauer, 2010).
Le principe de la fonction va être assez similaire à l’utilisation de
la fonction rma que nous avons déjà largement utilisée. Il
faut néanmoins expliquer quelques arguments supplémentaires :
mods: Nous avons vu que les analyses en sous-groupes représente des analyses de modération. On va donc préciser que la variable “type” est un modérateur de la manière suivante :mods=~type1.random: on doit préciser les effets aléatoires. Habituellement, dans une analyse de variance, l’effet aléatoire est les participants. Dans beaucoup de logiciels, c’est sous-entendu et il n’est pas nécessaire de le préciser. Étant donné que, dans cette fonction, il est possible de préciser les effets aléatoires (voir le chapitre sur les modèles linéaires mixtes pour plus d’informations), il est nécessaire de préciser l’effet aléatoire. On préciser la partie aléatoire en la faisant précéder par une barre vertical (un pipe). Ainsi,random = ~ type | studysignifie que les études représentent l’effet aléatoire et que ces études appartiennent à une structure supérieure qu’est le type de tutorat.struct: cet argument permet de spécifier la structure de la matrice de variance/covariance. Il ne semble pas opportun de rentrer dans les détails ici puisque ce qui nous intéresse est un effet entre les groupes. Dans une matrice de variance covariance, les variances correspondent à la diagonale2, on va donc préciser que la structure eststr=“DIAG”.tau2: cet argument permet de spécifier la manière dont le \(\tau^2\) doit être estimé. Si vous ne spécifier rien, la fonction fera une estimation de \(\tau^2\). Cependant, pour reproduire exactement les résultats que nous avons décrits précédemment, nous allons préciser ces valeurs. Ainsi, on crée un vecteur avec les \(tau^2\) des deux modèles que nou avons estimés précédemment.
dif <- rma.mv(g, var, mods = ~ type, random = ~ type | study, struct="DIAG",
data=dat, slab=study, tau2=c(resA$tau2, resB$tau2))
dif##
## Multivariate Meta-Analysis Model (k = 10; method: REML)
##
## Variance Components:
##
## outer factor: study (nlvls = 10)
## inner factor: type (nlvls = 2)
##
## estim sqrt k.lvl fixed level
## tau^2.1 0.0000 0.0000 5 yes A
## tau^2.2 0.0000 0.0000 5 yes B
##
## Test for Residual Heterogeneity:
## QE(df = 8) = 12.9745, p-val = 0.1127
##
## Test of Moderators (coefficient 2):
## QM(df = 1) = 13.4626, p-val = 0.0002
##
## Model Results:
##
## estimate se zval pval ci.lb ci.ub
## intrcpt 0.3241 0.0535 6.0633 <.0001 0.2193 0.4289 ***
## typeB 0.2870 0.0782 3.6691 0.0002 0.1337 0.4403 ***
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Exercice 3
Des chercheurs se sont intéressés à l’efficacité du vaccin BCG contre la tuberculose. Treize études ont été inclues. Les données peuvent être chargées dans la mémoire de R en utilisant la commande
data(dat.bcg)Les données sont présentées ci-dessous.
trial | author | year | tpos | tneg | cpos | cneg | ablat | alloc |
|---|---|---|---|---|---|---|---|---|
1 | Aronson | 1,948 | 4 | 119 | 11 | 128 | 44 | random |
2 | Ferguson & Simes | 1,949 | 6 | 300 | 29 | 274 | 55 | random |
3 | Rosenthal et al | 1,960 | 3 | 228 | 11 | 209 | 42 | random |
4 | Hart & Sutherland | 1,977 | 62 | 13,536 | 248 | 12,619 | 52 | random |
5 | Frimodt-Moller et al | 1,973 | 33 | 5,036 | 47 | 5,761 | 13 | alternate |
6 | Stein & Aronson | 1,953 | 180 | 1,361 | 372 | 1,079 | 44 | alternate |
7 | Vandiviere et al | 1,973 | 8 | 2,537 | 10 | 619 | 19 | random |
8 | TPT Madras | 1,980 | 505 | 87,886 | 499 | 87,892 | 13 | random |
9 | Coetzee & Berjak | 1,968 | 29 | 7,470 | 45 | 7,232 | 27 | random |
10 | Rosenthal et al | 1,961 | 17 | 1,699 | 65 | 1,600 | 42 | systematic |
11 | Comstock et al | 1,974 | 186 | 50,448 | 141 | 27,197 | 18 | systematic |
12 | Comstock & Webster | 1,969 | 5 | 2,493 | 3 | 2,338 | 33 | systematic |
13 | Comstock et al | 1,976 | 27 | 16,886 | 29 | 17,825 | 33 | systematic |
Dans ce tableau, les variables correspondent à :
- trial : le numéro de l’essai
- author : le auteurs de l’étude
- year : l’année de réalisation de l’étude
- tpos : le nombre de personnes positives à tuberculose chez les personnes vaccinées
- tneg : le nombre de personnes négatives à tuberculose chez les personnes vaccinées
- cpos : le nombre de personnes positives à tuberculose dans le groupe de contrôle
- cneg : le nombre de personnes négatives à tuberculose dans le groupe de contrôle
- alloc : variable qualitative indiquant si l’administration du vaccin a été fait de manière aléatoire (random), en alternance (alternate) ou si l’administration suivait un ordre systématique (systematic).
Faites une méta-analyse permettant d’identifier si le mode d’administration du vaccin influence la taille d’effet. Pour cette méta-analyse, utiliser un modèle à effet fixes.
Faites une représentation graphique de ces effets. Le graphique doit faire apparaître :
- le nom des études (et l’année)
- la moyenne pondérée de la taille d’effet pour l’ensemble des études mais également pour chaque sous-groupe
- les indices d’hétérogénéité
- le poids des études
Attention à vérifier vos résultats.
Cliquez pour la solution
On commence par calculer les tailles d’effet de sorte à ce que ce soit utilisable par la fonction <coderma :
dat <- escalc(measure="OR", ai=tpos, bi=tneg, ci=cpos, di=cneg, data=dat.bcg)Ensuite, et en prévision du forest plot, où on veut que les études
d’un même traitement se suivent dans le jeu de données, nous allons
trier les données. On peut le faire avec la fonction
arrange du package ‘dplyr’ (Wickham et al., 2023).
require(dplyr)
dat<-dat%>%arrange(alloc)On continue en réalisant l’analyse sur l’ensemble des données.
rma.out <- rma(yi = yi, # la colonne dans le jeu de données où on trouve la taille d'effet
vi = vi, # la colonne où on trouve l'erreur de mesure
data=dat, # le nom du jeu de données
method="EE", # utiliser la méthode de Paule-Mandel
digits=2) # le nombre de décimale à imprimer dans la sortie de résultats
summary(rma.out)##
## Equal-Effects Model (k = 13)
##
## logLik deviance AIC BIC AICc
## -76.03 163.16 154.06 154.62 154.42
##
## I^2 (total heterogeneity / total variability): 92.65%
## H^2 (total variability / sampling variability): 13.60
##
## Test for Heterogeneity:
## Q(df = 12) = 163.16, p-val < .01
##
## Model Results:
##
## estimate se zval pval ci.lb ci.ub
## -0.44 0.04 -10.32 <.01 -0.52 -0.35 ***
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1On calcule ensuite l’intervalle de prédiction et on obtient les résultats sous une forme plus interprétable que sous une forme logarithmique.
predict(object = rma.out, transf=exp) ##
## pred ci.lb ci.ub
## 0.65 0.60 0.70On répète l’opération pour chacun des sous-groupes.
rma.out.random <- rma(yi = yi, # la colonne dans le jeu de données où on trouve la taille d'effet
vi = vi, # la colonne où on trouve l'erreur de mesure
data=dat, # le nom du jeu de données
method="EE", # utiliser la méthode de Paule-Mandel
digits=2,# le nombre de décimale à imprimer dans la sortie de résultats
subset=alloc=="random")
summary(rma.out.random)##
## Equal-Effects Model (k = 7)
##
## logLik deviance AIC BIC AICc
## -53.49 111.31 108.98 108.92 109.78
##
## I^2 (total heterogeneity / total variability): 94.61%
## H^2 (total variability / sampling variability): 18.55
##
## Test for Heterogeneity:
## Q(df = 6) = 111.31, p-val < .01
##
## Model Results:
##
## estimate se zval pval ci.lb ci.ub
## -0.30 0.05 -5.39 <.01 -0.40 -0.19 ***
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1predict(object = rma.out.random, transf=exp) ##
## pred ci.lb ci.ub
## 0.74 0.67 0.83rma.out.alternate <- rma(yi = yi, # la colonne dans le jeu de données où on trouve la taille d'effet
vi = vi, # la colonne où on trouve l'erreur de mesure
data=dat, # le nom du jeu de données
method="EE", # utiliser la méthode de Paule-Mandel
digits=2,# le nombre de décimale à imprimer dans la sortie de résultats
subset=alloc=="alternate")
summary(rma.out.alternate)##
## Equal-Effects Model (k = 2)
##
## logLik deviance AIC BIC AICc
## -2.47 8.83 6.93 5.63 10.93
##
## I^2 (total heterogeneity / total variability): 88.67%
## H^2 (total variability / sampling variability): 8.83
##
## Test for Heterogeneity:
## Q(df = 1) = 8.83, p-val < .01
##
## Model Results:
##
## estimate se zval pval ci.lb ci.ub
## -0.84 0.09 -9.21 <.01 -1.02 -0.66 ***
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1predict(object = rma.out.alternate, transf=exp) ##
## pred ci.lb ci.ub
## 0.43 0.36 0.52rma.out.systematic<- rma(yi = yi, # la colonne dans le jeu de données où on trouve la taille d'effet
vi = vi, # la colonne où on trouve l'erreur de mesure
data=dat, # le nom du jeu de données
method="EE", # utiliser la méthode de Paule-Mandel
digits=2,# le nombre de décimale à imprimer dans la sortie de résultats
subset=alloc=="systematic")
summary(rma.out.systematic)##
## Equal-Effects Model (k = 4)
##
## logLik deviance AIC BIC AICc
## -7.03 16.93 16.06 15.44 18.06
##
## I^2 (total heterogeneity / total variability): 82.28%
## H^2 (total variability / sampling variability): 5.64
##
## Test for Heterogeneity:
## Q(df = 3) = 16.93, p-val < .01
##
## Model Results:
##
## estimate se zval pval ci.lb ci.ub
## -0.41 0.10 -4.33 <.01 -0.60 -0.23 ***
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1predict(object = rma.out.systematic, transf=exp) ##
## pred ci.lb ci.ub
## 0.66 0.55 0.80Nous allons commencer par tester si l’effet moyen pondéré dépend de la méthode d’administration du traitement.
slab<-paste0(dat.bcg$author," (",dat.bcg$year, ")")
dif <- rma.mv(yi, vi, mods = ~ alloc, random = ~ alloc | trial, struct="DIAG",
data=dat, slab= slab)
dif##
## Multivariate Meta-Analysis Model (k = 13; method: REML)
##
## Variance Components:
##
## outer factor: trial (nlvls = 13)
## inner factor: alloc (nlvls = 3)
##
## estim sqrt k.lvl fixed level
## tau^2.1 0.2421 0.4921 2 no alternate
## tau^2.2 0.4131 0.6427 7 no random
## tau^2.3 0.4200 0.6480 4 no systematic
##
## Test for Residual Heterogeneity:
## QE(df = 10) = 137.0728, p-val < .0001
##
## Test of Moderators (coefficients 2:3):
## QM(df = 2) = 1.6402, p-val = 0.4404
##
## Model Results:
##
## estimate se zval pval ci.lb ci.ub
## intrcpt -0.6171 0.3684 -1.6751 0.0939 -1.3391 0.1050 .
## allocrandom -0.3783 0.4646 -0.8143 0.4155 -1.2888 0.5322
## allocsystematic 0.1876 0.5202 0.3606 0.7184 -0.8320 1.2072
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Notez que nous devons pas mettre le \(\tau^2\) puisque nous avons modélisé la méta-analyse avec un effet fixe. Cette valeur n’a donc pas été calculée.
par(mar=c(4,4,2,2))
dat$weights <- paste0(round(weights(rma.out)), "%") # weights in % (rounded)
hetéro.text<-paste0("I² =",
round(rma.out$I2,3),
", Q =", round(rma.out$QE,1),
", p ",
ifelse(round(rma.out$QEp, 3)==0, "< .001", paste0("= ",round(rma.out$QEp, 3) )))
forest(rma.out, xlim=c(-4,3.5), cex=.5, refline=0.5,
header=c("Etude", "OR [IC à 95%]"),
xlab="Rapport de cote", mlab="Combiné",
ilab=dat$weights, ilab.xpos=2.6, ilab.pos=2.6,
#at=seq(-0.2, 1.2, by=0.2), efac=c(0,1,1.5),
rows=c(22:21, 17:11, 7:4), slab=slab,
ylim=c(-1.2,27))
text(-4, 23.5, substitute(paste(bold('Alternance'))), , pos=4, font=1, cex=.5)
text(-4, 18.5, substitute(paste(bold('Aléatoire'))), , pos=4, font=1, cex=.5)
text(-4, 8.5, substitute(paste(bold('Systématique'))), , pos=4, font=1, cex=.5)
text(2.6, -1, "100%", pos=2, cex = .5)
text(2.6, 26 , "Weight", pos=2, font=1, cex = .5)
text(-3, -1 ,hetéro.text , pos=4, font=1, cex =0.5)
addpoly(rma.out.random, row=20, mlab="Combiné")
addpoly(rma.out.alternate, row= 10, mlab="Combiné")
addpoly(rma.out.systematic, row= 2, mlab="Combiné")
Figure 3. Graphique en forêt de l’efficacité du vaccin BCG en fonction du type d’administration
Les résultats de la méta-analyse montrent que le vaccin BCG est globalement efficace pour réduire le risque de tuberculose. En effet, les personnes vaccinées ont 1.54 moins de chances de développer une tuberculose que les personnes qui ne le sont pas. Si les indices d’hétérogénéité, I² =92.645, Q =163.2, p < .001indique la préence d’une hétérogénétié, les analyses de modération ne permettent pas d’attribuer cette hétérogénéité à des différences méthodologiques relatives au mode d’administration du vaccin. En effet, bien les bénéfices du vaccin (\(OR_{aléatoire}=.74\),\(OR_{alternance}=.43\),\(OR_{systématique}=.66\)) varient sensiblement entre les différents modes d’administration, nous devons tolérer l’hypothèse nulle selon laquelle l’efficacité du vaccin dépend du type d’études qui a été réalisée, \(Q_M\) (2)=1.64, p =0.44.
Les analyses en sous-groupes pour les modèles à effets aléatoires.
Le raisonnement qui prévalait dans le chapitre consacré aux effets fixes et aléatoires vaut également lorsque nous réalisons des analyses de modération.
Concrètement, quand nous avons utilisé des modèles à effets fixes dans la section précédente, nous avons fait l’hypothèse que les tailles d’effet étaient homogènes à l’intérieur de chaque type de tutorat. Dans les modèles aléatoires, l’hypothèse est qu’à l’intérieur de chaque mode de tutorat persistera une certaine variabilité, ce qui pourrait par exemple se justifier par le fait que tout le monde n’a pas les mêmes compétences pédagogiques, et que les apprenants n’ont pas tous les mêmes prédispositions au départ.
On identifie que les modèles aléatoires restent le choix de prédilection dans la plupart des situations en psychologie.
Lorsque nous utilisons des modèles à effets aléatoires, la manière de calculer les statistiques ainsi que les conséquences liées à l’utilisation de modèles aléatoires sont les mêmes que celles qui ont été présentées dans les chapitres précédents. Ainsi, il va être nécessaire de calculer la valeur \(\tau2\). Pour rappel, le \(\tau^2\) renvoie à la variance des vraies taille d’effet pour un ensemble d’études. Le point qui est critique ici est la manière de définir l’ensemble d’études dont nous venons de parler.
Concrètement, deux choix s’offrent à nous : dans le premier cas, on calcule \(\tau^2\) sur l’ensemble des études de la méta-analyse, qu’importe le sous-groupe auquel appartient une étude ; dans le second cas, on caclule \(\tau^2\) pour chacun des sous-groupes. On comprend ici que Dans le premier cas, la dispersion sera plus large alors que dans le second cas, elle tendra à être plus petite, du moins si chacun des sets d’études représentent des clustes qu’on peut clairement identifier.
Dans l’exemple qui a été développé depuis le début de ce chapitre, on veut calculer la moyenne du sous-groupe A, celle du sous-groupe B ainsi que leur distribution respective. Il apparaît donc clairement que \(\tau^2\) doit être estimé pour chacun des sous-groupes. Ceci s’explique par le fait que nous voulons pouvoir séparer la part de variance qui explique la différence de tailles d’effet moyennes entre les sous-groupe. Le fait de calculer \(\tau^2\) à l’intérieur de chaque sous-groupe fait en sorte qu’on ne peut pas attribuer à la méthode utilisée une mesure plus large de la dispersion étant donné qu’il n’y a qu’une seule méthode utilisée à l’intérieur de chacun des sous-groupes.
Faut-il ou non combiner les valeurs de \(\tau^2\) ?
Pour rappel, la formule du \(\tau2\) correspond à l’équation 3 \[T^2=\frac{Q-ddl}{C} \text{ (3)}\]
où C vaut :
\[C=\sum W_i - \frac{\sum W^2_i}{\sum W_i} \text{ (4)}\] Pour estimer le \(\tau2\) de chacun des sous-groupes, nous devons commencer par calculer C en premier.
On commence par calculer le carré des pondérations.
dat$W2<-dat$pondérations^2study | type | g | var | pondérations | g_x_pond | W2 |
|---|---|---|---|---|---|---|
Thornhill | A | 0.110 | 0.010 | 100.00000 | 11.000000 | 10,000.000 |
Kendall | A | 0.224 | 0.030 | 33.33333 | 7.466667 | 1,111.111 |
Vandamm | A | 0.338 | 0.020 | 50.00000 | 16.900000 | 2,500.000 |
Leonard | A | 0.451 | 0.015 | 66.66667 | 30.066667 | 4,444.444 |
Professor | A | 0.480 | 0.010 | 100.00000 | 48.000000 | 10,000.000 |
Jeffries | B | 0.440 | 0.015 | 66.66667 | 29.333333 | 4,444.444 |
Fremont | B | 0.492 | 0.020 | 50.00000 | 24.600000 | 2,500.000 |
Doyle | B | 0.651 | 0.015 | 66.66667 | 43.400000 | 4,444.444 |
Stella | B | 0.710 | 0.025 | 40.00000 | 28.400000 | 1,600.000 |
Thorwald | B | 0.740 | 0.012 | 83.33333 | 61.666667 | 6,944.444 |
Pour le traitement A, on obtient \(\tau2\) de la manière suivante :
\[C=(100+33.33+50+66.67+100) - \frac{(10000.000 + 1111.111 + 2500.000 + 4444.444 + 10000.000)}{(100+33.33+50+66.67+100)}=269.8413 \]
Nous pouvons à présent calculer \(\tau^2\) :
\[\tau^2 = \frac{8.43-4}{269.8413} = 0.01642297\] En appliquant la même procédure pour le groupe B, nous obtenons 0.002. Il faut à présent choisir si nous combinons ou non ces estimations de \(\tau^2\). si nous ne le faisons pas, l’aternative est d’appliquer l’estimation de \(\tau^2\) de chaque sous-groupe à ce sous-groupe.
Il faut néanmoins noter que ce n’est pas les \(\tau2\) en tant que tels qui sont regroupés mais les estimations de Q, des degrés de liberté et des C.
La décision de regrouper ou non les valeurs de \(\tau^2\) va dépendre de la capacité que nous avons à pouvoirf faire l’hypothèse que les différences observées dans les tailles d’effet entre les études de chacun des sous-groupe est identique ou non. Si c’est le cas, alors, nous pouvons regrouper les valeurs de \(\tau^2\). pour obtenir une estimation commune, qui sera ensuite utilisée pour l’ensemble des sous-groupes. Dans l’exemple sur le tutorat, cela semble être une hypothèse raisonnable.
Á l’inverse, si nous faisons l’hypothèse que les variations à l’intérieur de chacun des sous-groupes ne sont pas les mêmes, alors il est préférable d’utiliser 2 estimations de \(\tau^2\). Cela peut survenir s’il existe une plus grande variabilité dans les populations ou les outils utilisés dans un groupe que dans un autre.
Le problème quand on veut faire une estimation de \(\tau^2\) pour chacun des sous-groupes est qu’il peut arriver que certains sous-groupes ne contiennent que peu d’études (<5). La conséquence dans ce cas est l’estimation de \(\tau^2\) sera très imprécise. Ainsi, quand on décide de faire des estimations séparées par sous-groupes, on doit pouvoir s’assurer qu’il y a assez d’études à l’intérieur de chaque sous-groupe pour atteindre une estimation raisonnablement précise de \(\tau^2\). Dans ce cas, il est plus raisonnable de faire une estimation combinée.
En résumé, si vous pouvez faire l’hypothèse que les études à l’intérieur de chaque sous-groupe partagent une même taille d’effet, on utilise un modèle à effet fixe, sinon un modèle à effet aléatoire. Dans ce dernier cas, on utilise un \(\tau^2\) combiné si on fait l’hypothèse que la variance à l’intérieur de chaque sous-groupe est identique, sinon on utilise des estimations de variance séparées, sauf si le nombre d’études ne le permet pas.
Les modèles à effets aléatoires avec une estimation séparée de \(\tau^2\)
La différence majeure avec les modèles à effets fixes est que les pondérations vont prendre en compte à la fois la variance de chaque étude et celle entre les études. Nous allons donc repartir du jeu de données initial mais en y ajoutant les informations relative à la variances entre les études (var_bet) à la variance spécifique à chaque étude (var_with) et ainsi que la variance totale, qui est la somme des deux précédentes.
Etude | Type | Taille d'effet | Variance intra (Vy) | Variance inter (tau²) | Variance totale |
|---|---|---|---|---|---|
Thornhill | A | 0.110 | 0.010 | 0.0164 | 0.0264 |
Kendall | A | 0.224 | 0.030 | 0.0164 | 0.0464 |
Vandamm | A | 0.338 | 0.020 | 0.0164 | 0.0364 |
Leonard | A | 0.451 | 0.015 | 0.0164 | 0.0314 |
Professor | A | 0.480 | 0.010 | 0.0164 | 0.0264 |
Jeffries | B | 0.440 | 0.015 | 0.0022 | 0.0172 |
Fremont | B | 0.492 | 0.020 | 0.0022 | 0.0222 |
Doyle | B | 0.651 | 0.015 | 0.0022 | 0.0172 |
Stella | B | 0.710 | 0.025 | 0.0022 | 0.0272 |
Thorwald | B | 0.740 | 0.012 | 0.0022 | 0.0142 |
Nous pouvons à présent calculer les pondérations, les pondérations au carré et le produit entre les tailles d’effet et les pondérations Notez que les données sont stockés dans un objet R qui n’a pas le même nom des variables que les variables dans le tableau pour des raisons purement pratique : R n’aime pas les caractères spéciaux dans les noms des variables. Néanmoins, leur nom est relativement transparent. Remarquez simplement que les codes sont strictement identiques à ce qu’on a utilisé précédemment pour les modèles à effets fixes.
Etude | Type | Taille d'effet | Variance intra (Vy) | Variance inter (tau²) | Variance totale | Pondérations (W) | Pondérations au carré | Produit entre la taille d'effet et la pondération | Taille d'effet au carré multiplié par la pondération |
|---|---|---|---|---|---|---|---|---|---|
Thornhill | A | 0.1100 | 0.0100 | 0.0164 | 0.0264 | 37.8788 | 1,434.8026 | 4.1667 | 0.4583 |
Kendall | A | 0.2240 | 0.0300 | 0.0164 | 0.0464 | 21.5517 | 464.4768 | 4.8276 | 1.0814 |
Vandamm | A | 0.3380 | 0.0200 | 0.0164 | 0.0364 | 27.4725 | 754.7398 | 9.2857 | 3.1386 |
Leonard | A | 0.4510 | 0.0150 | 0.0164 | 0.0314 | 31.8471 | 1,014.2399 | 14.3631 | 6.4777 |
Professor | A | 0.4800 | 0.0100 | 0.0164 | 0.0264 | 37.8788 | 1,434.8026 | 18.1818 | 8.7273 |
Jeffries | B | 0.4400 | 0.0150 | 0.0022 | 0.0172 | 58.1395 | 3,380.2055 | 25.5814 | 11.2558 |
Fremont | B | 0.4920 | 0.0200 | 0.0022 | 0.0222 | 45.0450 | 2,029.0561 | 22.1622 | 10.9038 |
Doyle | B | 0.6510 | 0.0150 | 0.0022 | 0.0172 | 58.1395 | 3,380.2055 | 37.8488 | 24.6396 |
Stella | B | 0.7100 | 0.0250 | 0.0022 | 0.0272 | 36.7647 | 1,351.6436 | 26.1029 | 18.5331 |
Thorwald | B | 0.7400 | 0.0120 | 0.0022 | 0.0142 | 70.4225 | 4,959.3335 | 52.1127 | 38.5634 |
On peut à présent calculer l’ensemble des indicateurs que nous avons calculé précédemment pour les modèles à effets fixes. Ainsi, on commence par la taille d’effet moyenne pondérée pour le groupe A et pour le groupe B
M_a. = sum(dat$WY[1:5])/sum(dat$W[1:5]) #le point indique qu'il s'agit de la taille d'effet moyenne pondérée pour un modèle aléatoire
M_a.## [1] 0.324492M_b. = sum(dat$WY[6:10])/sum(dat$W[6:10]) #le point indique qu'il s'agit de la taille d'effet moyenne pondérée pour un modèle aléatoire
M_b.## [1] 0.6100599Continue avec leur variance respective :
V_M_a. = 1/sum(dat$W[1:5]) #le point indique qu'il s'agit de la taille d'effet moyenne pondérée pour un modèle aléatoire
V_M_a.## [1] 0.006384515V_M_b. = 1/sum(dat$W[6:10]) #le point indique qu'il s'agit de la taille d'effet moyenne pondérée pour un modèle aléatoire
V_M_b.## [1] 0.003724237Leur erreur-type est la racine carrée des variances que nous venons de calculer :
se_Ma.<-sqrt(V_M_a.)
se_Ma.## [1] 0.07990316se_Mb.<-sqrt(V_M_b.)
se_Mb.## [1] 0.06102653On peut donc calculer les intervalles de confiance pour chaque taille d’effet.
M_a.-1.96*se_Ma.# limite inférieure de l'intervalle de confiance## [1] 0.1678818M_a.+1.96*se_Ma.# limite supérieure de l'intervalle de confiance## [1] 0.4811022M_b.-1.96*se_Mb.## [1] 0.4904479M_b.+1.96*se_Mb.## [1] 0.7296719On peut également calculer la statistique z en faisant le rapport entre la taille d’effet moyenne pondérée et son erreur-type :
Z_Ma.<-M_a./se_Ma.
Z_Ma.## [1] 4.061066Z_Mb.<-M_b./se_Mb.
Z_Mb.## [1] 9.996634Ces deux valeurs amènent à rejeter l’hypothèse nulle selon laquelle la taille d’effet vaut 0.
Il ne reste plus qu’à calculer la statistique Q pour chaque groupe, non pas pour tester l’hétérogénéité mais pour partionner la variance entre ces différentes composantes.
Q_a.<-sum(dat$WY2[1:5])-(sum(dat$WY[1:5])^2/sum(dat$W[1:5]))
Q_a.## [1] 3.391042Q_b.<-sum(dat$WY2[6:10])-(sum(dat$WY[6:10])^2/sum(dat$W[6:10]))
Q_b.## [1] 3.962959