加快正态分布概率质量分配

Question

我们有N个用户，每个用户有P个点，其中每个点都是0到1之间的单个值。我们需要使用已知密度为1的正态分布来分布每个点的质量。 0.05点有一定的不确定性。此外，我们需要将质量包裹在0和1周围，例如一个位于0.95的点也将在0附近分配质量。我在下面提供了一个工作示例，该示例将正态分布分为D = 50个bin。该示例使用Python输入模块，但如果需要，您可以忽略它。

from typing import List, Any
import numpy as np
import scipy.stats
import matplotlib.pyplot as plt

D = 50
BINS: List[float] = np.linspace(0, 1, D + 1).tolist()


def probability_mass(distribution: Any, x0: float, x1: float) -> float:
    """
    Computes the area under the distribution, wrapping at 1.
    The wrapping is done by adding the PDF at +- 1.
    """
    assert x1 > x0
    return (
        (distribution.cdf(x1) - distribution.cdf(x0))
        + (distribution.cdf(x1 + 1) - distribution.cdf(x0 + 1))
        + (distribution.cdf(x1 - 1) - distribution.cdf(x0 - 1))
    )


def point_density(x: float) -> List[float]:
    distribution: Any = scipy.stats.norm(loc=x, scale=0.05)
    density: List[float] = []
    for i in range(D):
        density.append(probability_mass(distribution, BINS[i], BINS[i + 1]))
    return density


def user_density(points: List[float]) -> Any:

    # Find the density of each point
    density: Any = np.array([point_density(p) for p in points])

    # Combine points and normalize
    combined = density.sum(axis=0)
    return combined / combined.sum()


if __name__ == "__main__":

    data: List[float] = [0.05, 0.3, 0.5, 0.5]
    density = user_density(data)

    # Plot the density as a bar chat
    middle: List[float] = []
    for i in range(D):
        middle.append((BINS[i] + BINS[i + 1]) / 2)
    plt.bar(x=middle, height=density, width=1.0 / D + 0.001)
    plt.xlim(0, 1)
    plt.xlabel("x")
    plt.ylabel("Density")
    plt.show()

在此示例中，N = 1，D = 50，P = 4。但是，我们希望将此方法扩展到N = 10000和P = 100，同时尽可能快。我不清楚我们如何将这种方法矢量化。我们如何最好地加快速度？

编辑

更快的解决方案可能会有稍微不同的结果。例如，它可以近似正态分布而不是使用精确的正态分布。

EDIT2

我们只关心使用density函数计算user_density()。该图仅用于帮助解释该方法。我们不在乎情节本身：）

Answer 1

这将是我的向量化方法：

data = np.array([0.05, 0.3, 0.5, 0.5])

np.random.seed(31415)
# random noise
randoms = np.random.normal(0,1,(len(data), int(1e5))) * 0.05

# samples with noise
samples = data[:,None] + randoms

# wrap [0,1]
samples = (samples % 1).ravel()


# histogram
hist, bins, patches = plt.hist(samples, bins=BINS, density=True)

输出：