Prueba de bondad de ajuste

Las pruebas de bondad de ajuste son un contraste de hipótesis para determinar el grado o nivel de ajuste de nuestros datos a una distribución teórica.

Estas pruebas se basan en la comparación de las frecuencias de ocurrencia observadas en una muestra empírica y las frecuencias esperadas de una distribución teórica. En donde, el objetivo será si existe o no discrepancia entre los valores observados y los valores esperados de la distribución de interés.

La hipótesis a probar de interés estará dada por \begin{align*} H_0: X \sim F_0(x) \quad \text{vs} \quad H_1: X \nsim F_0(x) \end{align*}

Siendo $F_0(x)$ la distribución de probabilidad hipótetica que se quiere probar.

Entre las pruebas de bondad de ajuste más usadas se tiene:

Prueba Librería Función
Kolmogorov–Smirnov stats ks.test()
Cramer–von Mises goftest cvm.test()
Anderson–Darling goftest ad.test()
QQPlot car qqPlot()

Donde éstas pruebas, requieren de los parámetros de la distribución que se quiere probar, y para encontrarlos, es posible emplear métodos de optimización que nos permitan observar cuales son los parámetros ajustados para un conjunto de datos determinado. Para emplear estos método de optimización es posible usar la función fitdistr() de la librería MASS.

Ejercicio

Suponga que el grupo de Economía de la Salud ha realizado un estudio sobre los contagios que se han presentado durante las últimas semanas por el COVID-19, ya que tienen la hipótesis de la tasa de contagios posee una distribución Weibull.

Para probar si dicha hipótesis se encuentra apoyada o no por la información empírica, el grupo de Economía de la Salud, ha decidido tomar como una muestra de los últimos 2 meses (61 días) sobre el número de contagios que se han presentado la región antioqueña, encontrado los siguientes registros (en miles)

714 754 679 735 684 759 759 741 713 722
739 710 737 708 711 744 733 739 681 741
723 689 699 680 690 678 694 710 690 694
722 725 747 693 701 724 713 720 696 682
713 657 714 717 695 724 735 700 713 739
737 721 744 668 728 727 695 725 705 678
713

Si se decide emplear un nivel de significancia del \(5\%\), pruebe si las hipótesis planteadas por el grupo de la economía de la salud se encuentran o no apoyadas por la información empírica.

Solución

En este caso estamos interesados en emplear un nivel de significancia \(\alpha=0.05\) para probar si la información obtenida en \(61\) días apoya o no la hipótesis sobre que los contagios que se han presentado durante las últimas semanas por el COVID-19 poseen una distribución Weibull, esto es \[\begin{align*} H_0:X \sim Wei(\alpha,\beta)\\ H_1:X \not\sim Wei(\alpha,\beta) \end{align*}\] y para ello debemos inicialmente ajustar los parámetros de una distribución Weibull que mejor ajustan los datos, mediante la función fitdistr de la librería MASS, tal que 

## Se cargan los datos
datos <- c(714, 754, 679, 735, 684, 759, 759, 741, 713, 722,
           739, 710, 737, 708, 711, 744, 733, 739, 681, 741,
           723, 689, 699, 680, 690, 678, 694, 710, 690, 694,
           722, 725, 747, 693, 701, 724, 713, 720, 696, 682,
           713, 657, 714, 717, 695, 724, 735, 700, 713, 739,
           737, 721, 744, 668, 728, 727, 695, 725, 705, 678,
           713)


## Se realiza el ajuste de los parámetros a una distribución Weibull
library(MASS)
parametros <- fitdistr(x = datos, densfun = "weibull")
parametros

     shape        scale   
   33.781908   724.649562 
 (  3.286809) (  2.906706)

Una vez conocidos el valor de los parámetros de la distribución Weibull que ofrecen un mejor ajuste a los datos muestrales, tal que \(\alpha = 33.781908\) y \(\beta = 724.649562\), se procede al ajuste de la distribución mediante el empleo de la prueba Cramer-Von Mises, mediante la función cvm.test de la librería goftest 

## Se hace el ajuste de los datos a la distribución Weibull
library(goftest)
cvm.test(x = datos, null = "pweibull", shape = parametros$estimate[1],
    scale = parametros$estimate[2])

 Cramer-von Mises test of goodness-of-fit
    Null hypothesis: Weibull distribution
    with parameters shape = 33.781908, scale =
    724.649562
    Parameters assumed to be fixed

data:  datos
omega2 = 0.067026, p-value = 0.7718

Donde se observa que la prueba de bondad de ajuste arrojó un P-valor\(=0.7718\), por lo cual al ser mayor que al nivel de significancia \(\alpha=0.05\), se concluirá que no hay evidencia en contra de rechazar la hipótesis nula de que los datos poseen una distribución Weibull con parámetro de forma \(\alpha = 33.781908\) y parámetro de escala \(\beta = 724.649562\).

Adicionalmente, podemos ver el ajuste que tiene la distribución Weibull de parámetros \(\alpha = 33.781908\) y \(\beta = 724.649562\), mediante un QQ-plot, tal que 

## Se muestra el ajuste de la distribución Weibull de forma gráfica
library(car)
qqPlot(datos, dist = "weibull", shape = parametros$estimate[1], scale = parametros$estimate[2],
    envelope = 0.95)

en donde se observa que todos los puntos asociados a la distirbución empírica caen dentro de los límites de confianza del \(95\%\) del QQ-plot, y por tanto se tendrá con un nivel de significancia \(\alpha=0.05\), que el número de contagios que se han presentado durante las últimas semanas por el COVID-19, posee una distribución Weibull con parámetro de forma \(\alpha = 33.781908\) y parámetro de escala \(\beta = 724.649562\).

Pruebas de normalidad

Un caso particular de las pruebas de bondad de ajuste son las pruebas específicas que permiten probar si un conjunto de datos se distribuyen o no normalmente, en donde el juego de hipótesis de interés estará dado por

\begin{align*} H_0: X \sim N(\mu,\sigma^2) \quad \text{vs} \quad H_1: X \nsim N(\mu,\sigma^2) \end{align*}

A saber, las funciones más usadas de bondad de ajuste para probar normalidad son:

Prueba Librería Función
Shapiro-Wilk stats shapiro.test()
Lilliefors nortest lillie.test()
Shapiro-Francia nortest sf.test()
Cramer Von-Mises nortest cvm.test()
Anderson-Darling nortest ad.test()
QQPlot car qqPlot()

en donde, a diferencia de las pruebas de bondad de ajuste anteriormente presentadas, éstas pruebas no requieren ajustar con anterioridad los parámetros de la distribución normal, ya que dicho procedimiento se realiza de forma interna dentro de las funciones de ajuste.

Ejercicio

Suponga que el grupo de Economía de la Salud ha realizado un estudio sobre los contagios que se han presentado durante las últimas semanas por el COVID-19, ya que tienen la hipótesis de la tasa de contagios posee una distribución Normal.

Para probar si dicha hipótesis se encuentra apoyada o no por la información empírica, el grupo de Economía de la Salud, ha decidido tomar como una muestra de los últimos 2 meses (61 días) sobre el número de contagios que se han presentado la región antioqueña, encontrado los siguientes registros (en miles)

714 754 679 735 684 759 759 741 713 722
739 710 737 708 711 744 733 739 681 741
723 689 699 680 690 678 694 710 690 694
722 725 747 693 701 724 713 720 696 682
713 657 714 717 695 724 735 700 713 739
737 721 744 668 728 727 695 725 705 678
713

Si se decide emplear un nivel de significancia del \(5\%\), pruebe si las hipótesis planteadas por el grupo de la economía de la salud se encuentran o no apoyadas por la información empírica. 

## Se cargan los datos
datos <- c(714, 754, 679, 735, 684, 759, 759, 741, 713, 722, 739, 710, 737, 708,
    711, 744, 733, 739, 681, 741, 723, 689, 699, 680, 690, 678, 694, 710, 690, 694,
    722, 725, 747, 693, 701, 724, 713, 720, 696, 682, 713, 657, 714, 717, 695, 724,
    735, 700, 713, 739, 737, 721, 744, 668, 728, 727, 695, 725, 705, 678, 713)

## Se hace el ajuste de los datos a la distribución Normal
library(nortest)
lillie.test(datos)

    Lilliefors (Kolmogorov-Smirnov) normality test

data:  datos
D = 0.066013, p-value = 0.7321

Observando que al obtener un P-valor\(=0.7321\) se tendrá que no hay evidencia en contra de la hipótesis nula, ya que el P-valor obtenido en la prueba de normalidad es superior al nivel de significancia \(\alpha=0.05\), concluyendo entonces que el número de contagios que se han presentado durante las últimas semanas por el COVID-19, posee una distribución Normal.

Además, si realizamos el QQ-plot para observar el ajuste de las observaciones a la distribución normal se tiene que 

## Se muestra el ajuste de la distribución Normal de forma gráfica
library(car)
qqPlot(datos, envelope = 0.95)

En donde se aprecia que la totalidad de los puntos asociados a los datos muestrales caen dentro de las bandas de confianza del \(95\%\) de la distribución normal, por lo cual se concluye con un nivel de significancia del \(5\%\) que el número de contagios que se han presentado durante las últimas semanas por el COVID-19, si posee una distribución Normal.

Tablas de contingencia

Las tablas de contingencia se emplean para probar si uno o más atributos de una población son o no independientes. Para ello se plantea el siguiente juego de hipótesis

\begin{align*} H_0&: \text{Los atributos de la población son independientes}\\ H_1&: \text{Los atributos de la población son dependientes} \end{align*}

Suponga que en una población se clasifica según los atributos $A$ y $B$, y además suponga que existen $r$ categorías en el atributo $A$ y $c$ categorías en el atributo $B$, de forma que hay un total de $r\times c$ posibles combinaciones entre los atributos, tal que

Donde $O_{ij}$: es el número de observaciones entre la i-ésima categoría de $A$ y la j-ésima categoría de $B$, $R_{i}$: es la suma de las observaciones por fila, y $C_{j}$: es la suma de las observaciones por columnas, tal que \begin{align*} R_{i}=\sum_{j=1}^{c}O_{ij} \qquad \qquad C_{j}=\sum_{i=1}^{r}O_{ij} \end{align*}

Ahora, con el fin de probar si hay o no independencia entre las categorías, es necesario calcular $E_{ij}$ el número de observaciones esperadas en la casilla $ij$, tal que \begin{align*} E_{ij}=np_{ij} \end{align*}

donde $pij$ es la probabilidad conjunta bajo el supuesto de independiencia, entonces \begin{align*} p_{ij} = \frac{R_{i}}{n}\times \frac{C_{j}}{n} =\frac{R_{i}C_{j}}{n^2} \end{align*}

y por tanto \begin{align*} E_{ij} = \frac{R_{i}C_{j}}{n} \end{align*}

Una vez calculados cada uno de los $E_{ij}$, es posible reescribir dichos valores de forma en forma de matríz de la forma

Finalmente, al comparar las frecuencias observadas $O_{ij}$ con las frecuencias esperadas $E_{ij}$, se tendrá que el estadístico de prueba para probar la independencia, será de la forma \begin{align*} \chi^2_{c}=\sum_{i=1}^{r}\sum_{j=1}^{c} \frac{(O_{ij}- E_{ij})^2}{E_{ij}} \sim \chi^2_{(r-1)(c-1)} \end{align*}

la región crítica \begin{align*} RC:\{\chi^2|\chi^2>\chi^2_{\alpha,(r-1)(c-1)}\} \end{align*}

y el P-valor \begin{align*} \text{P-valor}=\mathbb{P}(\chi^2_{(r-1)(c-1)}>\chi^2_c) \end{align*}

Ejercicio

Suponga que se desea clasificar los defectos de fabricación en una empresa manufacturera de acuerdo al tipo de defecto y la jornada en la cual fueron producidas.

Suponga que se obtienen los siguientes resultados de una muestra aleatoria de \(350\) productos defectuosos

Contrastar a un nivel de significancia del \(5%\), si existe independencia entre el tipo de defecto y el turno de producción en el cual fue manufacturado el producto.

Solución

Con el fin de probar si existe o no independencia entre el tipo de defecto y el turno de producción en el cual fue manufacturado el producto, es necesario calcular inicialmente las sumas marginales por filas y columnas de la varible observada, tal que

Una vez calculadas las sumas marginales por fila y por columna, se procede a calcular las frecuencias esperadas, mediante la ecuación \[\begin{align*} E_{ij} = \frac{R_{i}C_{j}}{n} \end{align*}\] para obtener la siguiente tabla de frecuencias esperadas, de la forma

Y ahora, al emplear los valores observados con los esperados, se realiza el cálculo del estadístico de prueba, tal que \[\begin{align*} \chi^2_{c}&=\sum_{i=1}^{r}\sum_{j=1}^{c} \frac{(O_{ij}- E_{ij})^2}{E_{ij}} \sim \chi^2_{(r-1)(c-1)}\\ &=\frac{(16 - 26.71429)^2}{26.71429} + \frac{(27 - 24.51429)^2}{24.51429} + \ldots + \frac{(23 - 17.45714)^2}{17.45714}\\ &=22.98415 \end{align*}\] Finalmente se plantea el criterio de decisión, en donde la región crítica asociada para un nivel de significancia del \(5\%\) está dada por \[\begin{align*} RC&:\{\chi^2|\chi^2>\chi^2_{\alpha,(r-1)(c-1)}\}\\ RC&:\{\chi^2|\chi^2>\chi^2_{0.05,(3-1)(4-1)}\}\\ RC&:\{\chi^2|\chi^2>\chi^2_{0.05,6}\}\\ RC&:\{\chi^2|\chi^2>12.5915872\}\\ \end{align*}\] mientras que el P-valor asociado al estadístico de prueba es igual a \[\begin{align*} \text{P-valor}&=\mathbb{P}(\chi^2_{(r-1)(c-1)} > \chi^2_c)\\ &=\mathbb{P}(\chi^2_{(3-1)(4-1)} > 22.98415)\\ &=\mathbb{P}(\chi^2_{6} > 22.98415)\\ &=0.0008018046 \end{align*}\] Encontrando con un nivel de significancia del \(5\%\) que se rechaza la hipótesis nula de que el tipo de defecto y el turno sean categorías independiente, lo cual quiere decir que el turno laboral si influye en el tipo de defecto que tienen los productos producidos por la empresa maufacturera, debido a que el estadístico de prueba cae dentro de la región crítica, y el P-valor es menor al nivel de significancia \(\alpha = 0.05\).