pymc 确定随机变量的总和

Question

我有两个独立的正态分布随机变量a, b。在 pymc 中是这样的：

from pymc import Normal


def model():
    a = Normal('a', tau=0.01)
    b = Normal('b', tau=0.1)

如果我们可以将其视为正态分布，我想知道a+b 是什么，即：

from pymc import Normal


def model():
    a = Normal('a', tau=0.01)
    b = Normal('b', tau=0.1)

    tau_c = Uniform("tau_c", lower=0.0, upper=1.0)
    c = Normal("a+b", tau=tau_c, observed=True, value=a+b)

然后我想估计 tau_c，但它不适用于 pymc，因为 a 和 b 是随机的（如果它们是数组，这是可能的，但我没有观察到 @987654329 @或b，我只知道他们的分布）。

我认为我可以做到的一种方法是使用每个 a 和 b 的分布生成随机值，然后执行以下操作：

def model(a, b):
    tau_c = Uniform("tau_c", lower=0.0, upper=1.0)
    c = Normal("a+b", tau=tau_c, observed=True, value=a+b)

但我认为使用 pymc 有更好的方法。

谢谢！

Answer 1

如果我正确理解了您的问题和代码，您应该做一些更简单的事情。如果要估计由 a 和 b 之和给出的分布参数，则仅使用以下示例中的第一个块。如果您还想独立于变量 b 的参数估计变量 a 的参数，则使用其他两个块

with pm.Model() as model:
    mu = pm.Normal('mu', mu=0, sd=10)
    sd = pm.HalfNormal('sd', 10)
    alpha = pm.Normal('alpha', mu=0, sd=10)
    ab = pm.SkewNormal('ab', mu=mu, sd=sd, alpha=alpha, observed=a+b)

    mu_a = pm.Normal('mu_a', mu=0, sd=10)
    sd_a = pm.HalfNormal('sd_a', 10)
    alpha_a = pm.Normal('alpha_a', mu=0, sd=10)
    a = pm.SkewNormal('a', mu=mu_a, sd=sd_a, alpha=alpha_a, observed=a)

    mu_b = pm.Normal('mu_b', mu=0, sd=10)
    sd_b = pm.HalfNormal('sd_b', 10)
    alpha_b = pm.Normal('alpha_b', mu=0, sd=10)
    b = pm.SkewNormal('b', mu=mu_b, sd=sd_b, alpha=alpha_b, observed=b)

    trace = pm.sample(1000)

请务必使用最新版本的 PyMC3，因为之前的版本不包含 SkewNormal 发行版。

更新：

鉴于你改变了你的问题：

如果 a 和 b 是独立的随机变量并且两者都服从正态分布，那么它们的总和将服从正态分布。

a ~ N(mu_a, sd_a²)

b ~ N(mu_b, sd_b²)

a+b ~ N(mu_a+mu_b, sd_a²+sd_b²)

也就是说，您将它们的均值相加，然后将它们的方差（而不是它们的标准差）相加。你不需要使用 PyMC3。

如果您仍想使用 PyMC3（可能是您的分布不是高斯分布，并且您不知道如何分析计算它们的总和）。您可以从您的 a 和 b 分布中生成合成数据，然后使用 PyMC3 来估计参数，类似于：

with pm.Model() as model:
    mu = pm.Normal('mu', mu=0, sd=10)
    sd = pm.HalfNormal('sd', 10)
    ab = pm.Normal('ab', mu=mu, sd=sd, observed=a+b)
    trace = pm.sample(1000)