在Haskell中实现Metropolis算法

Question

我正在尝试编写基于Haskell中的Metropolis算法的MCMC程序,我遇到了从概率分布(生成伪随机数)和构造程序的问题.现在,我很高兴使用带有硬编码种子的生成器,而不是处理处理IO的复杂性.

似乎我应该使用状态monad来跟踪随机生成器状态,先前的马尔可夫链状态,卡方值和算法的每个步骤之间的接受计数,然后最终收集所有马尔可夫链状态和最后接受计数.这是最好的/惯用的方式吗？如果是这样,程序的布局应该是什么(即提议函数的类型签名和大都市步骤函数等).

我已经看到一些处理随机数的示例程序,其中从某种概率monad生成特定长度的随机数列表,然后通过一些简单的函数来执行计算.如果可能的话,我真的想避免这种内部形式的程序.

编辑:暂时删除WIP代码.

Answer 1

以下是关于编写惯用Haskell的一些反馈.

• 除非你正在写一元代码,使用do纯函数(即constructMuTable,metropolis)是非常不地道的Haskell.

  代替

  foo = do
      let x = ...
          y = ...
          z = ...
      bar x y z
  

  写吧

  foo =
      let x = ...
          y = ...
          z = ...
      in bar x y z
  

  或使用where而不是let ... in ....

• ETA-减少.在一些地方(zVec,muVec,sigmaVec中main)你写(\x -> f x).这仅相当于f,模_|_,seq等等.

• 用Data.Vector.Unboxed.你有很多V.Vector Double,哪些商店装箱Doubles,效率低下.对于原始类型Double,使用未装箱的向量(可能)使用更少的内存更快的代码.

• (!!)尽可能避免索引列表.使用Data.Vector,而不是作为V.!是O(1)同时(!!)的O(n).

• 看起来您可以在State这里使用monad来清理代码.然而,就目前的势在必行形式而言,我很难看到转型.

也许您可以尝试应用我给出的一些建议并简化一些大型密集函数,然后对您的算法的更高级别反馈将变得更加明显.