viernes, 20 de abril de 2012

Pruebas estadísticas para los números pseudoaleatorios

En esta entrada hablaremos sobre los números pseudoaleatorios, generaremos un número pseudo-aleatorio por medio de un programa hecho en Ruby y de ahí nos pasaremos a hacer unas pruebas estadísticas sobre dichos datos.

Así es como se generan números aleatorios en Ruby:

archivo = File.new("rand.dat","w+")
for s in 1...100
      print s,"\n";
      archivo.print s," ",rand(),"\n";
end

Este genera un archivo llamado rand.dat en el cual se van guardando los datos, para que posteriormente sean gráficados.


Donde genera 100 números de manera pseudoaleatoria en un rango entre 0 y 1.

Para hacer dicha modificaré el código de Ruby para generar números que sean 0 ó 1:

archivo = File.new("rand.dat","w+")
for s in 1...100
      print s,"\n";
      archivo.print s," ",rand(2),"\n";
end

Ahora haremos una prueba llamada "Random Excursions Test" para comprobar que realmente son pseudoaleatorios, donde el propósito de está prueba es concentrarse en el número de ciclos que tienen exactamente K visitas en un paseo al azar la suma acumulada, para más información entra a está liga:


Donde se programó esta prueba en Python:

import math
import random
def tabla(x, k):
    val = 1.0 / (2.0 * math.fabs(x))
    return val * (1.0 - val)**(k - 1)

def lookup(lista):
    ciclo = list()
    c = None
    for e in lista:
        if e == 0:
            if c is not None:
                ciclo.append(c)
            c = list()
        else:
            c.append(e)
    return ciclo

file = open("rand.dat")
lista = list()
line = file.readline()
line = line.strip()
for elem in line.split():
    lista.append(int(elem)*2 -1)
n=len(lista)
s=list()

print lista
for e in range(0,n+1):
    s.append(sum(lista[:e]))

s.append(0)
print s
ciclo=lookup(s)
J = len(ciclo)

v=dict()
for x in range(-4, 5):
    if x != 0:
        for c in ciclo:
            oc = c.count(x)
            if oc > 5:
                oc = 5 # los mayores a cinco van con cinco
            if (x, oc) in v:
                v[(x, oc)] += 1
            else:
                v[(x, oc)] = 1
for x in range(-4, 5):
    if x != 0:
        ver = 0
        stat = 0.0
        for k in range(0, 6):
            if (x, k) in v:
                val = v[(x, k)] # frecuencia observada
            else:
                val = 0.0
            teor = J * tabla(x, k) # frecuencia esperada
            stat += (val - teor)**2 / teor
            ver += val
        print x, stat 
Donde se normalizan los números pseudoaleatorios generados entre 0 y 1 donde 0 se cambia a -1, después de ello se hacen sumatorias parciales comenzando de X1 así sucesivamente hasta Xn, para luego examinar la cantidad de repeticiónes de frecuencias y compararlos apartir de la estadística Chi-square, lo que se esperaría una probabilidad de una distribución aleatoria.

Como resultado obtenemos lo siguiente:


Después comparamos con esta applet si realmente tienen una alta probabilidad de distribución aleatoria:


Con esto concluimos que a mayor cantidad sean los valores, hay mayor probabilidad de que se rechazen los valores y que a menor cantidad sean hay mayor probabilidad de que se acepten.

Se comprobó que los datos generados son uniformes y que puede mejorarse la aleatoridad. Se llegó a esta conclusión gracias a que dividimos los -1 y 1 entre los valores esperados para así mismo comprobarlo con el applet y esta misma applet nos determinó en todas las pruebas que se hicieron que nos daba una probabilidad de 0,9999.



Link de referencia:

Información de Random Excursions Test
http://es.scribd.com/doc/11305936/133/Random-Excursions-Test
Applet
http://www.danielsoper.com/statcalc3/calc.aspx?id=11