Untitled

# Alias sampling

def alias_sample(probs):
    n = len(probs)
    F = probs * n

    alias_list = np.zeros(n) 
    prob_list = np.zeros(n)

    G = list(np.where(F>=1)[0])
    S = list(np.where(F<1)[0])

    while G and S:
        j = S.pop()
        i = G.pop()
        
        prob_list[j] = F[j]
        alias_list[j] = i
        F[i] = F[i] - (1.0 - F[j])

        if F[i] < 1:
            S.append(i)
        else:
            G.append(i)
    
    # rest should be assigned 1
    for i in G + S:
        prob_list[i] = 1.0

    return prob_list, alias_list

probs = np.array([7/48,5/48,1/8,1/16,1/4,5/16])

prob_list,alias_list = alias_sample(probs)


def alias_sample(prob, alias, n):
    K = len(prob)
    samples = np.zeros(n)
    for i in range(n):
        k = np.random.randint(K)
        if np.random.random(1) < prob[k]:
            samples[i] = k
        else:
            samples[i] = alias[k]
    return samples



plt.hist(alias_sample(prob_list,alias_list,10000)+1)