Desviación estándar de observaciones agrupadas

Tengo un conjunto de datos de observaciones de muestra, almacenadas como recuentos dentro de intervalos. p.ej:

min/max  count
40/44    1
45/49    2
50/54    3
55/59    4
70/74    1

Ahora, encontrar una estimación del promedio de esto es bastante sencillo. Simplemente use la media (o mediana) de cada intervalo de rango como observación y el recuento como un peso y encuentre el promedio ponderado:

$$\bar{x}^* = \frac{1}{\sum_{i=1}^N w_i} \sum_{i=1}^N w_ix_i$$

Para mi caso de prueba, esto me da 53.82.

Mi pregunta ahora es, ¿cuál es el método correcto para encontrar la desviación estándar (o varianza)?

A través de mi búsqueda, encontré varias respuestas, pero no estoy seguro de cuál, si alguna, es realmente apropiado para mi conjunto de datos. Pude encontrar la siguiente fórmula tanto en otra pregunta aquí como en un documento NIST aleatorio .

$$s^{2*} = \frac{ \sum_{i=1}^N w_i (x_i - \bar{x}^*)^2 }{ \frac{(M-1)}{M} \sum_{i=1}^N w_i }$$

Lo que da una desviación estándar de 8.35 para mi caso de prueba. Sin embargo, el artículo de Wikipedia sobre medios ponderados ofrece la fórmula:

$$s^{2*} = \frac{ \sum_{i=1}^N w_i}{(\sum_{i=1}^N w_i)^2 - \sum_{i=1}^N w_i^2} \sum_{i=1}^N w_i(x_i-\bar{x}^*)^2$$

y

$$s^{2*} = \frac{1}{(\sum_{i=1}^N w_i) - 1} \sum_{i=1}^N w_i(x_i-\bar{x}^*)^2$$

Que dan desviaciones estándar de 8.66 y 7.83, respectivamente, para mi caso de prueba.

Actualizar

Gracias a @whuber que sugirió investigar las correcciones de Sheppard y sus útiles comentarios relacionados con ellas. Desafortunadamente, me está costando entender los recursos que puedo encontrar al respecto (y no puedo encontrar ningún buen ejemplo). Sin embargo, para resumir, entiendo que lo siguiente es una estimación sesgada de la varianza:

$$s^{2*} = \frac{1}{\sum_{i=1}^N w_i} \sum_{i=1}^N w_i(x_i-\bar{x}^*)^2$$

También entiendo que la mayoría de las correcciones estándar para el sesgo son para muestras aleatorias directas de una distribución normal. Por lo tanto, veo dos problemas potenciales para mí:

1. Estas son muestras aleatorias agrupadas (que, estoy bastante seguro, es donde entran las correcciones de Sheppard).
2. No se sabe si los datos son o no para una distribución normal (por lo tanto, supongo que no, lo que, estoy bastante seguro, invalida las correcciones de Sheppard).

Entonces, mi pregunta actualizada es; ¿Cuál es el método apropiado para manejar el sesgo impuesto por la fórmula de desviación / varianza estándar ponderada "simple" en una distribución no normal? Más específicamente con respecto a los datos agrupados.

Nota: estoy usando los siguientes términos:

• es la varianza ponderada $s^{2*}$
• es el número de observaciones. (es decir, el número de contenedores)$N$
• es el número de pesos distintos de cero. (es decir, el número de contenedores con recuentos)$M$
• son los pesos (es decir, los recuentos)$w_i$
• son las observaciones. (es decir, el bin significa)$x_i$
• es la media ponderada.$\bar{x}^*$

Esta respuesta presenta dos soluciones: las correcciones de Sheppard y una estimación de máxima verosimilitud. Ambos coinciden estrechamente en una estimación de la desviación estándar: para el primero y 7,69 para el segundo (cuando se ajusta para ser comparable al estimador "imparcial" habitual).$7.70$$7.69$

---

Correcciones de Sheppard

Las "correcciones de Sheppard" son fórmulas que ajustan los momentos calculados a partir de datos agrupados (como estos) donde

• se supone que los datos se rigen por una distribución compatible con un intervalo finito $[a,b]$

• ese intervalo se divide secuencialmente en contenedores iguales de ancho común que es relativamente pequeño (ningún contenedor contiene una gran proporción de todos los datos)$h$

• La distribución tiene una función de densidad continua.

Se derivan de la fórmula de suma de Euler-Maclaurin, que se aproxima a integrales en términos de combinaciones lineales de valores del integrando en puntos regularmente espaciados, y por lo tanto generalmente aplicable (y no solo a distribuciones normales).

Aunque estrictamente hablando, una distribución Normal no se admite en un intervalo finito, sí lo es en una aproximación extremadamente cercana. Esencialmente, toda su probabilidad está contenida dentro de siete desviaciones estándar de la media. Por lo tanto, las correcciones de Sheppard son aplicables a los datos que se supone provienen de una distribución Normal.

Las dos primeras correcciones de Sheppard son

1. Utilice la media de los datos agrupados para la media de los datos (es decir, no se necesita corrección para la media).

2. $h^2/12$

$h^2/12$$h$$-h/2$$h/2$$h^2/12$

Hagamos los cálculos. Solía Rilustrarlos, comenzando especificando los recuentos y los contenedores:

counts <- c(1,2,3,4,1)
bin.lower <- c(40, 45, 50, 55, 70)
bin.upper <- c(45, 50, 55, 60, 75)

La fórmula adecuada para usar para los recuentos proviene de la replicación de los anchos de la papelera por las cantidades dadas por los recuentos; es decir, los datos agrupados son equivalentes a

42.5, 47.5, 47.5, 52.5, 52.5, 57.5, 57.5, 57.5, 57.5, 72.5

$x$$k$$kx^2$

bin.mid <- (bin.upper + bin.lower)/2
n <- sum(counts)
mu <- sum(bin.mid * counts) / n
sigma2 <- (sum(bin.mid^2 * counts) - n * mu^2) / (n-1)

mu$1195/22 \approx 54.32$sigma2$675/11 \approx 61.36$$7.83$$h=5$$h^2/12 = 25/12 \approx 2.08$$\sqrt{675/11 - 5^2/12} \approx 7.70$

---

Estimaciones de máxima verosimilitud

$F_\theta$$\theta$$(x_0, x_1]$$k$$F_\theta$

$$\log \prod_{i=1}^k \left(F_\theta(x_1) - F_\theta(x_0)\right) = k\log\left(F_\theta(x_1) - F_\theta(x_0)\right)$$

(Ver MLE / Probabilidad de intervalo lognormalmente distribuido ).

$\Lambda(\theta)$$\hat\theta$$-\Lambda(\theta)$$\theta$R

sigma <- sqrt(sigma2) # Crude starting estimate for the SD
likelihood.log <- function(theta, counts, bin.lower, bin.upper) {
  mu <- theta[1]; sigma <- theta[2]
  -sum(sapply(1:length(counts), function(i) {
    counts[i] * 
      log(pnorm(bin.upper[i], mu, sigma) - pnorm(bin.lower[i], mu, sigma))
  }))
}
coefficients <- optim(c(mu, sigma), function(theta) 
  likelihood.log(theta, counts, bin.lower, bin.upper))$par

$(\hat\mu, \hat\sigma) = (54.32, 7.33)$

$\sigma$$n/(n-1)$$\sigma$$\sqrt{n/(n-1)} \hat\sigma = \sqrt{11/10}\times 7.33 = 7.69$$7.70$

Verificando los supuestos

Para visualizar estos resultados, podemos trazar la densidad normal ajustada sobre un histograma:

hist(unlist(mapply(function(x,y) rep(x,y), bin.mid, counts)),
     breaks = breaks, xlab="Values", main="Data and Normal Fit")
curve(dnorm(x, coefficients[1], coefficients[2]), 
      from=min(bin.lower), to=max(bin.upper), 
      add=TRUE, col="Blue", lwd=2)

$11$

$\chi^2$$\chi^2$R

breaks <- sort(unique(c(bin.lower, bin.upper)))
fit <- mapply(function(l, u) exp(-likelihood.log(coefficients, 1, l, u)),
              c(-Inf, breaks), c(breaks, Inf))
observed <- sapply(breaks[-length(breaks)], function(x) sum((counts)[bin.lower <= x])) -
  sapply(breaks[-1], function(x) sum((counts)[bin.upper < x]))
chisq.test(c(0, observed, 0), p=fit, simulate.p.value=TRUE)

La salida es

Chi-squared test for given probabilities with simulated p-value (based on 2000 replicates)

data:  c(0, observed, 0) 
X-squared = 7.9581, df = NA, p-value = 0.2449

$0.245$