Comment calculer la valeur p d'un rapport de cotes dans R?

J'ai le tableau de valeurs suivant:

25  75
38  162

Le rapport de cotes est de 0,7037 et le log (OR) est de -0,3514. Pour un tableau de contingence avec les valeurs a, b, c et d, la variance de log (OR) est donnée par

(1/a + 1/b + 1/c + 1/d)

Comment puis-je calculer la valeur p de log (OR) à partir de ces données dans R (si elle est significativement différente de 0)?

r variance rnso
la source

Réponses:

Vous pouvez utiliser le test exact de Fisher, qui entre un tableau de contingence et génère une valeur de p, avec une hypothèse nulle que le rapport de cotes est 1 et une autre hypothèse selon laquelle le rapport de cotes n'est pas égal à 1.

(tab <- matrix(c(38, 25, 162, 75), nrow=2))
#      [,1] [,2]
# [1,]   38  162
# [2,]   25   75
fisher.test(tab)
# 
#   Fisher's Exact Test for Count Data
# 
# data:  tab
# p-value = 0.2329
# alternative hypothesis: true odds ratio is not equal to 1
# 95 percent confidence interval:
#  0.3827433 1.3116294
# sample estimates:
# odds ratio 
#  0.7045301

Dans ce cas, la valeur p est de 0,23.

josliber
la source

Merci pour une façon intelligente de déterminer la valeur p. Le test du chi carré peut également être utilisé d'une manière similaire.

rnso

@rnso bien sûr, bien que le test exact de Fisher soit préféré au chi carré lorsque vous avez de petites cellules dans votre tableau de contingence.

josliber

C'est un mythe de longue date, mais ce n'est malheureusement pas vrai. Le Pearson ordinaire fournit des valeurs plus précises que le soi-disant test "exact" de Fisher, même lorsque les fréquences attendues sont aussi basses que 1,0.

χ^{2}

$\chi^2$

P

$P$

Frank Harrell

pourriez-vous en dire un peu plus sur ce @FrankHarrell? Je sais que le serait un résultat asymptotique, alors que le test exact de Fisher repose sur la distribution exacte, comment la valeur est-elle plus "précise" en utilisant la méthode asymptotique?

χ^{2}

$\chi^2$

p

$p$

bdeonovic

Voir de nombreux commentaires à ce sujet sur le site. En bref, les valeurs P du test de Fisher sont trop grandes. L'erreur absolue moyenne dans les valeurs P du test de Pearson est plus petite. Fisher est seulement "exact" dans le sens où les valeurs P sont "garanties" de ne pas être trop petites.

Frank Harrell

Une autre façon de le faire (autre que le test exact de Fisher) est de mettre les valeurs dans un GLM binomial:

d <- data.frame(g=factor(1:2),
                s=c(25,75),
                f=c(38,162))
g <- glm(s/(s+f)~g,weights=s+f,data=d,
    family="binomial")
coef(summary(g))["g2",c("Estimate","Pr(>|z|)")]
##   Estimate   Pr(>|z|) 
## -0.3513979  0.2303337

Pour obtenir le test du rapport de vraisemblance (légèrement plus précis que la valeur Wald illustrée ci-dessus), procédez comme suit: $p$

anova(g,test="Chisq")

qui donne

##      Df Deviance Resid. Df Resid. Dev Pr(>Chi)
## NULL                     1     1.4178         
## g     1   1.4178         0     0.0000   0.2338

(LRT Wald Fisher dans ce cas parce que l'échantillon est assez grand) $p=0.2338 \approx$ $p=0.2303337 \approx$ $p=0.2329$

Ben Bolker
la source

Il vaut mieux généraliser la solution et utiliser le rapport de vraisemblance $\chi^2$ test à partir d'un modèle statistique tel que le modèle logistique. Le test LR fournit assez précis $P$ -valeurs. Cela gère également les cas où vous devez tester plus d'un paramètre, par exemple, des problèmes à 3 groupes, des effets continus non linéaires, etc. Le test LR pour le modèle global (qui est tout ce qui est nécessaire dans cet exemple car il n'y a pas d'ajustement variables) peuvent être facilement obtenus dans la base R ou en utilisant le rmspackage, par exemple

f <- lrm(y ~ groups, weights=freqs)
f  # prints LR chi-sq, d.f., P, many other quantities

Ici, les modèles imbriqués sont ce modèle et un modèle d'interception uniquement.

Frank Harrell
la source

J'ai pu constater que le test LR (lrtest) est utilisé pour comparer les modèles imbriqués. Comment pouvons-nous l'utiliser ici? Pourriez-vous écrire une ligne de code R pour cela?

rnso

car ce que ça vaut, c'est plus ou moins la même approche statistique (bien qu'avec une meilleure explication) que dans ma réponse ci-dessus. lrm()a des valeurs par défaut, des formats de sortie, etc. différents, mais le modèle statistique (IIUC) est le même queglm(...,family="binomial")

Ben Bolker