Seguir seleccionando un azar es una buena manera de eliminar el sesgo.

Actualización

Podríamos hacer el código rápido si buscamos una x en el rango divisible por n.

// Assumptions
// rand() in [0, RAND_MAX]
// n in (0, RAND_MAX]

int x; 

// Keep searching for an x in a range divisible by n 
do {
    x = rand();
} while (x >= RAND_MAX - (RAND_MAX % n)) 

x %= n;

El bucle anterior debe ser muy rápido, digamos 1 iteración en promedio.

@user1413793 es correcto sobre el problema. No voy a discutir eso más, excepto para hacer un punto: sí, para valores pequeños de n y valores grandes de RAND_MAX, el sesgo modulo puede ser muy pequeño. Pero el uso de un patrón que induce el sesgo significa que debe considerar el sesgo cada vez que calcule un número aleatorio y elija diferentes patrones para diferentes casos. Y si tomas la decisión equivocada, los errores que introduce son sutiles y casi imposibles de probar unitariamente. En comparación con solo usando la herramienta adecuada (como arc4random_uniform), eso es trabajo extra, no menos trabajo. Hacer más trabajo y obtener una solución peor es una ingeniería terrible, especialmente cuando hacerlo bien cada vez es fácil en la mayoría de las plataformas.

Desafortunadamente, las implementaciones de la solución son todas incorrectas o menos eficientes de lo que deberían ser. (Cada solución tiene varios comentarios que explican los problemas, pero ninguna de las soluciones se ha fijado para abordarlos.) Esto es probable que confunda el casual buscador de respuestas, así que estoy proporcionando una implementación bien conocida aquí.

Nuevamente, la mejor solución es solo usar arc4random_uniform en plataformas que lo proporcionan, o una solución de rango similar para su plataforma (como Random.nextInt en Java). Hará lo correcto sin costo de código para usted. Esta es casi siempre la llamada correcta a hacer.

Si no tiene arc4random_uniform, entonces puede usar el poder de opensource para ver exactamente cómo se implementa en la parte superior de un RNG de mayor rango (ar4random en este caso, pero un enfoque similar también podría funcionar sobre otros RNGS).

Aquí está la implementación de OpenBSD :

/*
 * Calculate a uniformly distributed random number less than upper_bound
 * avoiding "modulo bias".
 *
 * Uniformity is achieved by generating new random numbers until the one
 * returned is outside the range [0, 2**32 % upper_bound).  This
 * guarantees the selected random number will be inside
 * [2**32 % upper_bound, 2**32) which maps back to [0, upper_bound)
 * after reduction modulo upper_bound.
 */
u_int32_t
arc4random_uniform(u_int32_t upper_bound)
{
    u_int32_t r, min;

    if (upper_bound < 2)
        return 0;

    /* 2**32 % x == (2**32 - x) % x */
    min = -upper_bound % upper_bound;

    /*
     * This could theoretically loop forever but each retry has
     * p > 0.5 (worst case, usually far better) of selecting a
     * number inside the range we need, so it should rarely need
     * to re-roll.
     */
    for (;;) {
        r = arc4random();
        if (r >= min)
            break;
    }

    return r % upper_bound;
}

Vale la pena señalar el último comentario de commit sobre este código para aquellos que necesitan implementar cosas similares:

Cambie arc4random_uniform() para calcular 2**32 % upper_bound'' as -upper_bound % upper_bound". Simplifica el código y lo hace el lo mismo en ambas arquitecturas ILP32 y LP64, y también un poco más rápido en LP64 arquitecturas mediante el uso de un resto de 32 bits en lugar de un 64 bits resto.

Señalado por Jorden Verwer en tech@ ok deraadt; no hay objeciones de djm o otto

La implementación de Java también es fácilmente encontrable (ver enlace anterior):

public int nextInt(int n) {
   if (n <= 0)
     throw new IllegalArgumentException("n must be positive");

   if ((n & -n) == n)  // i.e., n is a power of 2
     return (int)((n * (long)next(31)) >> 31);

   int bits, val;
   do {
       bits = next(31);
       val = bits % n;
   } while (bits - val + (n-1) < 0);
   return val;
 }

Definición

Modulo Bias es el sesgo inherente en el uso de la aritmética modulo para reducir un conjunto de salida a un subconjunto del conjunto de entrada. En general, existe un sesgo cuando la asignación entre el conjunto de entrada y salida no está distribuida equitativamente, como en el caso de usar aritmética de módulo cuando el tamaño del conjunto de salida no es un divisor del tamaño del conjunto de entrada.

Este sesgo es particularmente difícil de evitar en la computación, donde los números se representan como cadenas de bits: 0s y 1s.Encontrar fuentes verdaderamente aleatorias de aleatoriedad también es extremadamente difícil, pero está más allá del alcance de esta discusión. Para el resto de esta respuesta, supongamos que existe una fuente ilimitada de bits verdaderamente aleatorios.

Ejemplo de problema

Consideremos simular un rollo de troquel (0 a 5) usando estos bits aleatorios. Hay 6 posibilidades, por lo que necesitamos suficientes bits para representar el número 6, que es 3 bits. Desafortunadamente, 3 bits aleatorios rinden 8 posibles resultados:

000 = 0, 001 = 1, 010 = 2, 011 = 3
100 = 4, 101 = 5, 110 = 6, 111 = 7

Podemos reducir el tamaño del conjunto de resultados a exactamente 6 tomando el valor módulo 6, sin embargo esto presenta el problema sesgo de módulo: 110 produce un 0, y 111 produce un 1. Este dado está cargado.

Soluciones potenciales

Enfoque 0:

En lugar de confiar en bits aleatorios, en teoría uno podría contratar a un pequeño ejército para tirar dados todo el día y registrar los resultados en una base de datos, y luego usar cada resultado solo una vez. Este es casi tan práctico como suena, y lo más probable es que no produzca resultados verdaderamente aleatorios de todos modos (juego de palabras).

Enfoque 1:

En lugar de usar el módulo, una solución ingenua pero matemáticamente correcta es descartar los resultados que producen 110 y 111 y simplemente intentarlo de nuevo con 3 bits nuevos. Desafortunadamente, esto significa que hay un 25% de probabilidad en cada tirada de que se requiera un re-roll, incluyendo cada uno de los re-rolls en sí. Esto es claramente poco práctico para todos los usos menos los más triviales.

Enfoque 2:

Use más bits: en lugar de 3 bits, use 4. Esto produce 16 posibles resultados. Por supuesto, volver a rodar cada vez que el resultado es mayor que 5 empeora las cosas (10/16 = 62.5%) por lo que solo no ayudará.

Observe que 2 * 6 = 12

Suena bien al principio, pero vamos a comprobar las matemáticas:{[17]]}

4 discarded results / 16 possibilities = 25%

En este caso, 1 bit extra no ayudó en absoluto!

Ese resultado es desafortunado, pero intentémoslo de nuevo con 5 bits:

32 % 6 = 2 discarded results; and
2 discarded results / 32 possibilities = 6.25%

Una mejora definitiva, pero no lo suficientemente buena en muchos casos prácticos. La buena noticia es, agregar más bits nunca aumentará las posibilidades de tener que descartar y volver a rodar. Esto vale no solo para los dados, sino en todos caso.

Como se demostró sin embargo, agregar 1 bit adicional puede no cambiar nada. De hecho, si aumentamos nuestro rollo a 6 bits, la probabilidad permanece 6.25%.

Esto plantea 2 preguntas adicionales:

1. Si agregamos suficientes bits, ¿hay una garantía de que la probabilidad de un descarte disminuirá?
2. ¿Cuántos bits son suficientes en el caso general?

Solución general

Afortunadamente la respuesta a la primera la pregunta es sí. El problema con 6 es que 2^x mod 6 voltea entre 2 y 4 que casualmente son un múltiplo de 2 entre sí, de modo que para un par x > 1,

[2^x mod 6] / 2^x == [2^(x+1) mod 6] / 2^(x+1)

Por lo tanto, 6 es una excepción más que la regla. Es posible encontrar módulos más grandes que producen potencias consecutivas de 2 de la misma manera, pero eventualmente esto debe envolver alrededor, y la probabilidad de un descarte se reducirá.

Sin ofrecer más pruebas, en general usando el doble de numero de bits requeridos proporcionará un menor, por lo general insignificante, posibilidad de un descarte.

Prueba de Concepto

Aquí hay un programa de ejemplo que usa libcrypo de OpenSSL para suministrar bytes aleatorios. Al compilar, asegúrese de enlazar a la biblioteca con -lcrypto que la mayoría de todos deberían tener disponible.

#include <iostream>
#include <assert.h>
#include <limits>
#include <openssl/rand.h>

volatile uint32_t dummy;
uint64_t discardCount;

uint32_t uniformRandomUint32(uint32_t upperBound)
{
    assert(RAND_status() == 1);
    uint64_t discard = (std::numeric_limits<uint64_t>::max() - upperBound) % upperBound;
    uint64_t randomPool = RAND_bytes((uint8_t*)(&randomPool), sizeof(randomPool));

    while(randomPool > (std::numeric_limits<uint64_t>::max() - discard)) {
        RAND_bytes((uint8_t*)(&randomPool), sizeof(randomPool));
        ++discardCount;
    }

    return randomPool % upperBound;
}

int main() {
    discardCount = 0;

    const uint32_t MODULUS = (1ul << 31)-1;
    const uint32_t ROLLS = 10000000;

    for(uint32_t i = 0; i < ROLLS; ++i) {
        dummy = uniformRandomUint32(MODULUS);
    }
    std::cout << "Discard count = " << discardCount << std::endl;
}

Animo a jugar con los valores MODULUS y ROLLS para ver cuántos re-rolls ocurren realmente bajo la mayoría de las condiciones. Una persona escéptica también puede desea guardar los valores calculados en el archivo y verificar que la distribución parezca normal.

Hay dos quejas habituales con el uso de modulo.

• Uno es válido para todos los generadores. Es más fácil de ver en un caso límite. Si su generador tiene un RAND_MAX que es 2 (que no es compatible con el estándar C) y desea solo 0 o 1 como valor, utilizando modulo generará 0 el doble de veces (cuando el generador genera 0 y 2) que generará 1 (cuando el generador genera 1). Tenga en cuenta que esto es cierto tan pronto como no suelte valores, cualquiera que sea la asignación que están utilizando desde los valores del generador hasta el deseado, uno ocurrirá dos veces más a menudo que el otro.

• Algún tipo de generador tiene sus bits menos significativos menos aleatorios que el otro, al menos para algunos de sus parámetros, pero tristemente esos parámetros tienen otra característica interesante (tal tiene ser capaz de tener RAND_MAX uno menos que una potencia de 2). El problema es bien conocido y durante mucho tiempo la implementación de la biblioteca probablemente evite el problema (por ejemplo, la muestra la implementación de rand() en el estándar C usa este tipo de generador, pero deja caer los 16 bits menos significativos), pero a algunos les gusta quejarse de eso y puede tener mala suerte

Usando algo como

int alea(int n){ 
 assert (0 < n && n <= RAND_MAX); 
 int partSize = 
      n == RAND_MAX ? 1 : 1 + (RAND_MAX-n)/(n+1); 
 int maxUsefull = partSize * n + (partSize-1); 
 int draw; 
 do { 
   draw = rand(); 
 } while (draw > maxUsefull); 
 return draw/partSize; 
}

Generar un número aleatorio entre 0 y n evitará ambos problemas (y evita el desbordamiento con RAND_MAX == INT_MAX)

POR cierto, C++11 introdujo formas estándar para la reducción y otro generador que rand().

La Solución de Mark (La solución aceptada) es Casi Perfecta.

int x;

do {
    x = rand();
} while (x >= (RAND_MAX - RAND_MAX % n));

x %= n;

Editado Mar 25 '16 a las 23:16

Mark Amery 39k21170211

Sin embargo, tiene una advertencia que descarta 1 conjunto válido de resultados descartados en todos los escenarios donde el valor del RAND_MAX (RM) es 1 menor que un múltiplo de N.

Es decir, cuando el número de valores que se descartarían como inválidos (I) es igual a N, entonces son en realidad un conjunto válido, no un inválido establecer.

EG:

RM = 255
N = 4

Discard X => RM - RM % N 

When X => 252, Discarded values = 252, 253, 254, 255

Number of discarded Values (I) = RM % N + 1

Como se puede ver en el ejemplo el recuento de los valores descartados = 4, cuando el recuento de valores descartados = N entonces el conjunto es válido para su uso.

Si describimos la diferencia entre los valores N y RM como D, es decir:

D = (RM - N)

Luego, a medida que el valor de D se vuelve más pequeño, el Porcentaje de vueltas innecesarias debido a este método aumenta en cada multiplicativo natural. (Así que cuando RAND_MAX NO es igual a un Número Primo esto es un válido preocupación)

EG:

RM=255 , N=2 Then: D = 253, Lost percentage = 0.78125%

RM=255 , N=4 Then: D = 251, Lost percentage = 1.5625%
RM=255 , N=8 Then: D = 247, Lost percentage = 3.125%
RM=255 , N=16 Then: D = 239, Lost percentage = 6.25%
RM=255 , N=32 Then: D = 223, Lost percentage = 12.5%
RM=255 , N=64 Then: D = 191, Lost percentage = 25%
RM=255 , N= 128 Then D = 127, Lost percentage = 50%

Para negar esto podemos hacer una simple enmienda Como se muestra aquí:

 int x;

 do {
     x = rand();
 } while (x > (RAND_MAX - ( ( ( RAND_MAX % n ) + 1 ) % n) );

 x %= n;

Esto proporciona una versión más general de la fórmula que tiene en cuenta las peculiaridades adicionales de usar módulo para definir sus valores máximos.

Ejemplos de usar un valor pequeño para RAND_MAX que es un multiplicativo de N.

Mark Versión original:

RAND_MAX = 3, n = 2, Values in RAND_MAX = 0,1,2,3, Valid Sets = 0,1 and 2,3.
When X >= (RAND_MAX - ( RAND_MAX % n ) )
When X >= 2 the value will be discarded, even though the set is valid.

Versión enmendada:

RAND_MAX = 3, n = 2, Values in RAND_MAX = 0,1,2,3, Valid Sets = 0,1 and 2,3.
When X > (RAND_MAX - ( ( RAND_MAX % n  ) + 1 ) % n )
When X > 3 the value would be discarded, but this is not a vlue in the set RAND_MAX so there will be no discard.

Además, en el caso en que N debe ser el número de valores en RAND_MAX; en este caso, podría establecer N = RAND_MAX + 1, a menos que RAND_MAX = INT_MAX.

En lo que respecta al bucle, puedes usar N = 1, y cualquier valor de X será aceptado, sin embargo, y poner una sentencia IF para tu multiplicador final. Pero tal vez usted tiene código que puede tener una razón válida para devolver un 1 cuando la función se llama con n = 1...

Así que puede ser mejor usar 0, que normalmente proporcionaría un error Div 0, cuando desea tener n = RAND_MAX + 1

IE:

int x;

if n != 0 {
    do {
        x = rand();
    } while (x > (RAND_MAX - ( ( ( RAND_MAX % n ) + 1 ) % n) );

    x %= n;
} else {
    x = rand();
}

Ambas soluciones resuelven el problema con resultados válidos innecesariamente descartados que ocurrirán cuando RM + 1 es un producto de n.

La segunda versión también cubre el escenario de caso de borde cuando se necesita n para igualar el conjunto total posible de valores contenidos en RAND_MAX.

El enfoque modificado en ambos es el mismo y permite una solución más general a la necesidad de proporcionar números aleatorios válidos y minimizar los descartados valor.

Para reiterar:

La Solución General Básica que extiende el ejemplo de mark:

 int x;

 do {
     x = rand();
 } while (x > (RAND_MAX - ( ( ( RAND_MAX % n ) + 1 ) % n) );

 x %= n;

La Solución General Extendida que Permite un escenario adicional de RAND_MAX + 1 = n:

int x;

if n != 0 {
    do {
        x = rand();
    } while (x > (RAND_MAX - ( ( ( RAND_MAX % n ) + 1 ) % n) );

    x %= n;
} else {
    x = rand();
}

Con un valor RAND_MAX de 3 (en realidad debería ser mucho más alto que eso, pero el sesgo seguiría existiendo) tiene sentido a partir de estos cálculos que hay un sesgo:

1 % 2 = 1 2 % 2 = 0 3 % 2 = 1 random_between(1, 3) % 2 = more likely a 1

En este caso, el % 2 es lo que no debe hacer cuando desea un número aleatorio entre 0 y 1. Sin embargo, puedes obtener un número aleatorio entre 0 y 2 haciendo % 3, porque en este caso: RAND_MAX es un múltiplo de 3.

Otro método

Hay mucho más simple, pero para agregar a otras respuestas, aquí está mi solución para obtener un número aleatorio entre 0 y n - 1, por lo que n diferentes posibilidades, sin sesgo.

• el número de bits (no bytes) necesarios para codificar el número de posibilidades es el número de bits de datos aleatorios que necesitará
• codifica el número desde bits aleatorios
• si este número es >= n, reinicie (sin módulo).

Los datos realmente aleatorios no son fáciles de obtener, entonces por qué usar más bits de los necesarios.

A continuación se muestra un ejemplo en Smalltalk, usando una caché de bits de un generador de números pseudo-aleatorios. No soy experto en seguridad, así que úsalo bajo tu propio riesgo.

next: n

    | bitSize r from to |
    n < 0 ifTrue: [^0 - (self next: 0 - n)].
    n = 0 ifTrue: [^nil].
    n = 1 ifTrue: [^0].
    cache isNil ifTrue: [cache := OrderedCollection new].
    cache size < (self randmax highBit) ifTrue: [
        Security.DSSRandom default next asByteArray do: [ :byte |
            (1 to: 8) do: [ :i |    cache add: (byte bitAt: i)]
        ]
    ].
    r := 0.
    bitSize := n highBit.
    to := cache size.
    from := to - bitSize + 1.
    (from to: to) do: [ :i |
        r := r bitAt: i - from + 1 put: (cache at: i)
    ].
    cache removeFrom: from to: to.
    r >= n ifTrue: [^self next: n].
    ^r

Como indica la respuesta aceptada, "sesgo modulo" tiene sus raíces en el bajo valor de RAND_MAX. Usa un valor extremadamente pequeño de RAND_MAX (10) para mostrar que si RAND_MAX fuera 10, entonces intentaste generar un número entre 0 y 2 usando %, los siguientes resultados resultarían:

rand() % 3   // if RAND_MAX were only 10, gives
output of rand()   |   rand()%3
0                  |   0
1                  |   1
2                  |   2
3                  |   0
4                  |   1
5                  |   2
6                  |   0
7                  |   1
8                  |   2
9                  |   0

Así que hay 4 salidas de 0 (4/10 de probabilidad) y solo 3 salidas de 1 y 2 (3/10 de probabilidad cada una).

Así que es parcial. Los números más bajos tienen una mejor oportunidad de venir fuera.

Pero eso solo aparece tan obviamente cuando RAND_MAX es pequeño. O más específicamente, cuando el número que está modificando es grande comparado con RAND_MAX.

Una solución mucho mejor que looping (que es increíblemente ineficiente y ni siquiera debería sugerirse) es usar un PRNG con un rango de salida mucho mayor. El algoritmo Mersenne Twister tiene una salida máxima de 4,294,967,295. Como tal hacer MersenneTwister::genrand_int32() % 10 para todos los efectos, será distribuido equitativamente y el efecto de sesgo modulo desaparecerá.

Acabo de escribir un código para el Método de Volteo de Moneda Imparcial de Von Neumann, que teóricamente debería eliminar cualquier sesgo en el proceso de generación de números aleatorios. Puede encontrar más información en ( http://en.wikipedia.org/wiki/Fair_coin )

int unbiased_random_bit() {    
    int x1, x2, prev;
    prev = 2;
    x1 = rand() % 2;
    x2 = rand() % 2;

    for (;; x1 = rand() % 2, x2 = rand() % 2)
    {
        if (x1 ^ x2)      // 01 -> 1, or 10 -> 0.
        {
            return x2;        
        }
        else if (x1 & x2)
        {
            if (!prev)    // 0011
                return 1;
            else
                prev = 1; // 1111 -> continue, bias unresolved
        }
        else
        {
            if (prev == 1)// 1100
                return 0;
            else          // 0000 -> continue, bias unresolved
                prev = 0;
        }
    }
}