Cuantiles de la combinación de distribuciones normales

Question

Tengo información sobre las distribuciones de las dimensiones antropométricas (como la envergadura de los hombros) para niños de diferentes edades. Para cada edad y dimensión, tengo la media y la desviación estándar. (También tengo ocho cuantiles, pero no creo que pueda obtener lo que quiero de ellos).

Para cada dimensión, me gustaría estimar determinados cuantiles de la distribución de longitudes. Si asumo que cada una de las dimensiones se distribuye normalmente, puedo hacerlo con las medias y las desviaciones estándar. ¿Existe una fórmula bonita que pueda utilizar para obtener el valor asociado a un cuantil particular de la distribución?

La operación inversa es bastante fácil: para un valor concreto, obtén el área a la derecha del valor para cada una de las distribuciones normales (edades). Suma los resultados y divídelos por el número de distribuciones.

Actualización : Aquí está la misma pregunta en forma de gráfico. Supongamos que cada una de las distribuciones de color se distribuye normalmente.

Además, obviamente puedo probar un montón de longitudes diferentes y seguir cambiándolas hasta que consiga una que se acerque lo suficiente al cuantil deseado para mi precisión. Me pregunto si hay una manera mejor que esta. Y si este es el enfoque correcto, ¿hay un nombre para él?

Answer 1

Desgraciadamente, la función cuantil normal (a partir de la cual se pueden determinar todas las demás, ya que la normal es una familia de escala de localización) no admite una forma cerrada (es decir, una "fórmula bonita"). Lo más parecido a una forma cerrada es que la función cuantil normal estándar es la función $w$ que satisface la ecuación diferencial

$$ \frac{d^2 w}{d p^2} = w \left(\frac{d w}{d p}\right)^2 $$

y las condiciones iniciales $w(1/2) = 0$ y $w'(1/2) = \sqrt{2 \pi}$ . En la mayoría de los entornos informáticos existe una función que calcula numéricamente la función cuantil normal. En R, se escribiría

qnorm(p, mean=mu, sd=sigma)

para obtener el $p$ 'th quantile of the $N(\mu, \sigma^2)$ distribución.

---

Editar: Con una comprensión modificada del problema, los datos se generan a partir de una mezcla de normales, de modo que la densidad de los datos observados es:

$$ p(x) = \sum_{i} w_{i} p_{i}(x) $$

donde $\sum_{i} w_{i} = 1$ y cada $p_{i}(x)$ es una densidad normal con media $\mu_{i}$ y la desviación estándar $\sigma_{i}$ . Se deduce que la FCD de los datos observados es

$$ F(y) = \int_{-\infty}^{y} \sum_{i} w_{i} p_{i}(x) dx = \sum_{i} w_{i} \int_{-\infty}^{y} p_{i}(x) = \sum_{i} w_{i} F_{i}(y) $$

donde $F_{i}(x)$ es la FCD normal con media $\mu_{i}$ y la desviación estándar $\sigma_{i}$ . La integración y la suma pueden intercambiarse porque estas integrales son finitas. Esta FCD es continua y bastante fácil de calcular en un ordenador, por lo que la FCD inversa, $F^{-1}$ , también conocida como función cuantílica, se puede calcular haciendo una búsqueda de líneas. Opto por esta opción porque no se me ocurre ninguna fórmula sencilla para la función cuantil de una mezcla de normales, como función de los cuantiles de las distribuciones constituyentes.

El siguiente código R calcula numéricamente $F^{-1}$ utilizando la bisección para la búsqueda de líneas. La función F_inv() es la función de cuantiles, es necesario suministrar el vector que contiene cada $w_{i}, \mu_{i}, \sigma_{i}$ y el cuantil a resolver, $p$ .

# evaluate the function at the point x, where the components 
# of the mixture have weights w, means stored in u, and std deviations
# stored in s - all must have the same length.
F = function(x,w,u,s) sum( w*pnorm(x,mean=u,sd=s) )

# provide an initial bracket for the quantile. default is c(-1000,1000). 
F_inv = function(p,w,u,s,br=c(-1000,1000))
{
   G = function(x) F(x,w,u,s) - p
   return( uniroot(G,br)$root ) 
}

#test 
# data is 50% N(0,1), 25% N(2,1), 20% N(5,1), 5% N(10,1)
X = c(rnorm(5000), rnorm(2500,mean=2,sd=1),rnorm(2000,mean=5,sd=1),rnorm(500,mean=10,sd=1))
quantile(X,.95)
    95% 
7.69205 
F_inv(.95,c(.5,.25,.2,.05),c(0,2,5,10),c(1,1,1,1))
[1] 7.745526

# data is 20% N(-5,1), 45% N(5,1), 30% N(10,1), 5% N(15,1)
X = c(rnorm(5000,mean=-5,sd=1), rnorm(2500,mean=5,sd=1),
      rnorm(2000,mean=10,sd=1), rnorm(500, mean=15,sd=1))
quantile(X,.95)
     95% 
12.69563 
F_inv(.95,c(.2,.45,.3,.05),c(-5,5,10,15),c(1,1,1,1))
[1] 12.81730