Max Hernandez: Tarea 3: Pruebas estadísticas para los números pseudoaleatorios

Las aplicaciones de los números aleatorios en la actualidad se ha vuelto tan diversa y tan importante en la actualidad que se necesita de cierto grado de calidad para que sea utilizable para las personas. Por eso esta semana quiero hablar sobre los generadores de números aleatorios y algunas pruebas para tener una mayor seguridad de que los generadores que utilizamos tienen la calidad ya mencionada.

Las aplicaciones de los números aleatorios son :

Criptografía.
Simulaciones.
Videojuegos.
Muestras en experimentos.
Obras de arte.

Generadores de numeros aleatorios

Los generadores de números aleatorios es un algoritmo determinista que producen una variable al azar especificada por una distribución..

Pruebas de aleatoriedad

PPara las pruebas de aleatoriedad elegí tres pruebas como ejemplo del conjunto de pruebas de NSIT los cuales son:

Frequency Test - Monobit: El propósito de esta prueba es determinar si la proporción de ceros y unos obtenidos.
Frequency Test - Block: El propósito de esta prueba es determinar la proporción de unos en bloques de M-bits (M=cantidad de bits por bloque)
Runs Test: El propósito de este es medir la cantidad de secuencias de bits idénticos interrumpidos.

import math
import random
from scipy.stats import chi2 

def frecuency_monobit(n ):
   sn = 0
   for i in range(n):
       bit = random.randint(0, 1) # numeros aleatorios 0 y 1
       if bit == 1:
           sn = sn + 1
       else:
           sn = sn - 1
   sobs = math.fabs(sn)/math.sqrt(n) 
   p_value = math.erfc(sobs/math.sqrt(2))
   
   return p_value

# M-> logitud de cada bloque
# n-> numero total de muestras
# N-> numero de bloques
def frecuency_block(n, M):
    N = int( n / M)
    bloques = []
    for i in range(N):
        bloques.append(0.0)
        for j in range(M):
            bit = random.randint(0, 1)
            if bit == 1:
                bloques[i] = bloques[i] + 1.0
        bloques[i] = bloques[i] / M

    x2obs = 0
    for i in bloques:
        x2obs = x2obs + (i-0.5)**2
    x2obs = 4*M*x2obs
    return (1.0 - chi2.cdf(x2obs, M-1))

#numero de muestras
#verificar el contador de cambios
def run_test(n):
    media = 0.0     vn_ops = 1.0

    for i in range(n):
        nuevo = random.randint(0,1) #int(algo[i])
        media = media + nuevo
        if i > 0 and anterior != nuevo:
            vn_ops = vn_ops + 1
        anterior = nuevo
    media = media / n

    if math.fabs(media-0.5) <= (2 / math.sqrt(n)):
        return math.erfc( math.fabs(vn_ops - 2*n*media*(1-media) )  / ( 2*math.sqrt(2*n)*media*(1-media) )  )
    else:
        return "run_test(): Prueba rechazada"

print frecuency_monobit(1000)
print frecuency_block(100,10)
print run_test(100)

Estas pruebas generan un numero el cuales la regla de decisión es que si los valores son menores a 0.01 entonces se considera que la muestra es aleatoria.

Para revisar esto genere 10 veces cada prueba y como se puede observar, la mayoría de estas pasan la prueba lo cuál nos entrega un poco mas de seguridad de el generador.

Generador de numeros Poisson

Una vez que tengo una cierta seguridad de que los números que son el producto de el generador tienen las características necesarias para considerarse les aleatorios, ahora puedo intentar generar números con una distribución para mis fines de simulación, en mi caso practique haciendo un ejemplo convirtiendo los números producto del generador de python a la distribución poisson.

import math
import random

def poisson(lamb):
    acumulador = exponencial( lamb)
#    acumulador = 0
    contador = 0
    while( acumulador < 1):
        acumulador = acumulador + exponencial( lamb)
        contador+=1
    return contador


def exponencial(lamb):
    return (-1*math.log( random.random() ))/lamb

def buscar(vector, buscado):
  for i in range(len(vector)):
      if vector[i][0] == buscado:
          return i
  return -1

def contar(fl, n, lamb):
    numeros = []
    for i in range(n):
         generado = poisson(lamb)
         pos = buscar(numeros ,generado )
         if pos != -1:
             numeros[pos][1] += 1
         else:
             numeros.append([generado ,0 ])
    try:
        fl = open(fl , "w")
    except:
        print "Error al arbir el archivo"

    for i in numeros:
        i[1]=float(i[1])/n
        fl.write(str(i[0])+" "+format(i[1], '.6f')+"\n")

    try:
        fl.close()
    except:
        print "Error al cerrar el archivo"

contar( "poisson.dat", 1000000, 20)
poisson_esperado("esperado.dat", 1, 50, 20)

-----------------------------
Este fue mi primer resultado, para comparar los números que esperaría que salieran con los números obtenidos genere 1000000 números con mi generador de números pPoisson y después conté el número de ocurrencias de cada uno y estas dividirlas entre el total para obtener una probabilidad práctica de que genere ese número mi programa. Después lo gráfique junto con la formula de la distribución Poisson dándole como parámetro una lambda semejante y el resultado fue el siguiente:

Después de algunas pruebas y modificaciones el resultado fue el siguiente donde los números tienen una probabilidad práctica muy parecida a la de la formula de Poisson, con lo cuál compruebo que este generador tiene cierta calidad para usarla en futuras simulaciones.