AJO_1415_Resumen_estimación_parámetros

1002 days ago by sevillad1415

Estimación de los parámetros estadísticos de una distribución normal

Supongamos que la variable $X$, que toma los valores resultantes de un cierto proceso de medición, sigue una distribución normal $N(\mu,\sigma)$. Casi siempre $\mu$ es desconocido (es el valor real o ideal a medir). En algunas situaciones conoceremos $\sigma$ (por las especificaciones del aparato de medida, por ejemplo), en otras no. Queremos estimar los parámetros de $X$ a partir de un conjunto finito de valores (una muestra). En esta hoja ilustraremos varios casos por medio de ejemplos, que serán generados artificialmente a modo de ilustración. En lugar de ejecutarlos sin más, estudia las fórmulas teóricas y su contrapartida en Sage.

En los casos en que hagamos estimaciones por intervalo de confianza, por defecto trabajamos con probabilidad 90%. En concreto, esto significa que de cada 10 ejecuciones del código, en 9 de ellas el valor real está contenido en el intervalo de estimación, y 1 de cada 10 ese intervalo no contiene el valor real. Puedes experimentar tú mismo/a.

 
       

Estimación puntual de la media

La estimación puntual de $\mu$ que usamos es la media muestral, $\overline{x}$. A continuación generamos varios datos aleatoriamente a partir de una distribución normal arbitraria. Puedes ejecutar el código varias veces para ver cómo cambia la estimación. También puedes cambiar los parámetros.

mu = 234 sigma = 2 n = 5 T = RealDistribution('gaussian',sigma) d = [mu + T.get_random_element() for i in (1..n)] d 
       
[237.192550675, 234.073122209, 235.623864404, 235.89263184,
233.351143727]
[237.192550675, 234.073122209, 235.623864404, 235.89263184, 233.351143727]
sum(d)/n 
       
235.226662571
235.226662571
 
       

Estimación puntual de la varianza

La estimación puntual de $\sigma^2$ que usamos es la cuasivarianza muestral, $S^2$. La varianza muestral no es un buen indicador: se queda corto, en general, por un factor $\frac{n-1}{n}$ (compárense las fórmulas de la varianza y cuasivarianza muestrales).

A continuación un ejemplo de estimación puntual de la varianza, hecho "a mano" y también con la función de Sage correspondiente.

mu = 234 sigma = 2 n = 5 T = RealDistribution('gaussian',sigma) d = [mu + T.get_random_element() for i in (1..n)] d 
       
[232.313391424, 234.063658988, 235.80399119, 236.21076118,
237.601227227]
[232.313391424, 234.063658988, 235.80399119, 236.21076118, 237.601227227]
xbarra = mean(d) sum([(x_i-xbarra)^2 for x_i in d])/(n-1); variance(d) 
       
4.19402649085
4.19402649085
4.19402649085
4.19402649085
 
       

Estimación por intervalo de la media (caso de varianza conocida)

Supongamos que conocemos $\sigma^2$. Hemos visto que la variable \[ \overline{X}=\frac{X_1+\cdots+X_n}{n} \] sigue una distribución $N(\mu,\frac{\sigma}{\sqrt{n}})$. Estandarizando, la nueva variable \[ Y=\frac{\overline{X}-\mu}{\sigma/\sqrt{n}} \] sigue una distribución $N(0,1)$. Usando el valor de $\overline{X}$ dado por la muestra, del intervalo de confianza que calculemos para $Y$ despejamos $\mu$ y sacamos su intervalo (nótese que en el último paso al despejar una desigualdad del tipo $a \leq -\mu \leq b$ se cambiar el orden y se obtiene $-b\leq \mu \leq -a$).

mu = 234 sigma = 2 n = 5 T = RealDistribution('gaussian',sigma) d = [mu + T.get_random_element() for i in (1..n)] d 
       
[228.519691401, 236.508351522, 234.444911781, 234.493869319,
232.817441619]
[228.519691401, 236.508351522, 234.444911781, 234.493869319, 232.817441619]
xbarra = mean(d) xbarra 
       
233.356853128
233.356853128
conf = 0.90 a_Y = r.qnorm((1-conf)/2)._sage_() b_Y = -a_Y a_Y, b_Y 
       
(-1.64485362695147, 1.64485362695147)
(-1.64485362695147, 1.64485362695147)
a_mu = -b_Y * (sigma/sqrt(n).n()) + xbarra b_mu = -a_Y * (sigma/sqrt(n).n()) + xbarra a_mu, b_mu 
       
(231.885651319329, 234.828054937649)
(231.885651319329, 234.828054937649)
 
       

Estimación por intervalo de la media (caso de varianza desconocida)

Si no conocemos $\sigma^2$, no podemos despejar $\mu$ en la fórmula de arriba. Así que lo hacemos de manera ligeramente distinta: usando \[ Y=\frac{\overline{X}-\mu}{S/\sqrt{n}} \] donde ahora sustituimos la desviación típica $\sigma$, que nos es desconocida, por la estimación $S$, la raíz cuadrada de la nueva variable \[ S^2 = \frac{1}{n-1}\sum (X_i-\overline{X}) . \] Lo que cambia es que ahora la distribución de $Y$ no es la normal $N(0,1)$ sino una t de Student con $n-1$ grados de libertad. A modo de ilustración aquí tienes una comparación animada entre $N(0,1)$ y las $t$ de Student para distintos grados de libertad.

pnormal = plot(RealDistribution('gaussian', 1),xmin=-4,xmax=4,color="black",figsize=6) frames = [plot(RealDistribution('t', i),xmin=-4,xmax=4,color="blue",legend_label=str(i)+" grados de libertad")+pnormal for i in [1..20]] animate(frames).show(delay=20) 
       

De nuevo, calcularemos un intervalo de confianza para $Y$ con los valores de $\overline{X}$ y $S$ de la muestra, y despejaremos $\mu$.

mu = 234 sigma = 2 n = 5 T = RealDistribution('gaussian',sigma) d = [mu + T.get_random_element() for i in (1..n)] xbarra = mean(d) s2 = variance(d) xbarra; s2 
       
233.808484212
2.53714735672
233.808484212
2.53714735672
conf = 0.90 a_Y = r.qt((1-conf)/2, n-1)._sage_() b_Y = -a_Y a_Y, b_Y 
       
(-2.13184678632665, 2.13184678632665)
(-2.13184678632665, 2.13184678632665)
a_mu = -b_Y * (sqrt(s2)/sqrt(n).n()) + xbarra b_mu = -a_Y * (sqrt(s2)/sqrt(n).n()) + xbarra a_mu, b_mu 
       
(232.289882683, 235.327085741)
(232.289882683, 235.327085741)
 
       

Estimación por intervalo de la varianza

Usamos que la variable \[ Y=\frac{(n-1)S^2}{\sigma^2} \] sigue una distribución chi cuadrado ($\chi^2$) con $n-1$ grados de libertad. A continuación las gráficas de varias distribuciones chi cuadrado.

frames = [plot(RealDistribution('chisquared', i),xmin=0,xmax=30,ymax=1,aspect_ratio=10,color="blue",legend_label=str(i)+" grados de libertad") for i in [1..20]] animate(frames).show(delay=10) 
       

Del intervalo de confianza para $Y$ y el valor de $S^2$ obtenido de la muestra despejamos $\sigma^2$ y obtenemos su intervalo de confianza. A la hora de despejar usamos que una desigualdad del tipo $a \leq \frac{1}{\sigma^2} \leq b$ se transforma en $\frac1b \leq \sigma^2 \leq \frac1a$.

mu = 234 sigma2 = 2 n = 5 T = RealDistribution('gaussian',sigma) d = [mu + T.get_random_element() for i in (1..n)] s2 = variance(d) s2 
       
3.72343328294
3.72343328294
conf = 0.90 a_Y = r.qchisq((1-conf)/2, n-1)._sage_() b_Y = r.qchisq(1-(1-conf)/2, n-1)._sage_() a_Y, b_Y 
       
(0.710723021397324, 9.48772903678115)
(0.710723021397324, 9.48772903678115)
a_sigma2 = (n-1)*s2/b_Y b_sigma2 = (n-1)*s2/a_Y a_sigma2, b_sigma2 
       
(1.56978904794, 20.9557488408)
(1.56978904794, 20.9557488408)
 
       

Estimación por intervalo del cociente de varianzas

Tenemos dos variables $X,Y$ que siguen distribuciones normales, cuyas varianzas son $\sigma_X^2$ y $\sigma_Y^2$ respectivamente. Usamos que la variable \[ Z=\frac{S_X^2/\sigma_X^2}{S_Y^2/\sigma_Y^2} \] sigue una distribución F de Snedecor con $n_X-1$ y $n_Y-1$ grados de libertad. Del intervalo de confianza para $Z$ y los valores muestrales de $S_X^2,S_Y^2$ despejamos $\sigma_X^2/\sigma_Y^2$ y obtenemos su intervalo de confianza de manera parecida al caso de arriba.

mu_X = 234 sigma_X = 2 n_X = 5 T_X = RealDistribution('gaussian',sigma_X) d_X = [mu_X + T_X.get_random_element() for i in (1..n_X)] mu_Y = 157 sigma_Y = 3 n_Y = 15 T_Y = RealDistribution('gaussian',sigma_Y) d_Y = [mu_Y + T_Y.get_random_element() for i in (1..n_Y)] s2_X = variance(d_X) s2_Y = variance(d_Y) sigma_X/sigma_Y.n(); n_X, sigma_X, s2_X; n_Y, sigma_Y, s2_Y 
       
0.666666666666667
(5, 2, 2.55138724205)
(15, 3, 11.0639007663)
0.666666666666667
(5, 2, 2.55138724205)
(15, 3, 11.0639007663)
conf = 0.90 a_Z = r.qf((1-conf)/2, n_X-1, n_Y-1)._sage_() b_Z = r.qf(1-(1-conf)/2, n_X-1, n_Y-1)._sage_() a_Z, b_Z 
       
(0.170260692120781, 3.11224984796139)
(0.170260692120781, 3.11224984796139)
a_cociente = s2_X/s2_Y * 1/b_Z b_cociente = s2_X/s2_Y * 1/a_Z a_cociente, b_cociente 
       
(0.0740958105582, 1.35442110725)
(0.0740958105582, 1.35442110725)
 
       

Estimación por intervalo de la diferencia de medias (caso de observaciones independientes y varianzas conocidas)

Tenemos dos variables $X,Y$ que siguen distribuciones normales, con medias $\mu_X$ y $\mu_Y$ y varianzas $\sigma_X^2$ y $\sigma_Y^2$ respectivamente. Como la variable $\overline{X} - \overline{Y}$ sigue una distribución normal de media $\mu_X - \mu_Y$ y varianza $\frac{\sigma_X^2}{n_X} + \frac{\sigma_Y^2}{n_Y}$, tenemos que \[ Z = \frac{(\overline{X} - \overline{Y}) - (\mu_X - \mu_Y)}{\sqrt{\sigma_X^2/n_X + \sigma_Y^2/n_Y}} \] sigue una distribución $N(0,1)$. Del intervalo de confianza para $Z$ y los distintos valores muestrales despejamos $\mu_X - \mu_Y$ y obtenemos su intervalo de confianza.

mu_X = 234 sigma_X = 2 n_X = 5 T_X = RealDistribution('gaussian',sigma_X) d_X = [mu_X + T_X.get_random_element() for i in (1..n_X)] mu_Y = 157 sigma_Y = 3 n_Y = 15 T_Y = RealDistribution('gaussian',sigma_Y) d_Y = [mu_Y + T_Y.get_random_element() for i in (1..n_Y)] xbarra = mean(d_X) ybarra = mean(d_Y) s2_X = variance(d_X) s2_Y = variance(d_Y) mu_X-mu_Y; n_X, xbarra, s2_X; n_Y, ybarra, s2_Y 
       
77
(5, 233.641330391, 6.26955215289)
(15, 157.043967167, 4.47165138349)
77
(5, 233.641330391, 6.26955215289)
(15, 157.043967167, 4.47165138349)
conf = 0.90 a_Z = r.qnorm((1-conf)/2)._sage_() b_Z = -a_Z a_Z, b_Z 
       
(-1.64485362695147, 1.64485362695147)
(-1.64485362695147, 1.64485362695147)
denom = sqrt(s2_X/n_X + s2_Y/n_Y) a_diferencia = -b_Z*denom + (xbarra-ybarra) b_diferencia = -a_Z*denom + (xbarra-ybarra) a_diferencia, b_diferencia 
       
(74.5482026717, 78.6465237748)
(74.5482026717, 78.6465237748)