Contrastes de normalidad

Traballo Fin de Grao Contrastes de normalidad Víctor Platas García 2020/2021 UNIVERSIDADE DE SANTIAGO DE COMPOSTELA GRAO DE MATEMÁTICAS Traballo Fin de Grao Contrastes de normalidad Víctor Platas García Julio 2021 UNIVERSIDADE DE SANTIAGO DE COMPOSTELA Trabajo propuesto Área de Coñecemento: Estadística e Investigación Operativa Título: Contrastes de normalidad Breve descrición do contido Muchos procedimientos estadísticos parten del supuesto de que los datos proceden de una distribución normal. Sin embargo, es muy habitual que los datos reales no cumplan esta suposición, lo cual repercutirá en la validez de tales procedimientos. Por este motivo son importantes los contrastes de normalidad. En este trabajo se revisarán los contrastes más conocidos, como el test de Kolmogorov-Smirnov- Lilliefors o el test de Shapiro-Wilk, y se completarán con otros tests, como los basados en la asimetría o la curtosis. Se estudiarán sus propiedades teóricas y se realizarán estudios de simulación para comprobar su funcionamiento en datos simulados. Recomendacións Outras observacións iii Índice general Resumen viii Introducción xi 1. La distribución normal 1 1.1. Contextohistórico ................................ 1 1.2. Denición..................................... 2 1.3. Teorema Central del Límite . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 1.4. Estadísticainferencial .............................. 6 1.4.1. Estimación de parámetros de máxima verosimilitud . . . . . . . . . . 7 2. El test Chi Cuadrado de Pearson 9 2.1. Contraste de una hipótesis nula simple . . . . . . . . . . . . . . . . . . . . . 9 2.2. Contraste de un modelo paramétrico . . . . . . . . . . . . . . . . . . . . . . 11 3. Test basados en la función de distribución 15 3.1. Test de Kolmogorov-Smirnov y versión de Lilliefors . . . . . . . . . . . . . . 15 3.1.1. Prueba de Kolmogorov-Smirnov . . . . . . . . . . . . . . . . . . . . . 15 3.1.2. TestdeLilliefors ............................. 17 3.2. TestdeCramer.................................. 17 3.3. Test de Anderson-Darling . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18 3.4. Relación Kolmogorov-Cramer-Anderson . . . . . . . . . . . . . . . . . . . . 19 4. Test basados en los momentos 21 4.1. Asimetría ..................................... 21 4.1.1. Denición................................. 21 4.1.2. Medidas de asimetría . . . . . . . . . . . . . . . . . . . . . . . . . . . 22 4.1.3. Asimetría en la distribución normal . . . . . . . . . . . . . . . . . . . 23 4.1.4. Test de normalidad basado en la asimetría . . . . . . . . . . . . . . . 24 v vi ÍNDICE GENERAL 4.2. Curtosis...................................... 25 4.2.1. Denición................................. 25 4.2.2. Medidas de curtosis . . . . . . . . . . . . . . . . . . . . . . . . . . . 25 4.2.3. Curtosis en la distribución normal . . . . . . . . . . . . . . . . . . . 26 4.2.4. Test de normalidad basado en la curtosis . . . . . . . . . . . . . . . . 27 4.3. TestdeJarque-Bera ............................... 28 5. Métodos grácos 31 5.1. Histograma y curva normal . . . . . . . . . . . . . . . . . . . . . . . . . . . 31 5.2. Gráco de cuantiles teóricos . . . . . . . . . . . . . . . . . . . . . . . . . . . 33 6. Test Shapiro-Wilk y Shapiro-Francia 37 6.1. TestdeShapiro-Wilk............................... 37 6.2. Test de Shapiro-Francia . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39 7. Estudio comparativo mediante simulaciones 41 7.1. Resultados de las simulaciones . . . . . . . . . . . . . . . . . . . . . . . . . . 41 7.2. Conclusiones ................................... 50 A. Códigos de R 51 Bibliografía 61 2 CAPÍTULO 1. LA DISTRIBUCIÓN NORMAL a mediados del siglo XVIII en distintos estudios hasta la actualidad, cuyo estudio y formalización ha hecho de esta distribución una de las bases de la estadística que conocemos hoy en día. 1.2. Denición La distribución normal, también conocida comoe distribución de Gauss o distribución gaussiana es una distribución de probabilidad de variable continua que suele aparecer frecuentemente tanto en teoría de probabilidades como en la propia estadística. Para denir formalmente la distribución normal veamos en primer lugar los parámetros que la componen: µ : es la media (aunque también puede ser la moda, la mediana o el valor esperado). σ : es la desviación típica. σ2 : es la varianza. ϕ : es la función de densidad de probabilidad. De los parámetros anteriores, destacamos la función de densidad de probabilidad, más comúnmente conocida como campana de Gauss que viene dada por: ϕ(x) = 1 σ√2πe−(x−µ)2/2σ2 Y cuya gráca, también conocida globalmente, es la siguiente: Figura 1.1: Función de densidad de probabilidad 1.2. DEFINICIÓN 3 Una vez denida la función de densidad de probabilidad, podemos denir la función de distribución de la distribución normal: [2] Φµ,σ2(x) = x Z −∞ ϕµ,σ2(u)du =1 σ√2π x Z −∞ e−(u−µ)2 2σ2du (1.1) La cual, del mismo modo que vimos para la función de densidad, presenta la siguiente gráca en función de los distintos valores de los parámetros µ y σ2 : Figura 1.2: Función de distribución Veamos a continuación las principales propiedades de la distribución normal. Propiedades 1. Es simétrica respecto de la media µ . 2. El valor de la media µ coincide con el valor de la moda y mediana. 3. Los puntos de inexión de la función de densidad se dan en x=µ−σ y en x=µ+σ . 4. Como podemos observar en la gura 1.1 en el intervalo [µ−2σ, µ + 2σ] encontramos cerca de el 95,5% de la distribución, de donde se obtienen los limites de las tablas utilizadas habitualmente en la normal estándar. 5. Dado X∼N(µ, σ2) y a , b∈R , entonces aX +b∼N(aµ +b, a2σ2) 6. Suponiendo X∼N(µX, σ2 X) y Y∼N(µY, σ2 Y) son dos variables aleatorias independientes, entonces: 4 CAPÍTULO 1. LA DISTRIBUCIÓN NORMAL La suma de ambas también seguirá una distribución normal, es decir: U=X+Y∼N(µX+µY, σ2 X+σ2 Y) Además, se tiene que dadas dos variables aleatorias independientes cuya suma sigue una distribución normal, podemos asegurar que ambas variables serán normales. Al igual que la suma, la diferencia de ambas variables también sigue una distribución normal, es decir: V=X−Y∼N(µX−µY, σ2 X+σ2 Y) Si las varianzas de las variables X e Y son iguales, entonces U y V son independientes entre sí. 7. Dadas n variables normales estándar independientes X1, ..., Xn , la suma de dichas variables X2 1+···+X2 n seguirá una distribución χ2 con n grados de libertad. 8. Dadas n variables normales estándar independientes X1, ..., Xn , la media muestral ¯ X= (X1+··· +Xn)/n y la varianza muestral X2= ((X1−¯ X)2+··· + (Xn− ¯ X)2)/(n−1) son independientes. Adicionalmente, debemos destacar el caso especial cuando µ= 0 y σ= 1 , ya que, en ese caso, la función de distribución recibe el nombre de función de distribución normal estándar, calculada de la siguiente forma sustituyendo dichos valores en la ecuación 1.1: Φ0,1(x) = x Z −∞ ϕ0,1(u)du =1 √2π x Z −∞ e−u2 2du Finalmente, una de las deniciones que también admite la función de distribución viene dada en términos de la función error: Φ(x) = 1 21 + erf x−µ σ√2 El complemento de la función de distribución Q(x) denida como 1−Φ(x) sirve para representar la cola de la probabilidad de la distribución gaussiana y puede venir dada en función del error Φ−1(p) = √2erf−1(2p−1), p ∈(0,1) o, en cambio, expresarse como: Φ−1 µ,σ2(p) = µ+σΦ−1(p) = µ+σ√2erf−1(2p−1), p ∈(0,1) 1.3. TEOREMA CENTRAL DEL LÍMITE 5 Estandarización Como consecuencia de las propiedades de la distribución normal previamente mencionadas, tenemos que es posible relacionar cualquier variable aleatoria normal con la distribución normal estándar de la siguiente forma. Dada una variable aleatoria X∼N(µ, σ2) , entonces Z=X−µ σ es una variable aleatoria normal estándar, es decir, Z∼N(0,1) . Esta transformación, que recibe también los nombres de normalización ó tipicación, sirve de gran ayuda dado que los valores de la función de distribución de la normal estándar están tabulados y de esta forma podemos usar dichos resultados para hallar valores de la función de distribución de variables no estandarizadas a través transformaciones simples. 1.3. Teorema Central del Límite Consideremos en primer lugar n variables aleatorias Xi, i = 1, . . . , n idénticamente distribuídas e independientes de forma que E[Xi] = µ y V ar[Xi] = σ2<∞∀i= 1, . . . , n . Podemos denir entonces la variable aleatoria Sn como la suma de esas n variables aleatorias, es decir: Sn:= X1+···+Xn= n X i=1 Xi Estandarizando Sn obtenemos una nueva variable denida como: Zn:= Sn−nµ σ√n=Pn i=1 Xi−nµ σ√n (1.2) cuya media es igual a 0 y cuya desviación estandar será igual a 1. Estamos pues en disposición de formalizar el enunciado de dicho teorema. Teorema 1.1 (Central del Límite) . Sean n variables aleatorias Xi, i = 1, . . . , n independientes e idénticamente distribuidas con E[Xi] = µ y V ar[Xi] = σ2<∞∀i= 1, . . . , n y sea Zn la distribución denida en 1.2. Entonces, la función de distribución de Zn converge hacia la función de distribución normal estándar cuando n→ ∞ , es decir, l´ım n→∞ P(Zn≤z) = Φ(z) = Zz −∞ 1 √2πe−x2 2dx 6 CAPÍTULO 1. LA DISTRIBUCIÓN NORMAL Otra versión que admite el enunciado de este teorema es el que adopta la versión estandarizada de Zn en función de la media muestral ¯ X , es decir Zn=¯ X−µ σ √n Propiedades del Teorema Central del Límite Cuando n es sucientemente grande, el Teorema Central del Límite nos garantiza una distribución aproximadamente normal. Para asegurar la convergencia, es suciente con que las variables que se suman sean independientes, idénticamente distribuidas, con valor esperado y varianza nitas. La aproximación entre las dos distribuciones es, en general, mayor en el centro de las mismas que en sus extremos o colas. 1.4. Estadística inferencial En numerosas ocasiones, precisamos de parámetros de la normal para realizar métodos inferencia estadística que nos permitan estimar intervalos de conanza para parámetros como la media o contrastes de hipótesis para valores propuestos para la misma. Dada una muestra de tamaño n , con media muestral ¯ X y desviación estándar muestral s , podemos calcular los límites del intervalo de conanza de la media poblacional para el nivel de conanza α elegido, y el p-valor del contraste de hipótesis sobre la hipótesis nula H0:µ=µ0 , siendo µ0 el valor propuesto para la media poblacional. El estadístico de contraste sobre la media seguirá una distribución t-Student con n−1 grados de libertad por lo que, en virtud del Teorema Central del Límite, la distribución de este estadístico tiende a la normal (es decir, será aproximadamente normal cuando el tamaño de la muestra sea sucientemente grande). De la misma forma, se puede realizar un proceso similar para estimar proporciones de un suceso que aconteció en la muestra. A la hora de comparar parámetros de dos muestras, también encontraremos que para realizar pruebas tipo la prueba t, precisaremos suponer que ambas muestras proceden de poblaciones normalmente distribuídas o, en el caso de la comparación de proporciones de un suceso en dos muestras, también aprecerá la normal en la distribución que sigue el estadístico empleado en dicha prueba. 1.4. ESTADÍSTICA INFERENCIAL 7 1.4.1. Estimación de parámetros de máxima verosimilitud La estimación por máxima verosimilitud, conocida también como EMV y, en ocasiones, MLE por sus siglas en inglés, es un método habitual para ajustar un modelo y estimar sus parámetros. Así pues, supongamos X1, . . . , Xn variables normalmente distribuidas e independientes con media µ y varianza σ2>0 . Tenemos por tanto que los valores observados de estas n variables aleatorias constituyen una muestra de tamaño n de una población normalmente distribuida. Intentaremos estimar la media poblacional µ y la desviación típica poblacional σ , basándose en las valores observados de esta muestra. La función de densidad conjunta de estas n variables aleatorias independientes es f(x1, ..., xn;µ, σ) = n Y i=1 ϕµ,σ2(xi) = 1 (σ√2π)n n Y i=1 exp −1 2xi−µ σn,(x1, ..., xn)∈Rn. La función de verosimilitud basada en las observaciones X1, . . . , Xn es L(µ, σ) = C σnexp −Pn i=1(Xi−µ)2 2σ2, µ ∈R, σ > 0, con alguna constante C > 0 . En el método de máxima verosimilitud, los valores de µ y σ que maximizan la función de verosimilitud se toman como estimadores de los parámetros poblacionales µ y σ . Normalmente, a la hora de maximizar una función de dos variables consideraríamos las derivadas parciales, pero en este caso, el valor de µ que maximiza la función de verosimilitud con σ jo no depende de σ . Sin embargo, tenemos que ese valor de µ se sustituye por µ en la función de verosimilitud con lo que nalmente encontramos el valor de σ que maximiza la expresión resultante. Veámoslo: En primer lugar observamos que la función de verosimilitud es una función decreciente de Pn i=1(Xi−µ)2 por lo que queremos obtener el valor de µ que minimice esta suma. Sea por tanto ¯ Xn= (X1+···+Xn)/n la media muestral basada en las n observaciones. Observamos que n X i=1 (Xi−µ)2= n X i=1 (Xi−¯ Xn)2+n(¯ Xn−µ)2, donde tan sólo el último elemento depende de µ y se minimiza por ˆµn=¯ Xn . Esta es, por lo tanto, la estimación de máxima verosimilitud de µ basada en las n observaciones X1, ..., Xn . Cuando sustituimos esta estimación por µ en la función de verosimilitud, obtenemos 8 CAPÍTULO 1. LA DISTRIBUCIÓN NORMAL L(¯ Xn, σ) = C σnexp −Pn i=1(Xi−¯ Xn)2 2σ2, σ > 0 Para obtener el estimador de máxima verosimilitud tanto de σ2 como de σ debemos calcular la log-función de verosimilitud, esto es, el logaritmo de la función de verosimilitud, que denotaremos como ` : `(¯ Xn, σ) = log C −n log σ −Pn i=1(Xi−¯ Xn)2 2σ2, σ > 0, de donde se obtiene ∂ ∂µ `(¯ Xn, σ) = −n σ3 σ2−1 n n X i=1 (Xi−¯ Xn)2!, σ > 0. Observamos que esta derivada tendrá valor positivo, cero o negativo según σ2 esté entre 0 y ˆσ2 n:= 1 nPn i=1(Xi−¯ Xn)2 , sea igual a dicha cantidad o sea mayor que esa cantidad respectivamente. Esta media de cuadrados de residuos es pues el estimador de máxima verosimilitud de σ2 , y su raíz cuadrada es el estimador de máxima verosimilitud de σ basado en las n observaciones. Observemos nalmente que si hay solamente una observación, es decir, que n=1 o X1=··· =Xn , lo cual solo ocurre con probabilidad cero, entonces ˆσ2 n= 0 lo que reeja el hecho de que en estos casos la función de verosimilitud es ilimitada cuando σ decrece hasta cero. Capítulo 2 El test Chi Cuadrado de Pearson La prueba Chi-cuadrado, también llamada Ji cuadrado ( χ2 ), fue publicada en el año 1900 por el matemático británico Karl Pearson y se considera una prueba no paramétrica que mide la discrepancia entre una distribución observada y otra teórica (bondad de ajuste), indicando en qué medida las diferencias existentes entre ambas, de haberlas, se deben al azar en el contraste de hipótesis. Una de las mayores ventajas de este test es que puede ser utilizado para trabajar tanto con distribuciones discretas como continuas como, por ejemplo, la distribución normal. Para una correcta realización de la prueba necesitaremos unos pequeños requerimientos tal y como veremos más adelante en profundidad: Necesitaremos que los datos se encuentren agrupados por categorías ó clases, es decir, que si nuestros datos no se encuentran agrupados originalmente, deberemos construir una tabla de frecuencias o histograma. También necesitaremos que el tamaño muestral sea sucientemente grande para que la prueba sea válida, lo cual supone una desventaja frente a otro tipo de test. 2.1. Contraste de una hipótesis nula simple Sea por lo tanto una muestra aleatoria simple X1, . . . , Xn donde todas las Xi siguen una misma distribución F y todas son observaciones independientes. Si la distribución F se corresponde a una variable aleatoria discreta, denotemos A1, . . . , Ak sus valores posibles. La idea del test consiste en realizar un recuento de las frecuencias de cada Ai en la muestra, y compararlas con las frecuencias que cabría esperar si fuera cierta la hipótesis nula, esto es, si la distribución fuera F0 . 9 10 CAPÍTULO 2. EL TEST CHI CUADRADO DE PEARSON En el caso de que F sea una distribución continua, se divide la recta real en k conjuntos que denotaremos de la misma forma como A1, . . . , Ak . De esta forma, podemos calcular las frecuencias de observaciones que pertenecen a cada Ai , y llevar a cabo el test igual que en el caso discreto, pues este test trabaja únicamente con k valores y sus frecuencias observadas y esperadas bajo F0 . Por tanto, consideremos Nj el número de observaciones muestrales que se encuentran en Aj . Suponiendo que la distribución fuese F0 , podemos calcular la probabilidad de cada Aj bajo esta distribución, la cual denotaremos como p0 j . De ser cierta la hipótesis nula H0:F=F0 , cada Nj tendría una distribución binomial con n intentos y probabilidad de éxito p0 j . Además, también bajo H0 , se tendría E(Nj) = np0 j , de modo que el número de observaciones que cabría esperar en el conjunto Aj sería np0 j por lo que denotaremos frecuencias observadas a los N1, . . . , Nk y frecuencias esperadas a np0 1, . . . , np0 k . Como estadístico de contraste usaremos el siguiente: D= k X j=1 (Nj−np0 j)2 np0 j que usaremos para rechazar H0 si D > c siendo c cierta constante, aceptando la hipótesis nula en caso contrario. Dicho estadístico se aproxima a una distribución Chi-cuadrado con (k−1) grados de libertad cuando el tamaño muestral n→ ∞ por lo que la constante c que considerábamos antes vendrá dada por c=χ2 k−1,α y, por ende, el test consistirá en: Rechazar H0:F=F0si D = k X j=1 (Nj−np0 j)2 np0 j > χ2 k−1,α Además, en el test chi-cuadrado, si la distribución es continua la forma en la que escogemos los conjuntos A1, . . . , Ak pueden suponer diferencias importantes por lo que a la hora de formar dichos conjuntos debemos tener en cuenta que: [3] Para que la aproximación de la distribución por la chi-cuadrado sea buena, es conveniente que las frecuencias esperadas sean grandes por lo que es recomendable que np0 j≥5,∀j∈1, . . . , k Para poder distinguir si la distribución verdadera es F0 , conviene que se divida en muchos conjuntos, por lo que se aconseja que k≥5 2.2. CONTRASTE DE UN MODELO PARAMÉTRICO 11 Observamos por lo tanto que ambas condiciones implican que el tamaño muestral n sea superior a 25, tal y como comentábamos al inicio del capítulo. Además, generalmente se respeta la primera condición, ya que sin ella la aproximación de la distribución no sería válida mientras que la segunda condición únicamente es utilizada para obtener resultados signicativos ya que no afecta a la validez de la prueba por lo que usualmente se toma k=√n 2.2. Contraste de un modelo paramétrico También se puede considerar una hipótesis nula del tipo H0:F∈Fθ/θ ∈Θ donde se plantea si la distribución pertenece al modelo paramétrico representado por Fθ donde θ representa los parámetros del modelo. Podemos considerar de la misma forma que para el caso de hipótesis nula simple las mismas observaciones independientes con la misma distribución F con Ai, i = 1, . . . , k la división de los conjuntos en el caso continuo y Nj la frecuencia de Aj . Considerando el estimador ˆ θ del parámetro θ podremos calcular las probabilidades pj(ˆ θ) de cada Aj bajo la distribución Fˆ θ y con frecuencias esperadas npi(ˆ θ), i = 1, . . . , k , por lo que podremos calcular el estadístico chi-cuadrado de la misma forma, con la diferencia de que ahora su aproximación será una distribución chi-cuadrado pero con menos grados de libertad como consecuencia de de la estimación de los parámetros. Por lo tanto, D= k X j=1 (Nj−npj(ˆ θ))2 npj(ˆ θ)≈χ2 k−q−1 donde q será el número de parámetros estimados que en el caso de la normal serán 2 si estimamos la media y la varianza. Con la intención de entender mejor el funcionamiento del estadístico chi-cuadrado, presentamos el siguiente ejemplo realizado en R , en el que se generan 30 números aleatorios que siguen una distribución normal y después se realiza la prueba, primero bajo hipótesis de que los datos pertenecen a una normal estándar, luego realizando una estimación de la media y varianza y, por último, empleando la función pearson.test() , cuyo funcionamiento se explica en profundidad en el anexo, pues también será la utilizada en las simulaciones del Capítulo 7. Por lo tanto, veamos en el siguiente código las diferencias de la implementación del estadístico y el resultado proporcionado por la función. 18 CAPÍTULO 3. TEST BASADOS EN LA FUNCIÓN DE DISTRIBUCIÓN acumulada Fn comparada con una función de distribución empírica F0 , o para comparar dos distribuciones empíricas. El estimador de Cramér von Mises acepta o rechaza la hipótesis nula H0 : Los datos provienen de una distribución conocida F0 , mediante el cálculo de la distancia entre las dos funciones, de modo que H0 debe aceptarse para valores de CV Mn pequeños: CV Mn=n ∞ Z −∞ [Fn(x)−F0(x)]2dF0(x) (3.1) donde n denota el número de datos. De forma alternativa, se puede demostrar que el estadístico se puede escribir de la siguiente forma: T=CV Mn=1 12n+ n X i=1 2i−1 2n−F0(Xi) Test de Watson Una pequeña modicación del anterior estadístico viene dada por este resultado, el cual emplea el estadístico: [6] U2=T−n(¯ F−1 2)2 , donde ¯ F=1 n n X i=1 F(Xi) 3.3. Test de Anderson-Darling El test de Anderson-Darling, nombrado así por Theodore Wilbur Anderson y Donald A. Darling , sirve para comprobar si una muestra procedente de una serie de datos proviene de un distribución determinada. Esta prueba se basa en la comparación de la distribución de probabilidades acumulada empírica (resultado de los datos) con la distribución de probabilidades acumulada teórica (denida en H0 ). Así pues, y de la misma forma que en los test anteriores, nuestras hipótesis serán que la muestra sigue una distribución de probabilidad especíca ( H0 ) o, si por el contrario, no lo hace ( H1 ). En nuestro caso, usaremos esta prueba para comprobar si la muestra sigue sigue una distribución normal. Este test, que es una variación de la prueba de Kolmogorov-Smirnov explicada anteriormente, le da más peso a las colas de lo que lo hace el dicho test K-S. 3.4. RELACIÓN KOLMOGOROV-CRAMER-ANDERSON 19 En el test de Kolmogorov-Smirno los valores críticos no dependen de la distribución especíca a la que estamos sometiendo el test, mientras que el de Anderson-Darling sí emplea dicha distribución a la hora de calcular los valores críticos. Como consecuencia, el test de Anderson-Darling tiene la ventaja de tener una mayor sensibilidad mientras que, por otro lado, los valores críticos variarán dependiendo de la distribución que hayamos escogido, por lo que deberán ser calculados para cada distribución, aunque hoy en día los valores críticos están tabulados para algunas distribuciones conocidas como la normal, la uniforme o la exponencial, entre otras. El estadístico de Anderson-Darling está denido como: A2=−n−S (3.2) donde S= n X i=1 (2i−1) n[ln(F(Yi)) + ln(1 −F(Yn+1−i))] y donde F es la función de distribución acumulativa de distribución especíca que estamos tratando de comprobar que siguen nuestros datos.[7] También debemos destacar que los datos de nuestra muestra ( Yi ) deben estar ordenados Por tanto, rechazaremos la hipótesis de que nuestra muestra sigue una determinada distribución si el valor del estadístico 3.2 es mayor que el valor crítico que, recordemos, depende de la distribución. 3.4. Relación Kolmogorov-Cramer-Anderson Puesto que los tres test presentados en este capítulo se basan en la función de distribución, prestaremos atención a las similitudes de ambos test, así como un breve repaso a sus principales diferencias. En primer lugar, debemos remontarnos a 1928 cuando Cramér propuso cuanticar la evidencia en contra de la hipótesis nula de que la función de distribución de nuestra muestra seguía una función de distribución determinada ( H0:F=F0 ). Para dichos cálculos, usó la siguiente cantidad: Z R [Fn(x)−F0(x)]2dK(x) donde K es una función no negativa. Esta cantidad le sirvió para medir la discrepancia entre Fn y F0 , ponderadas por K(x) . Cramér, correctamente, siguió el criterio de rechazar H0 cuando esta discrepancia era muy grande. 20 CAPÍTULO 3. TEST BASADOS EN LA FUNCIÓN DE DISTRIBUCIÓN A la vez que Cramér, von Mises llegó a resultados similares aunque de manera independiente. Más tarde, al rededor de 1937, Smirnov consiguió generalizar el criterio de Cramér y von Mises como: [8] W2 n=nZ R [Fn(x)−F0(x)]2ψ(x)dF0(x) (3.3) donde ψ era una función que ponderaba las diferencias entre las distribuciones Fn y F0 . La expresión anterior se trata de un estadístico, pues Fn está denida a partir de la muestra aleatoria y puede denotarse también como: W2 n=EF0[ψ(X)(Fn−F(X))2] donde EF0 signica que la esperanza se calcula suponiendo cierta H0 , es decir, que la distribución de X es F0 . W2 n cuantica el valor esperado de la diferencia cuadrática entre Fn y F0 ponderada por ψ(x) . Si escogiésemos ψ(x) = ψ(F0(x)) , es decir que la función que pondera las diferencias dependa de F0(x) entonces W2 n es libre de distribución bajo H0 . Así pues, sustituyendo ψ(x)=1 en 3.3 obtendríamos el estadístico de Cramér-von Mises 3.1. En cambio, Anderson y Darling propusieron utilizar para su estadístico la siguiente función: ψ(x) = 1 F0(x)(1 −F0(x)) Esta elección de ψ es lo que implica que las diferencias en las colas de la distribución tienen mayor peso. Además, y de la misma forma que en el caso anterior, basta con sustituir la denición de Anderson y Darling de ψ(x) en 3.3 para obtener el estadístico de su test. Pese a que los test de Cramér-von Mises y Anderson -Darling guardan una gran similitud con respecto al de Kolmogorov-Smirnov (K-S) y, por ende, al de Lilliefors, una de las principales diferencias que presentan es que en el test K-S los valores críticos son independientes de la distribución empírica escogida mientras que para las otras dos pruebas no lo son. Esta última diferencia implica que a la hora de realizar los tests, podamos tener tabulados los valores críticos del test K-S mientras que no podremos contar con esa ventaja en el resto de pruebas, pues dichos valores variarán según las distintas distribuciones. Capítulo 4 Test de normalidad basados en los momentos En este capítulo deniremos y analizaremos en profundidad los conceptos de asimetría y curtosis y su inuencia en la distribución normal. Brevemente repasaremos los test de normalidad basados en estos dos parámetros para nalmente estudiar el test de Jarque- Bera que combina ambos conceptos. 4.1. Asimetría Veamos qué es la asimetría estadística junto con algunos valores que la denen. 4.1.1. Denición Las medidas de asimetría permiten medir el grado de simetría (o asimetría) que presenta una distribución de probabilidad de una variable aleatoria sin la necesidad de realizar su representación gráca. Puesto que la simetría se mide en torno a un eje de simetría, consideramos como tal una recta paralela al eje de ordenadas que pasa por la media de la distribución. Por lo tanto, denimos a una distribución como simétrica si existe el mismo número de valores a la derecha que a la izquierda de la media y, por lo tanto, tiene el mismo número de desviaciones con signo positivo que con signo negativo. En cambio, diremos que hay asimetría positiva (o a la derecha) si la cola a la derecha de la media es más larga que la de la izquierda, es decir, si hay valores más separados de la media a la derecha. Del mismo modo, armaremos que hay asimetría negativa (o a la izquierda) si la cola a la izquierda de la media es más larga que la de la derecha, es decir, 21 22 CAPÍTULO 4. TEST BASADOS EN LOS MOMENTOS Figura 4.1: Representación gráca de ambos tipos de asimetría. si hay valores más separados de la media a la izquierda. 4.1.2. Medidas de asimetría A continuación presentamos algunas de las medidas de asimetría más conocidas y utilizadas. Coeciente de asimetría de Fisher Esta medida, que es una de las más empleadas, está basada en el tercer momento estándar dado que nos interesa mantener el signo de las desviaciones con respecto a la media, para obtener si son mayores las que ocurren a la derecha de la media que las de la izquierda. Aún así no se debe tomar el momento estándar con respecto a la media de orden 1 debido a que una simple suma de todas las desviaciones siempre es cero. Como ejemplo de esto tenemos que si los datos están agrupados en k clases obtenemos que: [9] k X i=1 fi(xi−µ) = k X i=1 fixi−µ k X i=1 fi=µ−µ= 0 donde xi representa la marca de la clase i -ésima y fi es la frecuencia relativa de dicha clase. Por lo tanto, se toman las desviaciones al cubo y el coeciente de asimetría de Fisher queda denido como: S=γ1=µ3 σ3= 1 nPn i=1(xi−¯x)3 (1 nPn i=1(xi−¯x)2)3/2 (4.1) donde µ3 es el tercer momento entorno a la media y σ la desviación típica. Así pues los valores de γ1 determinan los distintos tipos de asimetría que sufre nuestra distribución: 4.1. ASIMETRÍA 23 S > 0 : la distribución es asimétrica positiva o a la derecha. S < 0 : la distribución es asimétrica negativa o a la izquierda. Por otro lado, sabemos que si la distribución es simétrica entonces S= 0 pero el recíproco es falso, es decir, si S= 0 no podemos armar que la distribución es simétrica. Coeciente de asimetría de Pearson El coeciente de asimetría de Pearson está basado en que en aquellas distribuciones que son simétricas, la media de la distribución coincide con la moda. Dicho coeciente sólo se debe utilizar en distribuciones uniformes, unimodales y moderadamente asimétricas y se dene como: Ap=µ−moda µ donde µ es el momento ordinario de orden 1, es decir la media. Por lo tanto, si estamos ante una distribución simétrica la media coincidirá con la moda y, consecuentemente, Ap= 0 . En cambio si la distribución tiene asimetría positiva, el valor de la media estará por encima de la moda, lo que signica que el valor de Ap>0 . Coeciente de asimetría de Bowley-Yule Este coeciente se basa en la posición de los cuartiles y la mediana y está denido como: ABY =Q3/4+Q1/4−2Me Q3/4−Q1/4 Ahora bien, si tenemos en cuenta que el primer y el tercer cuartil están a la misma distancia de la media si la distribución es simétrica, obtenemos que, en caso de simetría, ABY = 0 , mientras que si la asimetría es positiva, ABY >0 . 4.1.3. Asimetría en la distribución normal Una vez denidos los distintos coecientes de asimetría veamos que el valor de éste cuando nos encontramos en una distribución es exactamente 0. Para simplicar los cálculos y pues cualquier tipo de distribución normal tendrá un coeciente de asimetría igual a 0 usaremos la distribución normal estándar. Proposición 4.1. El coeciente de asimetría de la distribución normal estándar es 0 24 CAPÍTULO 4. TEST BASADOS EN LOS MOMENTOS Demostración. En primer lugar, sea X∼N(0,1) . Tenemos que en el caso de la distribución normal estándar, la desviación típica tiene valor 1 y la media 0, por lo que el coeciente de asimetría quedará de la siguiente forma: S=µ3(X) = E[(X−E[X])3] = E[X3] = α3(X) Así pues, nos bastará con calcular el momento de tercer orden de la distribución normal en torno a la media, es decir, el valor de la tercera derivada de la Función Generatriz de Momentos de X para t= 0 . Por lo tanto y volviendo a que nos encontramos trabajando con la distribución estándar, tenemos que ϕX(t) = et2 2 . De donde la primera derivada será: ϕ0 X(t) = 2t 22tϕX(t) = tϕX(t) , y puesto que estamos centrándonos en t= 0 , ϕ0 X(t= 0) = 0 Para la segunda derivada tenemos que ϕ00 X(t) = 1ϕX(t) + tϕ0 X(t) = ϕX(t) + t2ϕX(t) = (1 + t2)ϕX(t) por lo que, haciendo t= 0 obtenemos que ϕ00 X(t= 0) = 1 para la tercera derivada: ϕ000 X(t)=2tϕX(t) + (1 + t2)ϕ0 X(t) = (3t+t3)ϕX(t) donde nuevamente haciendo t=0, ϕ000 X(t)=0 por lo tanto, µ3(X) = α3(X) = S(X)=0 4.1.4. Test de normalidad basado en la asimetría Una vez denido el concepto de asimetría, podremos usarlo para determinar si una variable aleatoria se ajusta a una distribución normal o no. Puesto que la distribución gaussiana es simétrica en torno a la media y que una de las formas más habituales de desviarse respecto de la normalidad es por falta de simetría, podemos construir un método de contraste basándonos en este concepto. 4.2. CURTOSIS 25 Por lo tanto, si la distribución de los datos es normal, el coeciente de asimetría seguirá una distribución asintótica normal de media igual a cero y varianza 6 n , por lo que emplearemos como estadístico de contraste el siguiente: rn 6S∼N(0,1) Por lo tanto, se rechazará la hipótesis nula de normalidad cuando el valor del estadístico anterior sea muy grande (positiva o negativamente), en comparación con los cuantiles de la normal estándar. 4.2. Curtosis 4.2.1. Denición Denimos la curtosis de una variable aleatoria como una característica de la forma que adquiere la función de distribución de dicha variable. Un valor grande en la curtosis signica que existe una gran acumulación de valores de la variable tanto cerca de la media como lejos de ella, a la vez que implica una menor frecuencia de los valores intermedios. Esta característica explica algunas formas de la función de distribución, con colas más gruesas, centro más apuntado y menor proporción de valores entre el pico y las colas. No debemos caer en el error de asumir que cuanto mayor sea la curtosis, mayor será la varianza o viceversa, pues no es cierto. 4.2.2. Medidas de curtosis Para calcular un valor de la curtosis empleamos el siguiente indicador: K=µ4 σ4= 1 nPn i=1(xi−¯x)4 (1 nPn i=1(xi−¯x)2)2 (4.2) Este coeciente viene denido por el cuarto momento centrado con respecto a la media ( µ4 ) y la desviación estándar ( σ ). Puesto que en la distribución normal se verica que µ4= 3σ4 , existe una segunda denición para el exceso de curtosis: g2=µ4 σ4−3 = K−3 de forma que este coeciente valdrá 0 cuando estemos trabajando con una distribución normal y se tome a esta como referencia para la curtosis. 26 CAPÍTULO 4. TEST BASADOS EN LOS MOMENTOS Así pues, si tomamos como referencia la distribución normal, una distribución podrá ser:[10] Leptocúrtica: se da cuando K > 3 y g2>0 y la distribución será más apuntada y con colas mas gruesas que la dist. Normal. Platicúrtica: se da cuando K < 3 y g2<0 y la distribución será menos apuntada y con colas menos gruesas que la dist. Normal. Mesocúrtica: se da cuando K= 3 y g2= 0 y la distribución seguirá una distribución normal. Figura 4.2: Representación gráca de los tipos de curtosis. Por último tenemos que existe otra forma de medir la curtosis, examinando dicho parámetro de la suma de variables aleatorias, es decir, suponiendo Y la suma de n variables aleatorias independientes con la misma distribución X , se tendrá: Kurt[Y] = Kurt[X] n 4.2.3. Curtosis en la distribución normal Del mismo modo que hicimos con la asimetría, centraremos nuestros cálculos en la distribución normal estándar, para comprobar que dicha distribución tiene un coeciente de curtosis K= 3 o, en su caso, g2= 0 . Proposición 4.2. El coeciente de curtosis K de la distribución normal estándar X∼ N(0,1) es 3, y, por lo tanto, g2= 0 . Demostración. El coeciente de curtosis K=µ4 σ4 es el momento central de cuarto orden en torno a la media, pero dado que la distribución normal estándar tiene media igual a 0 y desviación típica 1, coincidirá con la cuarta derivada de la Función Generatriz de Momentos para el valor t= 0 . 4.2. CURTOSIS 27 Aprovechando los cálculos de la demostración de la Proposición 4.1 tenemos que la cuarta derivada de la Función Generatriz de Momentos será: ϕIV X(t) = (3+3t2)ϕX(t)+(3t+t3)ϕ0 X(t) = (3+3t2)ϕX(t)+(3t+3t3)tϕX(t) = (3+6t2+t4)ϕX(t) donde, para t= 0 , K=µ4(X) = α4(X) = ϕIV X(t= 0) = 3 . De donde, por la propia denición de g2= 0 : g2=µ4 σ4−3 = K−3=3−3=0 Por lo tanto, y a la vista del resultado anterior, la distribución normal posee la forma tipo en cuanto a aplastamiento o apuntamiento y sirve de modelo para las demás distribuciones por lo que numerosos autores consideran que el valor referencia del coeciente de curtosis es 3 y no 0. 4.2.4. Test de normalidad basado en la curtosis Tras denir el concepto de curtosis, podremos usarlo para determinar si una variable aleatoria se ajusta a una distribución normal o no. Así pues, podremos rechazar la hipótesis de la normalidad si la curtosis de nuestra muestra es mucho mayor o mucho menor que cero, que es el valor que se le asigna a la distribución normal como vimos en la Proposición 4.2. Además, si la distribución de los datos muestrales es normal, entonces la curtosis tendrá distribución asintótica normal de media cero y varianza 24 n , por lo que se puede emplear como estadístico de contraste el siguiente: rn 24g2∼N(0,1) Por lo tanto, se rechazará la hipótesis nula de normalidad cuando el estadístico anterior sea muy grande (tanto en positivo como en negativo), en comparación con los cuantiles de N(0,1) . 34 CAPÍTULO 5. MÉTODOS GRÁFICOS Por último podemos relacionar los QQ-plots con los grácos de probabilidad ya que dichos grácos emplean la esperanza del k -ésimo orden de la distribución frente a los QQ- plots que usan los cuantiles, por lo que solo existe una diferencia importante entre ambos al trabajar con muestras de tamaño n pequeño. Fueron simuladas tres muestras de tamaño 100 que siguen una distribución normal, normal logarítmica y T-Student con 20 grados de libertad respectivamente y cuyos grácos son los siguientes: Figura 5.4: Ejemplo de QQ-Plot de una distribución normal. Como podemos observar, los datos se adecúan y ajustan a la gráca de y=x dado que, evidentemente, siguen una distribución normal. A continuación, veremos lo que ocurre con la distribución normal logarítmica, donde cabe esperar que los cuantiles de nuestros datos apenas coincidan con los de la normal, desviándose completamente de la trayectoria de la diagonal. 5.2. GRÁFICO DE CUANTILES TEÓRICOS 35 Figura 5.5: Ejemplo de QQ-Plot de una distribución log-normal. Pero no siempre se obtienen grácas con un desfase de los cuantiles de la muestra y los de la normal tan grande. Veamos un ejemplo de esta situación con el caso de la distribución T-Student con 20 grados de libertad: Figura 5.6: Ejemplo de QQ-Plot de una distribución T-Student. Observamos que en la gráca los puntos se adecúan bastante bien y siguen la diagonal, lo que puede dar lugar a pensar que nos encontramos ante una muestra que sigue una distribución normal del mismo modo que en la gura 5.4. En cambio, nuestra muestra sigue una distribución T-Student con 20 grados de libertad, con lo que comprobamos que los métodos grácos pueden ser confusos y dar lugar a más errores que los estudiados previamente. 36 CAPÍTULO 5. MÉTODOS GRÁFICOS Capítulo 6 Test de Shapiro-Wilk y variante de Shapiro-Francia En este capítulo, estudiaremos el test de Shapiro-Wilk junto con su variante de Shapiro- Francia, además de que observaremos las diferencias de resultados entre ambos test, realizando ambas pruebas a una misma muestra. 6.1. Test de Shapiro-Wilk En 1965 los estadísticos Samuel Shapiro y Martin Wilk publicaron el test que lleva hoy su nombre y que, además, es considerado uno de los más potentes. De la misma forma que los test anteriores, contrasta las hipótesis de que el conjunto de datos de estudio proviene de una distribución normal ( H0 ) o que, en cambio sigue cualquier otra distribución ( H1 ). Para realizar el test, en primer lugar debemos considerar los datos estandarizados: Zi=Xi−ˆ X S, i ∈1, . . . , n Y, una vez estandarizados podemos construir el estadístico de la siguiente manera:[13] W= [n/2] X i=1 ai,n Z(n−i+1):n−Zi:n donde Z1:n< . . . < Zn:n es la muestra ordenada de los datos estandarizados y ai , n constantes. Se trata de calcular las distancias entre los datos de la muestra ordenada, simétricos respecto de la mediana, es decir, la distancia entre el primero y el último, el 37 38 CAPÍTULO 6. TEST SHAPIRO-WILK Y SHAPIRO-FRANCIA segundo y el penúltimo, y así sucesivamente, de forma genérica, se comparan el Zi:n y el Z(n−i+1):n . La nalidad es comparar estas distancias con las que habría en una muestra de observaciones normales. Además, el coeciente de determinación guarda una estrecha relación con el QQ-Plot, ya que lo podemos interpretar como el coeciente de determinación de los puntos representados en el QQ-Plot. Por lo tanto, rechazaremos la normalidad cuando el estadístico de Shapiro-Wilk sea pequeño, ya que eso será debido a que los puntos del QQ-Plot se alejan de la recta esperada Para visualizar de mejor forma la aplicación del test de Shapiro-Wilk nos ayudamos del siguiente ejemplo. Ejemplo: Se han medido los tiempos de ejecución de cierta tarea por parte de seis trabajadores, con los siguientes resultados en horas: 6,55,5,58,4,33,2,53,5,49,8,64,3,87,4,82 Se quiere saber si las mediciones provienen de una población normal. Solución: Ciertos parámetros de estas distribuciones, como media, desviación típica u otros, que en las secciones anteriores se determinaron con el propósito de facilitar la comparación, ahora son irrelevantes, puesto que el test de Shapiro-Wilk comienza con una estandarización de los datos. Ayudándonos de R, utilizaremos la función shapiro.test: x = c(6.55, 5.58, 4.33, 2.53, 5.49, 8.64, 3.87, 4.82) shapiro.test(x) Ejecutando lo anterior obtenemos que: Shapiro-Wilk normality test data: x W = 0.9827, p-value = 0.9643 El test nos devuelve un p-valor muy grande, 0.9643, que es claramente mayor que 0.05 ya que estamos suponiendo un nivel de signicación del 5 por ciento, por lo tanto, no tenemos evidencias sucientes para rechazar la hipótesis de normalidad. 6.2. TEST DE SHAPIRO-FRANCIA 39 6.2. Test de Shapiro-Francia Introducido por Samuel Shapiro y R. S. Francia en 1972, es una simplicación del test de Shapiro-Wilk, por lo que su funcionamiento es muy similar y son pocas las diferencias entre ambos. Así pues, sea x(i) el i -ésimo valor ordenado de nuestra muestra (véase como ejemplo que si nuestra muestra fuese 5, -1, 2, 0 nuestro x2= 0 ), tendremos que el estadístico está denido como:[14] W0=Cor2(x(i), mi) = σ(x,m) σxσm =Pn i=1(x(i)−¯x)(mi−¯m) q(Pn i=1(x(i)−¯x)2)(Pn i=1(mi−¯m)2) Además, se tiene que dicho estadístico es el cuadrado del coeciente de correlación entre los estadísticos de orden muestrales (es decir, nuestros cuantiles en un gráco QQ-plot) x(1), . . . , x(n) y los valores esperados de los estadísticos de orden de la distribución N(0,1) m1, . . . , mn , es decir, nuestros cuantiles teóricos en el QQ-plot. Del mismo modo que hicimos con el test de Shapiro-Wilk, podemos seguir el mismo ejemplo que el realizado para dicho test y ver que ocurre para la variante de Shapiro- Francia. Ejemplo: Sean los mismos tiempos de ejecución de cierta tarea por parte de seis operarios, con los mismos resultados en horas: 6,55,5,58,4,33,2,53,5,49,8,64,3,87,4,82 Veamos si siguen una distribución normal. Solución: Ayudándonos de R, utilizaremos la función sf.test, de la librería nortest : library(nortest) x = c(6.55, 5.58, 4.33, 2.53, 5.49, 8.64, 3.87, 4.82) sf.test(x) Que, tras ejecutar, nos dará el siguiente resultado: Shapiro-Francia normality test data: x W = 0.96053, p-value = 0.826 40 CAPÍTULO 6. TEST SHAPIRO-WILK Y SHAPIRO-FRANCIA El test nos devuelve un p-valor grande, 0.826, que es claramente mayor que 0.05 ya que estamos suponiendo un nivel de signicación del 5 por ciento igual que para el test de Shapiro-Wilk, por lo tanto, no tenemos evidencias sucientes para rechazar la hipótesis de normalidad. Capítulo 7 Estudio comparativo mediante simulaciones En el siguiente capítulo se comparan los resultados obtenidos en las simulaciones realizadas en R , donde para diferentes distribuciones de probabilidad se realizan los distintos test estudiados previamente y se calcula la potencia de cada método aplicado a cada distribución para distintos tamaños muestrales. Para cada distribución fueron realizadas 1000 simulaciones al 5 por ciento de signicación por lo que cuanto más cercano sea el tanto por uno de rechazo a 0.05 más signicativo será el resultado de que se acepta la hipótesis nula. Dichos resultados han sido tabulados de forma que en cada tabla se observan las potencias de los distintos test para tamaños muestrales diferentes de muestras que siguen una distribución especíca. Adicionalmente, los datos de las tablas son representados mediante grácas donde las potencias de los test se superponen para tener una forma más visual de poder compararlos. Finalmente, se concluye el capítulo con las conclusiones extraídas de las simulaciones. 7.1. Resultados de las simulaciones Para comenzar y como ejemplo, podemos observar lo que ocurre al aplicar los distintos test a una distribución normal. 41 42 CAPÍTULO 7. ESTUDIO COMPARATIVO MEDIANTE SIMULACIONES Figura 7.1: Tantos por uno de rechazo para una distribución normal. Al tratarse de una distribución normal, la totalidad de los test dan resultados muy similares y cercanos al 0.05 de rechazo, pues todos detectan la normalidad de la distribución, pero veamos que ocurre cuando realizamos el test a una muestra que sigue otra distribución. Figura 7.2: Tantos por uno de rechazo para una distribución exponencial. Claramente el tanto por uno de rechazo se dispara con respecto al ejemplo anterior pues estamos trabajando con una distribución exponencial, pero podemos sacar algunas conclusiones de los datos anteriores ya que podemos observar que cuanto mayor es el 7.1. RESULTADOS DE LAS SIMULACIONES 43 tamaño muestral, mayor es el rechazo de la hipótesis nula de que nuestra distribución sigue una normal. Además, observamos que la prueba de Shapiro-Wilk es la que mejor y más rápido detecta que los datos no siguen una distribución normal, aunque no hay una gran diferencia con respecto al resto de pruebas dado que esta distribución solo tiene una cola, por lo que es más fácil de detectar para los test. Veamos que ocurre con otras distribuciones. A continuación, analizaremos que ocurre en el caso de que nuestra muestra siga una distribución T-Student. Dado que esta distribución se aproxima a una distribución normal cuanto mayor es el número de sus grados de libertad, por lo que la lógica nos hace pensar que el rechazo de la hipótesis nula irá siendo menor según aumenten los grados de libertad. Veamos que ocurre con dos distribuciones T-Student con 5 y 50 grados de libertad respectivamente y veremos que diferencias hay entre ambas. Figura 7.3: Tantos por uno de rechazo para una distribución T-Student con 5 grados de libertad. Como podemos observar, para tamaños muestrales pequeños, muchos de los test son incapaces de detectar que nos encontramos ante una distribución no normal, incluso algunos, como el Chi-Cuadrado de Pearson precisa de un tamaño muestral muy grande para obtener un resultado signicativo que permita rechazar la hipótesis nula. Esta situación se complica, naturalmente, cuando aumentan los grados de libertad: 50 CAPÍTULO 7. ESTUDIO COMPARATIVO MEDIANTE SIMULACIONES 7.2. Conclusiones Analizados los resultados arrojados del comportamiento de los test ante diferentes distribuciones, podemos extraer algunas conclusiones sobre cuales de ellos arrojan resultados mas ables y tienen una mayor potencia respecto a los demás. Observamos que, salvo para la distribución binomial, el test de Shapiro-Wilk es el que mejor porcentaje de rechazo arroja para cada tamaño muestral, con lo que podemos concluir que es el test que más potencia tiene de los estudiados. Además su variante de Shapiro-Francia también presenta muy buenos resultados para la mayoría de distribuciones. Respecto al test de Jarque-Bera, tal y como analizábamos previamente, su abilidad se ve trastocada cuando se enfrenta a distribuciones con parámetros de asimetría y curtosis similares a la normal, por lo que puede arrojar resultados erróneos que nos hagan aceptar la hipótesis de normalidad cuando esta no sea cierta. Por otro lado, el test Chi Cuadrado de Pearson se presenta como el más able en la distribución binomial, pero su rendimiento baja notoriamente en el resto de distribuciones, sobre todo en las T-Student y Chi Cuadrado, independientemente de sus grados de libertad que, como vimos, también afectan a los resultados obtenidos. El resto de test arrojaron resultados muy similares, por lo que el test de Shapiro-Wilk y su variante de Shapiro-Francia se presuponen como los test más ables a la hora de realizar un contraste de normalidad cuando desconozcamos la distribución de la muestra a la que le apliquemos el test. Por último, recalcar el papel fundamental que juega el tamaño muestral a la hora de realizar las pruebas, ya que es evidente, tanto en las tablas como en las grácas, que cuanto mayor es el volumen de la muestra, mejores y más ables son los resultados de los porcentajes de rechazo independientemente del test realizado. Apéndice A Códigos de R En el siguiente anexo se presentan todas las líneas de código utilizadas tanto en las simulaciones como en la obtención de ejemplos a lo largo de todo el trabajo, junto con las explicaciones de las funciones y paquetes utilizados. A continuación se presentan los códigos utilizados en R para obtener los resultados simulados. En primer lugar, destacaremos las funciones necesarias para realizar los diferentes test, ya que más tarde observaremos que el código es el mismo para las distintas distribuciones, con la única variación de los comandos necesarios para la generación de números aleatorios que siguen cada una de las distribuciones especícas, dichas funciones proceden de los paquetes normtest y nortest . Así pues: ad.test(x) # Para realizar el test de Anderson-Darling cvm.test(x) # Para realizar el test de Cramer-von-Mises sf.test(x) # Para realizar el test de Shapiro-Francia shapiro.test(x) # Para realizar el test de Shapiro-Wilk lillie.test(x) # Para realizar el test de Lilliefors pearson.test(x) # Para realizar el test Chi-Cuadrado de Pearson jb.norm.test(x) # Para realizar el test de Jarque-Bera A continuación se explica brevemente el funcionamiento de las funciones anteriores, junto con los parámetros que requieren. La función ad.test() emplea el estadístico de prueba A=−n−1 nPn i=1(2i−1)(ln(F(x(i)))+ ln(1 −F(x(n−i+1)))) , tal y como denimos en el capítulo 3 donde ¯x y s son la media y la desviación estándar de los valores de los datos. El p-valor se calcula a partir del estadístico modicado Z=A(1,0 + 0,75 n+2,25 n2) de acuerdo con la tabla 4.9 en Stephens (1986). [15] 51 52 APÉNDICE A. CÓDIGOS DE R La función cvm.test() requiere de varios parámetros para su correcto funcionamiento: x : el vector que contiene los valores de nuestra muestra. null : es la distribución especíca que queremos comprobar que sigue nuestra muestra estimated : es un valor lógico, de ser "TRUE"(que viene por defecto) asumirá que los valores de media y varianza son estimados, de ser "FALSE"los da por conocidos. Las funciones shapiro.test() y sf.test() tienen un funcionamiento idéntico cuya única diferencia es la misma que en el propio método ya que para realizar el test de Shapiro- Francia sí estimará media y varianza. En ambos casos el p-valor viene calculado por la fórmula dada por Royston (1993). [16] La función pearson.test() requiere de tres parámetros para su correcto funcionamiento: x : el vector que contiene los valores de nuestra muestra. n.classes : el número de intervalos en los que se dividirá la muestra (por defecto se dividirá en k intervalos con k= d(2(n2/5)) ) adjust : es un valor lógico, de ser "TRUE"(que viene por defecto) la función realizará el test a una distribución chi-cuadrado con n.classes-3 grados de libertad (es decir, como si estuviésemos estimando la media y la varianza) y de ser "FALSE"serán n.classes-1 grados de libertad Dicha función usará como estadístico P=P(Ci−Ei)2 Ei donde Ci son las frecuencias observadas y Ei las esperadas. La función jb.norm.test(x) requiere de los siguientes parámetros: x : el vector que contiene los valores de nuestra muestra. n.repl : el número de replicaciones aleatorias por el método de Monte Carlo. Esta función se basa en el siguiente estadístico: JB =n 6S2+1 4(K−3)2 tal y como está explicado en el capítulo 4, con S el coeciente de asimetría y de la muestra y K el de curtosis. Por último, el p-valor se calcula mediante simulación de Monte Carlo. A continuación, se muestra como fueron generadas las variables aleatorias para cada una de las simulaciones. Dichas muestras son números aleatorios que siguen una distribución especíca tal como sigue: 53 x=rnorm(muestra) # Dist. Normal y=runif(muestra) # Dist. Unif z=rexp(muestra) # Dist. Expo chi1=rchisq(muestra,5) #Chi 5 grados chi2=rchisq(muestra,50) #Chi 50 grados t1=rt(muestra,5) # TStudent 5 grados t2=rt(muestra,50) # TSTudent 50 grados bi=rbinom(muestra,muestra,0.5) # Binomial l=rlnorm(muestra) # Log-normal b=rbeta(muestra,1.5,1.5) # Beta g=rgamma(muestra,1.5,1.5) # Gamma Como podemos observar, el tamaño muestral viene determinado por la variable muestra cuyo valor se modica dentro de un bucle como veremos a continuación. Como ejemplo, veremos como sigue el código para la realización del test a la muestra x que sigue una distribución normal, ya que el código es exactamente el mismo y se sigue de la misma forma para el resto. set.seed(123456) n=c(10,30,50,100,500) # Tamaño muestral ns=1000 # Numero de muestras simuladas # Librerías necesarias para los distintos test library(normtest) library(nortest) # Creación de vectores donde irán los p-valores de cada test pv_anderson_norm=c() pv_cvm_norm=c() pv_sf_norm=c() pv_lil_norm=c() pv_pear_norm=c() pv_sw_norm=c() pv_jb_x=c() # Creación de vectores de potencias de cada test 54 APÉNDICE A. CÓDIGOS DE R pot_anderson_norm=c() pot_cvm_norm=c() pot_sf_norm=c() pot_lil_norm=c() pot_pear_norm=c() pot_sw_norm=c() pot_jb_x=c() # Bucle de la realización de los test a las muestras for (i in 1:length(n)){ muestra=n[i] for (is in 1:ns){ x=rnorm(muestra,0,1) anderson_norm=ad.test(x) cvm_norm=cvm.test(x) sf_norm=sf.test(x) lil_norm=lillie.test(x) pear_norm=pearson.test(x) sw_norm=shapiro.test(x) jb_norm=jb.norm.test(x) # Vectores con los p-valores pv_anderson_norm[is]= anderson_norm$p.value pv_cvm_norm[is]=cvm_norm$p.value pv_sf_norm[is]=sf_norm$p.value pv_lil_norm[is]=lil_norm$p.value pv_pear_norm[is]=pear_norm$p.value pv_sw_norm[is]=sw_norm$p.value pv_jb_x[is]=jb_norm$p.value } # Vectores con las potencias de cada test pot_anderson_norm[i]= sum(pv_anderson_norm<0.05)/ns pot_cvm_norm[i]=sum(pv_cvm_norm<0.05)/ns 55 pot_sf_norm[i]=sum(pv_sf_norm<0.05)/ns pot_lil_norm[i]=sum(pv_lil_norm<0.05)/ns pot_pear_norm[i]=sum(pv_pear_norm<0.05)/ns pot_sw_norm[i]=sum(pv_sw_norm<0.05)/ns pot_jb_x[i]=sum(pv_jb_x<0.05)/ns } # Creación de un data frame con los valores de las potencias Potencias_Normal= data.frame(rbind(pot_pear_norm,pot_lil_norm,pot_cvm_norm, pot_anderson_norm,pot_jb_x,pot_sw_norm,pot_sf_norm)) rownames(Potencias_Normal)<-c("Test Chi Cuadrado de Pearson","Test de Lilliefors", "Test de Cramer-von-Mises","Test de Anderson-Darling","Test de Jarque-Bera", "Test de Shapiro-Wilk","Test de Shapiro-Francia") colnames(Potencias_Normal)<-c("Tamaño 10", "Tamaño 30", "Tamaño 50", "Tamaño 100", "Tamaño 500") #Visualización del data frame View(Potencias_Normal) Una vez generado el data frame con los datos correspondientes, usamos el siguiente código para la creación de las tablas y las grácas para representar los resultados obtenidos. De la misma forma que para el código anterior, únicamente se presenta el código empleado para una distribución normal, pues la única variación que hay es el uso de los data frames correspondientes a las distintas distribuciones. # PARA LAS GRAFICAS # Division del área de trabajo para que salgan 4 gráficas en un mismo plot Conf2x2 = matrix(c(1:4), nrow=2, byrow=TRUE) layout(Conf2x2) layout.show(4) # GRAFICA EXPONENCIAL plot(1:5,xlab="Tamaños Muestrales",ylab="Porcentaje de Rechazo",main= "Dist. Exponencial", ylim=c(0,1),xaxt="n",las=2) 56 APÉNDICE A. CÓDIGOS DE R axis(1, at = 1:5, labels = n) lines(1:5,pot_pear_exp,lwd=2,col="orange") lines(1:5,pot_jb_exp,lwd=2,col="blue") lines(1:5,pot_sw_exp,lwd=2) lines(1:5,pot_sf_exp,lwd=2,col="yellow") lines(1:5,pot_lil_exp,lwd=2,col="green") lines(1:5,pot_cvm_exp,lwd=2,col="red") lines(1:5,pot_anderson_exp,lwd=2,col="light blue") # GRAFICA BINOMIAL plot(1:5,xlab="Tamaños Muestrales",ylab="Porcentaje de Rechazo",main= "Dist. Binomial", ylim=c(0,1),xaxt="n",las=2) axis(1, at = 1:5, labels = n) lines(1:5,pot_pear_binom,lwd=2,col="orange") lines(1:5,pot_jb_binom,lwd=2,col="blue") lines(1:5,pot_sw_binom,lwd=2) lines(1:5,pot_sf_binom,lwd=2,col="yellow") lines(1:5,pot_lil_binom,lwd=2,col="green") lines(1:5,pot_cvm_binom,lwd=2,col="red") lines(1:5,pot_anderson_binom,lwd=2,col="light blue") # GRAFICA UNIFORME plot(1:5,xlab="Tamaños Muestrales",ylab="Porcentaje de Rechazo",main= "Dist. Uniforme", ylim=c(0,1),xaxt="n",las=2) axis(1, at = 1:5, labels = n) lines(1:5,pot_pear_unif,lwd=2,col="orange") lines(1:5,pot_jb_unif,lwd=2,col="blue") lines(1:5,pot_sw_unif,lwd=2) lines(1:5,pot_sf_unif,lwd=2,col="yellow") lines(1:5,pot_lil_unif,lwd=2,col="green") lines(1:5,pot_cvm_unif,lwd=2,col="red") lines(1:5,pot_anderson_unif,lwd=2,col="light blue") # GRAFICA LOG plot(1:5,xlab="Tamaños Muestrales",ylab="Porcentaje de Rechazo",main= "Dist.Log-Normal", ylim=c(0,1),xaxt="n",las=2) 57 axis(1, at = 1:5, labels = n) lines(1:5,pot_pear_log,lwd=2,col="orange") lines(1:5,pot_jb_log,lwd=2,col="blue") lines(1:5,pot_sw_log,lwd=2) lines(1:5,pot_sf_log,lwd=2,col="yellow") lines(1:5,pot_lil_log,lwd=2,col="green") lines(1:5,pot_cvm_log,lwd=2,col="red") lines(1:5,pot_anderson_log,lwd=2,col="light blue") # GRAFICA TSTUDENT5 plot(1:5,xlab="Tamaños Muestrales",ylab="Porcentaje de Rechazo",main= "Dist. T-Student 5 Grados de Libertad", ylim=c(0,1),xaxt="n",las=2) axis(1, at = 1:5, labels = n) lines(1:5,pot_pear_t5,lwd=2,col="orange") lines(1:5,pot_jb_t5,lwd=2,col="blue") lines(1:5,pot_sw_t5,lwd=2) lines(1:5,pot_sf_t5,lwd=2,col="yellow") lines(1:5,pot_lil_t5,lwd=2,col="green") lines(1:5,pot_cvm_t5,lwd=2,col="red") lines(1:5,pot_anderson_t5,lwd=2,col="light blue") # GRAFICA TSTUDENT 50 plot(1:5,xlab="Tamaños Muestrales",ylab="Porcentaje de Rechazo",main= "Dist. T-Student 50 Grados de Libertad", ylim=c(0,1),xaxt="n",las=2) axis(1, at = 1:5, labels = n) lines(1:5,pot_pear_t50,lwd=2,col="orange") lines(1:5,pot_jb_t50,lwd=2,col="blue") lines(1:5,pot_sw_t50,lwd=2) lines(1:5,pot_sf_t50,lwd=2,col="yellow") lines(1:5,pot_lil_t50,lwd=2,col="green") lines(1:5,pot_cvm_t50,lwd=2,col="red") lines(1:5,pot_anderson_t50,lwd=2,col="light blue") # GRAFICA CHI CUADRADO 5 plot(1:5,xlab="Tamaños Muestrales",ylab="Porcentaje de Rechazo",main= "Dist. Chi Cuadrado 5 Grados de Libertad", ylim=c(0,1),xaxt="n",las=2) 58 APÉNDICE A. CÓDIGOS DE R axis(1, at = 1:5, labels = n) lines(1:5,pot_pear_c5,lwd=2,col="orange") lines(1:5,pot_jb_c5,lwd=2,col="blue") lines(1:5,pot_sw_c5,lwd=2) lines(1:5,pot_sf_c5,lwd=2,col="yellow") lines(1:5,pot_lil_c5,lwd=2,col="green") lines(1:5,pot_cvm_c5,lwd=2,col="red") lines(1:5,pot_anderson_c5,lwd=2,col="light blue") # GRAFICA CHI CUADRADO 50 plot(1:5,xlab="Tamaños Muestrales",ylab="Porcentaje de Rechazo",main= "Dist. Chi Cuadrado 50 Grados de Libertad", ylim=c(0,1),xaxt="n",las=2) axis(1, at = 1:5, labels = n) lines(1:5,pot_pear_c50,lwd=2,col="orange") lines(1:5,pot_jb_c50,lwd=2,col="blue") lines(1:5,pot_sw_c50,lwd=2) lines(1:5,pot_sf_c50,lwd=2,col="yellow") lines(1:5,pot_lil_c50,lwd=2,col="green") lines(1:5,pot_cvm_c50,lwd=2,col="red") lines(1:5,pot_anderson_c50,lwd=2,col="light blue") # PARA LAS TABLAS # Uso de la librería necesaria library(gt) # Implementación de los datos en las tablas tab_1_ad<-gt(Potencias_Normal, rowname_col = "Test",rownames_to_stub = TRUE) tab_2_ad<- tab_1_ad %>% tab_header( title = md("Tantos por uno de rechazo de los diferentes test de normalidad para muestras simuladas que siguen una **distribución Normal**"), subtitle = md("Test hechos al 5% de significación y 1000 simulaciones para cada tamaño muestral") ) %>% 59 tab_style( style = list( cell_fill(color= "light blue",alpha=0.3), cell_text(weight = "normal") ), locations = cells_body( columns = TRUE, rows = c(1,3,5,7))) %>% tab_style( style = list( cell_fill(color= "light blue",alpha=0.3), cell_text(weight = "normal") ), locations = cells_stub( rows = c(1,3,5,7))) tab_2_ad Para el apartado de métodos grácos, el siguiente código fue empleado para obtener los histogramas que aparecen en dicho capítulo: set.seed(123456) # Creación de muestras con distintos tamaños x= seq(-4, 4, length = 500) normal10=rnorm(10,0,1) normal50=rnorm(50,0,1) normal250=rnorm(250,0,1) normal500=rnorm(500,0,1) Conf2x2 = matrix(c(1:4), nrow=2, byrow=TRUE) layout(Conf2x2) layout.show(4) # Creación de los histogramas con la función de densidad de la normal superpuesta hist(normal10, freq = F,xlim = c(-4, 4),main="Tamaño muestral n=10",xlab="") lines(x, dnorm(x, mean = 0, sd = 1), col = "red", lty = 2, lwd = 2)