返回

算法 | 相等概率随机分布设计

后端

我们经常会在编写算法时，设计一个给定固定概率的随机分布，例如 [0, 1] 或者 [-1, 1]。通常情况下，编写一个产生 [0, 1] 或 [-1, 1] 均匀随机分布的函数很简单，但如果我们需要生成一个不同于 [0, 1] 的非均匀分布，事情就会变得复杂一些。

如何根据不相等概率[0,1]设计出相等概率的[0,1]？例如我们已知产生随机数 x1 的概率分布为 p1(x1) ， x1 ∈ [a,b]，如果 x1 满足某条件，比如 x1 > c，其产生的概率为 p2(x1|x1>c)，其中 a, b, c 是常数，求当 x1>c 时的 x1 的概率分布 p2(x1) 。
这个过程可以使用蒙特卡罗方法和拒绝采样法来完成。

第一步：使用蒙特卡罗方法生成一个随机数 x1，其概率分布为 p1(x1)。

第二步：判断 x1 是否满足条件 x1 > c。

第三步：如果 x1 满足条件 x1 > c，则接受 x1 并将 p2(x1) 设为 p2(x1) = p1(x1)/p2(x1|x1>c)。

第四步：否则，拒绝 x1 并重复步骤 1 到 3。

重复这些步骤，直到生成一定数量的 x1。x1 的概率分布将接近 p2(x1)。

代码示例

import numpy as np

def inverse_cdf(p, x, a, b):
  """
  Calculate the inverse CDF of a given distribution.

  Args:
    p: The probability value.
    x: The random variable.
    a: The lower bound of the distribution.
    b: The upper bound of the distribution.

  Returns:
    The value of x corresponding to the given probability.
  """

  # Calculate the cumulative distribution function.
  cdf = np.zeros_like(p)
  cdf[0] = p[0]
  for i in range(1, len(p)):
    cdf[i] = cdf[i - 1] + p[i]

  # Find the value of x corresponding to the given probability.
  x_inv = np.zeros_like(p)
  for i in range(len(p)):
    x_inv[i] = a + (b - a) * np.searchsorted(cdf, p[i])

  return x_inv


def generate_non_uniform_random_numbers(p, a, b, n):
  """
  Generate a list of n random numbers from a given non-uniform distribution.

  Args:
    p: The probability distribution.
    a: The lower bound of the distribution.
    b: The upper bound of the distribution.
    n: The number of random numbers to generate.

  Returns:
    A list of n random numbers.
  """

  # Generate a list of n random numbers from a uniform distribution.
  u = np.random.uniform(size=n)

  # Calculate the inverse CDF of the given distribution.
  x = inverse_cdf(u, p, a, b)

  return x


# Define the non-uniform distribution.
p = [0.2, 0.3, 0.4, 0.1]

# Generate 1000 random numbers from the non-uniform distribution.
x = generate_non_uniform_random_numbers(p, 0, 1, 1000)

# Plot the histogram of the generated random numbers.
plt.hist(x, bins=10)
plt.show()

以上给出了如何根据不相等概率 [0, 1] 设计出相等概率的 [0, 1] 的算法。这种算法在许多应用中都有用，例如模拟随机变量或生成随机样本。