Pruebas para normalidad

library(MASS)
library(nortest)
library(normtest)
library(DescTools)
library(moments)

1. Introducción

Muchos procedimientos estadísticos requieren partir del supuesto de normalidad. Al seleccionar una muestra aleatoria de una distribución es necesario validar si el supuesto de normalidad se cumple. Se presentan en este trabajo algunas de las pruebas que más se utilizan para probar la normalidad. Se presentan tanto los procedimientos manuales como utilizando los paquetes disponibles en R.

Para ilustrar las pruebas se utilizan 40 datos que corresponden al tiempo (minutos) que tardaron 40 estudiantes en resolver un problema de cálculo.

datos<-c(3.5,7.2,2.6,3.8,6.2,2.8,5.4,3.9,3.9,5.0,4.1,3.4,4.2,3.6,4.7,3.7,3.9,4.5,2.7,3.6,4.2,4.0,4.4,1.5,3.3,3.3,4.1,6.6,2.2,2.9,5.2,2.8,5.4,4.1,2.0,5.5,2.7,4.4,5.9,2.8)
n<-length(datos)

2. Métodos gráficos

Los gráficos son una forma visual de observar el comportamiento de una variable. Obviamente la observación no basta para establecer si una variable cumple con el supuesto de normalidad, pero en algunos casos puede ser una guía.

2.1 histograma

Un histograma presenta la forma general de una distribución así como su dispersión, dividiendo los datos en \(m\) clases mutuamente excluyentes. Es una representación gráfica en forma de barras, que permiten tener una visión general de la forma de la distribución de los datos muestrales. Cuando los datos muestrales provienen de una distribución normal se espera que las mayores frecuencias se encuentren en los valores centrales y que tienden a decrecer hacia los valores extremos. Se le suele agregar una curva de la distribución normal con media y desviación estándar de los datos bajo estudio.

hist(datos,col="green",border="red",freq=FALSE,main="")
curve(dnorm(x,mean(datos),sd(datos)),add=T,lwd=2)

2.2 Diagrama de tallos y hojas

El diagrama de tallos y hojas Stem and Leaf Diagram trata de mostrar la distribución de los datos y su representación gráfica, la cual es muy parecida a la de un histograma y que muestra una visión general de los datos.

stem(datos,scale=1)

## 
##   The decimal point is at the |
## 
##   1 | 5
##   2 | 026778889
##   3 | 33456678999
##   4 | 0111224457
##   5 | 024459
##   6 | 26
##   7 | 2

2.3 Diagrama de cajas

Los diagramas de cajas (boxplot) son una forma de resumen de los datos que muestra su distribución, centramiento, variabilidad, posibles valores atípicos y asimetría. Para su construcción utiliza cinco valores(tres cuartiles y valores máximo y mínimo).

boxplot(datos,col="yellow",border="red",horizontal=TRUE,notch=FALSE)

2.4 Gráfico de cuantil-cuantil

El gráfico cuantil-cuantil permite comparar la distribución empírica de un conjunto de datos con la distribución teórica normal.

Si los datos de una muestra provienen de una distribución normal, cabría esperarse que los cuantiles empíricos se correspondan linealmente con los correspondientes cuantiles teóricos de la distribución normal. Entonces, si se grafican los percentiles teóricos en el eje horizontal y los cuantiles empíricos en el eje vertical, el intercepto corresponde a la media y la pendiente corresponde a la desviación estándar de los datos de la muestra.

Si los datos observados en la muestra son \(x_{1},x_{2}, \cdots,x_{n}\) , el gráfico de cuantiles normales se construye con los puntos \((z_{i},x_{(i)})\), donde \(x_{(i)}\), \(i=1,2,3,\cdots,n\) son los datos ordenados en forma ascendente de magnitud, y \(z_{i}\) es el percentil correspondiente a la función de distribución normal estándar, obtenido mediante la transformación: \(z_{i}=\Phi^{-1}\left( \dfrac{i-0.375}{n+0.25}\right)\)

Donde \(\Phi\) es la función de distribución normal estándar.

Las distancias de los puntos del gráfico respecto de la línea recta muestran algunas características de los datos, tales como asimetría, apuntamiento y la presencia de datos atípicos.

x<-sort(datos)
i<-c(1:40)
fi<-round((i-0.375)/(n+0.25),4)
do<-qnorm(fi)

par(mfrow=c(1,2))
plot(do,x,main="qq calculado")
qqline(datos,lwd=2)
qqnorm(datos)
qqline(datos,lwd=2)

3. Métodos analíticos

Los estadísticos de asimetría y kurtosis se pueden utilizar para determinar si en los datos de la muestra se presenta desviación a la normalidad.

El coeficiente de asimetría, se define como:

\[\begin{equation*} g_{1}=\dfrac{m_{3}}{m_{2}^{\frac{3}{2}}}=\dfrac{\frac{1}{n}\sum (x_{i}-\bar{x})^{3}}{(\frac{1}{n}\sum (x_{i}-\bar{x})^{2})^{\frac{3}{2}}} \end{equation*}\]

Cuando el coeficiente de asimetría tiende a cero, la distribución de los datos tiende a ser simétrica respecto de la media. Si el coeficiente de asimetría se encuentra entre -1 y 1, se considera que existe una ligera desviación de la normalidad. En la medida en que el coeficiente de asimetría tienda a aumentar, se puede afirmar que los datos muestrales se alejan de la normalidad.

El coeficiente de kurtosis (o apuntamiento) se define como:

\[\begin{equation*} g_{2}=\dfrac{m_{4}}{m_{2}^{2}} =\dfrac{\frac{1}{n}\sum (x_{i}-\bar{x})^{4}}{(\frac{1}{n}\sum (x_{i}-\bar{x})^{2})^{2}} \end{equation*}\]

El coeficiente de kurtosis es una medida de qué tan apuntada o achatada es la distribución de la variable. Una variable que presente distribución normal tendrá un coeficiente de kurtosis de 3. En la medida en que se aleje de este valor se tendrá una mayor desviación respecto a la normalidad.

m<-mean(datos)
g1<-(sum(((datos-m)^3))/n)/(((sum((datos-m)^2))/n)^(3/2));g1

## [1] 0.4670349

g2<-(sum(((datos-m)^4))/n)/(((sum((datos-m)^2))/n)^(2));g2

## [1] 3.016572

4 Pruebas estadísticas

Las pruebas estadísticas son técnicas de estadística inferencial que se utilizan para obtener una conclusión obtenida a partir de un conjunto de datos muestrales puedan ser inferidas a la población de la cual se obtuvo la muestra.

\(H_{0}\): la muestra proviene de una población con distribución normal

\(H_{.1}\): la muestra no proviene de una población con distribución normal

4.1 Prueba chi-cuadrado de bondad de ajuste

La prueba chi-cuadrado de bondad de ajuste fue propuesta en 1900 por Karl Pearson . Dada una muestra aleatoria, se divide en \(m\) grupos mutuamente excluyentes y se obtiene la frecuencia observada en cada una de las clases. Se determinan las frecuencias esperadas bajo el supuesto de normalidad o de otra distribución que se suponga. Las frecuencias observadas y esperadas se comparan mediante el estadístico de prueba:

\[\begin{equation*} \chi_{p}^{2}=\sum \dfrac{(fo_{i}-fe_{i})^{2}}{fe_{i}} \end{equation*}\]

Donde : \(fo_{i}\) son las frecuencias observadas en la muestra en cada una de las $ m $ clases y$ fe_{i}$ son las frecuencias esperadas bajo la distribución propuesta.

Este estadístico presenta una distribución \(\chi^{2}\) con $ m-1 $ grados de libertad si se especifican de antemano los parámetros de la distribución propuesta o $ m-p-1 $ si los parámetros de la distribución propuesta se estiman a partir de la muestra ($ p $ es el número de parámetros que deben estimarse) que se compara con el valor teórico de una distribución $ ^{2} $ según el nivel de significancia deseada.

Procedimiento.

1. Utilizar la muestra para estimar los parámetros de la distribución propuesta, si no se suponen de antemano.

2. Construir con los datos una distribución de frecuencias con las \(m\) clases supuestas.

3. Obtener la probabilidad \(p_{i}\) para cada valor o intervalos de valores de acuerdo con la distribución supuesta.

4. Obtener las frecuencias esperadas para cada valor o intervalo de valores: \(fe_{i}=np_{i}\).

5. Calcular el estadístico de prueba.

6. De acuerdo con el estadístico de prueba y el nivel de significancia utilizado, rechazar o no rechazar la hipótesis nula.

A continuación se muestran los resultados obtenidos, utilizando los datos presentados inicialmente

min(datos)

## [1] 1.5

max(datos)

## [1] 7.2

c<-round(((max(datos)-min(datos))/6),2);c

## [1] 0.95

a<-seq(min(datos),6*c+min(datos),c);a

## [1] 1.50 2.45 3.40 4.35 5.30 6.25 7.20

T5<-table(cut(datos,breaks=c(a),include.lowest=TRUE,right=FALSE))
de<-sd(datos)
n1<-round(pnorm((2.45-m)/de)*n,0)
n2<-round(pnorm((3.4-m)/de)*n,0)-round(pnorm((2.45-m)/de)*n,0)
n3<-round(pnorm((4.35-m)/de)*n,0)-round(pnorm((3.4-m)/de)*n,0)
n4<-round(pnorm((5.3-m)/de)*n,0)-round(pnorm((4.35-m)/de)*n,0)
n5<-round(pnorm((6.25-m)/de)*n,0)-round(pnorm((5.3-m)/de)*n,0)
n6<-round(pnorm((7.2-m)/de)*n,0)-round(pnorm((6.25-m)/de)*n,0)
fe<-c(n1,n2,n3,n4,n5,n6)
fo<-c(3,9,15,6,5,1)
sum(((fo-fe)^2)/fe)

## [1] 3.304545

cbind(T5,fe)

##            T5 fe
## [1.5,2.45)  3  4
## [2.45,3.4)  9  9
## [3.4,4.35) 15 11
## [4.35,5.3)  6 10
## [5.3,6.25)  5  5
## [6.25,7.2]  1  1

pearson.test(datos,n.classes = 6)

## 
##  Pearson chi-square normality test
## 
## data:  datos
## P = 3.2, p-value = 0.3618

4.2 Prueba de Kolmogorov-Smirnov

La prueba de Kolmogorov-Smirnov fue propuesta en 1933 y es una de las pruebas más utilizadas para bondad de ajuste y a diferencia de la prueba chi-cuadrado no necesita agrupar en clases los datos y es aplicable a muestras pequeñas. Se basa en la comparación entre las funciones de distribución acumulada que se observa en la muestra ordenada y la distribución supuesta bajo la hipótesis nula. Si la comparación revela una diferencia significativa entre las distribuciones muestral y supuesta, entonces la hipótesis nula debe ser rechazada, es decir la distribución no es la supuesta.

La prueba de Kolmogorov-Smirnov tiene la desventaja de requerir que se tenga completamente especificada la distribución supuesta.

\[\begin{equation*} F_{n}(x)= \begin{cases} 0 &\text{si $x < x_{(1)}$}\\ \dfrac{i}{n} &\text{si $x_{(i)} \leq x_{(n-1)}$}\\ 1 &\text{si $x \geq x_{(n)}$} \end{cases} \end{equation*}\]

El estadístico de prueba es:

\[\begin{equation*} D_{n}=max _{x} \vert F_{n}-F_{0} \vert \end{equation*}\]

Donde \(F_{0} (x)\) es la distribución acumulada de la distribución supuesta.\ La estadística \(D_{n}\) tiene una distribución que es independiente del modelo supuesto bajo la hipótesis nula por lo cual puede evaluarse sólo en función del tamaño de muestra. 6 La tabla muestra los cálculos necesarios para el cálculo de la prueba de Kolmogorov-Smirnov. puede observarse que en la línea 26 ocurre el máximo.

Los valores críticos para la prueba se encuentran tabulados, que para esta situación es 0.215, para un nivel de significancia de 0.05.

y<-i/n
x2<-round(pnorm(x,m,de),5)
k<-abs(y-x2)
cbind(i,y,x,x2,k)

##        i     y   x      x2       k
##  [1,]  1 0.025 1.5 0.02339 0.00161
##  [2,]  2 0.050 2.0 0.05585 0.00585
##  [3,]  3 0.075 2.2 0.07614 0.00114
##  [4,]  4 0.100 2.6 0.13276 0.03276
##  [5,]  5 0.125 2.7 0.15059 0.02559
##  [6,]  6 0.150 2.7 0.15059 0.00059
##  [7,]  7 0.175 2.8 0.16995 0.00505
##  [8,]  8 0.200 2.8 0.16995 0.03005
##  [9,]  9 0.225 2.8 0.16995 0.05505
## [10,] 10 0.250 2.9 0.19083 0.05917
## [11,] 11 0.275 3.3 0.28886 0.01386
## [12,] 12 0.300 3.3 0.28886 0.01114
## [13,] 13 0.325 3.4 0.31662 0.00838
## [14,] 14 0.350 3.5 0.34544 0.00456
## [15,] 15 0.375 3.6 0.37520 0.00020
## [16,] 16 0.400 3.6 0.37520 0.02480
## [17,] 17 0.425 3.7 0.40571 0.01929
## [18,] 18 0.450 3.8 0.43681 0.01319
## [19,] 19 0.475 3.9 0.46831 0.00669
## [20,] 20 0.500 3.9 0.46831 0.03169
## [21,] 21 0.525 3.9 0.46831 0.05669
## [22,] 22 0.550 4.0 0.50000 0.05000
## [23,] 23 0.575 4.1 0.53169 0.04331
## [24,] 24 0.600 4.1 0.53169 0.06831
## [25,] 25 0.625 4.1 0.53169 0.09331
## [26,] 26 0.650 4.2 0.56319 0.08681
## [27,] 27 0.675 4.2 0.56319 0.11181
## [28,] 28 0.700 4.4 0.62480 0.07520
## [29,] 29 0.725 4.4 0.62480 0.10020
## [30,] 30 0.750 4.5 0.65456 0.09544
## [31,] 31 0.775 4.7 0.71114 0.06386
## [32,] 32 0.800 5.0 0.78678 0.01322
## [33,] 33 0.825 5.2 0.83005 0.00505
## [34,] 34 0.850 5.4 0.86724 0.01724
## [35,] 35 0.875 5.4 0.86724 0.00776
## [36,] 36 0.900 5.5 0.88356 0.01644
## [37,] 37 0.925 5.9 0.93461 0.00961
## [38,] 38 0.950 6.2 0.95991 0.00991
## [39,] 39 0.975 6.6 0.98067 0.00567
## [40,] 40 1.000 7.2 0.99454 0.00546

max(k)

## [1] 0.11181

ks.test(datos,"pnorm",m,de,alternative="two.sided")

## Warning in ks.test(datos, "pnorm", m, de, alternative = "two.sided"): ties
## should not be present for the Kolmogorov-Smirnov test

## 
##  One-sample Kolmogorov-Smirnov test
## 
## data:  datos
## D = 0.11181, p-value = 0.6993
## alternative hypothesis: two-sided

4.3 Prueba de Liliefors

Esta prueba fue introducida por Lilliefors (1967) y es una modificación de la prueba de Kolmogorov-Smirnov para normalidad cuando cuando la media y la varianza de la distribución no son conocidas y deben ser estimadas a partir de los datos.

Utiliza la diferencia máxima entre entre la función de distribución empírica y la propuesta.

\(D=max(D^{+},D^{-})\)

\(D^{+}=max[\dfrac{i}{n}-F(x_{i}]\)

\(D^{-}= max[F(x_{i}-\dfrac{i}{n}]\)

\(F(x_{i})= \Phi\left(\dfrac{x_{i}-\bar{x}}{s}\right)\)

El estadístico de prueba es el mismo que el de Kolmogorov-Smirnov, pero los valores críticos son diferentes, que se encuentran tabulados. Para este ejemplo el valor crítico es 0.1389.

x2<-round(pnorm(x,m,de),5)
k<-abs(y-x2)
cbind(i,y,x,x2,k)

##        i     y   x      x2       k
##  [1,]  1 0.025 1.5 0.02339 0.00161
##  [2,]  2 0.050 2.0 0.05585 0.00585
##  [3,]  3 0.075 2.2 0.07614 0.00114
##  [4,]  4 0.100 2.6 0.13276 0.03276
##  [5,]  5 0.125 2.7 0.15059 0.02559
##  [6,]  6 0.150 2.7 0.15059 0.00059
##  [7,]  7 0.175 2.8 0.16995 0.00505
##  [8,]  8 0.200 2.8 0.16995 0.03005
##  [9,]  9 0.225 2.8 0.16995 0.05505
## [10,] 10 0.250 2.9 0.19083 0.05917
## [11,] 11 0.275 3.3 0.28886 0.01386
## [12,] 12 0.300 3.3 0.28886 0.01114
## [13,] 13 0.325 3.4 0.31662 0.00838
## [14,] 14 0.350 3.5 0.34544 0.00456
## [15,] 15 0.375 3.6 0.37520 0.00020
## [16,] 16 0.400 3.6 0.37520 0.02480
## [17,] 17 0.425 3.7 0.40571 0.01929
## [18,] 18 0.450 3.8 0.43681 0.01319
## [19,] 19 0.475 3.9 0.46831 0.00669
## [20,] 20 0.500 3.9 0.46831 0.03169
## [21,] 21 0.525 3.9 0.46831 0.05669
## [22,] 22 0.550 4.0 0.50000 0.05000
## [23,] 23 0.575 4.1 0.53169 0.04331
## [24,] 24 0.600 4.1 0.53169 0.06831
## [25,] 25 0.625 4.1 0.53169 0.09331
## [26,] 26 0.650 4.2 0.56319 0.08681
## [27,] 27 0.675 4.2 0.56319 0.11181
## [28,] 28 0.700 4.4 0.62480 0.07520
## [29,] 29 0.725 4.4 0.62480 0.10020
## [30,] 30 0.750 4.5 0.65456 0.09544
## [31,] 31 0.775 4.7 0.71114 0.06386
## [32,] 32 0.800 5.0 0.78678 0.01322
## [33,] 33 0.825 5.2 0.83005 0.00505
## [34,] 34 0.850 5.4 0.86724 0.01724
## [35,] 35 0.875 5.4 0.86724 0.00776
## [36,] 36 0.900 5.5 0.88356 0.01644
## [37,] 37 0.925 5.9 0.93461 0.00961
## [38,] 38 0.950 6.2 0.95991 0.00991
## [39,] 39 0.975 6.6 0.98067 0.00567
## [40,] 40 1.000 7.2 0.99454 0.00546

max(k)

## [1] 0.11181

lillie.test(datos)

## 
##  Lilliefors (Kolmogorov-Smirnov) normality test
## 
## data:  datos
## D = 0.11181, p-value = 0.2364

4.4 Prueba de Shapiro-Wilks

La prueba de Shapiro-Wilk es una prueba de bondad de ajuste que se utiliza para comprobar si unos datos determinados \((x_{1}, x_{2}, \cdots ,x_{n})\) han sido extraídos de una población normal. Los parámetros de la distribución no tienen porqué ser conocidos y es adecuada para muestras pequeñas$ (n<50)$.

Se ordenan de menor a mayor los datos de la muestra. Se obtienen los coeficientes \(a_{1}, a_{2}, \cdots a_{k}\), donde \(k\) es aproximadamente \(\frac{n}{2}\). Los \(a_{i}\) se encuentran tabulados.

Sean

\(SC=\displaystyle\sum_{i=1}^{n}(x_{i}-\bar{x})^{2}\)

$D= ^{2} $

El estadístico de prueba es:

\[\begin{equation*} SW=\dfrac{D}{SC} \end{equation*}\]

a<-c(0.3964,0.2737,0.2368,0.2098,0.1878,0.1691,
     0.1526,0.1376,0.1237,0.1108,0.0986,0.0870,
     0.0759,0.0651,0.0546,0.0444,0.0343,0.0244,
     0.0146,0.0049)
me<-head(x,20)
ma<-tail(x,20)
mi<-sort(ma,decreasing=TRUE)
sc<-sum((datos-m)^2)
dif<-mi-me
ca<-a*dif
cbind(a,me,mi,ca)

##            a  me  mi      ca
##  [1,] 0.3964 1.5 7.2 2.25948
##  [2,] 0.2737 2.0 6.6 1.25902
##  [3,] 0.2368 2.2 6.2 0.94720
##  [4,] 0.2098 2.6 5.9 0.69234
##  [5,] 0.1878 2.7 5.5 0.52584
##  [6,] 0.1691 2.7 5.4 0.45657
##  [7,] 0.1526 2.8 5.4 0.39676
##  [8,] 0.1376 2.8 5.2 0.33024
##  [9,] 0.1237 2.8 5.0 0.27214
## [10,] 0.1108 2.9 4.7 0.19944
## [11,] 0.0986 3.3 4.5 0.11832
## [12,] 0.0870 3.3 4.4 0.09570
## [13,] 0.0759 3.4 4.4 0.07590
## [14,] 0.0651 3.5 4.2 0.04557
## [15,] 0.0546 3.6 4.2 0.03276
## [16,] 0.0444 3.6 4.1 0.02220
## [17,] 0.0343 3.7 4.1 0.01372
## [18,] 0.0244 3.8 4.1 0.00732
## [19,] 0.0146 3.9 4.0 0.00146
## [20,] 0.0049 3.9 3.9 0.00000

sum(ca)

## [1] 7.75198

sw<-(sum(ca))^2/sc;sw

## [1] 0.9745896

shapiro.test(datos)

## 
##  Shapiro-Wilk normality test
## 
## data:  datos
## W = 0.97661, p-value = 0.5655

4.5 Prueba de Cramér-von Mises

Desarrollada por Harald Cramér y Richard Edler von Mises, quienes la expusieron entre los años 1928-1930. Compara la función de distribución empírica con la función de distribución teórica.

\[\begin{equation*} T=nw^{2}=\dfrac{1}{12n}+\displaystyle\sum_{i=1}^{n}\left[\dfrac{2i-1}{2n}-F(x_{i})\right]^{2} \end{equation*}\]

z<-pnorm((x-m)/de)
cvm<-1/(12*n)+sum(((2*i-1)/(2*n)-z)^2)
cvm

## [1] 0.06536419

cvm.test(datos)

## 
##  Cramer-von Mises normality test
## 
## data:  datos
## W = 0.065364, p-value = 0.3138

4.6 Prueba de Watson

Es una versión modificada de la prueba de Cramér-von Mises.

\[\begin{equation*} U^{2}= T-n(\bar{F}-\dfrac{1}{2})^{2} \end{equation*}\]

donde:

\(\bar{F}=\frac{1}{n}\displaystyle\sum_{i=1}^{n}F(x_{i})\)

mf<-1/n*sum(z)
U2<-cvm-n*(mf-1/2)^2
U2

## [1] 0.05831482

Los valores críticos se encuentran tabulados. El valor critico para la prueba es 0.117

4.7 Prueba de Shapiro-Francia

Es una modificación de la prueba de Shapiro-Wilk propuesta por S. S. Shapiro y R. S. Francia (1972).

Está dada como el cuadrado de la correlación entre los valores de la muestra ordenados (de menor a mayor) y los cuantiles esperados de la distribución normal estándar.

Los cuantiles ordenados esperados bajo la distribución normal se aproximan mediante qnorm(ppoints(x, a = 3/8)).

\[\begin{equation*} W_{SF}=\dfrac{\left[\displaystyle\sum_{i=1}^{n}m_{i}x_{i}\right]^{2}}{\displaystyle\sum_{i=1}^{n}(x_{i}-\bar{x})^{2}\displaystyle\sum_{i=1}^{n}m_{i^{2}}} \end{equation*}\]

El valor de esta prueba cae entre 0 y 1. Los valores cercanos a 1 indican la normalidad de los datos y los valores cercanos a 0 indican no normalidad en los datos.

po<-qnorm(ppoints(x, a = 3/8))
sf<-(cov(x,po)/(sd(x)*sd(po)))^2;sf

## [1] 0.9762366

ShapiroFranciaTest(datos)

## 
##  Shapiro-Francia normality test
## 
## data:  datos
## W = 0.97624, p-value = 0.4698

4.8 Prueba de Geary

Desarrollada en la década de 1930, esta prueba es muy sencilla en su aplicación y recientemente se ha estado utilizando bastante. Se basa en un estadístico muy sencillo que es un cociente entre dos estimadores de la desviación estándar poblacional \((\sigma)\).

\[\begin{equation*} d=\dfrac{1}{n\cdot s}\sum _{i=1}^{n}|x_{i}-\bar{x}| \end{equation*}\]

donde:

\[\begin{equation*} s^{2}=\dfrac{1}{n}\sum _{i=1}^{n}(x_{i}-\bar{x})^{2} \end{equation*}\]

El denominador es un estimador (no insesgado) de la desviación estándar poblacional , sin importar si la población es normal o no y el numerador (la desviación media) es un buen estimador de la desviación estándar si la distribución es normal, sin embargo podría sobreestimar o subestimar a la desviación estándar cuando se presenten grandes desviaciones de la normalidad, de tal manera que los valores de \(d\) que difieran considerablemente de 1 representan una probable desviación de la normalidad.

Si la hipótesis nula de normalidad es verdadera, el valor esperado de \(d\) es aproximadamente \(\sqrt{\dfrac{2}{\pi}}\approx 0.7979\). Entonces se rechaza la hipótesis nula si \(d\) es muy grande ( \(d > 0.7979)\) o muy pequeño (\(d < 0.7979)\)

d2<-sum((datos-m)^2)/n
d3<-sqrt(d2)
d<-1/(n*d3)*sum(abs(datos-m));d

## [1] 0.7691864

geary.norm.test(datos)

## 
##  Geary test for normality
## 
## data:  datos
## d = 0.76919, p-value = 0.861

4.9 Prueba de Frosini

El estadístico de prueba es:

\[\begin{equation*} B= \dfrac{1}{\sqrt{n}}\displaystyle\sum_{i=1}^{n}\vert\Phi(Y_{i})-\dfrac{i-0.5}{n}\vert \end{equation*}\]

donde

\(Y_{i}=\dfrac{x_{(i)-\bar{x}}}{s}\)

\(s=\sqrt{\dfrac{1}{n}\displaystyle\sum_{i=1}^{n}(x_{i}-\bar{x})^{2}}\)

El valor \(p\) se obtiene por simulación Montecarlo.

y1<-pnorm((x-m)/d3)
b1<-abs(y1-(i-0.5)/n)
B=1/sqrt(n)*sum(b1);B

## [1] 0.1857081

frosini.norm.test(datos)

## 
##  Frosini test for normality
## 
## data:  datos
## B = 0.18571, p-value = 0.4315

4.10 Prueba de Hegazy-Green

Se basa en el estadístico

\[\begin{equation*} T_{1}=\dfrac{1}{n}\displaystyle\sum_{i=1}^{n}\mid Y_{i}-\Phi^{-1}\left(\dfrac{i}{n+1}\right)\mid \end{equation*}\]

Donde

\(Y_{i}= \dfrac{x_{(i)}-\bar{x}}{s}\)

\(s^{2}=\dfrac{1}{n}\displaystyle\sum_{i=1}^{n}(x_{i}-\bar{x})^{2}\)

El valor se obtiene por simulación Montecarlo.

y2<-((x-m)/d3)
c1<-qnorm(i/(n+1))
T1<-1/n*(sum(abs(y2-c1)));T1

## [1] 0.123191

hegazy1.norm.test(datos)

## 
##  Hegazy-Green test for normality
## 
## data:  datos
## T = 0.12319, p-value = 0.4655

4.11 Prueba de Hegazy-Green 2

Es una modificación de la prueba original.

\[\begin{equation*} T_{1}=\dfrac{1}{n}\displaystyle\sum_{i=1}^{n}\left( Y_{i}-\Phi^{-1}\left(\dfrac{i}{n+1}\right)\right)^{2} \end{equation*}\]

Donde

\(Y_{i}= \dfrac{x_{(i)}-\bar{x}}{s}\)

\(s^{2}=\dfrac{1}{n}\displaystyle\sum_{i=1}^{n}(x_{i}-\bar{x})^{2}\)

El valor se obtiene por simulación Montecarlo.

y2<-((x-m)/d3)
c1<-qnorm(i/(n+1))
T2<-1/n*(sum((y2-c1)^2));T2

## [1] 0.0305858

hegazy2.norm.test(datos)

## 
##  Hegazy-Green test for normality
## 
## data:  datos
## T = 0.030586, p-value = 0.3785

4.12 Prueba de Ryan-Joiner

El gráfico cuantil - cuantil, como se vio anteriormente, es un método informal para evaluar la forma de la distribución de la población. Cuánto más cercano se encuentren los puntos a la recta, es más probable suponer la distribución en la cual se basa el gráfico.

Sea \(x_{(i)}\) las observaciones muestrales ordenadas. El gráfico para comprobar la normalidad sería el representado por las parejas \((x_{(i), y_{i} })\) donde \(y_{i}=\Phi^{ - 1}\left( \dfrac{i-0.5}{n+0.25}\right)\)

Una medida de la magnitud con la que los puntos se agrupan alrededor de la recta es el coeficiente de correlación muestral \(r\). En la medida en que los puntos se alejen de la recta el coeficiente de correlación tenderá a ser más pequeño. Esta idea se puede utilizar para realizar un procedimiento de prueba, tal como rechazar la hipótesis de normalidad si el coeficiente de correlación es menor que algún valor crítico, a un nivel de significancia especificado.

\[\begin{equation*} R_{J}=\dfrac{\sum(y_{i}-\bar{y})b_{i}}{\sqrt{s^{2}(n-1)\sum b_{i}^{2}}} \end{equation*}\]

Se define como el coeficiente de correlación entre los datos y las puntuaciones normales de los datos. Si el coeficiente de correlación está cercano a 1 los datos se encuentran cerca de la gráfica de probabilidad normal.

Los valores críticos se encuentran tabulados. Para este caso el valor crítico es 0.9715 con un nivel de significancia de 0.05.

Rj<-cor(x,do);Rj

## [1] 0.9880512

Otra forma

va<-var(datos)
faz<-sum(x*do)
ni<-sum(do^2)
Rj2<-faz/(sqrt((n-1)*va*ni));Rj2

## [1] 0.9880512

4.13 Prueba de Spiegelhalter

El estadístico de prueba está dado por

\[\begin{equation*} T=\left( (c_{n}u)^{-(n-1)}+g^{-n-1} \right)^{\frac{1}{n-1}} \end{equation*}\]

Donde

\(u=\dfrac{x_{(n)}-x_{(1)}}{s}\)

\(g=\dfrac{\sum|x_{i}-\bar{x}}{s\sqrt{n(n-1)}}\)

\(c_{n}=\dfrac{(n!)^{\frac{1}{n-1}}}{2n}\)

\(s^{2}=\dfrac{1}{n-1}\displaystyle\sum(x_{i}-\bar{x})^{2}\)

u<-max(datos)-min(datos);u

## [1] 5.7

g<-sum(abs(datos-m))/(de*sqrt(n*(n-1)));g

## [1] 0.7691864

cn<-(factorial(n))^(1/(n-1))/(2*n);cn

## [1] 0.2115503

T<-((cn*u)^(-(n-1))+g^(-(n-1)))^(1/(n-1));T

## [1] 1.300075

spiegelhalter.norm.test(datos)

## 
##  Spiegelhalter test for normality
## 
## data:  datos
## T = 1.3001, p-value = 0.194

4.14 Prueba de Weisberg - Bingham

El estadístico de prueba está dado por

\[\begin{equation*} WB=\dfrac{\dfrac{\left(\sum m_{i}x_{(i)} \right)^{2}}{\sum m_{i}^{2}}}{\sum(x_{i}-\bar{x})^{2}} \end{equation*}\]

\(m_{i}=\Phi^{-1}\left(\dfrac{i-\frac{3}{8}}{n+\frac{1}{4}}\right)\)

El valor se obtiene por simulación Montecarlo.

vi<-(sum(do*x)^2)/(sum(do^2))
vo<-sum((datos-m)^2)
wb<-vi/vo;wb

## [1] 0.9762452

wb.norm.test(datos)

## 
##  Weisberg-Bingham test for normality
## 
## data:  datos
## WB = 0.97624, p-value = 0.4695

4.15 Prueba de Jarque-Bera

Esta prueba se basa en los coeficientes de asimetría y kurtosis. bajo el supuesto de que la muestra proviene de una población normal, el estadístico JB tiene una distribución chi-cuadrado con dos grados de libertad.

\[\begin{equation*} JB=n \left[ \dfrac{g_{1}^{2}}{6}+\dfrac{(g_{2}-3)^{2}}{24}\right] \end{equation*}\]

jb<-n*(g1^2/6+(g2-3)^2/24);jb

## [1] 1.454602

jb.norm.test(datos)

## 
##  Jarque-Bera test for normality
## 
## data:  datos
## JB = 1.4546, p-value = 0.3415

4.16 Prueba de Jarque- Bera ajustada

\[\begin{equation*} AJB=\dfrac{g_{1}^{2}}{Var(g_{1})}+\dfrac{(g_{2}-E(g_{2}))^{2}}{Var(g_{2})} \end{equation*}\]

\(Var(g_{1})= \dfrac{6(n-2)}{(n+1)(n+3)}\)

\(E(g_{2})= \dfrac{3(n-1)}{n+1}\)

\(Var(g_{2})=\dfrac{24n(n-2)(n-3)}{(n+1)^{2})(n+3)(n+5)}\)

Vb1<-(6*(n-2))/((n+1)*(n+3))
Eb2<-((3*(n-1)))/((n+1))
Vb2<-((24*n*(n-2)*(n-3)))/((n+1)^2*(n+3)*(n+5))
ajb<-(g1^2/Vb1)+((g2-Eb2)^2)/Vb2;ajb

## [1] 1.750575

ajb.norm.test(datos)

## 
##  Adjusted Jarque-Bera test for normality
## 
## data:  datos
## AJB = 1.7506, p-value = 0.3065

4.17 Prueba de Anderson-Darling

Esta prueba data de 1954, pero en realidad se trata de una modificación de la prueba de Kolmogorov-Smirnov. Esta prueba presenta una gran sensibilidad hacia los extremos de la distribución, aunque no funciona bien en distribuciones de tipo discreto. \\(x_{i}\), \(i=1,2,\cdots,n\)

Procedimiento.

1. Utilizar la muestra para estimar los parámetros de la distribución.

2. Organizar los datos en forma ascendente:\(x_{i}\), \(i=1,2,\cdots,n\)

3. Organizar los datos en forma descendente:\(x_{n+1-i}\), \(i=1,2,\cdots,n\)

4. Calcular, según la distribución supuesta, la probabilidad para valor \(x_{i}\) y para cada valor \(x_{n+1-i}\)

5. Calcular el estadístico de prueba.

6. De acuerdo con el estadístico de prueba y el nivel de significancia utilizado, rechazar o no rechazar la hipótesis nula.

El estadístico de Anderson-Darling está dado por

\[\begin{equation*} A_{n}^{2}=-\left[n+\dfrac{1}{n}\displaystyle\sum_{i=1}^{n}(2i-1)[lnF(x_{(i)})+ln(1-Fx_{(n+1-i)})] \right] \end{equation*}\]

t2<-sort(datos,decreasing = TRUE)#3
h1<-2*i-1;h1#4

##  [1]  1  3  5  7  9 11 13 15 17 19 21 23 25 27 29 31 33 35 37 39 41 43 45 47 49
## [26] 51 53 55 57 59 61 63 65 67 69 71 73 75 77 79

z1<-round((x-m)/de,4)#5
z2<-round((t2-m)/de,4)#6
fi1<-round(pnorm(z1),5)#7
fi2<-round(1-pnorm(z2),5)#8
fi3<-round(log(fi1),5)#9
fi4<-round(log(fi2),5)#10
u1<-h1*(fi3+fi4)
cbind(i,x,t2,h1,z1,z2,fi1,fi2,fi3,fi4,u1)

##        i   x  t2 h1      z1      z2     fi1     fi2      fi3      fi4        u1
##  [1,]  1 1.5 7.2  1 -1.9882  2.5450 0.02339 0.00546 -3.75545 -5.21031  -8.96576
##  [2,]  2 2.0 6.6  3 -1.5906  2.0678 0.05585 0.01933 -2.88509 -3.94610 -20.49357
##  [3,]  3 2.2 6.2  5 -1.4315  1.7497 0.07614 0.04009 -2.57518 -3.21663 -28.95905
##  [4,]  4 2.6 5.9  7 -1.1134  1.5111 0.13277 0.06538 -2.01914 -2.72754 -33.22676
##  [5,]  5 2.7 5.5  9 -1.0339  1.1929 0.15059 0.11645 -1.89319 -2.15029 -36.39132
##  [6,]  6 2.7 5.4 11 -1.0339  1.1134 0.15059 0.13277 -1.89319 -2.01914 -43.03563
##  [7,]  7 2.8 5.4 13 -0.9544  1.1134 0.16994 0.13277 -1.77231 -2.01914 -49.28885
##  [8,]  8 2.8 5.2 15 -0.9544  0.9544 0.16994 0.16994 -1.77231 -1.77231 -53.16930
##  [9,]  9 2.8 5.0 17 -0.9544  0.7953 0.16994 0.21322 -1.77231 -1.54543 -56.40158
## [10,] 10 2.9 4.7 19 -0.8748  0.5567 0.19084 0.28887 -1.65632 -1.24178 -55.06390
## [11,] 11 3.3 4.5 21 -0.5567  0.3976 0.28887 0.34546 -1.24178 -1.06288 -48.39786
## [12,] 12 3.3 4.4 23 -0.5567  0.3181 0.28887 0.37520 -1.24178 -0.98030 -51.10784
## [13,] 13 3.4 4.4 25 -0.4772  0.3181 0.31661 0.37520 -1.15008 -0.98030 -53.25950
## [14,] 14 3.5 4.2 27 -0.3976  0.1591 0.34546 0.43680 -1.06288 -0.82828 -51.06132
## [15,] 15 3.6 4.2 29 -0.3181  0.1591 0.37520 0.43680 -0.98030 -0.82828 -52.44882
## [16,] 16 3.6 4.1 31 -0.3181  0.0795 0.37520 0.46832 -0.98030 -0.75860 -53.90590
## [17,] 17 3.7 4.1 33 -0.2386  0.0795 0.40571 0.46832 -0.90212 -0.75860 -54.80376
## [18,] 18 3.8 4.1 35 -0.1591  0.0795 0.43680 0.46832 -0.82828 -0.75860 -55.54080
## [19,] 19 3.9 4.0 37 -0.0795  0.0000 0.46832 0.50000 -0.75860 -0.69315 -53.71475
## [20,] 20 3.9 3.9 39 -0.0795 -0.0795 0.46832 0.53168 -0.75860 -0.63171 -54.22209
## [21,] 21 3.9 3.9 41 -0.0795 -0.0795 0.46832 0.53168 -0.75860 -0.63171 -57.00271
## [22,] 22 4.0 3.9 43  0.0000 -0.0795 0.50000 0.53168 -0.69315 -0.63171 -56.96898
## [23,] 23 4.1 3.8 45  0.0795 -0.1591 0.53168 0.56320 -0.63171 -0.57412 -54.26235
## [24,] 24 4.1 3.7 47  0.0795 -0.2386 0.53168 0.59429 -0.63171 -0.52039 -54.14870
## [25,] 25 4.1 3.6 49  0.0795 -0.3181 0.53168 0.62480 -0.63171 -0.47032 -53.99947
## [26,] 26 4.2 3.6 51  0.1591 -0.3181 0.56320 0.62480 -0.57412 -0.47032 -53.26644
## [27,] 27 4.2 3.5 53  0.1591 -0.3976 0.56320 0.65454 -0.57412 -0.42382 -52.89082
## [28,] 28 4.4 3.4 55  0.3181 -0.4772 0.62480 0.68339 -0.47032 -0.38069 -46.80555
## [29,] 29 4.4 3.3 57  0.3181 -0.5567 0.62480 0.71113 -0.47032 -0.34090 -46.23954
## [30,] 30 4.5 3.3 59  0.3976 -0.5567 0.65454 0.71113 -0.42382 -0.34090 -45.11848
## [31,] 31 4.7 2.9 61  0.5567 -0.8748 0.71113 0.80916 -0.34090 -0.21176 -33.71226
## [32,] 32 5.0 2.8 63  0.7953 -0.9544 0.78678 0.83006 -0.23981 -0.18626 -26.84241
## [33,] 33 5.2 2.8 65  0.9544 -0.9544 0.83006 0.83006 -0.18626 -0.18626 -24.21380
## [34,] 34 5.4 2.8 67  1.1134 -0.9544 0.86723 0.83006 -0.14245 -0.18626 -22.02357
## [35,] 35 5.4 2.7 69  1.1134 -1.0339 0.86723 0.84941 -0.14245 -0.16321 -21.09054
## [36,] 36 5.5 2.7 71  1.1929 -1.0339 0.88355 0.84941 -0.12381 -0.16321 -20.37842
## [37,] 37 5.9 2.6 73  1.5111 -1.1134 0.93462 0.86723 -0.06762 -0.14245 -15.33511
## [38,] 38 6.2 2.2 75  1.7497 -1.4315 0.95991 0.92386 -0.04092 -0.07919  -9.00825
## [39,] 39 6.6 2.0 77  2.0678 -1.5906 0.98067 0.94415 -0.01952 -0.05747  -5.92823
## [40,] 40 7.2 1.5 79  2.5450 -1.9882 0.99454 0.97661 -0.00547 -0.02367  -2.30206

An<--(sum(u1)/n+n);An

## [1] 0.3749013

ad.test(datos)

## 
##  Anderson-Darling normality test
## 
## data:  datos
## A = 0.37489, p-value = 0.3983

Los valores críticos para la prueba A-D se encuentran tabulados. Para un nivel de significancia de 0.05 el valor crítico es 0.751.

4.18 Prueba de Bonett-Seier

Esta prueba es una modificación de la prueba de Geary.

\[\begin{equation*} Z_{w}=\dfrac{\sqrt{n+2}(\hat{w}-3)}{3.54} \end{equation*}\]

Donde

\(\hat{w}=13.29\cdot(ln\sqrt{m_{2}}-ln(\frac{1}{n}\sum|x_{i}-\bar{x}|))\)

\(Z_{w}\) sigue aproximadamente una distribución normal estándar.

s<-1/n*sum((datos-m)^2);s

## [1] 1.5415

dam<-1/n*sum(abs(datos-m)) 
ww<-13.29*(log(sqrt(s))-log(dam));ww

## [1] 3.487587

TW<-sqrt(n+2)*(ww-3)/3.54;TW

## [1] 0.8926348

bonett.test(datos,alternative="two.sided")

## 
##  Bonett-Seier test for Geary kurtosis
## 
## data:  datos
## tau = 0.95500, z = 0.89263, p-value = 0.3721
## alternative hypothesis: kurtosis is not equal to sqrt(2/pi)

4.19 Prueba de kurtosis

Se basa en el estadístico

\[\begin{equation*} g_{2}=\dfrac{m_{4}}{m_{2}^{2}} =\dfrac{\frac{1}{n}\sum (x_{i}-\bar{x})^{4}}{(\frac{1}{n}\sum (x_{i}-\bar{x})^{2})^{2}} \end{equation*}\]

Que tiene asintóticamente una distribución normal con media 3 y varianza \(24/n\).

pl<-sqrt(24/n)
z1<-(g2-3)/pl;z1

## [1] 0.02139409

2*(1-pnorm(z1))

## [1] 0.9829313

kurtosis.norm.test(datos)

## 
##  Kurtosis test for normality
## 
## data:  datos
## T = 3.0166, p-value = 0.9805

4.20 Prueba de asimetría

Se basa en el estadístico

\[\begin{equation*} g_{1}=\dfrac{m_{3}}{m_{2}^{\frac{3}{2}}}=\dfrac{\frac{1}{n}\sum (x_{i}-\bar{x})^{3}}{(\frac{1}{n}\sum (x_{i}-\bar{x})^{2})^{\frac{3}{2}}} \end{equation*}\]

Que tiene asintóticamente una distribución normal con media 0 y varianza \(6/n\).

pm<-sqrt(6/n)
z20<-(g1-0)/pm;z20

## [1] 1.205879

2*(1-pnorm(z20))

## [1] 0.2278642

skewness.norm.test(datos)

## 
##  Skewness test for normality
## 
## data:  datos
## T = 0.46703, p-value = 0.195

4.21 Prueba de D’Agostino

Basada en una transformación sobre la definición del coeficiente de asimetría.

\[\begin{equation*} Zs=b\cdot ln(u+\sqrt{u^{2}+1}) \end{equation*}\]

Que tiene una distribución normal estándar.

Donde

\(c=\dfrac{3(n^{2}+27n-70)(n+1)(n+3)}{(n-2)(n+5)(n+7)(n+9)}\)

\(w^{2}=\sqrt{2(c-1)}-1\)

\(a=\sqrt{\dfrac{W^{2}-1}{2}}\)

\(b=\dfrac{1}{\sqrt{ln(W)}}\)

\(u=a\cdot g_{1}\cdot\sqrt{\dfrac{(n+1)(n+3)}{6(n-2)}}\)

c<-(3*(n^2+27*n-70)*(n+1)*(n+3))/((n-2)*(n+5)*(n+7)*(n+9))
W2<--1+sqrt(2*(c-1))
W1<-sqrt(W2)
a<-sqrt((W2-1)/2)
b<-1/sqrt(log(W1))  
u<-a*g1*sqrt(((n+1)*(n+3))/((6*(n-2))))
Zs<-b*log(u+sqrt(u^2+1));Zs

## [1] 1.329087

agostino.test(datos,alternative="two.sided")

## 
##  D'Agostino skewness test
## 
## data:  datos
## skew = 0.46703, z = 1.32909, p-value = 0.1838
## alternative hypothesis: data have a skewness

4.22 Prueba de kurtosis (otra versión)

Se basa en una transformación del coeficiente de kurtosis.

\[\begin{equation*} Z_{k}=\dfrac{1-r-v^{\frac{1}{3}}}{\sqrt{r}} \end{equation*}\]

Que tiene una distribución normal estándar.

Donde

\(d=\sqrt{\dfrac{(n+1)^{2}(n+3)(n+5)}{24n(n-2)(n-3)}}\)

\(e=\dfrac{6(n^{2}-5n+2)}{(n+7)(n+9)}\sqrt{\dfrac{6(n+3)(n+5)}{n(n-2)(n-3)}}\)

\(f=6+\dfrac{8}{e}\left(\dfrac{2}{e} +\sqrt{1+\dfrac{4}{e^{2}}}\right)\)

\(g=d\left( g_{2}-\dfrac{3(n-1)}{n+1}\right) \sqrt{\dfrac{2}{f-4}}\)

\(v=\dfrac{1-\dfrac{2}{f}}{1+g}\)

\(r=\dfrac{2}{9f}\)

d<-sqrt(((n+1)^2*(n+3)*(n+5))/((24*n*(n-2)*(n-3))))
e<-(6*(n^2-5*n+2))/((n+7)*(n+9))*sqrt((6*(n+3)*(n+5))/(n*(n-2)*(n-3)))
f<-6+8/e*(2/e+sqrt(1+4/e^2))
g<-d*(g2-3*(n-1)/(n+1))*sqrt(2/(f-4))
v<-(1-2/f)/(1+g)
r<-2/(9*f)
Zk<-(1-r-v^(1/3))/sqrt(r);Zk

## [1] 0.4842748

4.23 Prueba omnibus de D’Agostino-Pearson

Se basa en el hecho de que los estadísticos \(Z_{s}\) y \(Z_{k}\) tienen distribución normal estándar, entonces el estadístico

\[\begin{equation*} DP=Z_{s}^{2}+Z_{k}^{2} \end{equation*}\]

tendrá una distribución \(\chi^{2}\) con 2 grados de libertad.

DP<-Zs^2+Zk^2; DP

## [1] 2.000995