随机器构建

n随机器：取值有n种情况的随机器
均匀n随机器：各种取值概率相等的n随机器

问题：由n随机器构造m随机器
具体问题：

• 均匀到不均匀
• 不均匀到均匀
• 不均匀到均匀
• 均匀到均匀
• 大均匀到小均匀
• 小均匀到大均匀

本文研究不均匀到均匀随机器。
重复试验k次，得到若干个均匀n随机器，然后利用这些均匀随机器构建想要的均匀随机器。

下面算法证明这个思路是正确的。

import random

import numpy as np

x = 0.3
y = 0.7


def rg():
    if random.random() < x:
        return 0
    else:
        return 1


k = 10
a = dict()


def build_table():
    cnt_map = dict()
    for i in range(2 ** k):
        s = np.array(list(map(int, bin(i)[2:])))
        cnt0 = np.count_nonzero(s == 1)
        if cnt0 not in cnt_map:
            cnt_map[cnt0] = []
        cnt_map[cnt0].append(i)
    for group in cnt_map.values():
        while len(group) > 1:
            bit_cnt = int(np.log2(len(group)))
            for ind in range(2 ** bit_cnt):
                now = np.zeros(bit_cnt, dtype=np.int32)
                for j in range(bit_cnt):
                    now[j] = 1 if ind & (1 << j) else 0
                if len(now):
                    a[group[ind]] = now
            group = group[2 ** bit_cnt:]


cache = []
called = 0


def my_rg():
    global called
    if len(cache):
        return cache.pop()
    while 1:
        called += k
        v = [rg() for i in range(k)]
        v = int(''.join(map(str, v)), 2)
        if v in a:
            cache.extend(a[v])
            return my_rg()


build_table()
ans = np.array([my_rg() for i in range(10000)])
print(np.count_nonzero(ans))
print(called / len(ans))

如何寻找k？怎样找到最合适的试验次数，使得获得每个随机数调用的随机器次数尽量少。

import numpy as np
import pylab as plt
from scipy.misc import comb

x = 0.3
y = 0.7


def try_k(k):
    # 尝试k次
    p = np.array([x ** i * y ** (k - i) for i in range(k + 1)])  # 每种类型发生的概率
    cnt = np.array([comb(k, i) for i in range(k + 1)])  # 每种类型的个数
    break_prob_list = []
    can_get_list = []
    for group_cnt, prob in zip(cnt, p):
        while group_cnt > 1:
            can_get = int(np.log2(group_cnt))
            if not can_get: break
            break_prob = 2 ** can_get * prob
            break_prob_list.append(break_prob)
            can_get_list.append(can_get)
            group_cnt -= 2 ** can_get
    break_prob = np.array(break_prob_list)
    can_get = np.array(can_get_list)
    return k / np.sum(break_prob * can_get)


min_index = None
min_value = 0
ans = []
max_try = 100
for i in range(2, max_try):
    now = try_k(i)
    print(i, now)
    ans.append(now)
    if min_index is None or min_value > now:
        min_value = now
        min_index = i
print('最佳答案', min_index, min_value)
plt.plot(list(range(2, max_try)), ans)
plt.show()