Comment résumer des intervalles crédibles pour un public médical

21

Avec des forfaits Stan et frontend rstanarmou brmsje peux facilement analyser les données de la bayésien comme je l'ai fait avant avec-modèles mixtes tels que lme. Bien que j'ai la plupart du livre et des articles de Kruschke-Gelman-Wagenmakers-etc. sur mon bureau, ceux-ci ne me disent pas comment résumer les résultats pour un public médical, partagé entre la colère de la Skylla de Bayesian et le Charybdis des examinateurs médicaux ( "nous voulons des significations, pas des choses diffuses").

Un exemple: la fréquence gastrique (1 / min) est mesurée en trois groupes; des contrôles sains sont la référence. Il y a plusieurs mesures pour chaque participant, donc à la fréquentiste j'ai utilisé le modèle mixte suivant lme:

summary(lme(freq_min~ group, random = ~1|study_id, data = mo))

Résultats légèrement modifiés:

Fixed effects: freq_min ~ group 
                   Value Std.Error DF t-value p-value
(Intercept)        2.712    0.0804 70    33.7  0.0000
groupno_symptoms   0.353    0.1180 27     3.0  0.0058
groupwith_symptoms 0.195    0.1174 27     1.7  0.1086

Pour plus de simplicité, j'utiliserai une erreur 2 * std en tant qu'IC à 95%.

Dans un contexte fréquentiste, j'aurais résumé cela comme suit:

  • Dans le groupe témoin, la fréquence estimée était de 2,7 / min (peut-être ajouter CI ici, mais j'évite parfois cela à cause de la confusion créée par CI absolu et différence).
  • Dans le groupe no_symptoms, la fréquence était plus élevée de 0,4 / min, IC (0,11 à 0,59) / min, p = 0,006 que le témoin.
  • Dans le groupe with_symptoms, la fréquence était plus élevée de 0,2 / min, IC (-0,04 à 0,4) / min, p = 0,11 que le contrôle.

Il s'agit de la complexité maximale acceptable pour une publication médicale, le critique me demandera probablement d'ajouter "non significatif" dans le deuxième cas.

Voici la même chose avec les stan_lmerpriors par défaut.

freq_stan = stan_lmer(freq_min~ group + (1|study_id), data = mo)


           contrast lower_CredI frequency upper_CredI
        (Intercept)     2.58322     2.714       2.846
   groupno_symptoms     0.15579     0.346       0.535
 groupwith_symptoms    -0.00382     0.188       0.384

où CredI sont des intervalles crédibles à 90% (voir la vignette rstanarm pourquoi 90% est utilisé par défaut.)

Des questions:

  • Comment traduire le résumé ci-dessus dans le monde bayésien?
  • Dans quelle mesure une discussion préalable est-elle nécessaire? Je suis tout à fait sûr que le papier reviendra avec l '"hypothèse subjective" habituelle lorsque je mentionnerai les prieurs; ou du moins sans "discussion technique, s'il vous plaît". Mais toutes les autorités bayésiennes demandent que l'interprétation ne soit valable que dans le contexte des prieurs.
  • Comment puis-je fournir un substitut de «signification» dans la formulation, sans trahir les concepts bayésiens? Quelque chose comme "crédiblement différent" (uuuh ...) ou presque crédiblement différent (buoha ..., sonne comme "au bord de la signification).

Jonah Gabry et Ben Goodrich (2016). rstanarm: Modélisation de régression appliquée bayésienne via Stan. Version du package R 2.9.0-3. https://CRAN.R-project.org/package=rstanarm

Équipe de développement Stan (2015). Stan: une bibliothèque C ++ pour la probabilité et l'échantillonnage, version 2.8.0. URL http://mc-stan.org/ .

Paul-Christian Buerkner (2016). brms: modèles de régression bayésienne utilisant Stan. Version du package R 0.8.0. https://CRAN.R-project.org/package=brms

Pinheiro J, Bates D, DebRoy S, Sarkar D et R Core Team (2016). nlme: Modèles d'effets mixtes linéaires et non linéaires . Version du package R 3.1-124, http://CRAN.R-project.org/package=nlme>.

Dieter Menne
la source
1
Je n'ai aucune expérience avec les examinateurs / rédacteurs en chef de revues médicales, mais vous pourriez peut-être essayer de dire qu'il n'y a aucune probabilité que l'ordonnée à l'origine soit négative, zéro probabilité que le coefficient de la variable muette "aucun symptôme" soit négatif et une probabilité d'environ 5% que le coefficient de la variable muette "avec symptômes" est négatif. Vous pouvez quantifier environ 5% plus précisément en faisant mean(as.matrix(freq_stan)[,"groupwith_symptoms"] < 0).
Ben Goodrich
Nous y avons pensé et les 5% sonnaient bien; les chercheurs le traduiront par «signification», mais comme ils comprennent mal la signification, ils auront raison par double négation. La «probabilité zéro», en revanche, est un tueur: l'accepteriez-vous? Peut-être que <1 / Reff (p <0,001) serait une approximation? Mais encore une fois: quand j'écris p <xxx, je suis dans le monde de la signification.
Dieter Menne
Corrigez Reff à n_eff ci-dessus.
Dieter Menne
1
Personnellement, je ne parlerais pas d'une probabilité de queue comme ayant "moins de 1 chance sur n_eff" parce que n_eff se rapporte à la précision avec laquelle la moyenne est estimée. Vous pourriez peut-être faire fonctionner vos chaînes assez longtemps pour obtenir 1 tirage négatif pour le coefficient group_nosymptoms, puis dire que la probabilité qu'il soit négatif est 1 / draws. Mais pour l'interception, la chaîne ne va jamais errer en territoire négatif pour ces données, donc je suppose que vous pourriez dire que la probabilité est inférieure à 1 / draws.
Ben Goodrich
J'ai obtenu quelques bons conseils sur l'inclusion de valeurs p pour un expert en domaine mais pas un réviseur expert en statistiques ici: stats.stackexchange.com/questions/148649/… . Nous avons utilisé p <minimum (n_eff de tous les paramètres) comme limite supérieure conservatrice lorsque p = 0.
stijn

Réponses:

16

Réflexions rapides:

1) Le problème clé est la question appliquée à laquelle vous essayez de répondre pour votre public, car cela détermine les informations que vous souhaitez obtenir de votre analyse statistique. Dans ce cas, il me semble que vous souhaitez estimer l'ampleur des différences entre les groupes (ou peut-être l'ampleur des ratios des groupes si c'est la mesure la plus familière à votre public). L'ampleur des différences n'est pas directement fournie par les analyses que vous avez présentées dans la question. Mais il est simple d'obtenir ce que vous voulez de l'analyse bayésienne: vous voulez la distribution postérieure des différences (ou ratios). Ensuite, à partir de la distribution postérieure des différences (ou ratios), vous pouvez faire une déclaration de probabilité directe comme celle-ci:

"Les différences les plus crédibles à 95% se situent entre [limite HDI basse à 95%] et [limite HDI élevée à 95%]" (ici j'utilise l'intervalle de densité la plus élevée à 95% [HDI] comme intervalle crédible, et parce que ce sont par définition des valeurs des paramètres de densité les plus élevées, ils sont considérés comme «les plus crédibles»)

Un public de revues médicales comprendrait intuitivement et correctement cette déclaration, car c'est ce que le public pense généralement que c'est la signification d'un intervalle de confiance fréquentiste (même si ce n'est pas la signification d'un intervalle de confiance fréquentiste).

Comment obtenez-vous les différences (ou ratios) de Stan ou JAGS? Simplement par le post-traitement de la chaîne MCMC terminée. À chaque étape de la chaîne, calculez les différences (ou ratios) pertinents, puis examinez la distribution postérieure des différences (ou ratios). Des exemples sont donnés dans DBDA2E https://sites.google.com/site/doingbayesiandataanalysis/ pour MCMC généralement dans la figure 7.9 (p. 177), pour JAGS dans la figure 8.6 (p. 211) et pour Stan dans la section 16.3 (p . 468), etc.!

2) Si vous êtes obligé par la tradition de déclarer si une différence de zéro est rejetée ou non, vous avez deux options bayésiennes.

2A) Une option consiste à faire des énoncés de probabilité concernant les intervalles proches de zéro et leur relation avec l'IDH. Pour cela, vous définissez une région d'équivalence pratique (CORDE) autour de zéro, qui n'est qu'un seuil de décision approprié pour votre domaine appliqué --- quelle est la différence d'une importance triviale? La définition de telles limites se fait systématiquement dans les tests cliniques de non-infériorité, par exemple. Si vous avez une mesure de «taille d'effet» dans votre champ, il peut y avoir des conventions pour une «petite» taille d'effet, et les limites de CORDE peuvent être, disons, la moitié d'un petit effet. Ensuite, vous pouvez faire des énoncés de probabilité directs tels que ceux-ci:

"Seulement 1,2% de la distribution postérieure des différences est pratiquement équivalent à zéro"

et

"Les différences les plus crédibles à 95% ne sont pas pratiquement toutes équivalentes à zéro (c'est-à-dire que l'IDH et la CORDE à 95% ne se chevauchent pas) et nous rejetons donc zéro." (remarquez la distinction entre l'énoncé de probabilité de la distribution postérieure et la décision ultérieure basée sur cet énoncé)

Vous pouvez également accepter une différence de zéro, à des fins pratiques, si les valeurs les plus crédibles à 95% sont toutes pratiquement équivalentes à zéro.

2B) Une deuxième option bayésienne est le test d'hypothèse nulle bayésienne. (Notez que la méthode ci-dessus n'était pasappelé "test d'hypothèse"!) Le test d'hypothèse nulle bayésienne fait une comparaison du modèle bayésien d'une distribution antérieure qui suppose que la différence ne peut être nulle qu'avec une distribution antérieure alternative qui suppose que la différence pourrait être une gamme diffuse de possibilités. Le résultat d'une telle comparaison de modèle (généralement) dépend très fortement du choix particulier de la distribution alternative, et donc une justification minutieuse doit être faite pour le choix de l'alternative préalable. Il est préférable d'utiliser des a priori au moins modérément informés pour les valeurs nulle et alternative afin que la comparaison des modèles soit réellement significative. Notez que la comparaison du modèle fournit des informations différentes de l'estimation des différences entre les groupes car la comparaison du modèle aborde une question différente. Ainsi, même avec une comparaison de modèles,

Il pourrait y avoir des moyens de faire un test d'hypothèse nulle bayésienne à partir de la sortie Stan / JAGS / MCMC, mais je ne sais pas dans ce cas. Par exemple, on pourrait essayer une approximation de Savage-Dickey à un facteur de Bayes, mais cela reposerait sur la connaissance de la densité antérieure sur les différences, ce qui nécessiterait une analyse mathématique ou une approximation MCMC supplémentaire de l'a priori.

Les deux méthodes pour décider des valeurs nulles sont discutées dans Ch. 12 de DBDA2E https://sites.google.com/site/doingbayesiandataanalysis/ . Mais je ne veux vraiment pas que cette discussion soit mise de côté par un débat sur la "bonne" façon d'évaluer les valeurs nulles; ils sont juste différents et ils fournissent des informations différentes. Le point principal de ma réponse est le point 1 ci-dessus: regardez la distribution postérieure des différences entre les groupes.

John K. Kruschke
la source
3
Bienvenue sur notre site! C'est formidable que vous deveniez membre de notre communauté!
Tim
Si vous souhaitez fusionner votre compte avec celui-ci stats.stackexchange.com/users/16592 (qui semble être le vôtre aussi), vous pouvez le faire automatiquement via stats.stackexchange.com/contact .
Amoeba dit Reinstate Monica
Vous pouvez effectuer le test d'hypothèse décrit ici à l'aide de brms. Voir: github.com/paul-buerkner/brms
bjw
3

Après l'étiquette SO, cela aurait dû être écrit comme un commentaire à @John K. Kruschke, mais les commentaires plus longs sont difficiles à structurer. Pardon.

  • @John K. Kruschke écrit: Simplement par le post-traitement de la chaîne MCMC terminée ...

lower_CredIet upper_CredIdans le post d'origine ont été calculés comme vous l'avez mentionné à partir des chaînes MCMC complètes et ne sont que légèrement reformatés pour une meilleure comparaison avec la lmesortie. Bien que vous préfériez l'IDH, ce sont de simples quantiles; avec le postérieur symétrique dans cet exemple, cela ne fait pas une grande différence.

  • CORDE et taille d'effet

J'ai vu des candidatures aux comités d'éthique où le pouvoir statistique a été calculé sans énoncer l'hypothèse de la taille de l'effet. Même dans le cas où il n'existe aucun moyen de définir un «effet cliniquement pertinent», il est difficile d'expliquer le concept aux chercheurs en médecine. C'est un peu plus facile pour les essais de non-infériorité, mais ceux-ci ne font pas aussi souvent l'objet d'une étude.

Je suis donc tout à fait sûr que l'introduction de CORDES ne sera pas acceptable - une autre hypothèse, les gens ne peuvent pas garder plus d'un numéro à l'esprit. Les facteurs de Bayes pourraient fonctionner, car il n'y a qu'un seul nombre à ramener à la maison comme les valeurs p auparavant.

  • Prieurs

Je suis surpris que ni @John K. Kruschke ni @Ben Goodrich de l'équipe Stan ne mentionnent de prieurs; la plupart des articles sur le sujet demandent une discussion détaillée de la sensibilité préalable lors de la présentation des résultats.

Ce serait bien si dans la prochaine édition de votre livre - avec un peu de chance avec Stan - vous pouviez ajouter des cases "Comment publier ceci (dans un document non statistique) avec 100 mots" pour des exemples sélectionnés. Quand je prendrais votre chapitre 23.1 par mot, un document de recherche médicale typique aurait 100 pages et des chiffres ...

Dieter Menne
la source
* Le point principal était de regarder la distribution postérieure des différences (entre les groupes, entre les combinaisons de groupes). C'est ce qui nécessite le post-traitement de la chaîne MCMC.
John K. Kruschke
* CORDE: Vous "êtes tout à fait sûr que les CORDES ne seront pas acceptables" et "il est difficile d'expliquer le concept aux chercheurs médicaux". Je ne vois pas alors comment les facteurs de Bayes seront plus faciles à expliquer ou à accepter, car un facteur de Bayes nécessite une explication et une justification encore plus élaborées d'un certain seuil BF particulier pour la décision !! Il me semble que vous avez supposé que votre public est ossifié en permanence dans un cadre fréquentiste; si tel est le cas, utilisez simplement les statistiques fréquentistes ou soumettez votre travail à un journal plus éclairé.
John K. Kruschke
* Vous exagérez sévèrement les recommandations du chapitre 23.1, qui peuvent en fait être traitées de manière concise dans une petite quantité de texte, en particulier pour les modèles simples tels que ceux que vous utilisez ici. Suite dans le commentaire suivant ...
John K. Kruschke
1
(i) Motiver l'utilisation du bayésien - il vous donne des distributions postérieures riches en informations. (ii) Expliquez le modèle et ses paramètres, ce qui est facile dans ce cas. (iii) Justifiez l'avant - encore une fois trivial dans ce cas juste pour dire que vous avez utilisé des prieurs diffus qui n'ont essentiellement aucun impact sur le postérieur. (Mais PAS si vous utilisez des facteurs Bayes, pour lesquels la priorité est cruciale.) (Iv) Signalez la fluidité de la chaîne MCMC - trivial de dire que l'ESS était d'environ 10 000 pour tous les paramètres et différences. Suite dans le commentaire suivant ...
John K. Kruschke
1
(v) Interpréter le postérieur: Indiquez simplement la tendance centrale (par exemple le mode) du postérieur et son 95% HDI, pour chaque différence d'intérêt. Ce n'est pas aussi court qu'un tweet, mais ce ne sont que quelques paragraphes.
John K. Kruschke