r - 似然积太小 - R 只给出 0

Question

我试图在 R 中计算这个后验分布。问题是分子，它是一堆 dbern(p_i, y_i) < 1 的乘积，太小了。 (我的 n 大约是 1500)。因此，R 吐出 0，并且所有\theta 的后验值也是 0。

澄清一下，每个 y_i 都有自己的 p_i，这些 p_i 一起构成了一个包含 n 个元素的向量，用于 n y。每个 theta 都有自己的 p_i 的 n 元素向量。



一个可重复的例子(分子的)

p <- sample(seq(0.001,0.999,by=0.01), 1500, replace=T)
y <- sample(c(0,1), 1500, replace=T)
dbern(y, p) # 1500-element vector, each element < 1
prod(dbern(y, p)) # produces 0
exp(sum(log(dbern(y, p)))) # produces 0

编辑(上下文):我正在做贝叶斯变化点分析(jstor.org/stable/25791783 - Western and Kleykamp 2004)。与论文中的连续 y 不同，我的 y 是二进制的，因此我使用了 Albert 和 Chib (1993) 中的数据增强方法。使用该方法，y 的似然为伯努利，其中 p = cdf-normal(x'B)。

那么 p 如何取决于 theta 呢？这是因为 theta 是变化点。其中一个 x 是时间虚拟变量——例如，如果 theta=10，那么对于第 10 天之后的所有观察值，时间虚拟变量 = 1，对于第 10 天之前的所有观察值，时间虚拟变量 = 0。

因此，p 取决于 x，x 取决于 theta——因此，p 取决于 theta。

我需要上述数量，因为是 Gibbs 采样中 theta 的完整条件。

Answer 1

解决此类精度问题的一种方法是在日志空间中工作。然而，这引入了分母的对数和积，这通常可能很痛苦。

如果您出于优化目的计算后验值，请注意您可能可以完全删除分母:您不需要归一化即可找到 argmax .