contraste Jarque Bera » G. R. Serrano

GEBR 7.2: Simulación de la potencia del contraste Jarque-Bera

Docencia, Econometría Grados No Responses »

may 302012

Otro sencillo ejemplo, pero relacionado con el anterior, nos permitirá estudiar la potencia del contraste de Jarque-Bera para distintos tamaños muestrales. La potencia es la probabilidad de rechazar la hipótesis nula cuando es falsa, así que vamos a generar 1000 muestras tamaño 100 de una variable aleatoria U(0,1) y calcularemos si la nula de normalidad se rechaza para cada una de ellas; es decir, vamos a contar el número de veces que el estadístico JB es mayor que 5.99 (el valor crítico de una χ²₂ para un nivel de significación del α = 0.05).

> set.seed(301)  ## Para hacer el ejemplo reproducible
> nrep  <- 1000
> n     <- 100
> Z     <- matrix(runif(n*nrep, 0, 1), nrow=n)
> jbrechaza <- function(x) {
+     return(jarque.bera.test(x)$statistic > qchisq(0.95, 2))
+ }
> rechaza  <- apply(Z, 2, jbrechaza)
> (potencia <- sum(rechaza)/length(rechaza))

[1] 0.549

Incluso con muestras relativamente grandes para lo que se maneja a menudo en econometría, el estadístico de Jarque-Bera sólo es capaz de rechazar la nula de normalidad en una de cada dos ocasiones en que la muestra es uniforme. Como cabe esperar, con muestras menores la potencia se reduce apreciablemente, acercándose a cero con n=20. Puede ir cambiando la variable n para trazar la curva de potencia del contraste para una alternativa U(0,1).

GEBR 7.1: Simulación del teorema central del límite

Docencia, Econometría Grados No Responses »

may 292012

En R resulta especialmente cómodo programar simulaciones para estudiar propiedades estadísticas de contrastes y estimadores. En este post y los dos siguientes se presentan algunos ejemplos.

Vamos a ilustrar el teorema central del límite, que en una versión sencilla podemos formular del siguiente modo. Sean X_i variables aleatorias independientes y no muy diferentes:

$W_n = \sum _{i=1}^n X_i \qquad \Rightarrow \qquad \frac{W_n - \mathrm{E}(W_n)}{\sqrt{\mathrm{V}(W_n)}} \underset{n \to \infty}{\longrightarrow} N(0,1).$

Para una aplicación empírica, se generan 1000 muestras, tamaño 2000, de una población log-normal (supuesto muy frecuente en econometría, donde se toman logaritmos para inducir normalidad), que son las X_i, y las almacenamos en una matriz 2000×1000. Si sumamos los primeros 10 elementos de cada columna, obtenemos una muestra tamaño 1000 de W₁₀, que representamos en un histograma, no se parece a una normal (Gráfico 1(a)).

(a) Histograma de W10

(b) Histograma de W500

Gráfico 1: Convergencia de la suma de log-normales a la normal.

Después tomamos los primeros 500 elementos de cada muestra de log normales y los sumamos, tendríamos una muestra de W₅₀₀. El histograma de W₅₀₀ empieza a mostrar un comportamiento más simétrico (Gráfico 1(b)). Finalmente, cogemos los 2000 elementos de cada muestra y los sumamos. Ahora disponemos de una muestra tamaño 1000 de la variable aleatoria W₂₀₀₀ (Gráfico 1(c)). El script de R para hacer esta simulación es:

> set.seed(301) ## Para hacer el ejemplo reproducible
> n    <- 2000
> nrep <- 1000
> X    <- matrix(exp(rnorm(n*nrep)), nrow=n)
> W10   <- apply(X[1:10, ], 2, sum)
> W500  <- apply(X[1:500, ], 2, sum)
> W2000 <- apply(X[1:2000, ], 2, sum)
> sapply(data.frame(W10, W500, W2000), jarque.bera.test)[c(1,3),]

          W10      W500         W2000
statistic 7842.408 32.48236     1.754108
p.value   0        8.841903e-08 0.4160067

La función sapply() aplica la función que aparece como segundo argumento a las columnas del data.frame en su primer argumento. De todos los resultados que devuelve la función jarque.bera.test(), en este punto sólo nos interesan dos de ellos: el valor del estadístico y el p-valor, que aparecen en la primera y tercera posición.

Utilizando el contraste de Jarque y Bera, se rechaza la nula para W₁₀ y W₅₀₀, pero no se puede rechazar normalidad para W₂₀₀₀. Para comparar la tasa de convergencia, puede cambiar la generación de una log-normal (exp(rnorm(n*nrep))) por una uniforme: X <- matrix(runif(n*nrep), nrow=n), para comprobar que a partir de 5 ó 10 uniformes, la suma es aproximadamente normal.