为采样创建概率分布的混合

Question

是否有一种通用的方法来连接SciPy(或NumPy)概率分布来创建混合概率分布,然后可以从中进行采样？

我有这样的分布显示使用类似的东西:

mixture_gaussian = (norm.pdf(x_axis, -3, 1) + norm.pdf(x_axis, 3, 1)) / 2

如果然后绘制如下:

但是,我不能从这个生成的模型中采样,因为它只是一个点的列表,将绘制为曲线.

注意,这个特定的分布只是一个简单的例子.我希望能够生成几种分布(包括"子"分布,这不仅仅是正态分布).理想情况下,我希望有一些功能可以自动规范化(即不必/ 2像上面的代码那样明确地执行.

SciPy/NumPy是否提供了一些轻松实现此目的的方法？

这个答案提供了一种方法,可以从多个分布中进行这样的采样,但它确实需要对给定的混合分布进行一些手工制作,特别是当想要不同地加权不同的"子"分布时.这是可用的,但如果可能的话,我希望方法更清洁,更直接.谢谢!

Answer 1

根据评论中@PaulPanzer 的指示，我创建了以下子类，以便从 SciPy 分布轻松创建混合模型。请注意，pdf我的问题不需要 ，但对我来说很好。

class MixtureModel(rv_continuous):
    def __init__(self, submodels, *args, **kwargs):
        super().__init__(*args, **kwargs)
        self.submodels = submodels

    def _pdf(self, x):
        pdf = self.submodels[0].pdf(x)
        for submodel in self.submodels[1:]:
            pdf += submodel.pdf(x)
        pdf /= len(self.submodels)
        return pdf

    def rvs(self, size):
        submodel_choices = np.random.randint(len(self.submodels), size=size)
        submodel_samples = [submodel.rvs(size=size) for submodel in self.submodels]
        rvs = np.choose(submodel_choices, submodel_samples)
        return rvs

mixture_gaussian_model = MixtureModel([norm(-3, 1), norm(3, 1)])
x_axis = np.arange(-6, 6, 0.001)
mixture_pdf = mixture_gaussian_model.pdf(x_axis)
mixture_rvs = mixture_gaussian_model.rvs(10)

Answer 2

从分布混合中采样(其中PDF与一些系数c_1,c_2,... c_n相加)等效于每个独立采样,然后,对于每个索引,从概率c_k中挑选第k个样本的值.

后者,混合,步骤可以有效地完成numpy.random.choice.这是一个混合了三种分布的例子.分布列在中distributions,其系数在coefficients.存在脂肪正态分布,均匀分布和窄的正态分布,系数为0.5,0.2,0.3.根据给定的系数生成在data[np.arange(sample_size), random_idx]之后发生的混合random_idx.

import numpy as np
import matplotlib.pyplot as plt

distributions = [
    {"type": np.random.normal, "kwargs": {"loc": -3, "scale": 2}},
    {"type": np.random.uniform, "kwargs": {"low": 4, "high": 6}},
    {"type": np.random.normal, "kwargs": {"loc": 2, "scale": 1}},
]
coefficients = np.array([0.5, 0.2, 0.3])
coefficients /= coefficients.sum()      # in case these did not add up to 1
sample_size = 100000

num_distr = len(distributions)
data = np.zeros((sample_size, num_distr))
for idx, distr in enumerate(distributions):
    data[:, idx] = distr["type"](size=(sample_size,), **distr["kwargs"])
random_idx = np.random.choice(np.arange(num_distr), size=(sample_size,), p=coefficients)
sample = data[np.arange(sample_size), random_idx]
plt.hist(sample, bins=100, density=True)
plt.show()