Distribución binomial: ¿estimación del intervalo de confianza sin la media?

Question

Esta pregunta es probablemente fácil, pero no pude encontrar la respuesta, ni recordar mis clases de estadística.

Tengo una bolsa (infinita) de fichas rojas (A) y azules (B), es decir

$P(A) = p = 1 - P(B)$

Quiero estimar el mínimo de fichas rojas $k$ que debería esperar en una muestra aleatoria de tamaño $m$ para un nivel de confianza dado $\alpha$ es decir $P(k<=x<=m) = \alpha$ donde $x$ es el número de fichas rojas de la muestra total $m$

Conozco la respuesta al problema para un conocido $p$ valor (utilizado esta conferencia ), salvo que en este caso no lo sé. Lo he estimado a partir de una muestra grande $n$ donde tomé $j$ fichas rojas, es decir $\hat{p}=\frac{j}{n}$

¿Cuál es la fórmula para vincular $k = f(m, \alpha, j, n)$ o $\alpha = f(m, k, j, n)$ considerando que la probabilidad subyacente es una estimación con incertidumbre?

Answer 1

Su situación es que $X \sim \mathcal{Binomial}(n,p)$ se observa, y en base a eso se quiere predecir de alguna manera el futuro $Y \sim \mathcal{Binomial}(m,p)$ . $X, Y$ son estocásticamente independientes, pero están vinculadas por compartir el mismo valor de $p$ . Se ha propuesto encontrar una densidad predictiva para $Y$ sustituyendo la MLE de $p$ basado en $X$ , $x/n$ pero señala que eso olvida la incertidumbre de estimación de $p$ . ¿Qué hacer?

En primer lugar, una solución bayesiana. Supongamos una prior conjugada (por simplicidad, el principio será el mismo para cualquier prior), es decir, una distribución beta. Lo hacemos en dos etapas, primero hallamos la distribución posterior de $p$ dado que $X=x$ y lo utilizamos como una probabilidad a priori que combinamos con la probabilidad basada en $Y$ . $$ \pi(p \mid x)=\mathcal{beta}(\alpha+x,\beta+n-x) $$ A continuación hallamos la densidad predictiva posterior (bueno, en este caso masa de probabilidad) de $Y=y$ integrando $p$ : $$ f_\text{pred}(y \mid x)=\frac{\binom{m}{y}p^y(1-p)^{m-y}\mathcal{beta}(p,\alpha+x-1,n+\beta-x-1)}{\int_ 0^1 \binom{m}{y}p^y(1-p)^{m-y}\mathcal{beta}(p,\alpha+x-1,n+\beta-x-1)\; dp} \\ =\frac{\binom{m}{y}B(\alpha+x+y,n+m+\beta-x-y)}{B(\alpha+x,n+\beta-x)} $$ donde $B$ es la función beta. He omitido el álgebra intermedia ya que se trata de manipulaciones estándar.

Pero, ¿qué hacer si no queremos recurrir a Bayes? Existe el concepto de probabilidad predictiva se ha realizado una revisión Probabilidad predictiva: Una revisión de Jan F. Bjørnstad, en Statist. Sci. 5(2): 242-254 (mayo de 1990). DOI: 10.1214/ss/1177012175 . Esto es algo más complicado y controvertido, ya que existen múltiples versiones y no una unicidad esencial como en el caso de la verosimilitud paramétrica. Empezaremos con la verosimilitud conjunta de $y,p$ basado en $X=x$ , pero entonces debemos encontrar una manera de eliminar $p$ . Veremos tres formas (hay más):

1. sustituyendo $p$ su MLE, estimador de máxima verosimilitud. Esto es lo que has hecho.

2. Perfilar, es decir, eliminar $p$ llevándolo al máximo.

3. Condicionado a una estadística suficiente mínima

A continuación, las soluciones.

1. $$f_0(y\mid x)=\binom{m}{y}\hat{p}_x^y (1-\hat{p}_x)^{m-y}; \quad \hat{p}_x=x/n $$

2. El mle de $p$ basado en $X, Y$ es $\hat{p}_y=\frac{x+y}{n+m}$ . Sustituyendo esto se obtiene la probabilidad predictiva del perfil $$ f_\text{prof}(y \mid x) = \binom{m}{y} \hat{p}_y^y (1-\hat{p}_y)^{m-y} $$ Obsérvese que esta probabilidad de perfil no suma necesariamente 1, por lo que para utilizarla como densidad de predicción es habitual renormalizar.

3. La estadística mínima suficiente para $p$ basada en la muestra conjunta $X,Y$ es $X+Y$ . Calculando la probabilidad condicional $\DeclareMathOperator{\P}{\mathbb{P}} \P(Y=y \mid X+Y=x+y)$ por suficiencia el parámetro desconocido $p$ se eliminará, y encontraremos una distribución hipergeométrica $$ f_\text{cond}(y \mid x) = \frac{\binom{m}{y}\binom{n}{x}}{\binom{n+m}{x+y}}$$ Tenga en cuenta que esto también necesitará renormalización para ser utilizado directamente como una densidad. Obsérvese que esto coincide con la solución conjugada bayesiana para el caso $\alpha=1, \beta=1$ .

Veamos algunos ejemplos numéricos (en los gráficos las densidades están renormalizadas, la solución bayesiana mostrada es para la prior Jeffrey's $\alpha=1/2, \beta=1/2$ ):

Primero con $n=m=10; x=3$ :

A continuación, un ejemplo con mayor $n$ por lo que se estima con mayor precisión $p$ :

A continuación se muestra el código del último gráfico:

p_0 <- function(n, x, m) function(y, log=FALSE) {
    phat <- x/n
    dbinom(y, m, phat, log=log)
}
p_prof <- function(n, x, m) function(y, log=FALSE) {
    phat <- (x + y)/(n + m)
    dbinom(y, m, phat, log=log)
}
p_cond <- function(n, x, m) function(y, log=FALSE) {
    dhyper(y, m, n, x + y, log=log)
}

n <- 100; m <- 20; x <- 30

plot(0:m,  p_0(n, x, m)(0:m), col="red", main="Predictive densities",
     ylab="density", xlab="y", type="b",
     sub="n=100, m=20;  x=30")
points(0:m,  p_prof(n, x, m)(0:m)/ sum( p_prof(n, x, m)(0:m)), col="blue", type="b")
points(0:m,  p_cond(n, x, m)(0:m)/ sum( p_cond(n, x, m)(0:m)), col="orange", type="b")
legend("topright", c("naive", "profile", "cond"),
       col=c("red", "blue", "orange"),
       text.col=c("red", "blue", "orange"), lwd=2)

A continuación, pueden construirse intervalos de predicción basados en estas probabilidades de predicción.