Intervalles de confiance pour la médiane

J'ai une distribution d'échantillons avec un petit nombre de valeurs dans chacun (moins de ). J'ai calculé la médiane de chaque échantillon, que je veux comparer avec un modèle et obtenir la différence entre le modèle et la médiane de chaque échantillon. Pour avoir un résultat cohérent, j'ai besoin d'une erreur sur cette différence.$10$

Il en résulte que trouver l'écart type dans un tel cas peut être assez difficile, au moins pour un non-pro comme moi (voir par exemple ici ).

J'ai trouvé ce site Web qui explique comment calculer les intervalles de confiance pour la médiane, même s'il n'y a aucune référence officielle citée.

Cela me semble raisonnable, mais je ne peux pas vraiment juger, alors je voudrais savoir:

1. ces formules sont-elles correctes?
2. Il y a une référence pour ça?
3. Et si je veux trouver un IC différent de ?$95\%$

Merci d'avance

EDIT: J'ai également trouvé cet exemple de bootstrap pour des données non gaussiennes . Maintenant, je ne sais pas grand chose sur le bootstrap, mais il serait bon d'avoir une adresse sur sa validité.

Sommaire

Lorsque vous pouvez supposer peu ou rien de la vraie loi de probabilité et en déduire peu - ce qui est le cas pour de petits échantillons de observations - alors une paire de statistiques d'ordre convenablement choisie constituera un intervalle de confiance pour la médiane. Les statistiques d'ordre à choisir peuvent être facilement trouvées avec une analyse rapide de la distribution binomiale . Il y a quelques choix à faire dans la pratique: ils sont discutés et illustrés à la fin de ce post.$n$$(n, 1/2)$

Par ailleurs, la même analyse peut être utilisée pour construire des intervalles de confiance pour tout quantile (dont la médiane, correspondant à , est un exemple). La distribution binomiale régit la solution dans ce cas.$q$$q=50\%$$(n, q)$

introduction

Rappelez-vous ce que signifie un intervalle de confiance (IC). Le cadre est un échantillon aléatoire indépendant avec chaque régie par la même distribution . On suppose seulement que est un élément d'un ensemble de distributions possibles. Chacun d'eux a une médiane . Pour tout fixe entre et , un CI de niveau est une paire de fonctions (aka "statistiques"), et , telles que$X = (X_1, X_2, \ldots, X_n)$$X_i$$F$$F$$\Omega$$F_{1/2}$$\alpha$$0$$1$$\alpha$$L$$U$

$${\Pr}_F(L(X)\le F_{1/2} \le U(X)) \ge 1 - \alpha.$$

Le côté droit est la couverture de la CI pour la distribution .$F$

En plus: pour que cela soit utile, nous préférons également que (1) l'infimum des couvertures sur soit aussi petit que possible et (2) la longueur attendue de l'intervalle, , devrait avoir tendance à être court pour tous ou "la plupart" .$F\in\Omega$$\mathbb{E}_F(U(X)-L(X))$$F\in\Omega$

Une analyse

Supposons que nous n'assumions rien sur . $\Omega$ Dans cette situation, nous pouvons toujours exploiter les statistiques de commande . Ce sont les valeurs spécifiques de l'échantillon trié. Pour simplifier la notation, trions l'échantillon une fois pour toutes afin que

$$X_1 \le X_2 \le \cdots \le X_n.$$

La valeur est la statistique d'ordre de l'échantillon. Puisque nous ne supposons rien sur , nous ne savons rien sur au début, donc nous ne pouvons pas en déduire beaucoup sur les intervalles probables entre chaque et son voisin . Cependant, nous pouvons encore raisonner quantitativement sur les valeurs individuelles: quelle est la chance que ne dépasse pas la médiane de ? Pour comprendre cela, soit une variable aléatoire régie par , et soit$X_i$$i^\text{th}$$\Omega$$F$$X_i$$X_{i+1}$$X_i$$F$$Y$$F$

$$\pi_F = {\Pr}_F(Y \le F_{1/2})$$

la chance que ne dépasse pas la médiane de . Ensuite, lorsque nous savons (depuis ) que notre échantillon non ordonné d'origine de valeurs doit avoir contenu au moins valeurs ne dépassant pas .$Y$$F$$X_i \le F_{1/2}$$X_1\le \cdots \le X_i \le F_{1/2}$$n$$i$$F_{1/2}$

Il s'agit d'un problème binomial. Formellement, si nous définissons la variable aléatoire à lorsque et sinon, ce qui précède montre que a une distribution de Bernoulli avec le paramètre . Un «succès» consiste à observer une valeur égale ou inférieure à la médiane. Par conséquent, est donné par la probabilité binomiale associée à moins de succès:$Z$$1$$Y \le F_{1/2}$$0$$Z$$\pi_F$$\Pr(X_i \gt F_{1/2})$$i$

$$\Pr(X_i \gt F_{1/2}) = \sum_{j=0}^{i-1} \binom{n}{j} \pi_F^j(1-\pi_F)^{n-j}.$$

Vous avez probablement remarqué que . En fait, pour de nombreuses distributions, les deux valeurs sont égales: elles ne diffèrent que lorsque attribue une probabilité positive à la médiane . Pour analyser la différence, écrivez pour . Pour cela implique$\pi_F \ge 1/2$$F$$F_{1/2}$$\pi_F = 1/2 + \varepsilon$$\varepsilon \ge 0$$2(j-1) \le n$

$$\eqalign{ \pi_F^j(1-\pi_F)^{n-j} &= (1/2+\varepsilon)^j(1/2-\varepsilon)^{n-j} = (1/2+\varepsilon)^j[(1/2-\varepsilon)^j(1/2-\varepsilon)^{n-2j}]\\ &=(1/4-\varepsilon^2)^j(1/2-\varepsilon)^{n-2j} \le (1/4)^j(1/2)^{n-2j}=2^{-n}. }$$

Par conséquent, lorsque , on peut se débarrasser de la dépendance de la somme sur , au prix de remplacer l'égalité par une inégalité:$2(i-1) \le n$$F$

$$\Pr(X_i \gt F_{1/2}) \le 2^{-n}\sum_{j=0}^{i-1} \binom{n}{j}.$$

Exactement le même argument (appliqué en inversant les statistiques de commande) montre que lorsque ,$2(i+1) \ge n$

$$\Pr(X_i \lt F_{1/2}) \le 2^{-n}\sum_{j=i+1}^n \binom{n}{j}.$$

Le côté droit se réduit à zéro chaque fois que (dans le premier cas) ou (dans le second). Par conséquent, il est toujours possible de trouver des index pour lesquels$i \le 0$$i \ge n$$l \le u$

$$\eqalign{ \Pr(X_l \gt F_{1/2} \text{ or } X_u \lt F_{1/2}) &= \Pr(X_l \gt F_{1/2}) + \Pr( X_u \lt F_{1/2}) \\ &\le 2^{-n}\left(\sum_{j=0}^{l-1} \binom{n}{j} + \sum_{j=u+1}^n \binom{n}{j}\right). }$$

Solution

C'est le complément de la condition définissant un intervalle de confiance, et donc équivalent à lui:

$$\Pr(X_l \le F_{1/2}\le X_u ) \ge 2^{-n}\sum_{j=l}^u \binom{n}{j}.$$

En sélectionnant pour faire du côté droit au moins , nous aurons trouvé une procédure d'intervalle de confiance dont le niveau est au moins .$l \le u$$1-\alpha$ $1-\alpha$

En d'autres termes, lors du choix de ces indices et , en définissant et , l'intervalle sera un CI pour la médiane ayant une couverture au moins . Vous pouvez calculer sa couverture réelle en termes de probabilités binomiales. Cette couverture sera atteinte pour toute distribution qui attribue une probabilité nulle à (qui inclut toutes les distributions continues). Il sera dépassé par tout qui attribue une probabilité non nulle à .$l$$u$$L(X) = X_l$$U(X) = X_u$$[L(X), U(X)]$$F_{1/2}$$1-\alpha$$F$$F_{1/2}$$F$$F_{1/2}$

Discussion

À ce stade, nous avons quelques choix. Le plus courant est de rendre les limites symétriques en fixant raisonnablement proche de . En fait, en stipulant , les limites de confiance peuvent être trouvées pour tout avec une recherche rapide ou en appliquant la fonction quantile binomiale.$u$$n+1-l$$u=n+1-l$$n$

Par exemple, soit et (pour illustrer une procédure CI ). Comptons la partie inférieure de la distribution binomiale cumulative avec les paramètres et :$n=10$$\alpha=10\%$$1-\alpha=90\%$$10$$1/2$

> i <- 0:5; names(i) <- i; print(pbinom(i, 10, 1/2), digits=1)
    0     1     2     3     4     5   
0.001 0.011 0.055 0.172 0.377 0.623 

(Ceci est une Rcommande et sa réponse.) Parce que la valeur à , égale à , est proche de , il est tentant de prendre et , pour alors la couverture sera de ce qui est proche de l'objectif de . Si vous devez atteindre la couverture souhaitée, vous devez prendre et ou et , les deux avec une couverture .$2$$5.5\%$$\alpha/2$$l=3$$u=10+1-3=8$$1 - 0.055 - 0.055 = 0.89$$90\%$$l=2$$u=8$$l=3$$u=9$$1 - 0.011 - .055 = 0.935$

Pour vérifier, simulons un grand nombre d'ensembles de données à partir de n'importe quelle distribution, calculons ces IC pour les ensembles de données et calculons la proportion d'IC ​​qui couvrent la vraie médiane. Cet Rexemple utilise une distribution normale:

n <- 10
n.sim <- 1e4
x <- apply(matrix(rnorm(n*n.sim), nrow=n), 2, sort)
covers <- function(x, l, u) mean(x[l, ] <= 0 & x[u, ] >= 0)
c(l3.u8=covers(x,3,8), l2.u8=covers(x,2,8), l3.u9=covers(x,3,9))

La sortie est

 l3.u8  l2.u8  l3.u9 
 0.8904 0.9357 0.9319 

Les couvertures concordent étroitement avec les valeurs théoriques.

Comme autre exemple, tirons des échantillons d'une distribution discrète, comme un Poisson:

lambda <- 2
x <- apply(matrix(rpois(n*n.sim, 2), nrow=n), 2, sort)
med <- round(lambda + 1/3 - 0.02/lambda)
c(l3.u8=covers(x,3,8), l2.u8=covers(x,2,8), l3.u9=covers(x,3,9))

 l3.u8  l2.u8  l3.u9 
0.9830 0.9845 0.9964 

Cette fois, les couvertures sont beaucoup plus élevées que prévu. La raison en est qu'il y a chances qu'une valeur aléatoire soit égale à la médiane. Cela augmente considérablement les chances que l'IC couvre la médiane. Ce n'est ni un problème ni un paradoxe. Par définition, la couverture doit être d'au moins quelle que soit la distribution - mais il est possible (comme dans ce cas) que la couverture pour des distributions particulières soit sensiblement supérieure à .$27\%$$1-\alpha$$F$$1-\alpha$

C'est là que réside le compromis: lorsque vous ne présumez rien de , l'IC basé sur les statistiques de commande est le seul que vous pouvez construire. Sa couverture pour votre vrai (mais inconnu) peut être un peu plus élevée que ce à quoi vous vous attendez. Cela signifie que votre CI sera plus large que si vous aviez fait des hypothèses fortes sur en limitant les possibilités de .$F$$F$$\Omega$$F$

Si vous souhaitez utiliser des méthodes numériques, vous pouvez générer une estimation de la distribution d'échantillonnage des médianes en utilisant le bootstrap. Rééchantillonnez votre échantillon à plusieurs reprises et calculez de nombreuses médianes. Le stdev de ces médianes sert d'estimation du stdev de la distribution d'échantillonnage des médianes. J'ai utilisé une méthode similaire pour calculer l'incertitude des résultats des parties d'échecs dans mon article sur les gambits d'échecs qui peut être trouvé ici https://sonoma.academia.edu/JamalMunshi/papers