flat-mcmc

commit 46840992f48ad098199b914da1a67838b1bf1a25
parent d38dbddfabddbdbe8d915fde445c45e4f409ef8f
Author: Jared Tobin <jared@jtobin.ca>
Date:   Sun,  3 Apr 2016 19:56:20 +0700

Bug fixes.

Diffstat:
M
lib/Numeric/MCMC/Flat.hs
 | 
12
++++++++++--
M
test/Rosenbrock.hs
 | 
19
+++++++++----------

2 files changed, 19 insertions(+), 12 deletions(-)
diff --git a/lib/Numeric/MCMC/Flat.hs b/lib/Numeric/MCMC/Flat.hs
@@ -4,6 +4,8 @@
 module Numeric.MCMC.Flat (
     mcmc
   , flat
+  , Particle
+  , Ensemble
 
   , module Sampling.Types
   , Chain
@@ -33,7 +35,13 @@ data Chain = Chain {
   }
 
 instance Show Chain where
-  show Chain {..} = show chainPosition -- FIXME better?
+  show Chain {..} =
+      init
+    . filter (`notElem` "[]")
+    . unlines
+    . V.toList
+    . V.map show
+    $ chainPosition
 
 type Particle = Vector Double
 
@@ -56,7 +64,7 @@ move :: Target Particle -> Particle -> Particle -> Double -> Double -> Particle
 move target p0 p1 z zc =
   let proposal = stretch p0 p1 z
       pAccept  = acceptProb target p0 proposal z
-  in  if   zc <= min 0 pAccept
+  in  if   zc <= min 1 (exp pAccept)
       then proposal
       else p0
 
diff --git a/test/Rosenbrock.hs b/test/Rosenbrock.hs
@@ -10,15 +10,14 @@ rosenbrock :: Vector Double -> Double
 rosenbrock xs = negate (5  *(x1 - x0 ^ 2) ^ 2 + 0.05 * (1 - x0) ^ 2) where
   [x0, x1] = V.toList xs
 
+ensemble :: Ensemble
+ensemble = V.fromList [
+    V.fromList [negate 1.0, negate 1.0]
+  , V.fromList [negate 1.0, 1.0]
+  , V.fromList [1.0, negate 1.0]
+  , V.fromList [1.0, 1.0]
+  ]
+
 main :: IO ()
-main = withSystemRandom . asGenIO $ mcmc 10 ensemble rosenbrock where
-  ensemble = V.fromList [
-      V.fromList [negate 1.0 :: Double, negate 1.0]
-    , V.fromList [negate 1.0, 1.0]
-    , V.fromList [1.0, negate 1.0]
-    , V.fromList [1.0, 1.0]
-    ]
+main = withSystemRandom . asGenIO $ mcmc 25000 ensemble rosenbrock
 
--- ensemble = V.fromList [ V.fromList [negate 1.0 :: Double, negate 1.0], V.fromList [negate 1.0, 1.0], V.fromList [1.0, negate 1.0], V.fromList [1.0, 1.0]]
---
--- rosenbrock xs = negate (5  *(x1 - x0 ^ 2) ^ 2 + 0.05 * (1 - x0) ^ 2) where [x0, x1] = V.toList xs