library(kableExtra) # zhu2021
library(scales)El material que aquí se incluye está elaborado a partir de apuntes de clase tomados en la materia “Simulación”, impartida en el ITAM por el Dr. Juan José Fernández Durán, de enero a mayo del 2001, como parte de mis estudios en la licenciatura en actuaría. Estos apuntes han sido, desde luego, posteriormente modificados por lecturas y revisión de material diverso.
El material de estos apuntes está elaborado utilizando R.
Una gran variedad de situaciones de la vida real requieren que contemos con un mecanismo que nos permita considerar o introducir la aleatoriedad de un evento. En el más estricto sentido es prácticamente imposible que un proceso humano genere aleatoriedad real, por eso hablamos de eventos (números en nuestro caso) pseudo-aleatorios. En estas notas haré uso de los términos aleatorio y pseudo-aleatorio de manera intercambiable.
La habilidad de poder generar números pseudo-aleatorios que cumplan con ciertas propiedades es, por lo tanto, muy relevante para una gran variedad de aplicaciones hoy en día. En el pasado se diseñaron medios físicos para la generación de estos números pseudo-aleatorios, sin embargo hoy en día recurrimos principalmente a métodos computacionales (algoritmos) que presentan ventajas significativas sobre los métodos físicos (Rubinstein and Kroese (2016)):
En la mayoría de los algoritmos para generar números pseudo-aleatorios un punto de partida fundamental es la generación de números pseudo-aleatorios uniformes, por lo que un primer objetivo fundamental consiste en construir procedimientos que permitan generar una serie de números aleatorios que se aproximen a una secuencia de números aleatorios uniformes \([0,1]\).
Uno de los métodos más sencillos para la generación de números pseudo-aleatorios es mediante Generadores Lineales Congruenciales (GLC). El procedimiento es el siguiente:
Selecciona \(X_0 = x_0\)
Para \(t \in \{ 1, \dots , N\}\)
2.1. \(X_t = (aX_{t-1} + c) \mod m\)
2.2. \(U_t = \frac{X_t}{m}\)
2.3. \(t = t+1\)
Fin.
A \(x_0\) se le conoce como la semilla del proceso. Es importante notar que, si se conoce al semilla es posible generar toda la cadena de valores pseudo-aleatorios (lo cual puede ser tanto deseable como no). Todos los lenguajes de programación hoy en día cuentan con alguna manera de establecer la semilla para un proceso, lo cual permite garantizar su reproducibilidad. En R, por ejemplo, se utiliza la función seed().
Las propiedades de la serie de números aleatorios generados por el proceso arriba descrito dependerán en gran medida de los valores de \(a\) (multiplicador), \(c\) (incremento) y \(m\) (módulo), mismos que son fijados por el analista.
temp <- data.frame(t = 0:10, X = NA, U = NA)
temp$X[1] <- 1
a <- 2
c <- 3
m <- 5
for (i in 2:11){
temp$X[i] <- (a*temp$X[i-1] + c) %% m
temp$U[i] <- temp$X[i]/m
}
kable(x = temp, digits = 4, align = c('c', 'c', 'c')) |> kable_classic()| t | X | U |
|---|---|---|
| 0 | 1 | NA |
| 1 | 0 | 0.0 |
| 2 | 3 | 0.6 |
| 3 | 4 | 0.8 |
| 4 | 1 | 0.2 |
| 5 | 0 | 0.0 |
| 6 | 3 | 0.6 |
| 7 | 4 | 0.8 |
| 8 | 1 | 0.2 |
| 9 | 0 | 0.0 |
| 10 | 3 | 0.6 |
Observemos algunas propiedades relevantes de la serie \(U_t\) generada por este procedimiento:
\(U_t \in [0, 1]\).
la secuencia \(X_0, X_1, X_2, \dots\) se repetirá, a lo más, después de \(m\) pasos, por lo que será periódica. Por lo tanto, una adecuada selección de los parámetros será fundamental para que la serie generada pase (o no) las pruebas de aleatoriedad.
El método de generadores lineales congruenciales, dadas sus características, no es muy empleado hoy en día.
Un generador recursivo múltiple (MRG, por sus siglas en inglés) de orden k es un generador que satisface la siguiente fórmula recursiva:
\[X_t = (a_1 X_{t-1} + \dots + a_k X_{t-k}) \mod m\]
dada una semilla \(X_0 = (X_{-k+1}, \dots, X_0)\). El pseudocódigo para la generación de números aleatorios utilizando este tipo de generadores es exactamente el mismo que para el GLC.
temp <- data.frame(t = 0:20, X = NA, U = NA)
temp$X[1:3] <- 1:3
a <- c(2:4)
c <- 3
m <- 5
for (i in 4:21){
temp$X[i] <- (a[1]*temp$X[i-3] + a[2]*temp$X[i-2] + a[3]*temp$X[i-1] + c) %% m
temp$U[i] <- temp$X[i]/m
}
kable(x = temp, digits = 4, align = c('c', 'c', 'c')) |> kable_classic()| t | X | U |
|---|---|---|
| 0 | 1 | NA |
| 1 | 2 | NA |
| 2 | 3 | NA |
| 3 | 3 | 0.6 |
| 4 | 3 | 0.6 |
| 5 | 0 | 0.0 |
| 6 | 3 | 0.6 |
| 7 | 1 | 0.2 |
| 8 | 1 | 0.2 |
| 9 | 1 | 0.2 |
| 10 | 2 | 0.4 |
| 11 | 1 | 0.2 |
| 12 | 0 | 0.0 |
| 13 | 0 | 0.0 |
| 14 | 0 | 0.0 |
| 15 | 3 | 0.6 |
| 16 | 0 | 0.0 |
| 17 | 2 | 0.4 |
| 18 | 2 | 0.4 |
| 19 | 2 | 0.4 |
| 20 | 1 | 0.2 |
Que sea de ciclo completo. Se dice que un generador es de ciclo completo si genera todos los números en el conjunto \(\{0,1,\dots,m-1\}\) antes de volverse a repetir la secuencia.
Los LCG puros (c = 0) alcanzan ciclo completo únicamente cuando \(m\) es primo y \(a^{m-1}-1\) es divisible entre \(m\). Si \(m\) no es primo, a lo más tendrán ciclo de longitud \(2^{n-2}\), donde \(n\) es el tamaño de la sucesión.
Los MRG alcanzan ciclo completo cuando: 1) c y m son primos relativos; 2) \(a-1\) es divisible por todos los elementos que dividen a \(m\); 3) \(a-1\) es divisible entre 4 si \(m\) es divisible entre \(4\).
\(m\) grande. Si \(m\) es grande se garantiza una densidad grande en el intervalo \((0,1)\), algunos valores comúnmente utilizados, por ejemplo, son \(2^{31}\) o \(2^{31}-1\).
Que parezca una sucesión de observaciones independientes de una v.a. \(U(0,1)\).
(PENDIENTE. Rubinstein and Kroese 2016 pág. 52.)
Normalmente deseamos contar con indicadores que nos permitan determinar qué tan bueno es nuestro GNPA. El objetivo de estos indicadores/pruebas es verificar que nuestro GNPA arroja resultados que, en efecto, se comportan como una variable aleatoria. Es muy importante recordar, sin embargo, que (casi) todos los GNPA son procedimientos determinísticos por lo que generalmente habrá pruebas que no pasarán. Dependiendo de la finalidad del GNPA y de las pruebas que pase (o no pase) será el diagnóstico de si es un buen GNPA o no.
Una prueba muy recurrida en otros contextos, busca comparar la frecuencia de aparición de ciertos valores de la variable aleatoria contra la frecuencia esperada de dichos valores suponiendo que la variable aleatoria sigue una determinada distribución.
La prueba \(\chi^2\) lleva el apellido de Carl Pearson, matemático británico que propuso la prueba en 1900. Hoy en día Pearson forma parte de la lista de investigadores en temas de probabilidad y estadística que han caído en desgracia por sus puntos de vista en temas de igualdad racial.
La prueba consiste en el siguiente procedimiento. Sea \(\underline{X}_n = (x_1, x_2, \dots, x_n)\) \(n\) observaciones de una muestra que se desea verificar si es consistente con una muestra proveniente de una distribución \(\mathcal{F}\), entonces:
La prueba \(\chi^2\) es una prueba de una cola y rechaza la hipótesis nula cuando el valor de la estadística de prueba es mayor al valor crítico determinado para un valor de significancia \(\alpha\).
x.n <- c(8,3,6,7,2,0,4,5,6,9,8,3,2,0,4,5,5,6,9,6,7,1,0,3,2,0,0,5,9,1,8,3,9,2,7,6,3,2,0,1,1,4,7,9,8,3,2,0,1,2,4,3,5,6,9,0,1,0,3,9,4,6,9,3,2,8,7,7,5,3,2,4,0,1,3,3,5,9,8,6,3,4,2,0,1,1,3,8,7,5,9,9,0,8,8,7,3,1,2,0)¿Existe evidencia de que los datos provengan de una distribución discreta \(U(0,9)\)? Utiliza un nivel de significancia del 5%.
conteo <- table(x.n)
print(conteo)## x.n
## 0 1 2 3 4 5 6 7 8 9
## 13 10 11 15 7 8 8 8 9 11est.prueba <- sum(((conteo - 10)^2)/10)
print(est.prueba)## [1] 5.8vc <- qchisq(p = 0.95, df = 10 - 1)Dado que el valor de la estadística de prueba 5.80 es menor que el valor crítico de una \(\chi^2_{9,0.95} = 16.92\) concluimos que no existe evidencia para rechazar la hipótesis de que los datos provenga de una distribución \(U(0,9)\).
Recomendaciones al realizar una prueba \(\chi^2\):
Utilizar particiones tales que todos los intervalos contengan al menos un mínimo de 5 observaciones esperadas bajo la distribución teórica.
Verificar que los intervalos que contienen a las colas de la distribución acumulen la mayoría de los valores. De ser necesario, para estos casos es posible no cumplir con la recomendación anterior.
Un número guía de intervalos es \(4 \times n^{2/5}\).
\[ H_0: \text{Dist. Uniforme en} \ [0,1]^k \]
Los generadores congruenciales para \(k=1\) no tienen problemas si es un buen generador.
Definir celdas en \([0,1]^k\).
Obtener \(D_i\) y \(l_i\) para cada celda.
Calcular la estadística \(\chi^2\).
Dado que el número de celdas crece exponencialmente en \(k\), en la práctica se limita generalmente hasta \(k=3\) por las exigencias computacionales de valores más grandes de k.
\[ H_0: vv.aa.ii. \]
Generar \(R\) sucesiones de tamaño \(r\).
Para cada sucesión obtener la permutación de índices.
Bajo la hipótesis nula, cada permutación de índices tiene probabilidad \(\frac{1}{r!}\).
Aplicar la prueba \(\chi^2\) con \(e_i = \frac{R}{r!}\).
Cuando se rechaza la hipótesis es porque existe un patrón de orden entre las sucesiones, lo que las hace dependientes. # TÉCNICAS DE REDUCCIÓN DE VARIANZA