Alias method

参考链接：

1. python实现：https://lips.cs.princeton.edu/the-alias-method-efficient-sampling-with-many-discrete-outcomes/

import numpy        as np
import numpy.random as npr

def alias_setup(probs):
    K       = len(probs)
    q       = np.zeros(K)
    J       = np.zeros(K, dtype=np.int)

    # Sort the data into the outcomes with probabilities
    # that are larger and smaller than 1/K.
    smaller = []
    larger  = []
    for kk, prob in enumerate(probs):
        q[kk] = K*prob
        if q[kk] < 1.0:
            smaller.append(kk)
        else:
            larger.append(kk)

    # Loop though and create little binary mixtures that
    # appropriately allocate the larger outcomes over the
    # overall uniform mixture.
    while len(smaller) > 0 and len(larger) > 0:
        small = smaller.pop()
        large = larger.pop()

        J[small] = large
        q[large] = q[large] - (1.0 - q[small])

        if q[large] < 1.0:
            smaller.append(large)
        else:
            larger.append(large)

    return J, q

def alias_draw(J, q):
    K  = len(J)

    # Draw from the overall uniform mixture.
    kk = int(np.floor(npr.rand()*K))

    # Draw from the binary mixture, either keeping the
    # small one, or choosing the associated larger one.
    if npr.rand() < q[kk]:
        return kk
    else:
        return J[kk]

K = 5
N = 1000

# Get a random probability vector.
probs = npr.dirichlet(np.ones(K), 1).ravel()

# Construct the table.
J, q = alias_setup(probs)

# Generate variates.
X = np.zeros(N)
for nn in xrange(N):
    X[nn] = alias_draw(J, q)

1. 算法思想：https://www.keithschwarz.com/darts-dice-coins/

Algorithm: Alias Method

• Initialization:
  1. Create arrays AliasAlias and Prob , each of size n .
  2. Create a balanced binary search tree T .
  3. Insert n⋅pi into T for each probability i .
  4. For j=1 to n−1 :
     1. Find and remove the smallest value in T ; call it pl.
     2. Find and remove the largest value in T ; call it pg .
     3. Set Prob[l]=pl .
     4. Set Alias[l]=g .
     5. Set pg:=pg−(1−pl) .
     6. Add pg to T.
  5. Let i be the last probability remaining, which must have weight 1.
  6. Set Prob[i]=1 .
• Generation:
  1. Generate a fair die roll from an n-sided die; call the side i .
  2. Flip a biased coin that comes up heads with probability Prob[i] .
  3. If the coin comes up "heads," return i .
  4. Otherwise, return Alias[i] .

1. 思想简单介绍：https://blog.csdn.net/qq_33765907/article/details/79182355
2. 各版本代码实现：https://blog.csdn.net/z69183787/article/details/53782322