Comment les Bayésiens vérifient-ils leurs méthodes en utilisant des méthodes de simulation Monte Carlo?

Contexte : J'ai un doctorat en psychologie sociale, où les statistiques théoriques et les mathématiques étaient à peine couvertes dans mes cours quantitatifs. Au cours des études de premier cycle et des cycles supérieurs, on m'a enseigné (un peu comme beaucoup d'entre vous aussi dans les sciences sociales, probablement) à travers le cadre fréquentiste "classique". Maintenant, j'aime aussi R et utiliser des méthodes de simulation pour vérifier que les méthodes fonctionnent bienplus de sens pour moi que les preuves mathématiques (encore une fois: fond dans une science sociale quantitative, pas de statistiques théoriques). Les méthodes fréquentistes et les méthodes de simulation ensemble ont beaucoup de sens pour moi. Parce que les fréquentistes voient la probabilité comme une chance à long terme (par exemple, si je le fais un nombre arbitrairement élevé de fois, et cela arrive 50% du temps, alors il y a une probabilité de 50%). Nous pouvons simuler ce long terme avec les méthodes de Monte Carlo!

Complications : depuis le premier cycle, je suis très au courant des méthodes bayésiennes, et il y a toujours eu des gens dans ma vie qui m'appelaient du côté bayésien, disant que les résultats étaient plus faciles à interpréter, que nous obtenons une probabilité pour une hypothèse au lieu des données donné une hypothèse, etc. J'étais vraiment dans ceci et ai pris une classe Bayesian, ai lu quelques livres et papiers Bayesian et suis maintenant assez familier avec Stan et ses paquets R associés.

Enter Mayo : Après avoir pensé "le bayésien est probablement la voie de l'avenir" pendant un certain temps, j'ai lu l' inférence statistique de Deborah Mayo comme un test sévère . Elle dit qu'elle ne choisit pas de camp au début du livre, mais elle le fait: elle est fréquentiste, et une grande partie du livre défend les méthodologies fréquentistes. Je ne veux pas nécessairement entrer dans une discussion pour savoir si nous pensons ou non que la façon dont elle voit les preuves est valide, mais cela m'a fait penser: Bayes est-elle vraiment tout ce qui est annoncé? Je veux dire, la foule des Bayes est tellement fracturée qu'elle-même que je ne connais même pas la "bonne" façon d'analyser les données dans un cadre bayésien. Habituellement, j'utiliserais simplementrstanarmet présentent des estimations ponctuelles et des intervalles crédibles ... qui correspondent souvent étroitement aux estimations fréquentistes et aux intervalles de confiance. Je pourrais faire des comparaisons de modèles, mais j'ai toujours peur de décrire les facteurs de Bayes comme des comparaisons de probabilités postérieures, etc.

Plus de réflexion : Ce que je n'arrêtais pas de penser à travers le livre de Mayo était: Il y a un moyen d'utiliser les ordinateurs pour s'assurer que nos méthodes fréquentistes fonctionnent, car la probabilité est ce que nous voyons à long terme et nous pouvons simuler cela. Les bayésiens ne peuvent même pas s'entendre sur la probabilité réelle, semble-t-il, selon l'école bayésienne (par défaut, subjective, etc.). Ce qui m'amène à ma question:

Question : Comment les bayésiens vérifient-ils que leurs méthodes définissent correctement l'incertitude (c.-à-d., Calculent des intervalles crédibles valides et des distributions postérieures) en utilisant des méthodes de simulation de Monte Carlo, si la probabilité n'est pas définie comme taux à long terme?

Exemple : je crée un générateur de données. Cela va simplement simuler à partir d'une distribution de Bernoulli avec une probabilité de 0,5:

set.seed(1839)
p <- .50
n <- 100
gen_dat <- function(n, p) {
  rbinom(n, 1, p)
}

Maintenant, disons que je veux m'assurer que les intervalles de confiance dans une régression logistique sont réellement valides. Je peux simuler une régression un grand nombre de fois et m'assurer que la valeur réelle de la population tombe dans l'intervalle de confiance à 95% 95% du temps. C'est un modèle d'interception uniquement, donc je veux juste m'assurer que j'évalue pcorrectement:

set.seed(1839)
iter <- 10000
results <- sapply(seq_len(iter), function(zzz) {
  mod <- glm(gen_dat(n, p) ~ 1, binomial)
  conf <- suppressMessages(confint(mod))
  log(p / (1 - p)) < max(conf) & log(p / (1 - p)) > min(conf)
})
mean(results)

Cela prend quelques minutes à courir, mais nous nous retrouvons avec l' mean(results)appel qui nous donne 0.9416. C'est environ 95%, et je suis sûr de dire que la glmcommande décrit l'incertitude d'une manière valide. Je suis sûr qu'il se serait rapproché de 95% si j'avais augmenté iteret voulu attendre ici plus longtemps sur mon ordinateur portable.

D'un autre côté, ajustons un modèle bayésien pour la même chose:

library(rstanarm)
set.seed(1839)
dat <- data.frame(y = gen_dat(n, p))
stan_mod <- stan_glm(y ~ 1, binomial, dat)
summary(stan_mod)

Cela me donne en partie:

Estimates:
                mean   sd    2.5%   25%   50%   75%   97.5%
(Intercept)    -0.1    0.2  -0.5   -0.2  -0.1   0.0   0.3  
mean_PPD        0.5    0.1   0.3    0.4   0.5   0.5   0.6  
log-posterior -73.0    0.7 -75.1  -73.1 -72.7 -72.5 -72.5  

Étant donné que les bayésiens ne définissent pas la probabilité comme ce que nous voyons à long terme, comment puis-je utiliser des méthodes de simulation pour vérifier que stan_glmcapturer avec précision l'incertitude? Autrement dit, comment puis-je avoir confiance que ces intervalles crédibles sont valides, en utilisant des méthodes de simulation? Et pour le moment, je ne suis même pas en train de définir un a priori - comment l'inclusion des prieurs intervient-elle ici, car cela affectera nos mesures d'incertitude?

Lorsque j'essayais d'écrire une régression bêta avec un composant de modèle d'obstacle dans Stan à partir de zéro une fois, quelqu'un m'a recommandé: "Simulez les données. Faites-le plusieurs fois, et les vraies estimations devraient être dans l'intervalle crédible d'environ 95 % du temps." Mais pour moi, cela va à l'encontre de la croyance même des Bayésiens! Cela repose sur une compréhension fréquentiste des probabilités! Alors, comment un bayésien pourrait-il me convaincre que l'intervalle crédible que j'obtiens de l' summary()appel à mon modèle décrit avec précision l'incertitude, en utilisant des méthodes de simulation?

But de la question : Ceci est un exemple trivial, mais souvent les clients me posent des problèmes difficiles. Et j'essaie des choses que je ne connais pas, alors je lance souvent une étude de simulation pour m'assurer que ce que je fais est valide. Si je devais écrire un modèle personnalisé dans Stan, comment saurais-je que ce que je fais est légitime? Comment pourrais-je utiliser des méthodes de simulation pour vérifier que ce que je fais dans Stan va réellement me dire ce que je veux savoir?

Je pense que vous voyez le problème logique dans votre question. Dans le paradigme fréquentiste, il est bon de présumer une vérité démographique, de générer des données et de voir si les estimations ont une bonne couverture, car c'est ce qu'elles sont censées faire. Dans le paradigme bayésien, cependant, il n'y a pas de vérité fondamentale pour générer des données! Les bayésiens demandent la probabilité de telles vérités à partir de données, donc dans la simulation, nous avons besoin de vérités différentes qui donnent naissance à des données, puis conditionnent les données. En pratique, on finit par simuler la loi de probabilité conditionnelle, qui, heureusement, tient toujours par définition. Je reprends ce numéro exact dans Rouder, 2014, Psychonomic Bulletin and Review. https://dx.doi.org/10.3758/s13423-014-0595-4

Comment les bayésiens vérifient-ils que leurs méthodes définissent correctement l'incertitude (c.-à-d., Calculent des intervalles crédibles valides et des distributions postérieures) en utilisant des méthodes de simulation de Monte Carlo, si la probabilité n'est pas définie comme taux à long terme?

Je crois que la confusion ici concerne le but des méthodes de simulation dans les statistiques bayésiennes. Le seul but des méthodes de Markov Chain Monte Carlo telles que Gibbs Sampling ou Hamiltonian Monte Carlo est de calculer le dénominateur de la règle de Bayes.

Bien sûr, il existe souvent d'autres méthodes qui rendraient MCMC inutile. Certains modèles peuvent être exprimés en utilisant la conjugaison, d'autres en appliquant une grille fine sur l'espace des paramètres, mais d'autres peuvent être résolus avec des tests d'acceptation-rejet. Lorsque MCMC est utile, c'est lorsque l'intégrale est mal conduite.

Bien que j'aimerais éviter les mathématiques, cela ne peut vraiment pas être évité. En examinant la règle de Bayes le numérateur est composé de et . est une vraisemblance et non une probabilité, donc il ne se résume à un que par hasard. Le dénominateur nous assure que égal à un. Le but de MCMC est de déterminer le nombre inférieur. Notez que le nombre inférieur est une constante. C'est la probabilité attendue.

La précision de ce nombre détermine certaines estimations de paramètres, mais pas toutes. Si vous utilisiez l'estimateur maximal a posteriori, alors MCMC est une étape inutile. Vous devriez plutôt construire un algorithme d'escalade. En revanche, il est nécessaire de déterminer la moyenne postérieure ou un intervalle. C'est parce que l'intervalle de 95% doit être 95% de quelque chose et le dénominateur détermine quelle est l'échelle de ce quelque chose.

Le but du MCMC dans les méthodologies bayésiennes est de faire converger les chaînes de Markov vers la densité postérieure. C'est ça. Il ne teste la validité de rien. Il s'agit simplement d'une tentative de déterminer une valeur à virgule fixe. C'est une forme d'intégration numérique. Comme il n'y a aucun moyen de savoir sans laisser l'algorithme s'exécuter à l'infini si toutes les régions denses ont été couvertes, il y a un certain jugement humain. L'algorithme aura une coupure quand il pense que c'est fait, mais cela ne signifie pas qu'il est réellement fait.

Dans les méthodologies fréquentistes, MCMC est souvent utilisé pour tester le caractère raisonnable d'un modèle ou pour approximer numériquement une solution lorsqu'une solution analytique n'est pas disponible. Il ne sert aucun objectif similaire ici.

Si je devais écrire un modèle personnalisé dans Stan, comment saurais-je que ce que je fais est légitime? Comment pourrais-je utiliser des méthodes de simulation pour vérifier que ce que je fais dans Stan va réellement me dire ce que je veux savoir?

Cette question est beaucoup plus difficile. Stan est un algorithme rapide, ce qui signifie qu'il échange la vitesse contre un risque supplémentaire d'inexactitude. Par construction, Stan sera plus souvent correct que incorrect. Il existe d'autres algorithmes conçus pour rechercher largement dans l'espace des paramètres des maximums locaux qui peuvent être plus précis, mais qui seront très lents.

Ce que vous devez faire, avant d'utiliser un algorithme particulier, est de lire la littérature sur cet algorithme et de regarder ses limites fonctionnelles. Malheureusement, il s'agit généralement d'un travail mathématique car le seul véritable objectif de toute méthode non conjuguée est d'estimer

La deuxième chose que vous pouvez faire est de le valider avec un algorithme alternatif. Les chiffres ne correspondront jamais, mais si vous les jugez suffisamment proches, alors tout va bien.

Troisièmement, la plupart des packages prédéfinis préviennent que quelque chose ne va pas. Si un avertissement apparaît, utilisez autre chose après avoir recherché la source du problème, afin de ne pas le recréer dans un autre algorithme.

Quatrièmement, regardez votre densité antérieure. Imaginez que vous aviez une densité antérieure de avec connue juste pour le simplifier et une probabilité de . À tout le moins, vous devriez aller wow, soit je me suis trompé, l'échantillon était mauvais, soit il y a quelque chose d'autre en cours que je devrais étudier.$\Pr(\mu)=\mathcal{N}(7,2^2)$$\sigma^2$$\mathcal{N}(25,.1^2)$

Cinquièmement, et vous devriez le faire avant de commencer Stan en premier lieu, représentez graphiquement vos probabilités marginales en une ou deux dimensions. Y a-t-il des surprises quelque part qui peuvent interférer avec l'algorithme?

Étant donné que les bayésiens ne définissent pas la probabilité comme ce que nous voyons à long terme, comment puis-je utiliser des méthodes de simulation pour vérifier que stan_glm capture avec précision l'incertitude? Autrement dit, comment puis-je avoir confiance que ces intervalles crédibles sont valides, en utilisant des méthodes de simulation? Et pour le moment, je ne suis même pas en train de définir un a priori - comment l'inclusion des prieurs intervient-elle ici, car cela affectera nos mesures d'incertitude?

Si vous ne définissez pas de priorité, votre modèle n'est pas valide. Si vous ne définissez pas une densité a priori raisonnable, pourquoi utiliseriez-vous un modèle bayésien? Les modèles Frequentist minimisent le risque de perte maximale qui pourrait se produire lors de la collecte d'un mauvais échantillon. Ils sont très pessimistes et il faut souvent plus d'informations pour produire le même résultat qu'une méthode bayésienne.

Néanmoins, cela ne sert à rien sans utiliser une bonne densité préalable. La densité antérieure permet à la méthode bayésienne de minimiser la perte moyenne due au choix d'un mauvais échantillon. Les informations contenues dans les antérieurs agissent comme un schéma de pondération de sorte que si un échantillon extrême est choisi par hasard, le prieur affaiblit le rôle que jouent les données.

EDIT J'ai réalisé que je n'avais pas fourni de réponse spécifique. C'était à la question

Comment pourrais-je utiliser des méthodes de simulation pour vérifier que ce que je fais dans Stan va réellement me dire ce que je veux savoir?

Ce qui rend cette question difficile, c'est que dans le paradigme bayésien, les points fixes sont, , l'échantillon. Dans les méthodes Frequentist, les paramètres sont fixes et des milliers d'échantillons invisibles sont créés. Du côté bayésien de la pièce, c'est l'échantillon qui est fixe. Vous devez simuler des milliers d'univers parallèles.$X$

Pour voir à quoi cela peut ressembler, imaginez toutes les fonctions de densité d'un tirage au sort avec une probabilité inconnue d'être des têtes et d'être des queues. Vous observez six têtes et deux queues. Imaginez un petit espace de paramètres où . Votre simulation considérerait tous les cas où six têtes pourraient être obtenues sur les trois distributions binomiales objectives. La valeur postérieure serait la moyenne pondérée de chaque paramètre étant la vraie valeur. Votre distribution prédictive serait la somme des distributions binomiales pondérées.$p$$1-p$$p\in\{1/3,1/2,2/3\}$

Pour vous, il est impossible que la prédiction bayésienne soit la vraie distribution. L'une des trois distributions est la vraie distribution. Les méthodes bayésiennes pondèrent leur probabilité en fonction de la valeur observée et de l'a priori. Le postérieur ne peut jamais être la vraie distribution, ni la densité prédictive.

Il demande "quelle est la probabilité de voir six têtes et deux queues sur l'ensemble de toutes les explications possibles (paramètres, modèles, etc.)".

Le Frequentist affirmerait que l'un des trois choix était la vraie valeur en le rendant nul. Six têtes et deux queues fausseraient mais pas les autres. Si, par hasard, vous avez choisi la bonne des trois distributions, alors vous avez parfaitement raison. Sinon, vous vous tromperez.$H_0:p=1/3,$

Si vous utilisiez des simulations pour maintenir un échantillon fixe, vous constateriez que Stan fonctionnerait admirablement car le théorème de Bayes est un théorème mathématique. Elle est optimale a posteriori. Tout ce que vous trouveriez, c'est que l'algorithme a correctement mis en œuvre le théorème de Bayes jusqu'au niveau d'erreur naturelle dans l'estimation du dénominateur.

Vous pouvez faire trois choses. Tout d'abord, vous pouvez utiliser des méthodes de notation de modèle pour les données hors échantillon. Deuxièmement, vous pouvez utiliser un processus de sélection de modèle bayésien ou de moyennage de modèle. Troisièmement, vous pouvez le traiter comme un problème fréquentiste et construire la distribution d'échantillonnage des estimateurs.

Pour les premiers, les méthodes de notation sont une littérature à part entière. Vous devriez les rechercher. La sélection des modèles bayésiens et la moyenne des modèles traitent les modèles comme des paramètres. Pour la sélection des modèles, la probabilité que les modèles soient vrais est calculée. Pour la moyenne d'un modèle, la probabilité que chaque modèle soit vrai est calculée et cela sert de pondération sur l'espace du modèle. Enfin, vous pouvez le traiter comme un modèle Frequentist.

Le dernier sera un problème dans de nombreux cas standard à cause du précédent. Pour les modèles avec trois dimensions ou plus et une distribution normale, la densité postérieure ne s'intégrera pas à l'unité si la densité antérieure n'est pas une densité appropriée. En d'autres termes, vous devez mordre la balle et choisir un prior pour n'importe quel modèle avec une réelle complexité.

La présence d'un a priori propre correctement centré force le cas où la méthode bayésienne sera supérieure à la méthode Frequentist correspondante en raison de l'amélioration de l'information. La méthode bayésienne gagnera sous n'importe quelle norme raisonnable. Ce n'est pas dû à un défaut de la méthode fréquentiste, mais la méthode bayésienne suppose des informations extérieures. La méthode Frequentist, en ne considérant que les informations contenues dans l'échantillon, contiendra moins d'informations si vous avez un réel préalable.

Encore une fois, si vous n'avez pas de véritable a priori, alors pourquoi utilisez-vous une méthode bayésienne?