intervalle de confiance d'une probabilité

Bonjour,
Pour simplifier, supposons un échantillonnage gaussien $x_i \sim_{\text{iid}} {\cal N}(\mu, \sigma^2)$. Comment avoir un intervalle de confiance de, par exemple, la probabilité $P(X_i < 1)$ ?

Réponses

  • Salut Steven,

    Je réponds de manière un peu naïve pour te faire préciser un peu la question (mon petit doigt me dit que tu as quelque chose derrière la tête :) ). Déjà tes $x_i$ te donnent des échantillons de Bernouilli en prenant $b_i=1_{x_i<1}$, et tu as un intervalle de confiance classique autour de l'estimateur $\hat{p}=\frac{1}{n}(b_1+\cdots+b_n)$, donné par l'inégalité de Hoeffding (ou, pour $n$ grand, par l'approximation gaussienne de la loi binomiale).

    Autre idée, tu cherches un intervalle de confiance sur $\Phi\left( \frac{1+\mu}{\sigma} \right)$, et ça doit pouvoir se bidouiller à partir d'intervalles de confiance sur $\mu$ et $\sigma$ : si $[\underline{\mu},\overline{\mu}]$ et $[\underline{\sigma},\overline{\sigma}] \subset \R_+$ sont des IC au niveau $1-\alpha/2$, alors $\underline{p}=\Phi\left( \frac{1+\underline{\mu}}{\overline{\sigma}} \right)$ et $\overline{p}=\Phi\left( \frac{1+\overline{\mu}}{\underline{\sigma}} \right)$ définissent un IC sur $p$, au niveau $1-\alpha$.
  • Thanx egoroffski, je n'avais même pas pensé la loi binomiale 8-)
    Rien de spécial derrière la tête, sauf que mon truc est un peu plus compliqué qu'un échantillonnage gaussien.
  • Petite correction si quelqu'un nous lit : c'est $\Phi\left( \frac{1-\mu}{\sigma} \right)$ pour ta 2ème solution.
    Et pour ête plus précis, le niveau de l'intervalle que tu proposes n'est pas exactement $1-\alpha$ mais il est au moins de $1-\alpha$.
  • Ah oui exact, merci pour la correction Steven. D'autre part mon truc avec les $\Phi(\frac{truc}{machin})$ est faux aussi tel quel, puisqu'il faut faire gaffe au signe de $1-\underline{\mu}$ et $1-\overline{\mu}$ dans les encadrements.

    Bon, j'ai fait des petits tests, avec R pour te faire plaisir, en simulant un $n$-échantillon de $\mathcal{N}(\mu,\sigma^2)$, et en calculant des intervalles de confiance classiques pour $\mu$ et $\sigma^2$ à base de loi de Student et du $\chi^2$ (ce qui avec $\Phi$ me donne un $IC_1$), et un autre pour les binomiales à base d'inégalité de Hoeffding (ce qui me donne un $IC_2$). Je répète 10000 fois l'expérience pour estimer (par Monte Carlo :)-D) la largeur de chaque intervalle, et le taux de couverture pour $p$.

    Avec les paramètres $\mu=3$, $\sigma=2$, $\alpha=5\%$ et une taille d'échantillon de $n=30$, ça marche, mais les intervalles sont vraiment trop larges (en moyenne $30\%$ pour $IC_1$ et $40\%$ pour $IC_2$) et les taux de couverture bien trop élevés, conformément à ta remarque ($99,75\%$ pour $IC_1$ et $99,99\%$ pour $IC_2$).

    Ci-dessous un graphe représentant les 1000 premières réplications ($IC_1$ en bleu, $IC_2$ en rouge, $p$ en vert).
    21482
  • Hola amigos,

    Un petit coup de delta méthode ça ne resserrerait pas le truc ?

    Amicalement
  • Hello Kuja,

    Oui c'est sans doute une excellente idée, j'essaierai ça demain. Mais sauf erreur ça n'est valable qu'asymptotiquement, or jusqu'à maintenant on a considéré des intervalles exacts ?
  • Thanx egoroffski, je regarderai ça quand j'aurai le temps. En ce moment je n'ai plus le ouebbe chez moi, c'est aussi pourquoi je fréquente moins le forum.
    Dans mon cas ce n'était pas grave que la couverture effective soit bien plus grande que la couverture nominale, car ça n'aurait pas changé la conclusion générale du problème.

    Par ailleurs ce serait intéressant d'avoir des intervalles de confiances de $\mu+k\sigma$.
    En bayésien c'est bien simple, modulo les problèmes du bayésien :D
  • est faux aussi tel quel, puisqu'il faut faire gaffe au signe de $1-\underline{\mu}$ et $1-\overline{\mu}$
    $1-x$ est décroissant en $x$, n'est-ce pas la seule chose qui importe ?
  • Mince tu as raison egoroffski : le problème c'est que $a/\sigma$ est croissant en $\sigma$ quand $a<0$.
  • Par ailleurs ce serait intéressant d'avoir des intervalles de confiances de $\mu+k\sigma$.
    $\rightarrow$ régression quantile
  • Voilà des bornes de confiance de $P(X_i > q)$ :
    $$
    \Phi\left[\frac{q-\hat\mu}{\hat\sigma}\left(1-\Phi^{-1}(1-\alpha/2)\sqrt{\frac{1}{n\frac{q-\hat\mu}{\hat\sigma}}+\frac{1}{2(n-1)}}\right) \right]
    $$
    et
    $$
    \Phi\left[\frac{q-\hat\mu}{\hat\sigma}\left(1+\Phi^{-1}(1-\alpha/2)\sqrt{\frac{1}{n\frac{q-\hat\mu}{\hat\sigma}}+\frac{1}{2(n-1)}}\right) \right]
    $$
  • Hello Steven,

    Je te crois sur parole, et j'ai essayé de tester ton expression, mais le problème est que $\frac{1}{n\frac{q-\hat\mu}{\hat\sigma}}+\frac{1}{2(n-1)}$ n'est pas toujours positif, du coup j'ai du mal à prendre la racine carrée.
  • Hello,
    Alors tu as tort de me croire :D
    C'est ça les formules :
    $$ \Phi\left[\frac{q-\hat\mu}{\hat\sigma}\left(1-\Phi^{-1}(1-\alpha/2)\sqrt{\frac{1}{n{\left(\frac{q-\hat\mu}{\hat\sigma}\right)}^2}+\frac{1}{2(n-1)}}\right) \right]
    $$ et $$ \Phi\left[\frac{q-\hat\mu}{\hat\sigma}\left(1+\Phi^{-1}(1-\alpha/2)\sqrt{\frac{1}{n{\left(\frac{q-\hat\mu}{\hat\sigma}\right)}^2}+\frac{1}{2(n-1)}}\right) \right]$$
    Et je précise que ce sont des bornes de confiance approximatives.
  • Les références :
    Bissell, A. F. (1990), "How Reliable Is Your Capability Index?" Applied Statistics, 30, 331 - 340.
    Kushler, R. H. and Hurley, P. (1992), "Confidence Bounds for Capability Indices," Journal of Quality Technology, 24, 188 - 195.
  • Merci pour la correction, et merci pour les références. Je suppose que ça découle d'un intervalle de confiance asymptotique sur le $q$ renormalisé $\frac{q-\mu}{\sigma}$ (ça ressemble à un "$z$-score" ?). Et ça un rapport avec ton histoire de "quantile regression" ?

    En tous cas ça marche très bien : pour $n=30$ (taille d'échantillon), mêmes paramètres $\mu,\sigma,\alpha$ que précédemment, et avec 10000 réplications, j'ai un taux de couverture de $95,06\%$, et une largeur moyenne d'intervalle égale à $20\%$.

    Ci-dessous les 1000 premières réplications :
    21527
  • Cool (tu)
    Oui à la base c'est un IC sur "l'indice de capabilité" $C_{pk} := \frac{q-\mu}{3\sigma}$ qu'on utilise en contrôle de qualité.
    Faudrait que j'essaye un bayésien un jour pour comparer.

    La régression quantile, non c'est autre chose, c'est pour estimer les quantiles, mais je n'ai trouvé des choses là-dessus que dans un cadre non-paramétrique.

    Merci beaucoup pour tes investigations :)-D
  • Pas de quoi mon cher Steven, je me suis bien amusé :)-D

    Quand tu dis "en bayésien", tu parles d'intervalles de crédibilité c'est ça ?
  • Ouaip, mais un intervalle de crédibilité avec une {\it bonne} loi a priori non-informative est aussi un intervalle de confiance avec niveau de confiance $\approx$ niveau de crédibilité.
  • Hum, je te crois sur parole, le bayésien c'est un peu nébuleux pour moi.
  • Tu peux me croire :D Nébuleux sur quel aspect ?
  • Bonne question. Pas la partie math, ça j'aime bien ; la partie numérique aussi. C'est plutôt sur l'interprétation des résultats. Bon il faut aussi dire que je n'ai jamais étudié tout ça sérieusement.

    Mais par exemple pour cette histoire d'intervalle de crédibilité, même dans un cas plus simple : modèle $\mathcal{N}(\mu,1)$ avec $\mu$ inconnu. Je vois mal comment on peut se débarasser de l'influence du choix d'un loi à priori pour $\mu$, même si je comprends que pour un grand nombre d'échantillons la loi a posteriori se met à ressembler de plus en plus à la vraisemblance "classique" (ou du moins le maximum a posteriori se rapproche du max de vraisemblance).
  • Me revoilà après une alarme incendie :D
    Du coup mon boulot n'a pas avancé, je te répondrai plus tard.

    [C'est parce que tu cherchais à bruler les étapes ? :D AD]
  • OK, pas de souci (content que tu sois sain et sauf :) )
  • Re,
    Très bon exemple le modèle $y \sim \mathcal{N}(\mu,1)$.
    Une bonne loi a priori non-informative n'est pas une loi en fait ; c'est une fonction qui, quand tu l'injectes dans la formule de Bayes comme si c'était une loi, mène à une loi a posteriori qui reflète au mieux l'information apportée par les données.
    Pour ce modèle, la bonne loi a priori non-informative est celle qu'utilisait Laplace : $\pi(\mu)=1$, la mesure de Lebesgue sur $\mathbb{R}$, façon de parler.
    La loi a posteriori est alors $(\mu | y) \sim \mathcal{N}(\bar y,1/n)$ ; tu peux vérifier qu'ici, si tu prends l'intervalle de crédibilité unilatéral à $95\%$, son niveau de confiance est exactement $95\%$.
    La suite de la leçon sur demande.
  • Merci Steven pour ta réponse pas si tardive que ça.

    Il me semble avoir lu ça quelque part en effet. Dans cet exemple, prendre la loi a priori $\pi \equiv 1$ fait que la loi a posteriori est proportionnelle à la vraisemblance de l'échantillon, qui est en $e^{-\frac{1}{2}(\bar{y}-\mu)^2}$, et vu la propriété de symétrie entre moyenne et variable de la densité gaussienne, ça revient à une loi normale pour le paramètre $\mu$ (d'ailleurs il me semble que l'astuce marche encore si $\pi(d\mu)=\mathcal{N}(m_{\mu},\sigma_{\mu}^2)$, avec une moyenne et une variance a posteriori différentes).

    Comme la variance de la loi a posteriori est $1/n$, on trouve l'intervalle de crédibilité bayésien $[\bar{y} \pm z_{1-\alpha/2} /n]$, qui est la même chose que l'intervalle de confiance bilatéral classique.

    En attendant la suite, qui viendra quand tu auras le temps, j'ai deux questions :

    1) Lorsqu'on prend une "loi a priori" de masse totale infinie, qu'est-ce qui garantit que la loi a posteriori sera quand même de masse totale finie ? Si la loi a posteriori est de masse infinie, ça devient un peu emm..dant pour définir un intervalle de crédibilité.

    3) La question qui tue (mais qui est peut-être idiote) : si on s'arrange pour choisir $\pi$ de sorte que la l'intervalle de crédibilité coïncide avec l'intervalle de confiance "classique" (à supposer que ce soit toujours possible), quel est finalement l'interêt de faire du bayésien ?
  • 1) Rien ne le garantit, et la loi a posteriori doit toujours être une vraie loi, si ce n'est pas le cas la loi a priori n'est pas "admissible".

    3) (y'a pas de 2 ?). Par exemple en bayésien tu as à ta disposition la loi prédictive a posteriori, pas en fréquentiste. Autre avantage : reprends le problème de départ de ce fil, tu simules $\mu$ et $\sigma$ a posteriori, tu peux en tirer des simulations de $f(\mu, \sigma)$ pour toute fonction $f$, et un intervalle de confiance pour $f(\mu, \sigma)$.
  • NB : ma 2ème remarque dans 3) n'est pas précise, la discussion se complique quand il y a plusieurs paramètres
  • Merci (encore).

    1) OK, ça me rassure un peu. Par exemple si j'ai le modèle $\mathcal{U}([-2,-1] \cup [\theta,\theta+1])$, et que toutes les observations sont dans $[-2,-1]$, alors la loi a priori "uniforme sur $\R_+$" pour $\theta$ n'est pas admissible ? Par contre elle le devient dès que j'ai une observation positive (bon, il y a sans doute des exemples moins tordus).

    2) Le 3 au lieu du 2 c'était pour voir si tu suivais :)

    3) OK, ça me convainc déjà un peu. La loi de prédiction à posteriori, c'est l'intégrale de $p(dx|\theta)$ contre la loi a priori j'imagine ? Et pour la simulation, là je comprends bien l'interêt, si je sais simuler la loi jointe a posteriori de $\theta=(\mu,\sigma)$, alors je sais simuler n'importe quelle fonction $f(\theta)$ et trouver des bornes de confiance crédibilité$^*$ empiriques sans hypothèses de monotonie sur $f$.

    [Edit, $*$ : si la loi a priori est "bien choisie" au sens que tu as expliqué plus haut, alors est-ce les intervalles de crédibilité deviennent des intervalles de confiance fréquentistes pour n'importe quelle quantité $f(\theta)$ ?]
  • 1) Oui normalement il faut une loi a posteriori "propre" pour toute issue possible.

    3) Oui, la prédictive a priori est la loi des observations intégrée sur la loi a priori, idem pour posteriori.

    Pour ton edit, ce n'est pas garanti. Mais ça devrait pas mal marcher en général, je pense.
  • Retour sur ton edit : il y a un célèbre exemple de Stein 1959, où on utilise une "bonne" loi a priori non-informative pour $\theta$ (multidimensionnel), mais avec celle-ci on obtient un intervalle de crédibilité sur un $f(\theta)$ qui a un niveau de confiance asymptotiquement nul :D
    C'est pourquoi le développement de la statistique bayésienne objective (= lois a priori non-informatives) en est venu à définir des lois a priori pour un paramètre d'intérêt $f(\theta)$ qu'on se fixe.
  • Par ailleurs ce serait intéressant d'avoir des intervalles de confiances de $\mu+k\sigma$.
    En fait ça c'est dans le domaine des intervalles de tolérance, je pense. Mais j'ai jamais trouvé un bon document là-dessus.
  • Voilà un intervalle de confiance obayésien ("o" pour "objectif") de P(X>q) :
    bounds <- function(y, q, alpha=.05, nsims=10000){
            n <- length(y)
    	ssq <- (n-1)*var(y)
    	rho <- rchisq(nsims,n)/ssq
    	mu <- rnorm(nsims, mean(y), 1/sqrt(n*rho))
    	p <- 1-pnorm(q,mu,1/sqrt(rho))
    return(quantile(p, c(alpha/2,1-alpha/2)))
    }
    
    n <- 100
    y <- rnorm(n,0,1)
    bounds(y,q=1.96)
    
           2.5%       97.5% 
    0.009790207 0.049672354
    
    Egoroffski, si ça t'amuse encore, tu pourras l'essayer avec ton code ?
  • Bon comme ça n'a plus l'air d'amuser egoroffski j'ai fait un essai.
    n <- 30
    mu <- 0
    sigma <- 1
    nsims <- 10000
    q <- 1.96
    p <- 1-pnorm(q,mu,sigma)
    
    test <- rep(NA, nsims)
    longueur <- rep(NA, nsims)
    
    for(i in 1:nsims){
    	y <- rnorm(nsims, mu, sigma)
    	bornes <- bounds(y,q)
    	test[i] <- (p > bornes[1]) & (p < bornes[2])
    	longueur[i] <- bornes[2]-bornes[1]
    }
    
    > mean(test)
    [1] 0.9515
    > mean(longueur)
    [1] 0.00391341
    
  • Hello Steven,

    Non, non ça m'amuse encore et surtout ça m'intéresse beaucoup, justement je viens de tester ta fonction, mais chez moi ça ne marche pas très bien. Du coup j'étais en train d'essayer de comprendre ce qu'elle fait pas-à-pas. Si je comprends on simule des couples $(\rho_j,\mu_j)$ qui ont la loi a posteriori de $(1/\sigma^2,m)$ sachant les $y_k$, qui sont des observations d'une v.a. $Y$ de loi $p(y|m,\sigma^2)=\mathcal{N}(m,\sigma^2)$ ?

    Je vais essayer de trouver le courage de me lancer dans les calculs pour retrouver quelle loi a priori tu utilises, mais en attendant j'obtiens, pour une réalisation de mon vecteur $X=(y_1,...,y_{30})$, simulé suivant $\mathcal{N}(3,2^2)$, l'intervalle :
    > bounds(X,q=1)
         2.5%     97.5% 
    0.7363267 0.9359803 
    
    qui est loin de contenir la vraie valeur de $p=\mathbb{P}(Y \leq 1) \simeq 0,1587$.
  • C'est $\mathbb{P}(Y \geq 1)$ :D
  • J'utilise la loi a priori de Jeffreys, qui est à la limite de la famille conjuguée. Ceci se généralise aux modèles linéaires. Si tu veux mes docs à ce sujet, donne-moi un mail. Voici la partie que j'utilise ici.
    21539
  • OK, à cause du $1-\Phi(\cdots)$.. quel gros malin je fais 8-). Bon je modifie ta fonction et je re-teste.. oui, ça marche mieux :) Maintenant en répliquant le test 10000 fois, j'obtiens une longueur moyenne de $20,5\%$ et un taux de couverture de $94,3\%$, sensiblement la même chose que pour l'IC fréquentiste asymptotique (d'ailleurs sur le graphique il est assez difficile de les distinguer). Par contre on sent la différence en termes de temps de calcul, mais bon le code est sans doute loin d'être optimal.

    Merci pour le doc, oui ça m'intéresse beaucoup en effet ! Il y a une adresse mail attachée à mon profil.

    Je retourne essayer de comprendre ta fonction, à toute.
  • Forcément que c'est bien plus long, car la fonction bounds() génère 10000 simulations de la posteriori.
    Ton profil dit :
    Adresse électronique:
    Caché
  • Ah mince, désolé. Je l'ai rendu visible maintenant.

    Oui, forcément, 10000 réplications d'un truc contenant 10000 simulations à chaque fois, ça commence à faire long.
  • Ouais et imagine quand tu simules avec du Gibbs, et t'es même pas sûr que ça a convergé :D
    Je t'ai envoyé les docs.
  • Hey j'ai reçu une notification d'échec d'envoi, t'es sûr de ton mail ?
  • Arf, mince, bon je t'ai envoyé une adresse e-mail qui marche par MP, désolé pour la perte de temps.

    Oui pour du MCMC c'est pire j'imagine.. en même temps, le fait de répliquer un grand nombre de fois les simulations peut justement constituer un bon diagnostic de convergence non ?
  • Ppfff.. diagnostic de convergence, on en parlera une autre fois :D
    Je continue à faire ton éducation obayésienne. Mon "slide" ci-dessous dit que le niveau de confiance des intervalles de crédibilité est optimal pour la loi a priori de Jeffreys, dans le cas d'un modèle à un seul paramètre. Optimal au sens de la "coïncidence fréquentiste". Je pense que tu sauras décoder, surtout que nous sommes dans le vif du sujet :)-D
    Pour les intervalles bilatéraux, la question est plus complexe.

    Mais c'est encore bien plus compliqué dans le cas d'un modèle à plusieurs paramètres. L'exemple de Stein que j'ai précedemment mentionné, concerne précisément la loi a priori de Jeffreys.

    Il existe une meilleure généralisation de la loi a priori de Jeffreys au cas de plusieurs paramètres. C'est la loi a priori dite "de référence", ou de Bernardo, ou de Berger & Bernardo. Elle est définie relativement au choix d'un paramètre d'intérêt, et semble passer au travers de plusieurs problèmes ou paradoxes des lois a priori candidates à l'analyse bayésienne "objective". Il y a des "paradoxes" délirissimes qui découlent de l'utilisation de la formule de Bayes avec une mesure de masse infinie.

    Enfin il est important que je mentionne que ni la loi a priori de Jeffreys ni la loi a priori de référence ne sont fondées sur un espoir de "coïncidence fréquentiste".
    21544
  • Merci encore Steven, c'est de plus en plus intéressant.

    Dans ce goût-là, je connaissais l'existence du théorème de Bernstein-Von MIses, qui si je ne m'abuse dit qu'asymptotiquement la loi a posteriori et la loi de l'EMV se "confondent" ; j'imagine que le cas $\pi \neq \pi^J$ est plus ou moins une conséquence de ce théorème (si la convergence qu'il implique est assez forte pour passer aux quantiles), et effectivement le fait qu'on puisse passer à du $O(1/n)$ est assez impressionnant, surtout si ce n'était pas "prévu au départ" comme tu le suggères.

    Je vais regarder ce que je trouve sur Berger-Bernardo.
  • Sur le site de José Bernardo il y a pas mal de ses articles qui sont téléchargeables.
    D'ailleurs ce matin j'ai trouvé un article de l'ami Bernardo qui traite le problème de ce fil apparemment(***) : http://www.uv.es/~bernardo/Capability.pdf. Voilà une excellente lecture pour ce week-end.
    Tu peux aussi regarder la partie sur les lois a priori non-informatives dans l'excellent bouquin de Christian Robert, "Le choix bayésien", ou sa version british, "The Bayesian choice".

    (***) EDIT : non (après lecture) le sujet de l'article n'est pas celui de ce fil
Connectez-vous ou Inscrivez-vous pour répondre.