choisir au hasard est une bonne façon d'éliminer le biais.

mise à Jour

nous pourrions faire le code rapidement si nous cherchons un x dans la gamme divisible par 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;

la boucle ci-dessus devrait être très rapide, disons 1 itération en moyenne.

@user1413793 est correct sur le problème. Je ne vais pas en discuter plus avant, sauf pour faire une remarque: Oui, pour les petites valeurs de n et les grandes valeurs de RAND_MAX , le biais modulo peut être très faible. Cependant, l'utilisation d'un modèle induisant un biais signifie que vous devez tenir compte du biais chaque fois que vous calculez un nombre aléatoire et de choisir des modèles différents pour différents cas. Et si vous faites le mauvais choix, les bugs qu'il introduit sont subtils et presque impossible à unit test. Comparé à juste utiliser l'outil approprié (comme arc4random_uniform ), c'est du travail supplémentaire, Pas Moins de travail. Faire plus de travail et obtenir une solution pire est une ingénierie terrible, surtout quand le faire correctement chaque fois est facile sur la plupart des plates-formes.

malheureusement, les implémentations de la solution sont toutes incorrectes ou moins efficaces qu'elles ne devraient l'être. (Chaque solution comporte divers commentaires expliquant les problèmes, mais aucune des solutions n'a été ils.) Ceci est susceptible de confondre le chercheur de réponse occasionnel, donc je fournis une mise en œuvre connue-bonne ici.

encore une fois, la meilleure solution est juste d'utiliser arc4random_uniform sur les plates-formes qui le fournissent, ou une solution à distance similaire pour votre plate-forme (comme "1519170920 Random.nextInt sur Java). Il fera la bonne chose sans coût de code pour vous. C'est presque toujours à l'appel correct de faire.

si vous ne pas avoir arc4random_uniform , alors vous pouvez utiliser la puissance d'opensource pour voir exactement comment il est mis en œuvre au-dessus d'une portée plus large RNG ( ar4random dans ce cas, mais une approche similaire pourrait également travailler sur le dessus d'autres RNGs).

voici le OpenBSD implementation :

/*
 * 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;
}

Il est intéressant de noter le dernier commit commentaire sur ce code pour ceux qui en ont besoin pour mettre en œuvre des choses semblables:

modifier arc4random_uniform () pour calculer 2**32 % upper_bound'' as -upper_bound % upper_bound". Simplifie le code et le rend même chose sur les architectures ILP32 et LP64, et légèrement plus rapide sur Architectures LP64 en utilisant un reste 32 bits au lieu d'un 64 bits reste.

souligné par Jorden Verwer sur tech@ ok ajacoutot; absence d'objection de la djm ou otto

L'implémentation Java est également facile à trouver (voir le lien précédent):

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;
 }

définition

le biais de Modulo est le biais inhérent à l'utilisation de l'arithmétique de modulo pour réduire un ensemble de sortie à un sous-ensemble de l'ensemble d'entrée. En général, un biais existe lorsque la correspondance entre l'ensemble d'entrées et de sorties n'est pas également répartie, comme dans le cas de l'utilisation de l'arithmétique modulo lorsque la taille de l'ensemble de sorties n'est pas un diviseur de la taille de l'ensemble d'entrées.

Ce biais est particulièrement difficile à éviter dans informatique, où les nombres sont représentés sous forme de chaînes de bits: 0s et 1s. Il est également extrêmement difficile de trouver des sources vraiment aléatoires d'aléatoire, mais cela dépasse le cadre de la présente discussion. pour le reste de cette réponse, supposons qu'il existe une source illimitée de bits vraiment aléatoires.

Exemple De Problème

considérons la simulation d'un rouleau de matrice (0 à 5) en utilisant ces bits aléatoires. Il y a 6 possibilités, donc nous besoin de suffisamment de bits pour représenter le nombre 6, qui est de 3 bits. Malheureusement, 3 bits aléatoires donne 8 résultats possibles:

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

nous pouvons réduire la taille du résultat fixé à exactement 6 en prenant la valeur modulo 6, mais cela présente le biais modulo problème: 110 donne un 0, et 111 donne un 1. Ce dé est chargé.

Solutions Possibles

approche 0:

plutôt que de compter sur des bits aléatoires, en théorie on pourrait engager une petite armée pour lancer des dés toute la journée et enregistrer les résultats dans une base de données, puis utiliser chaque résultat une seule fois. C'est à peu près aussi pratique qu'il semble, et plus que probable ne produirait pas vraiment des résultats aléatoires de toute façon (jeu de mots prévu).

approche 1:

au lieu d'utiliser le module, une solution naïve mais mathématiquement correcte est de jeter les résultats qui donnent 110 et 111 et simplement essayer à nouveau avec 3 nouveaux bits. Malheureusement, cela signifie qu'il y a une 25% de chance sur chaque rouleau qu'un reliquage sera nécessaire, y compris chacun des reliquats eux-mêmes. Il est clair que cela n'est pas pratique pour toutes les utilisations, sauf la plus insignifiante.

approche 2:

utilisez plus de bits: au lieu de 3 bits, Utilisez 4. Cela donne 16 résultats possibles. Bien sûr, re-rouler chaque fois que le résultat est plus grand que 5 rend les choses pires (10/16 = 62,5%) de sorte que seul ne sera pas utile.

noter que 2 * 6 = 12 < 16, de sorte que nous pouvons prendre en toute sécurité tout résultat inférieur à 12 et de réduire ce modulo 6 pour répartir également les résultats. Les quatre autres résultats doivent être écartés, puis reportés comme dans l'approche précédente.

sonne bien au début, mais vérifions le calcul:

4 discarded results / 16 possibilities = 25%

dans ce cas, 1 morceau supplémentaire n'a pas aidé du tout!

ce résultat est malheureux, mais essayons à nouveau avec 5 bits:

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

une amélioration certaine, mais pas assez bonne dans de nombreux cas pratiques. La bonne nouvelle est, ajouter plus de bits n'augmentera jamais les chances d'avoir besoin de jeter et re-roll . Cela vaut non seulement pour les dés, mais dans tous les cas.

As démontré cependant, l'ajout d'un 1 bit supplémentaire ne peut rien changer. en fait si nous augmentons notre rouleau à 6 bits, la probabilité reste de 6,25%.

2 questions supplémentaires:

1. si nous ajoutons suffisamment de bits, y a-t-il une garantie que la probabilité d'un écart diminuera?
2. combien de bits suffisent dans le cas général?

Solution Générale

Heureusement, la réponse à la première question est oui. Le problème avec 6 est que 2^x mod 6 flips entre 2 et 4 qui par coïncidence sont un multiple de 2 de l'autre, de sorte que pour un X Pair > 1,

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

ainsi 6 est une exception plutôt que la règle. Il est possible de trouver de plus grands modules qui produisent des puissances consécutives de 2 de la même manière, mais finalement ceci doit s'enrouler autour, et le la probabilité d'un écart sera réduit.

sans autre preuve, en général en utilisant le double du nombre des bits nécessaires fournira un plus petit, généralement insignifiant, chance de les jeter.

preuve de Concept

voici un exemple de programme qui utilise libcrypo D'OpenSSL pour fournir des octets aléatoires. Lors de la compilation, assurez-vous de lien à la bibliothèque avec -lcrypto que tout le monde devrait avoir 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;
}

j'encourage à jouer avec les valeurs MODULUS et ROLLS pour voir combien de rouleaux se produisent réellement dans la plupart des conditions. Une personne sceptique peut également vouloir enregistrer les valeurs calculées pour classer et vérifier la distribution semble normale.

Il ya deux plaintes habituelles avec l'utilisation de modulo.

• un est valable pour tous les générateurs. Il est facile de voir dans un cas limite. Si votre générateur a un RAND_MAX qui est 2 (qui n'est pas compatible avec le standard C) et que vous voulez seulement 0 ou 1 comme valeur, l'utilisation de modulo générera 0 deux fois plus souvent (quand le générateur génère 0 et 2) qu'il générera 1 (Quand le générateur génère 1). Notez que cela est vrai dès que vous n'avez pas les valeurs de goutte, quelle que soit la correspondance que vous utilisez entre les valeurs du générateur et la valeur désirée, l'une se produit deux fois plus souvent que l'autre.

• une sorte de générateur ont leurs bits moins significatifs moins aléatoire que l'autre, au moins pour certains de leurs paramètres, mais malheureusement ces paramètres ont d'autres caractéristiques intéressantes (tel a la possibilité D'avoir RAND_MAX un moins qu'une puissance de 2). Le problème est bien connu et depuis longtemps bibliothèque l'implémentation évite probablement le problème (par exemple l'implémentation de RAND() dans le standard C utilise ce genre de générateur, mais laisse tomber les 16 bits moins significatifs), mais certains aiment se plaindre à ce sujet et vous pourriez avoir de la malchance

utilisant quelque chose comme

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; 
}

pour générer un nombre aléatoire entre 0 et n évitera les deux problèmes (et il évite le débordement avec RAND_MAX == INT_MAX)

BTW, C++11 a introduit des moyens standard à la réduction et autre générateur que rand ().

La Solution de Mark

(la solution acceptée) est presque parfaite.

int x;

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

x %= n;

édité Mar 25 '16 à 23:16

Mark Amery 39k21170211

Cependant, il a une mise en garde qui ignore les 1 valable ensemble de résultats jetés dans tous les scénarios où la valeur de LA RAND_MAX (RM) est inférieur de 1 un multiple de N.

ie, quand le nombre de valeurs qui seraient mis au rebut comme invalide (I) est égal à N, Alors ils sont en fait un ensemble valide, pas un ensemble invalide.

par exemple:

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

comme vous pouvez voir dans l'exemple le nombre de valeurs rejetées = 4, Quand le nombre de valeurs rejetées = N alors l'ensemble est valide pour l'utilisation.

si nous décrivons la différence entre les valeurs N et RM comme D, i.e.:

D = (RM - N)

alors que la valeur de D devient plus petite, le Le pourcentage de rouleaux réenroulés non nécessaires à cause de cette méthode augmente à chaque multiplicateur naturel. (Donc quand RAND_MAX N'est pas égal à un nombre premier c'est une préoccupation valide)

par exemple:

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%

pour annuler ceci nous pouvons faire un simple amendement comme montré ici:

 int x;

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

 x %= n;

cela fournit une version plus générale de la formule qui tient compte des particularités supplémentaires de l'utilisation du module pour définissez vos valeurs max.

Exemples de l'aide d'une petite valeur pour RAND_MAX qui est un multiplicateur de N.

version originale:

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.

Version Modifiée:

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.

de plus, dans le cas où N devrait être le nombre de valeurs dans RAND_MAX; dans ce cas, vous pouvez définir N = RAND_MAX +1, sauf si RAND_MAX = INT_MAX.

en boucle, vous pouvez utiliser N = 1, et toute valeur de X sera accepté, cependant, et de mettre une instruction if dans le final de votre multiplicateur. Mais peut-être avez-vous un code qui peut avoir une raison valable pour retourner un 1 lorsque la fonction est appelée avec n = 1...

donc il peut être préférable d'utiliser 0, qui fournirait normalement une erreur Div 0, lorsque vous souhaitez avoir n = RAND_MAX+1

c'est à dire:

int x;

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

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

ces deux solutions résolvent le problème avec inutilement résultats valides rejetés qui se produiront lorsque RM+1 est un produit de N.

la deuxième version couvre également le scénario edge case lorsque vous avez besoin de n Pour égaler l'ensemble des valeurs possibles contenues dans RAND_MAX.

l'approche modifiée dans les deux cas est la même et permet de trouver une solution plus générale au besoin de fournir des nombres aléatoires valides et de minimiser les valeurs rejetées.

pour réitérer:

la Solution générale de base qui étend l'exemple de mark:

 int x;

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

 x %= n;

la Solution générale étendue qui permet un scénario supplémentaire 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();
}

avec une valeur de RAND_MAX de 3 (en réalité elle devrait être beaucoup plus élevée que cela, mais le biais existerait toujours) il est logique de ces calculs qu'il y ait un biais:

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

dans ce cas, le % 2 est ce que vous ne devriez pas faire quand vous voulez un nombre aléatoire entre 0 et 1 . Vous pourriez obtenir un nombre aléatoire entre 0 et 2 en faisant % 3 cependant, parce que dans ce cas: RAND_MAX est un multiple de 3 .

Autre méthode

il est beaucoup plus simple, mais pour ajouter à d'autres réponses , voici ma solution pour obtenir un nombre aléatoire entre 0 et n - 1 , donc n différentes possibilités, sans parti pris.

• le nombre de bits (et non d'octets) nécessaires pour encoder le nombre de possibilités est le nombre de bits de données aléatoires dont vous aurez besoin
• encoder le nombre de bits aléatoires
• si ce numéro est >= n , redémarrez (pas de modulo).

les données vraiment aléatoires ne sont pas faciles à obtenir, alors pourquoi utiliser plus de bits que nécessaire.

ci-dessous est un exemple dans Smalltalk, en utilisant un cache de bits à partir d'un générateur de nombres pseudo-aléatoires. Je ne suis pas un expert en matière de sécurité afin d'utiliser au vos propres risques.

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

comme l'indique réponse acceptée , le" biais modulo "a ses racines dans la faible valeur de RAND_MAX . Il utilise une très petite valeur de RAND_MAX (10) pour montrer que si RAND_MAX était 10, alors vous avez essayé de générer un nombre entre 0 et 2 en utilisant %, les résultats suivants résulteraient:

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

il y a donc 4 sorties de 0 (4/10 de chance) et seulement 3 sorties de 1 et 2 (3/10 de chance chacune).

donc c'est partial. Les plus faibles ont de meilleures chances de sortir.

mais qui apparaît de façon évidente seulement quand RAND_MAX est petit . Ou plus précisément, quand le nombre de votre modding par est grand par rapport à RAND_MAX .

une bien meilleure solution que looping (qui est incroyablement inefficace et ne devrait même pas être suggéré) est d'utiliser un PRNG avec un beaucoup plus grand plage de sortie. L'algorithme Mersenne Twister a une sortie maximale de 4.294.967.295. En tant que tel faisant MersenneTwister::genrand_int32() % 10 à toutes fins pratiques, sera également distribué et l'effet de biais modulo sera tout sauf disparaître.

je viens d'écrire un code pour la méthode non biaisée Coin Flip De Von Neumann, qui devrait théoriquement éliminer tout biais dans le processus de génération de nombres aléatoires. De plus amples informations sont disponibles au ( 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;
        }
    }
}