机器算法验证 - 计算对数概率。高维狄利克雷分布 - 吾爱随笔录

计算对数概率。高维狄利克雷分布

机器算法验证狄利克雷分布

2022-04-08 05:51:43

我有兴趣计算从 Dirichlet 分布中提取的数据的对数概率。特别是，我对在高维（可能是 1000 维或更多维）中稳定的计算感兴趣。

我尝试使用 Tom Minka 的 fastfit 工具箱中的“dirichlet_logProb.m”函数：http ://research.microsoft.com/en-us/um/people/minka/software/fastfit/ 。然而，在高维度上，这似乎给我提供的每个数据向量提供了意想不到的大正密度值。

当然，当分布急剧达到峰值时，我们可以在支持的某些点处具有任意大（大于 1）的密度值。例如，考虑一个方差非常小的零均值正态分布......它会给 x 提供非常大的密度 (>>1) 接近于零，但对于 x 的密度非常小 (<<1) 多于几个远离标准差。

任何分布都应该如此：我们应该总是能够找到具有大密度和小密度的点，因为密度需要在支持上积分为一个。但是，在高维中，我无法以这种方式计算狄利克雷对数概率计算……我尝试的每个点的密度都大于 1。

这是我的实验草图：我考虑一个浓度参数<1的对称狄利克雷。这应该更喜欢稀疏的数据，其中只有一些维度具有显着的概率质量，而所有其他维度都接近于零。然后我看三个可能的绘图，我在这里使用 K=4 维度来说明。

极其稀疏：x = [1 0 0 0]
中途制服：x = [0.5 0.5 0 0]
均匀：x = [0.25 0.25 0.25 0.25]

给定浓度参数 lambda < 1，我们应该期望 1 的密度比 2 大得多，而 2 的密度比 3 大得多。由于均匀性基本上与稀疏性相反，我们也应该期望这个结果可以忽略不计（小于比一）密度......它对应于负对数概率。

这是计算这三个结果的对数概率的结果，跨维度 K 的许多值，浓度 = 0.01。请记住，这些是对数概率，因此正值意味着密度 > 1，负值意味着密度 < 1。

        [1 0...0] [统一 0..0] [统一]
K= 4 9.18e+01 5.59e+01 -9.71e+00
K= 10 2.77e+02 1.43e+02 -2.09e+01
K= 50 1.52e+03 7.85e+02 -3.58e+01
K= 100 3.07e+03 1.64e+03 -4.04e+00
K= 500 1.55e+04 8.98e+03 7.80e+02
K=1000 3.11e+04 1.87e+04 2.25e+03

您可以使用以下代码重现此表：http: //gist.github.com/3699842

当然，问题在于，当 K=500 和 K=1000时，所有结果 1、2、3 的对数概率都非常正（密度 >> 1）。这对我来说似乎很麻烦，因为我预计最不可能的结果（制服）的密度小于一……我是不是弄错了？

有人可以提出答案吗？这是一个数字问题吗？还是我对狄利克雷有误解？我的主要怀疑是支持区域有点有趣（单纯形，而不是实数）。

1个回答

Dirichlet 分布的 pdf 定义为

f (θ; α) = B^{- 1} \prod_{i = 1}^{K} θ_{i}^{α_{i} - 1}

$f(\theta; \alpha) = B^{-1} \prod_{i=1}^K \theta_i^{\alpha_i - 1}$ 在哪里

B (α)

$B(\alpha)$ 是广义的 Beta 函数。请注意，如果有

θ_{i}

$\theta_i$ 为 0，则整个乘积为零。换句话说，狄利克雷分布的支持度超过了向量

θ

$\theta$ 其中每个

θ_{i} \in (0, 1)

$\theta_i \in (0, 1)$ 和

\sum_{i = 1}^{K} θ_{i} = 1

$\sum_{i=1}^K \theta_i = 1$ . 我不熟悉 Minka 的工具包，但是包含 0 的数据肯定会出现问题。

至于统一列，我相信这些值是正确的。这是我用来测试的python代码：

import math

def lbeta(alpha):
    return sum(math.lgamma(a) for a in alpha) - math.lgamma(sum(alpha))

def ldirichlet_pdf(alpha, theta):
    kernel = sum((a - 1) * math.log(t) for a, t in zip(alpha, theta))
    return kernel - lbeta(alpha)

for k in [4, 10, 50, 100, 500, 1000]:
    print ldirichlet_pdf([.01] * k, [1.0 / k] * k)

运行此脚本会产生输出：

4 -9.71111566837
10 -20.946493708
50 -35.7564901905
100 -4.03613939138
500 779.669123528
1000 2251.99967563

现在让我们从我们的狄利克雷分布中生成一些更可能的向量，K=1000。代码非常简单：

def sample_dirichlet(alpha):
    gammas = [random.gammavariate(a, 1) for a in alpha]
    norm = sum(gammas)
    return [g / norm for g in gammas]

现在如果我们将这个函数与我们之前的ldirichlet_pdf函数结合使用，我们会看到对于 K=1000，2e+3 是一个相对较小的密度。例如以下代码的结果：

alpha = [.01] * 1000
ldirichlet_pdf(alpha, sample_dirichlet(alpha))

产生 9.4e+4 和 1e+5 之间的值。

这里的关键见解是要意识到您不需要小于 1 的值才能使积分计算为 1。一个简单的例子是 $\int_0^1 1 dx = 1$ . 碰巧的是，对于 K=1000 和浓度 0.01 的对称 Dirichlet，pdf 到处都大于 1，但整个支持的积分仍然是 1。在更高的维度中，你需要有很多更小的浓度以获得均匀的负对数 pdf 例如，浓度为 0.0001 且 K=1000，均匀向量的对数 pdf 约为 -2.3e+3。

其它你可能感兴趣的问题

上一篇K-means Mahalanobis 与欧几里得距离下一篇为什么只考虑 1 个随机截距参数的 k-1 截距是有效的？