Julia 中似然的高效网格近似

Question

在给定一些数据的情况下近似模型最大似然的一种简单方法是网格近似。例如，在 R 中，我们可以生成一个参数值网格，然后在给定一些数据的情况下评估每个值的可能性（示例来自 McElreath 的 Statistical Rethinking）：

p_grid <- seq(from=0, to=1, length.out=1000)
likelihood <- dbinom(6, size=9, prob=p_grid)

此处，likelihood 是一个包含 1000 个值的数组，我认为这是获取此类数组的有效方法。

我是 Julia 的新手（不太擅长 R），所以我的方法与上述相同，依赖于理解语法：

using Distributions

p_grid = collect(LinRange(0, 1, 1000))
likelihood = [pdf(Binomial(n9=, p=p), 6) for p in p_grid]

这不仅笨重而且似乎效率低下，因为新的二项式被构造了 1000 次。是否有更好的（也许是矢量化的）方法来完成相同的任务？

Answer 1

在 R 或 Python 等语言中，人们经常使用术语“矢量化”来表示“避免语言中的 for 循环”。我说“在语言中”是因为仍然有 for 循环，只是它们现在在 C 中而不是 R/Python.

在朱莉娅，没有什么可担心的。有时您仍然会听到“矢量化”，但 Julia 人员倾向于在硬件矢量化的原始意义上使用它。更多关于 here.

至于你的代码，我觉得还可以。可以肯定的是，让我们进行基准测试！

julia> using BenchmarkTools

julia> @btime [pdf(Binomial(9, p), 6) for p in $p_grid]
  111.352 μs (1 allocation: 7.94 KiB)

另一种写法是使用 map:

julia> @btime map($p_grid) do p
           pdf(Binomial(9, p), 6)
           end;
  111.623 μs (1 allocation: 7.94 KiB)

要检查构建开销，您可以 lower-level 调用 StatsFuns，就像这样

julia> using StatsFuns

julia> @btime map($p_grid) do p
           binompdf(9, p, 6)
           end;
  109.809 μs (1 allocation: 7.94 KiB)

看起来有一些差异，但差异很小，可能占总成本的 2% 左右。

Efficient grid approximation of likelihood in Julia