Comment obtenir un intervalle de confiance pour un centile?

J'ai un tas de valeurs de données brutes qui sont des montants en dollars et je veux trouver un intervalle de confiance pour un centile de ces données. Existe-t-il une formule pour un tel intervalle de confiance?

Cette question, qui couvre une situation courante, mérite une réponse simple et non approximative. Heureusement, il y en a un.

Supposons que sont des valeurs indépendantes d'une distribution inconnue dont le quantile j'écrirai . Cela signifie que chaque a une chance (au moins) d'être inférieur ou égal à . Par conséquent, le nombre de inférieur ou égal à a une distribution binomiale . F q th F - 1 ( q ) X i q F - 1 ( q ) X i F - 1 ( q ) ( n , q $X_1, \ldots, X_n$$F$$q^\text{th}$$F^{-1}(q)$$X_i$$q$$F^{-1}(q)$$X_i$$F^{-1}(q)$$(n,q)$

Motivés par cette simple considération, Gerald Hahn et William Meeker dans leur manuel Statistical Intervals (Wiley 1991) écrivent

Un intervalle de confiance bilatéral sans distribution conservateur de pour est obtenu ... commeF - 1 ( q ) [ X ( l ) , X ( u )$100(1-\alpha)\%$$F^{-1}(q)$$[X_{(l)}, X_{(u)}]$

où sont les statistiques d'ordre de l'échantillon. Ils disent$X_{(1)}\le X_{(2)}\le \cdots \le X_{(n)}$

On peut choisir des entiers symétriquement (ou presque symétriquement) autour de et le plus près possible sous réserve des exigences queq ( n + 1 ) B ( u - 1 ; n , q ) - B ( l - 1 ; n , q ) ≥ 1 - α $0 \le l \le u \le n$$q(n+1)$

$$B(u-1;n,q) - B(l-1;n,q) \ge 1-\alpha.\tag{1}$$

L'expression à gauche est la chance qu'une variable binomiale ait l'une des valeurs . Evidemment, c'est la chance que le nombre de valeurs de données comprises dans les inférieurs de la distribution ne soit ni trop petit (inférieur à ) ni trop grand ( ou supérieur).{ l , l + 1 , … , u - 1 } X i 100 q % l $(n,q)$$\{l, l+1, \ldots, u-1\}$$X_i$$100q\%$$l$$u$

Hahn et Meeker suivent avec quelques remarques utiles, que je citerai.

L'intervalle précédent est conservateur car le niveau de confiance réel, donné par le côté gauche de l'équation , est supérieur à la valeur spécifiée . ...1 - $(1)$$1-\alpha$

Il est parfois impossible de construire un intervalle statistique sans distribution ayant au moins le niveau de confiance souhaité. Ce problème est particulièrement aigu lors de l'estimation des centiles dans la queue d'une distribution à partir d'un petit échantillon. ... Dans certains cas, l'analyste peut faire face à ce problème en choisissant et non symétrique. Une autre alternative peut être d'utiliser un niveau de confiance réduit.$l$$u$

Examinons un exemple (également fourni par Hahn & Meeker). Ils fournissent un ensemble ordonné de «mesures d'un composé issu d'un processus chimique» et demandent un intervalle de confiance de pour le percentile. Ils affirment que et fonctionneront.100 ( 1 - α ) = 95 % q = 0,90 l = 85 u = $n=100$$100(1-\alpha)=95\%$$q=0.90$$l=85$$u=97$

La probabilité totale de cet intervalle, comme le montrent les barres bleues de la figure, est de : c'est aussi proche que possible de , mais toujours au-dessus, en choisissant deux seuils et en éliminant toutes les chances dans le la queue gauche et la queue droite qui sont au-delà de ces seuils.95 $95.3\%$$95\%$

Voici les données, présentées dans l'ordre, en laissant de côté des valeurs du milieu:$81$

Le plus grand est et le plus grand est . L'intervalle est donc . 24,33 97 e 33,24 [ 24,33 , 33,24 $85^\text{th}$$24.33$$97^\text{th}$$33.24$$[24.33, 33.24]$

Réinterprétons cela. Cette procédure était censée avoir au moins chances de couvrir le percentile. Si ce centile dépasse réellement , cela signifie que nous aurons observé ou plus sur valeurs dans notre échantillon qui sont inférieures au centile. C'est trop. Si ce centile est inférieur à , cela signifie que nous aurons observé ou moins de valeurs dans notre échantillon qui sont inférieures au centile. C'est trop peu.90 e 33,24 97 100 90 e 24,33 84 90 $95\%$$90^\text{th}$$33.24$$97$$100$$90^\text{th}$$24.33$$84$$90^\text{th}$ Dans les deux cas - exactement comme indiqué par les barres rouges sur la figure - ce serait une preuve contre le centile se trouvant dans cet intervalle.$90^\text{th}$

Une façon de trouver de bons choix de et est de rechercher en fonction de vos besoins. Voici une méthode qui commence par un intervalle approximatif symétrique puis recherche en faisant varier à la fois et jusqu'à afin de trouver un intervalle avec une bonne couverture (si possible). Il est illustré avec du code. Il est configuré pour vérifier la couverture dans l'exemple précédent pour une distribution normale. Sa sortie est$l$$u$$l$$u$$2$R

La couverture moyenne de simulation était de 0,9503; la couverture prévue est de 0,9523

L'accord entre simulation et attente est excellent.

#
# Near-symmetric distribution-free confidence interval for a quantile `q`.
# Returns indexes into the order statistics.
#
quantile.CI <- function(n, q, alpha=0.05) {
  #
  # Search over a small range of upper and lower order statistics for the 
  # closest coverage to 1-alpha (but not less than it, if possible).
  #
  u <- qbinom(1-alpha/2, n, q) + (-2:2) + 1
  l <- qbinom(alpha/2, n, q) + (-2:2)
  u[u > n] <- Inf
  l[l < 0] <- -Inf
  coverage <- outer(l, u, function(a,b) pbinom(b-1,n,q) - pbinom(a-1,n,q))
  if (max(coverage) < 1-alpha) i <- which(coverage==max(coverage)) else
    i <- which(coverage == min(coverage[coverage >= 1-alpha]))
  i <- i[1]
  #
  # Return the order statistics and the actual coverage.
  #
  u <- rep(u, each=5)[i]
  l <- rep(l, 5)[i]
  return(list(Interval=c(l,u), Coverage=coverage[i]))
}
#
# Example: test coverage via simulation.
#
n <- 100      # Sample size
q <- 0.90     # Percentile
#
# You only have to compute the order statistics once for any given (n,q).
#
lu <- quantile.CI(n, q)$Interval
#
# Generate many random samples from a known distribution and compute 
# CIs from those samples.
#
set.seed(17)
n.sim <- 1e4
index <- function(x, i) ifelse(i==Inf, Inf, ifelse(i==-Inf, -Inf, x[i]))
sim <- replicate(n.sim, index(sort(rnorm(n)), lu))
#
# Compute the proportion of those intervals that cover the percentile.
#
F.q <- qnorm(q)
covers <- sim[1, ] <= F.q & F.q <= sim[2, ]
#
# Report the result.
#
message("Simulation mean coverage was ", signif(mean(covers), 4), 
        "; expected coverage is ", signif(quantile.CI(n,q)$Coverage, 4))

Dérivation

Le -quantile (c'est le concept plus général que le centile) d'une variable aléatoire est donné par . L'échantillon homologue peut être écrit comme - ce n'est que le quantile échantillon. Nous sommes intéressés par la distribution de:$\tau$$q_\tau$$X$$F_X^{-1}(\tau)$$\hat{q}_\tau = \hat{F}^{-1}(\tau)$

$\sqrt{n}(\hat{q}_\tau - q_\tau)$

Premièrement, nous avons besoin de la distribution asymptotique du cdf empirique.

Puisque , vous pouvez utiliser le théorème de la limite centrale. est une variable aléatoire de Bernoulli, donc la moyenne est et la variance est .$\hat{F}(x) = \frac{1}{n} \sum 1\{X_i < x\}$$1\{X_i < x\}$$P(X_i < x) = F(x)$$F(x)(1-F(x))$

$\sqrt{n}(\hat{F}(x) - F(x)) \rightarrow N(0, F(x)(1-F(x))) \qquad (1)$

Maintenant, comme l'inverse est une fonction continue, nous pouvons utiliser la méthode delta.

[** La méthode delta dit que si , et est une fonction continue, alors **]$\sqrt{n}(\overline{y} - \mu_y) \rightarrow N(0,\sigma^2)$$g(\cdot)$$\sqrt{n}(g(\overline{y}) - g(\mu_y)) \rightarrow N(0, \sigma^2 (g'(\mu_y))^2)$

Dans la partie gauche de (1), prenez , et$x=q_\tau$$g(\cdot) = F^{-1}(\cdot)$

$\sqrt{n}(F^{-1}(\hat{F}(q_\tau)) - F^{-1}(F(q_\tau))) = \sqrt{n}(\hat{q}_\tau - q_\tau)$

[** notez qu'il y a un peu de main dans la dernière étape car , mais ils sont asymptotiquement égaux s'ils sont fastidieux à montrer **]$F^{-1}(\hat{F}(q_\tau)) \neq \hat{F}^{-1}(\hat{F}(q_\tau)) = \hat{q}_\tau$

Maintenant, appliquez la méthode delta mentionnée ci-dessus.

Puisque (fonction inverse théorème)$\frac{\textrm{d}}{\textrm{d}x} F^{-1}(x) = \frac{1}{f(F^{-1}(x))}$

$\sqrt{n}(\hat{q}_\tau - q_\tau) \rightarrow N\left(0, \frac{F(q_\tau)(1-F(q_\tau))}{f(F^{-1}(F(q_\tau)))^2}\right) = N\left(0, \frac{F(q_\tau)(1-F(q_\tau))}{f(q_\tau)^2}\right)$

Ensuite, pour construire l'intervalle de confiance, nous devons calculer l'erreur standard en branchant des homologues échantillons de chacun des termes de la variance ci-dessus:

Résultat

Donc$se(\hat{q}_\tau) = \sqrt{\frac{\hat{F}(\hat{q}_\tau)(1-\hat{F}(\hat{q}_\tau))}{n \hat{f}(\hat{q}_\tau)^2}} =$ $\sqrt{\frac{\tau (1 - \tau)}{n \hat{f}(\hat{q}_\tau)^2}}$

Et$CI_{0.95}(\hat{q}_\tau) = \hat{q}_\tau \pm 1.96 se(\hat{q}_\tau)$

Cela vous obligera à estimer la densité de , mais cela devrait être assez simple. Alternativement, vous pouvez également démarrer le CI assez facilement.$X$