4  Intervalos de confianza

… la verdad, a veces, solo aporta una profunda soledad.
- Entonces, no quiere saber la verdad, sino contemplar su retrato a diario y pensar en una posibilidad. ¿Se trata solo de eso?
- Exacto -asintió Menshiki-. En lugar de una verdad inamovible elijo una posibilidad con margen de variación. Elijo encomendarme a esa variación.

Haruki Murakami, La muerte del comendador (libro 1)

4.1 ¿Qué es un intervalo de confianza?

El objetivo de este tema es utilizar una muestra para obtener un intervalo tal que podemos confiar en que contiene al verdadero valor de un parámetro desconocido. Veremos cómo interpretar esa confianza un poco más adelante.

Comenzamos con una definición formal de lo que buscamos: consideramos una muestra X1,,Xn de una variable aleatoria con función de distribución G(;θ), siendo θΘR un parámetro desconocido.

Definición 4.1 Sean dos estadísticos a(X1,,Xn) y b(X1,,Xn) con a(X1,,Xn)<b(X1,,Xn) c.s. y un valor α(0,1). Supongamos que se verifica Pθ{a(X1,,Xn)<θ<b(X1,,Xn)}=1α, para todo  θΘ. Entonces para una realización concreta de la muestra, x1,,xn, se dice que [a(x1,,xn),b(x1,,xn)] es un intervalo de confianza para θ con nivel de confianza 1α. Lo denotaremos IC1α(θ).

4.2 Intervalo de confianza para la media

Veamos cómo construir un intervalo de confianza para la media poblacional. Con el fin de comprender mejor las ideas tratamos primero una situación sencilla aunque no muy realista.

Población normal con varianza conocida

Supongamos que X1,,Xn son v.a.i.i.d. N(μ,σ2), donde μ es un parámetro desconocido y σ2 es conocida (por simplificar, un poco más adelante trataremos el caso en el que la varianza es desconocida). Sabemos que X¯N(μ,σ2n),  y, estandarizando,  X¯μσ/nN(0,1). Dado α(0,1), zα denota el cuantil 1α en la normal estándar (es decir, Φ(zα)=1α, donde Φ es la función de distribución de la N(0,1)). Entonces, Pμ{zα/2<X¯μσ/n<zα/2}=1α y, despejando μ, Pμ{X¯zα/2σn<μ<X¯+zα/2σn}=1α. Se concluye que IC1α(μ)=(x¯zα/2σn,x¯+zα/2σn) es un intervalo de confianza de nivel 1α para μ.

Interpretación frecuentista

Fijamos, por ejemplo, 1α=0.95 y extraemos muchas muestras de tamaño n de una distribución N(μ,σ2). Para cada una de las muestras aplicamos la fórmula anterior y calculamos el correspondiente intervalo. Entonces, aproximadamente el 95% de estos intervalos de confianza contiene al verdadero valor del parámetro, μ. En la práctica, tendremos una única muestra y no podremos saber si contiene o no a μ. Sin embargo, confiamos en que esa única muestra está dentro del 95% de las muestras buenas que producen intervalos que sí contienen al parámetro.

En la siguiente simulación, extraemos m=100 muestras de tamaño n=30 de una población normal estándar y determinamos el número de ellas en las que el correspondiente intervalo de nivel 0.95 (calculado con la fórmula que hemos visto) contiene al verdadero valor del parámetro. Se observa que hay tres muestras de las 100 que dan lugar a intervalos que no contienen al verdadero valor del parámetro μ=0, lo que no coincide pero se aproxima al valor nominal del 5 % de veces en las que el intervalo no contiene al parámetro.

Margen de error

Al radio del intervalo de confianza anterior se le suele llamar margen de error, E=zα/2σ/n. El margen de error depende de:

  • El tamaño muestral. Para un nivel de confianza fijo, a mayor tamaño muestral menor margen de error.
  • El nivel de confianza. Si fijamos n, a mayor nivel de confianza mayor es también el margen de error. Como queremos estar más seguros de no fallar y no tenemos más información, la única posibilidad es hacer más grande el intervalo.
  • La homogeneidad de la población. Cuanto menor sea σ, menor es el margen de error para valores de n y α dados. Esto es bastante intuitivo, si todos los individuos de una población son similares para cierta variable, es más fácil estimar la media poblacional de la variable.

Población normal con varianza desconocida

En la práctica σ no es conocida. En este caso, se sustituye σ por S en la estandarización de la media. Por el lema de Fisher, X¯μS/ntn1, lo que lleva de inmediato al siguiente intervalo de confianza de nivel 1α: IC1α(μ)=(x¯tn1;α/2sn,  x¯+tn1;α/2sn). La notación tn1;α/2 representa el valor que deja a su derecha una probabilidad de α/2 en la distribución t de Student con n1 grados de libertad (o sea, el cuantil 1α/2 de la distribución).

Ejemplo

Se mide el tiempo de duración (en segundos) de un proceso químico realizado 20 veces en condiciones similares, obteniéndose los siguientes resultados:

resultados <- c(93, 90, 97, 90, 93, 91, 96, 94, 91, 91, 88, 93, 95, 91, 89, 92,
87, 88, 90, 86)
mean(resultados)
#> [1] 91.25

Suponemos que los datos proceden de una distribución normal. Para calcular un intervalo de confianza de nivel 0.95 para la duración media del proceso usamos el comando t.test:

t.test(resultados)$conf.int   # 95% es el nivel por defecto
#> [1] 89.87604 92.62396
#> attr(,"conf.level")
#> [1] 0.95
t.test(resultados, conf.level = 0.9)$conf.int   # nivel 90%
#> [1] 90.11492 92.38508
#> attr(,"conf.level")
#> [1] 0.9

Cuando la población no es normal y la varianza es desconocida

Aunque la población no sea normal siempre que tenga varianza finita podemos aplicar el TCL. Se cumple X¯μσ/ndN(0,1). Sustituyendo σ por un estimador consistente σ^ y usando el lema de Slutsky tenemos X¯μσ^/ndN(0,1). De la propiedad anterior se obtiene el siguiente intervalo de confianza para μ con nivel aproximado 1α, IC1α(μ)(x¯zα/2σ^n,x¯+zα/2σ^n) Este intervalo es (aproximadamente) válido para cualquier distribución con varianza finita, siempre que n sea lo bastante grande.

Caso particular importante: intervalo para una proporción

Sean X1,,Xn v.a.i.i.d. B(1,p). Si denotamos p^=X¯, por el TCL p^pp(1p)ndN(0,1), y reemplazando p por su estimador natural p^, obtenemos que el intervalo de confianza aproximado para p es, IC1α(p)(p^zα/2p^(1p^)n,p^+zα/2p^(1p^)n).

Ejemplo

Se estima la proporción p de piezas defectuosas en la producción de una fábrica con una muestra de 200 piezas de las cuales 8 resultan ser defectuosas. Calcula un intervalo de confianza de nivel 0.95 para p.

Sustituyendo en la fórmula obtenemos IC0.95(p)=(8200±1.960.040.96200)=(0.04±0.02716)=(0.01284,0.06716).

Supongamos que este margen de error (la mitad de la longitud del IC) se considera insatisfactorio y se desea obtener un intervalo con un error de, como mucho, 0.01. ¿Qué tamaño muestral habría que elegir?

Nuestro objetivo es conseguir E=1.96p^(1p^)n0.01.

Como valor de p^ podemos tomar a modo de aproximación el obtenido en la muestra anterior (que a veces se llama muestra piloto). Entonces, E1.960.040.96n0.01. Despejando, obtenemos n1.962(0.040.960.012)=1475.17. Por tanto, habría que tomar n1476.

Cuando se quiere determinar el tamaño muestral necesario para obtener un error E y no se tiene ninguna información previa sobre el valor de p (no hay muestra piloto disponible) podemos actuar poniéndonos en el peor de los casos (es decir, el que da un intervalo de confianza más amplio), que es p^=1/2.

En el ejemplo anterior se tendría n=1.962(0.50.50.012)=9604.

Ejemplo: una ficha técnica

Los cálculos del ejemplo anterior son los que subyacen a la información de las fichas técnicas de las encuestas. Veamos un ejemplo típico:

En la ficha técnica anterior se nos dice que el tamaño muestral es n=1100 y el nivel de confianza es 1α=0.9545. Para este nivel resulta zα/2=2:

alpha <- 1-0.9545
qnorm(1-alpha/2)
#> [1] 2.000002

Por ello el intervalo es de la forma [p^2p^(1p^)/n]. Finalmente el cálculo del margen de error que da la ficha técnica corresponde al caso más desfavorable en el que p^=1/2:

n <- 1100
error <- 2*sqrt(1/(4*n))
error
#> [1] 0.03015113

4.3 El método de la cantidad pivotal

En los ejemplos anteriores los resultados que nos han permitido calcular los intervalos eran de la forma: X¯μσ/nN(0,1),  o  X¯μS/ntn1.

En general, hemos usado funciones que dependían de la muestra y del parámetro de interés θ, T(X1,,Xn;θ) pero tales que su distribución era totalmente conocida para cualquier valor del parámetro. Estas funciones se llaman cantidades pivotales o pivotes. Si identificamos una cantidad pivotal, podemos encontrar (en el caso discreto, al menos aproximadamente) dos cuantiles q1(α) y q2(α) tales que 1α=Pθ(q1(α)<T(X1,,Xn;θ)<q2(α)). Si T no es una función demasiado complicada, podremos despejar θ para obtener una región de confianza. En el caso en que T sea una función monótona de θ esta región será un intervalo.


Ejemplo

Si X1,,Xn son v.a.i.i.d. de una distribución U(0,θ), entonces T=X(n)/θ es una cantidad pivotal para el parámetro ya que se distribuye como el máximo de n observaciones i.i.d. uniformes en el intervalo (0,1). De hecho, la función de distribución de T en el intervalo (0,1) viene dada por F(x)=xn.

  • Determina dos valores a y b tales que P(aTb)=1α, y tales que además el intervalo (a,b) tenga la menor longitud posible.
  • Determina un intervalo de confianza para θ de nivel 1α a partir de los valores a y b.

4.3.1 Ejemplos de cantidades pivotales en poblaciones normales

Además de las cantidades pivotales que ya hemos visto para la media, a continuación enumeramos otras cantidades pivotales que se usan en poblaciones normales.

Diferencia de dos medias (datos independientes)

Un problema relevante desde el punto de vista aplicado es el de la comparación de las medias de dos poblaciones. Supongamos que tenemos dos muestras independientes de v.a.i.i.d.. La primera, X1,,Xn1, procede de una distribución N(μ1,σ2) mientras que la segunda Y1,,Yn2 procede de una distribución N(μ2,σ2). El objetivo es encontrar una cantidad pivotal para μ1μ2. Obsérvese que estamos suponiendo que las varianzas de ambas poblaciones son iguales. Esta es una hipótesis de simplificación habitual que recibe el nombre de homocedasticidad.

La estrategia habitual es pensar en un estimador natural para el parámetro de interés y luego modificarlo (usualmente mediante algún tipo de estandarización) de forma que la distribución del estimador estandarizado sea conocida.

El estimador natural para μ1μ2 es X¯Y¯. Ahora, como las muestras son independientes, X¯Y¯N(μ1μ2, σ2(1n1+1n2)). A partir de aquí, (X¯Y¯)(μ1μ2)σ1n1+1n2N(0,1). Si conociéramos σ2 ya tendríamos la cantidad pivotal. Como en la práctica no la conocemos es necesario estimarla a partir de las dos muestras. Las varianzas muestrales de cada de ellas S12 y S22 son estimadores distintos del mismo parámetro (debido a la hipótesis de homocedasticidad). Ambos estimadores se combinan mediante una media que da pesos proporcionales a los tamaños muestrales: Sp2=(n11)S12+(n21)S22n1+n22. Se proponen los ejercicios siguientes:

  1. Demostrar que E(Sp2)=σ2, es decir, que Sp2 es insesgado para σ2.
  2. Demostrar que (n1+n22)Sp2/σ2χn1+n222.
  3. Demostrar que, como consecuencia, (X¯Y¯)(μ1μ2)Sp1n1+1n2tn1+n22.
  4. De la cantidad pivotal para μ1μ2 anterior deducir el intervalo de confianza: IC1α(μ1μ2)=[X¯Y¯tn1+n22;α/2Sp1n1+1n2].
Ejemplo

A un grupo de 20 pollos se les suministró pienso con harina de maíz de una nueva variedad transgénica. A otro grupo de 20 pollos (grupo de control) se le alimentó con un pienso que no contenía la variedad mejorada. Ganancias de peso de los pollos (en gramos) al cabo de 21 días de alimentación:

maiz.normal <- c(380, 321, 366, 356, 283, 349, 402, 462, 356, 410, 329, 399, 
    350, 384, 316, 272, 345, 455, 360, 431)
maiz.transgenico <- c(361, 447, 401, 375, 434, 403, 393, 426, 406, 318, 467, 
    407, 427, 420, 477, 392, 430, 339, 410, 326)

Tal y como tenemos los datos aún no están ordenados del todo. Recordamos que debemos tener:

  • Todos los datos del mismo análisis en el mismo data.frame
  • Cada fila del fichero corresponde a un individuo
  • Cada columna del fichero corresponde a una variable

Para ordenar los datos usamos el código siguiente:

peso <- c(maiz.normal, maiz.transgenico)
tipo <- gl(2, 20, labels = c('normal', 'transgénico'))   # genera un factor con dos niveles
datos_maiz <- data.frame(peso, tipo)
glimpse(datos_maiz)
#> Rows: 40
#> Columns: 2
#> $ peso <dbl> 380, 321, 366, 356, 283, 349, 402, 462, 356, 410, 329, 399, 350, ~
#> $ tipo <fct> normal, normal, normal, normal, normal, normal, normal, normal, n~

Podemos usar diagramas de cajas para comparar gráficamente las dos muestras:

ggplot(datos_maiz) +
  geom_boxplot(aes(x=tipo, y=peso), fill='olivedrab4')

Parece que los pollos alimentados con un pienso normal ganaron menos peso que los alimentados con el pienso transgénico. Finalmente, calculamos el intervalo de confianza para μ1μ2, donde μ1 es la ganancia media de peso con pienso normal y μ2 es la ganancia media de peso con pienso transgénico. Usamos de nuevo el comando t.test:

t.test(peso ~ tipo, data=datos_maiz, var.equal=TRUE)$conf.int
#> [1] -66.700161  -6.599839
#> attr(,"conf.level")
#> [1] 0.95

El hecho de que todos los valores del intervalo sean negativos apoya la afirmación de que μ1<μ2.

Diferencia de dos medias (datos emparejados)

Se observa una muestra (X1,Y1),,(Xn,Yn) de datos normales bidimensionales y no es posible suponer que las variables X e Y son independientes. En estos casos se trabaja con las diferencias D1,,Dn, donde Di=XiYi. Siempre se tiene que μ=E(Di)=E(Xi)E(Yi)=μ1μ2. Por lo tanto, el intervalo de confianza para la diferencia de medias se puede construir a partir del intervalo para μ.

Ejemplo

Tomamos medidas de la concentración de zinc en la superficie y en el fondo (en mg/l) de seis puntos de un río. ¿Es la concentración media igual en la superficie y en el fondo? Como cada medida en profundidad y fondo corresponde al mismo punto del río no podemos suponer independencia. Esto se indica usando el argumento paired = TRUEen el comando t.test:

fondo <- c(0.41, 0.24, 0.39, 0.41, 0.60, 0.61)
superficie <- c(0.43, 0.27, 0.57, 0.53, 0.71, 0.72)
t.test(fondo, superficie, paired = TRUE)$conf.int
#> [1] -0.15822795 -0.03177205
#> attr(,"conf.level")
#> [1] 0.95

Varianza

Si X1,,Xn son v.a.i.i.d. con distribución N(μ,σ2), entonces uno de los apartados del lema de Fisher garantiza: (n1)S2σ2χn12. Por lo tanto (n1)S2/σ2 es una cantidad pivotal para σ2. Si χn1,α2 es el cuantil 1α de la la distribución χn12, entonces 1α=P(χn1,1α/22<(n1)S2σ2<χn1,α/22). De aquí se deduce fácilmente un intervalo de confianza para σ2 (ejercicio).

La diferencia respecto a los intervalos para medias es que ahora la distribución no es simétrica por lo que los cuantiles necesarios no son iguales salvo el signo. De hecho hay infinitos pares de cuantiles que podríamos usar, todos los pares que encierran entre sí una probabilidad 1α). Los que se usan habitualmente (los de la expresión anterior) de hecho no proporcionan el intervalo de confianza más corto.

Cociente de varianzas

Volvemos a la situación de dos muestras independientes de v.a.i.i.d.. La primera, X1,,Xn1, procede de una distribución N(μ1,σ12) mientras que la segunda, Y1,,Yn2, procede de una distribución N(μ2,σ22). En este caso, queremos encontrar un intervalo de confianza para el cociente entre las varianzas σ12/σ22. La cantidad pivotal en este caso depende de una nueva distribución:

Definición 4.2 Sean Y1 e Y2 dos v.a. independientes con distribuciones χn12 y χn22, respectivamente. Se dice que la variable Y=Y1/n1Y2/n2 tiene distribución F con n1 y n2 grados de libertad. (Notación: YFn1,n2).

La forma de la función de densidad de una v.a. con distribución F es parecida a la de una χ2.

datos <- data.frame(x = c(0, 6))
ggplot(datos, aes(x)) +
  geom_function(fun = df, args = list(df1 = 5, df2= 20), linewidth = 1.05)  +
  geom_function(fun = df, args = list(df1 = 20, df2= 20), col = 'blue', linetype = 2, linewidth = 1.05) +
  geom_function(fun = df, args = list(df1 = 20, df2= 5), col = 'red', linetype = 3, linewidth = 1.05)

Cuestiones
  • ¿Qué relación hay entre la distribución F y la t de Student?
  • Demuestra que Fn1,n2;α=1/Fn2,n1;1α.

Sean S12 y S22 las varianzas muestrales de cada una de las dos muestras. Usando el lema de Fisher y la independencia de las muestras tenemos que F=S12/σ12S22/σ22=S12/S22σ12/σ22Fn11,n21. Por lo tanto, F es una cantidad pivotal para el cociente de las varianzas de la que se puede obtener el intervalo de confianza correspondiente: IC1α(σ12/σ22)=[S12/S22Fn11,n21;α/2,  S12/S22Fn11,n21;1α/2]=[S12S22Fn21,n11;1α/2,  S12S22Fn21,n11;α/2].

Para calcular este intervalo con R se usa el comando var.test:

# Generación de  datos
set.seed(123)
n <- 100
sigma1 <- sqrt(2)
sigma2 <- 1
x <- rnorm(n, sd = sigma1)   
y <- rnorm(n, sd = sigma2)

# Calculo de intervalos
var.test(x, y)$conf.int  # Por defecto 95% es el nivel de confianza
#> [1] 1.199136 2.648760
#> attr(,"conf.level")
#> [1] 0.95
var.test(x, y, conf.level = 0.99)$conf.int   # Nivel 99%
#> [1] 1.057455 3.003648
#> attr(,"conf.level")
#> [1] 0.99

4.4 Cantidades pivotales asintóticas basadas en el estimador de máxima verosimilitud

Si X1,,Xn son v.a.i.i.d. de una distribución F{Fθ:θΘR}, sabemos por el tema anterior que, bajo las condiciones de regularidad, el estimador de máxima verosimilitud θ^ verifica n(θ^θ)dN(0,1I(θ)), donde I(θ) es la cantidad de información de Fisher de cada observación. Multiplicando por I(θ) tenemos que nI(θ)(θ^θ)dN(0,1), y entonces nI(θ)(θ^θ) es una cantidad pivotal asintótica. Si es posible despejar el parámetro en la expresión 1αP(zα/2nI(θ)(θ^θ)zα/2), tendremos un intervalo de confianza para θ de nivel aproximado 1α.

A veces puede ser difícil o imposible despejar el parámetro. En ese caso aún podemos construir un intervalo de confianza asintótico si encontramos un estimador consistente de la información de Fisher, I(θ)^. En este caso, también se cumple (lema de Slutsky): nI(θ)^(θ^θ)dN(0,1), de donde se obtiene el siguiente intervalo de confianza de nivel aproximado 1α: IC1α(θ)[θ^(nI(θ)^)1/2zα/2, θ^+(nI(θ)^)1/2zα/2].

Ejemplo

Sean X1,,Xn v.a.i.i.d. con distribución exponencial de parámetro θ (f(x;θ)=θeθx, si x>0). Es muy fácil ver que el EMV del parámetro es θ^=1/X¯ y que la información de Fisher de cada observación vale I(θ)=1/θ2. Como consecuencia, n(θ^θ)dN(0,θ2)n(θ^θ))θdN(0,1). Despejando el parámetro en 1αP(zα/2n(θ^θ)/θzα/2) se deduce el intervalo de confianza aproximado IC1α(θ)[nθ^n+zα/2, nθ^nzα/2].

Hay al menos otros dos métodos para obtener intervalos de confianza en este ejemplo (se propone como ejercicio deducir los intervalos correspondientes):

  • Estimar la información de Fisher consistentemente. Por la LDGN y el teorema de la aplicación continua, θ^pθ. Por lo tanto, n(θ^θ)θ^dN(0,1).
  • Determinar una función estabilizadora de la varianza g tal que n(g(θ^)g(θ))dN(0,1). Aplicando el método delta, esta función tiene que verificar g(θ)2θ2=1, por lo que es muy fácil encontrarla.

Puede demostrarse (ejercicio) que los tres intervalos coinciden si despreciamos todos los términos O(1/n) o menores.

4.5 Intervalos desde el punto de vista bayesiano

En un problema de inferencia con un enfoque bayesiano el elemento fundamental para realizar la inferencia es la distribución a posteriori π(θ|x1,,xn). A partir de esta distribución se define una región creíble de nivel 1ϵ, con ϵ(0,1) como un subconjunto AΘ tal que P(θA|x1,,xn)=Aπ(θ|x1,,xn)dθ=1ϵ. Bajo el punto de vista bayesiano el parámetro es una v.a. y, por tanto, para una muestra fija puede hablarse propiamente de la probabilidad de que el parámetro esté dentro del intervalo. Por el contrario, en el enfoque frecuentista, si ya hemos obtenido la muestra y tenemos un intervalo concreto I=I(x1,,xn), no se puede decir estrictamente que la probabilidad de que el parámetro esté en I es 1ϵ porque en I ya no hay nada aleatorio y el valor verdadero θ (desconocido) del parámetro cumplirá θI o θI pero no tiene sentido asignar probabilidad a que el parámetro esté o no en el intervalo. La interpretación correcta es la frecuentista que hemos comentado anteriormente.

Ejemplo

Se desea obtener un intervalo creíble para el parámetro λ de una distribución de Poisson a partir de una muestra x1,,xn, suponiendo que λγ(α,β), siendo αN conocido. Estamos usando la siguiente parametrización de la distribución gamma: f(x;α,β)=βαΓ(α)xα1eβx,   x>0. La función característica de esta v.a. es φY(t)=(1it/β)α. Si c>0, φcY(t)=(1ict/β)α, que corresponde a una γ(α,β/c).

Es fácil comprobar que la distribución a posteriori de λ es γ(α=nx¯+α,β=n+β). Por tanto, usando la propiedad anterior y la relación de la distribución gamma con la χ2, la distribución a posteriori de 2(n+β)λ es γ(nx¯+α,1/2)χ2(nx¯+α)2. Así pues, P{χ2(xi+p);1ϵ/222(n+a)λχ2(xi+p);ϵ/22}=1ϵ, y un intervalo creíble de nivel 1ϵ es A=(χ2(nx¯+α);1ϵ/222(n+β),χ2(nx¯+α);ϵ/222(n+β)).

4.6 Referencias

Todo el material de este capítulo es muy estándar por lo que hay muchos libros de texto que se pueden leer para ampliar información o trabajar más ejemplos. Como referencias apropiadas citamos los capítulos 24 y 25 de Dekking et al (2005), que incluyen una sección sobre el uso del bootstrap para obtener intervalos, y el capítulo 7 de Thijssen (2016), que resume sintéticamente las principales ideas. Una visión más teórica se puede encontrar en Casella y Berger (2001) y en Knight (1999).

  • Casella, G., y Berger, R. L. (2001). Statistical inference, second edition. Cengage Learning.
  • Dekking, F. M., Kraaikamp, C., Lopuhaä, H. P. y Meester, L. E. (2005). A Modern Introduction to Probability and Statistics: Understanding why and how. Springer.
  • Knight, K. (1999). Mathematical statistics. CRC Press.
  • Thijssen, J. (2016). A Concise Introduction to Statistical Inference. Chapman and Hall/CRC.