{-# LANGUAGE OverloadedStrings #-}

-- |
-- Module      :  Mcmc.Mcmc
-- Description :  Mcmc helpers
-- Copyright   :  (c) Dominik Schrempf, 2020
-- License     :  GPL-3.0-or-later
--
-- Maintainer  :  dominik.schrempf@gmail.com
-- Stability   :  unstable
-- Portability :  portable
--
-- Creation date: Fri May 29 10:19:45 2020.
--
-- Functions to work with the 'Mcmc' state transformer.
module Mcmc.Mcmc
  ( Mcmc,
    mcmcOutB,
    mcmcOutS,
    mcmcWarnB,
    mcmcWarnS,
    mcmcInfoB,
    mcmcInfoS,
    mcmcDebugB,
    mcmcDebugS,
    mcmcAutotune,
    mcmcClean,
    mcmcResetA,
    mcmcSummarizeCycle,
    mcmcReport,
    mcmcMonitorExec,
    mcmcRun,
  )
where

import Control.Monad
import Control.Monad.IO.Class
import Control.Monad.Trans.State
import Data.Aeson
import qualified Data.ByteString.Lazy.Char8 as BL
import Data.Maybe
import Data.Time.Clock
import Data.Time.Format
import Mcmc.Item
import Mcmc.Monitor
import Mcmc.Monitor.Time
import Mcmc.Proposal
import Mcmc.Save
import Mcmc.Status hiding (debug)
import Mcmc.Verbosity
import Numeric.Log
import System.Directory
import System.IO
import Prelude hiding (cycle)

-- | An Mcmc state transformer; usually fiddling around with this type is not
-- required, but it is used by the different inference algorithms.
type Mcmc a = StateT (Status a) IO

msgPrepare :: Char -> BL.ByteString -> BL.ByteString
msgPrepare :: Char -> ByteString -> ByteString
msgPrepare Char
c ByteString
t = Char -> ByteString -> ByteString
BL.cons Char
c (ByteString -> ByteString) -> ByteString -> ByteString
forall a b. (a -> b) -> a -> b
$ ByteString
": " ByteString -> ByteString -> ByteString
forall a. Semigroup a => a -> a -> a
<> ByteString
t

-- | Write to standard output and log file.
mcmcOutB :: BL.ByteString -> Mcmc a ()
mcmcOutB :: ByteString -> Mcmc a ()
mcmcOutB ByteString
msg = do
  Handle
h <- Handle -> Maybe Handle -> Handle
forall a. a -> Maybe a -> a
fromMaybe ([Char] -> Handle
forall a. HasCallStack => [Char] -> a
error [Char]
"mcmcOut: Log handle is missing.") (Maybe Handle -> Handle)
-> StateT (Status a) IO (Maybe Handle)
-> StateT (Status a) IO Handle
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> (Status a -> Maybe Handle) -> StateT (Status a) IO (Maybe Handle)
forall (m :: * -> *) s a. Monad m => (s -> a) -> StateT s m a
gets Status a -> Maybe Handle
forall a. Status a -> Maybe Handle
logHandle
  IO () -> Mcmc a ()
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO () -> Mcmc a ()) -> IO () -> Mcmc a ()
forall a b. (a -> b) -> a -> b
$ ByteString -> IO ()
BL.putStrLn ByteString
msg IO () -> IO () -> IO ()
forall (m :: * -> *) a b. Monad m => m a -> m b -> m b
>> Handle -> ByteString -> IO ()
BL.hPutStrLn Handle
h ByteString
msg

-- | Write to standard output and log file.
mcmcOutS :: String -> Mcmc a ()
mcmcOutS :: [Char] -> Mcmc a ()
mcmcOutS = ByteString -> Mcmc a ()
forall a. ByteString -> Mcmc a ()
mcmcOutB (ByteString -> Mcmc a ())
-> ([Char] -> ByteString) -> [Char] -> Mcmc a ()
forall b c a. (b -> c) -> (a -> b) -> a -> c
. [Char] -> ByteString
BL.pack

-- Perform warning action.
mcmcWarnA :: Mcmc a () -> Mcmc a ()
mcmcWarnA :: Mcmc a () -> Mcmc a ()
mcmcWarnA Mcmc a ()
a = (Status a -> Verbosity) -> StateT (Status a) IO Verbosity
forall (m :: * -> *) s a. Monad m => (s -> a) -> StateT s m a
gets Status a -> Verbosity
forall a. Status a -> Verbosity
verbosity StateT (Status a) IO Verbosity
-> (Verbosity -> Mcmc a ()) -> Mcmc a ()
forall (m :: * -> *) a b. Monad m => m a -> (a -> m b) -> m b
>>= \Verbosity
v -> Verbosity -> Mcmc a () -> Mcmc a ()
forall (m :: * -> *). Applicative m => Verbosity -> m () -> m ()
info Verbosity
v Mcmc a ()
a

-- | Print warning message.
mcmcWarnB :: BL.ByteString -> Mcmc a ()
mcmcWarnB :: ByteString -> Mcmc a ()
mcmcWarnB = Mcmc a () -> Mcmc a ()
forall a. Mcmc a () -> Mcmc a ()
mcmcWarnA (Mcmc a () -> Mcmc a ())
-> (ByteString -> Mcmc a ()) -> ByteString -> Mcmc a ()
forall b c a. (b -> c) -> (a -> b) -> a -> c
. ByteString -> Mcmc a ()
forall a. ByteString -> Mcmc a ()
mcmcOutB (ByteString -> Mcmc a ())
-> (ByteString -> ByteString) -> ByteString -> Mcmc a ()
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Char -> ByteString -> ByteString
msgPrepare Char
'W'

-- | Print warning message.
mcmcWarnS :: String -> Mcmc a ()
mcmcWarnS :: [Char] -> Mcmc a ()
mcmcWarnS = ByteString -> Mcmc a ()
forall a. ByteString -> Mcmc a ()
mcmcWarnB (ByteString -> Mcmc a ())
-> ([Char] -> ByteString) -> [Char] -> Mcmc a ()
forall b c a. (b -> c) -> (a -> b) -> a -> c
. [Char] -> ByteString
BL.pack

-- Perform info action.
mcmcInfoA :: Mcmc a () -> Mcmc a ()
mcmcInfoA :: Mcmc a () -> Mcmc a ()
mcmcInfoA Mcmc a ()
a = (Status a -> Verbosity) -> StateT (Status a) IO Verbosity
forall (m :: * -> *) s a. Monad m => (s -> a) -> StateT s m a
gets Status a -> Verbosity
forall a. Status a -> Verbosity
verbosity StateT (Status a) IO Verbosity
-> (Verbosity -> Mcmc a ()) -> Mcmc a ()
forall (m :: * -> *) a b. Monad m => m a -> (a -> m b) -> m b
>>= \Verbosity
v -> Verbosity -> Mcmc a () -> Mcmc a ()
forall (m :: * -> *). Applicative m => Verbosity -> m () -> m ()
info Verbosity
v Mcmc a ()
a

-- | Print info message.
mcmcInfoB :: BL.ByteString -> Mcmc a ()
mcmcInfoB :: ByteString -> Mcmc a ()
mcmcInfoB = Mcmc a () -> Mcmc a ()
forall a. Mcmc a () -> Mcmc a ()
mcmcInfoA (Mcmc a () -> Mcmc a ())
-> (ByteString -> Mcmc a ()) -> ByteString -> Mcmc a ()
forall b c a. (b -> c) -> (a -> b) -> a -> c
. ByteString -> Mcmc a ()
forall a. ByteString -> Mcmc a ()
mcmcOutB (ByteString -> Mcmc a ())
-> (ByteString -> ByteString) -> ByteString -> Mcmc a ()
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Char -> ByteString -> ByteString
msgPrepare Char
'I'

-- | Print info message.
mcmcInfoS :: String -> Mcmc a ()
mcmcInfoS :: [Char] -> Mcmc a ()
mcmcInfoS = ByteString -> Mcmc a ()
forall a. ByteString -> Mcmc a ()
mcmcInfoB (ByteString -> Mcmc a ())
-> ([Char] -> ByteString) -> [Char] -> Mcmc a ()
forall b c a. (b -> c) -> (a -> b) -> a -> c
. [Char] -> ByteString
BL.pack

-- Perform debug action.
mcmcDebugA :: Mcmc a () -> Mcmc a ()
mcmcDebugA :: Mcmc a () -> Mcmc a ()
mcmcDebugA Mcmc a ()
a = (Status a -> Verbosity) -> StateT (Status a) IO Verbosity
forall (m :: * -> *) s a. Monad m => (s -> a) -> StateT s m a
gets Status a -> Verbosity
forall a. Status a -> Verbosity
verbosity StateT (Status a) IO Verbosity
-> (Verbosity -> Mcmc a ()) -> Mcmc a ()
forall (m :: * -> *) a b. Monad m => m a -> (a -> m b) -> m b
>>= \Verbosity
v -> Verbosity -> Mcmc a () -> Mcmc a ()
forall (m :: * -> *). Applicative m => Verbosity -> m () -> m ()
debug Verbosity
v Mcmc a ()
a

-- | Print debug message.
mcmcDebugB :: BL.ByteString -> Mcmc a ()
mcmcDebugB :: ByteString -> Mcmc a ()
mcmcDebugB = Mcmc a () -> Mcmc a ()
forall a. Mcmc a () -> Mcmc a ()
mcmcDebugA (Mcmc a () -> Mcmc a ())
-> (ByteString -> Mcmc a ()) -> ByteString -> Mcmc a ()
forall b c a. (b -> c) -> (a -> b) -> a -> c
. ByteString -> Mcmc a ()
forall a. ByteString -> Mcmc a ()
mcmcOutB (ByteString -> Mcmc a ())
-> (ByteString -> ByteString) -> ByteString -> Mcmc a ()
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Char -> ByteString -> ByteString
msgPrepare Char
'D'

-- | Print debug message.
mcmcDebugS :: String -> Mcmc a ()
mcmcDebugS :: [Char] -> Mcmc a ()
mcmcDebugS = ByteString -> Mcmc a ()
forall a. ByteString -> Mcmc a ()
mcmcDebugB (ByteString -> Mcmc a ())
-> ([Char] -> ByteString) -> [Char] -> Mcmc a ()
forall b c a. (b -> c) -> (a -> b) -> a -> c
. [Char] -> ByteString
BL.pack

-- | Auto tune the 'Proposal's in the 'Cycle' of the chain. Reset acceptance counts.
-- See 'autotuneCycle'.
mcmcAutotune :: Mcmc a ()
mcmcAutotune :: Mcmc a ()
mcmcAutotune = do
  ByteString -> Mcmc a ()
forall a. ByteString -> Mcmc a ()
mcmcDebugB ByteString
"Auto tune."
  Status a
s <- StateT (Status a) IO (Status a)
forall (m :: * -> *) s. Monad m => StateT s m s
get
  let a :: Acceptance (Proposal a)
a = Status a -> Acceptance (Proposal a)
forall a. Status a -> Acceptance (Proposal a)
acceptance Status a
s
      c :: Cycle a
c = Status a -> Cycle a
forall a. Status a -> Cycle a
cycle Status a
s
      c' :: Cycle a
c' = Acceptance (Proposal a) -> Cycle a -> Cycle a
forall a. Acceptance (Proposal a) -> Cycle a -> Cycle a
autotuneCycle Acceptance (Proposal a)
a Cycle a
c
  Status a -> Mcmc a ()
forall (m :: * -> *) s. Monad m => s -> StateT s m ()
put (Status a -> Mcmc a ()) -> Status a -> Mcmc a ()
forall a b. (a -> b) -> a -> b
$ Status a
s {cycle :: Cycle a
cycle = Cycle a
c'}

-- | Clean the state.
mcmcClean :: Mcmc a ()
mcmcClean :: Mcmc a ()
mcmcClean = do
  Status a
s <- StateT (Status a) IO (Status a)
forall (m :: * -> *) s. Monad m => StateT s m s
get
  let cl :: Maybe (Cleaner a)
cl = Status a -> Maybe (Cleaner a)
forall a. Status a -> Maybe (Cleaner a)
cleaner Status a
s
      i :: Int
i = Status a -> Int
forall a. Status a -> Int
iteration Status a
s
  case Maybe (Cleaner a)
cl of
    Just (Cleaner Int
n a -> a
f) | Int
i Int -> Int -> Int
forall a. Integral a => a -> a -> a
`mod` Int
n Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== Int
0 -> do
      ByteString -> Mcmc a ()
forall a. ByteString -> Mcmc a ()
mcmcDebugB ByteString
"Clean state."
      let (Item a
st Log Double
pr Log Double
lh) = Status a -> Item a
forall a. Status a -> Item a
item Status a
s
      [Char] -> Mcmc a ()
forall a. [Char] -> Mcmc a ()
mcmcDebugS ([Char] -> Mcmc a ()) -> [Char] -> Mcmc a ()
forall a b. (a -> b) -> a -> b
$
        [Char]
"Old log prior and log likelihood: " [Char] -> [Char] -> [Char]
forall a. [a] -> [a] -> [a]
++ Double -> [Char]
forall a. Show a => a -> [Char]
show (Log Double -> Double
forall a. Log a -> a
ln Log Double
pr) [Char] -> [Char] -> [Char]
forall a. [a] -> [a] -> [a]
++ [Char]
", " [Char] -> [Char] -> [Char]
forall a. [a] -> [a] -> [a]
++ Double -> [Char]
forall a. Show a => a -> [Char]
show (Log Double -> Double
forall a. Log a -> a
ln Log Double
lh) [Char] -> [Char] -> [Char]
forall a. [a] -> [a] -> [a]
++ [Char]
"."
      let prF :: a -> Log Double
prF = Status a -> a -> Log Double
forall a. Status a -> a -> Log Double
priorF Status a
s
          lhF :: a -> Log Double
lhF = Status a -> a -> Log Double
forall a. Status a -> a -> Log Double
likelihoodF Status a
s
          st' :: a
st' = a -> a
f a
st
          pr' :: Log Double
pr' = a -> Log Double
prF a
st'
          lh' :: Log Double
lh' = a -> Log Double
lhF a
st'
      [Char] -> Mcmc a ()
forall a. [Char] -> Mcmc a ()
mcmcDebugS ([Char] -> Mcmc a ()) -> [Char] -> Mcmc a ()
forall a b. (a -> b) -> a -> b
$
        [Char]
"New log prior and log likelihood: " [Char] -> [Char] -> [Char]
forall a. [a] -> [a] -> [a]
++ Double -> [Char]
forall a. Show a => a -> [Char]
show (Log Double -> Double
forall a. Log a -> a
ln Log Double
pr') [Char] -> [Char] -> [Char]
forall a. [a] -> [a] -> [a]
++ [Char]
", " [Char] -> [Char] -> [Char]
forall a. [a] -> [a] -> [a]
++ Double -> [Char]
forall a. Show a => a -> [Char]
show (Log Double -> Double
forall a. Log a -> a
ln Log Double
lh') [Char] -> [Char] -> [Char]
forall a. [a] -> [a] -> [a]
++ [Char]
"."
      let dLogPr :: Double
dLogPr = Double -> Double
forall a. Num a => a -> a
abs (Double -> Double) -> Double -> Double
forall a b. (a -> b) -> a -> b
$ Log Double -> Double
forall a. Log a -> a
ln Log Double
pr Double -> Double -> Double
forall a. Num a => a -> a -> a
- Log Double -> Double
forall a. Log a -> a
ln Log Double
pr'
          dLogLh :: Double
dLogLh = Double -> Double
forall a. Num a => a -> a
abs (Double -> Double) -> Double -> Double
forall a b. (a -> b) -> a -> b
$ Log Double -> Double
forall a. Log a -> a
ln Log Double
lh Double -> Double -> Double
forall a. Num a => a -> a -> a
- Log Double -> Double
forall a. Log a -> a
ln Log Double
lh'
      Bool -> Mcmc a () -> Mcmc a ()
forall (f :: * -> *). Applicative f => Bool -> f () -> f ()
when
        (Double
dLogPr Double -> Double -> Bool
forall a. Ord a => a -> a -> Bool
> Double
0.01)
        ([Char] -> Mcmc a ()
forall a. [Char] -> Mcmc a ()
mcmcWarnS ([Char] -> Mcmc a ()) -> [Char] -> Mcmc a ()
forall a b. (a -> b) -> a -> b
$ [Char]
"Log of old and new prior differ by " [Char] -> [Char] -> [Char]
forall a. [a] -> [a] -> [a]
++ Double -> [Char]
forall a. Show a => a -> [Char]
show Double
dLogPr [Char] -> [Char] -> [Char]
forall a. [a] -> [a] -> [a]
++ [Char]
".")
      Bool -> Mcmc a () -> Mcmc a ()
forall (f :: * -> *). Applicative f => Bool -> f () -> f ()
when
        (Double
dLogPr Double -> Double -> Bool
forall a. Ord a => a -> a -> Bool
> Double
0.01)
        ([Char] -> Mcmc a ()
forall a. [Char] -> Mcmc a ()
mcmcWarnS ([Char] -> Mcmc a ()) -> [Char] -> Mcmc a ()
forall a b. (a -> b) -> a -> b
$ [Char]
"Log of old and new likelihood differ by " [Char] -> [Char] -> [Char]
forall a. [a] -> [a] -> [a]
++ Double -> [Char]
forall a. Show a => a -> [Char]
show Double
dLogLh [Char] -> [Char] -> [Char]
forall a. [a] -> [a] -> [a]
++ [Char]
".")
      Status a -> Mcmc a ()
forall (m :: * -> *) s. Monad m => s -> StateT s m ()
put (Status a -> Mcmc a ()) -> Status a -> Mcmc a ()
forall a b. (a -> b) -> a -> b
$ Status a
s {item :: Item a
item = a -> Log Double -> Log Double -> Item a
forall a. a -> Log Double -> Log Double -> Item a
Item a
st' Log Double
pr' Log Double
lh'}
    Maybe (Cleaner a)
_ -> () -> Mcmc a ()
forall (m :: * -> *) a. Monad m => a -> m a
return ()

-- | Reset acceptance counts.
mcmcResetA :: Mcmc a ()
mcmcResetA :: Mcmc a ()
mcmcResetA = do
  ByteString -> Mcmc a ()
forall a. ByteString -> Mcmc a ()
mcmcDebugB ByteString
"Reset acceptance ratios."
  Status a
s <- StateT (Status a) IO (Status a)
forall (m :: * -> *) s. Monad m => StateT s m s
get
  let a :: Acceptance (Proposal a)
a = Status a -> Acceptance (Proposal a)
forall a. Status a -> Acceptance (Proposal a)
acceptance Status a
s
  Status a -> Mcmc a ()
forall (m :: * -> *) s. Monad m => s -> StateT s m ()
put (Status a -> Mcmc a ()) -> Status a -> Mcmc a ()
forall a b. (a -> b) -> a -> b
$ Status a
s {acceptance :: Acceptance (Proposal a)
acceptance = Acceptance (Proposal a) -> Acceptance (Proposal a)
forall k. Ord k => Acceptance k -> Acceptance k
resetA Acceptance (Proposal a)
a}

-- | Print short summary of 'Proposal's in 'Cycle'. See 'summarizeCycle'.
mcmcSummarizeCycle :: Mcmc a BL.ByteString
mcmcSummarizeCycle :: Mcmc a ByteString
mcmcSummarizeCycle = do
  Acceptance (Proposal a)
a <- (Status a -> Acceptance (Proposal a))
-> StateT (Status a) IO (Acceptance (Proposal a))
forall (m :: * -> *) s a. Monad m => (s -> a) -> StateT s m a
gets Status a -> Acceptance (Proposal a)
forall a. Status a -> Acceptance (Proposal a)
acceptance
  Cycle a
c <- (Status a -> Cycle a) -> StateT (Status a) IO (Cycle a)
forall (m :: * -> *) s a. Monad m => (s -> a) -> StateT s m a
gets Status a -> Cycle a
forall a. Status a -> Cycle a
cycle
  ByteString -> Mcmc a ByteString
forall (m :: * -> *) a. Monad m => a -> m a
return (ByteString -> Mcmc a ByteString)
-> ByteString -> Mcmc a ByteString
forall a b. (a -> b) -> a -> b
$ Acceptance (Proposal a) -> Cycle a -> ByteString
forall a. Acceptance (Proposal a) -> Cycle a -> ByteString
summarizeCycle Acceptance (Proposal a)
a Cycle a
c

fTime :: FormatTime t => t -> String
fTime :: t -> [Char]
fTime = TimeLocale -> [Char] -> t -> [Char]
forall t. FormatTime t => TimeLocale -> [Char] -> t -> [Char]
formatTime TimeLocale
defaultTimeLocale [Char]
"%B %-e, %Y, at %H:%M %P, %Z."

-- Open log file.
mcmcOpenLog :: Mcmc a ()
mcmcOpenLog :: Mcmc a ()
mcmcOpenLog = do
  Status a
s <- StateT (Status a) IO (Status a)
forall (m :: * -> *) s. Monad m => StateT s m s
get
  let lfn :: [Char]
lfn = Status a -> [Char]
forall a. Status a -> [Char]
name Status a
s [Char] -> [Char] -> [Char]
forall a. [a] -> [a] -> [a]
++ [Char]
".log"
      n :: Int
n = Status a -> Int
forall a. Status a -> Int
iteration Status a
s
      frc :: Bool
frc = Status a -> Bool
forall a. Status a -> Bool
forceOverwrite Status a
s
  Bool
fe <- IO Bool -> StateT (Status a) IO Bool
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO Bool -> StateT (Status a) IO Bool)
-> IO Bool -> StateT (Status a) IO Bool
forall a b. (a -> b) -> a -> b
$ [Char] -> IO Bool
doesFileExist [Char]
lfn
  Maybe Handle
mh <- IO (Maybe Handle) -> StateT (Status a) IO (Maybe Handle)
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO (Maybe Handle) -> StateT (Status a) IO (Maybe Handle))
-> IO (Maybe Handle) -> StateT (Status a) IO (Maybe Handle)
forall a b. (a -> b) -> a -> b
$ case Status a -> Verbosity
forall a. Status a -> Verbosity
verbosity Status a
s of
    Verbosity
Quiet -> Maybe Handle -> IO (Maybe Handle)
forall (m :: * -> *) a. Monad m => a -> m a
return Maybe Handle
forall a. Maybe a
Nothing
    Verbosity
_ ->
      Handle -> Maybe Handle
forall a. a -> Maybe a
Just (Handle -> Maybe Handle) -> IO Handle -> IO (Maybe Handle)
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> case (Bool
fe, Int
n, Bool
frc) of
        (Bool
False, Int
_, Bool
_) -> [Char] -> IOMode -> IO Handle
openFile [Char]
lfn IOMode
WriteMode
        (Bool
True, Int
0, Bool
True) -> [Char] -> IOMode -> IO Handle
openFile [Char]
lfn IOMode
WriteMode
        (Bool
True, Int
0, Bool
False) -> [Char] -> IO Handle
forall a. HasCallStack => [Char] -> a
error [Char]
"mcmcInit: Log file exists; use 'force' to overwrite output files."
        (Bool
True, Int
_, Bool
_) -> [Char] -> IOMode -> IO Handle
openFile [Char]
lfn IOMode
AppendMode
  Status a -> Mcmc a ()
forall (m :: * -> *) s. Monad m => s -> StateT s m ()
put Status a
s {logHandle :: Maybe Handle
logHandle = Maybe Handle
mh}
  [Char] -> Mcmc a ()
forall a. [Char] -> Mcmc a ()
mcmcDebugS ([Char] -> Mcmc a ()) -> [Char] -> Mcmc a ()
forall a b. (a -> b) -> a -> b
$ [Char]
"Log file name: " [Char] -> [Char] -> [Char]
forall a. [a] -> [a] -> [a]
++ [Char]
lfn [Char] -> [Char] -> [Char]
forall a. [a] -> [a] -> [a]
++ [Char]
"."
  ByteString -> Mcmc a ()
forall a. ByteString -> Mcmc a ()
mcmcDebugB ByteString
"Log file opened."

-- Set the total number of iterations, the current time and open the 'Monitor's
-- of the chain. See 'mOpen'.
mcmcInit :: Mcmc a ()
mcmcInit :: Mcmc a ()
mcmcInit = do
  Mcmc a ()
forall a. Mcmc a ()
mcmcOpenLog
  Status a
s <- StateT (Status a) IO (Status a)
forall (m :: * -> *) s. Monad m => StateT s m s
get
  -- Start time.
  UTCTime
t <- IO UTCTime -> StateT (Status a) IO UTCTime
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO IO UTCTime
getCurrentTime
  [Char] -> Mcmc a ()
forall a. [Char] -> Mcmc a ()
mcmcInfoS ([Char] -> Mcmc a ()) -> [Char] -> Mcmc a ()
forall a b. (a -> b) -> a -> b
$ [Char]
"Start time: " [Char] -> [Char] -> [Char]
forall a. Semigroup a => a -> a -> a
<> UTCTime -> [Char]
forall t. FormatTime t => t -> [Char]
fTime UTCTime
t
  -- Monitor.
  let m :: Monitor a
m = Status a -> Monitor a
forall a. Status a -> Monitor a
monitor Status a
s
      n :: Int
n = Status a -> Int
forall a. Status a -> Int
iteration Status a
s
      nm :: [Char]
nm = Status a -> [Char]
forall a. Status a -> [Char]
name Status a
s
      frc :: Bool
frc = Status a -> Bool
forall a. Status a -> Bool
forceOverwrite Status a
s
  Monitor a
m' <- if Int
n Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== Int
0 then IO (Monitor a) -> StateT (Status a) IO (Monitor a)
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO (Monitor a) -> StateT (Status a) IO (Monitor a))
-> IO (Monitor a) -> StateT (Status a) IO (Monitor a)
forall a b. (a -> b) -> a -> b
$ [Char] -> Bool -> Monitor a -> IO (Monitor a)
forall a. [Char] -> Bool -> Monitor a -> IO (Monitor a)
mOpen [Char]
nm Bool
frc Monitor a
m else IO (Monitor a) -> StateT (Status a) IO (Monitor a)
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO (Monitor a) -> StateT (Status a) IO (Monitor a))
-> IO (Monitor a) -> StateT (Status a) IO (Monitor a)
forall a b. (a -> b) -> a -> b
$ [Char] -> Monitor a -> IO (Monitor a)
forall a. [Char] -> Monitor a -> IO (Monitor a)
mAppend [Char]
nm Monitor a
m
  Status a -> Mcmc a ()
forall (m :: * -> *) s. Monad m => s -> StateT s m ()
put (Status a -> Mcmc a ()) -> Status a -> Mcmc a ()
forall a b. (a -> b) -> a -> b
$ Status a
s {monitor :: Monitor a
monitor = Monitor a
m', start :: Maybe (Int, UTCTime)
start = (Int, UTCTime) -> Maybe (Int, UTCTime)
forall a. a -> Maybe a
Just (Int
n, UTCTime
t)}

-- | Report what is going to be done.
mcmcReport :: ToJSON a => Mcmc a ()
mcmcReport :: Mcmc a ()
mcmcReport = do
  Status a
s <- StateT (Status a) IO (Status a)
forall (m :: * -> *) s. Monad m => StateT s m s
get
  let b :: Maybe Int
b = Status a -> Maybe Int
forall a. Status a -> Maybe Int
burnInIterations Status a
s
      t :: Maybe Int
t = Status a -> Maybe Int
forall a. Status a -> Maybe Int
autoTuningPeriod Status a
s
      n :: Int
n = Status a -> Int
forall a. Status a -> Int
iterations Status a
s
      c :: Maybe (Cleaner a)
c = Status a -> Maybe (Cleaner a)
forall a. Status a -> Maybe (Cleaner a)
cleaner Status a
s
  case Maybe Int
b of
    Just Int
b' -> [Char] -> Mcmc a ()
forall a. [Char] -> Mcmc a ()
mcmcInfoS ([Char] -> Mcmc a ()) -> [Char] -> Mcmc a ()
forall a b. (a -> b) -> a -> b
$ [Char]
"Burn in for " [Char] -> [Char] -> [Char]
forall a. Semigroup a => a -> a -> a
<> Int -> [Char]
forall a. Show a => a -> [Char]
show Int
b' [Char] -> [Char] -> [Char]
forall a. Semigroup a => a -> a -> a
<> [Char]
" iterations."
    Maybe Int
Nothing -> () -> Mcmc a ()
forall (m :: * -> *) a. Monad m => a -> m a
return ()
  case Maybe Int
t of
    Just Int
t' -> [Char] -> Mcmc a ()
forall a. [Char] -> Mcmc a ()
mcmcInfoS ([Char] -> Mcmc a ()) -> [Char] -> Mcmc a ()
forall a b. (a -> b) -> a -> b
$ [Char]
"Auto tune every " [Char] -> [Char] -> [Char]
forall a. Semigroup a => a -> a -> a
<> Int -> [Char]
forall a. Show a => a -> [Char]
show Int
t' [Char] -> [Char] -> [Char]
forall a. Semigroup a => a -> a -> a
<> [Char]
" iterations (during burn in only)."
    Maybe Int
Nothing -> () -> Mcmc a ()
forall (m :: * -> *) a. Monad m => a -> m a
return ()
  case Maybe (Cleaner a)
c of
    Just (Cleaner Int
c' a -> a
_) -> [Char] -> Mcmc a ()
forall a. [Char] -> Mcmc a ()
mcmcInfoS ([Char] -> Mcmc a ()) -> [Char] -> Mcmc a ()
forall a b. (a -> b) -> a -> b
$ [Char]
"Clean state every " [Char] -> [Char] -> [Char]
forall a. Semigroup a => a -> a -> a
<> Int -> [Char]
forall a. Show a => a -> [Char]
show Int
c' [Char] -> [Char] -> [Char]
forall a. Semigroup a => a -> a -> a
<> [Char]
" iterations."
    Maybe (Cleaner a)
Nothing -> () -> Mcmc a ()
forall (m :: * -> *) a. Monad m => a -> m a
return ()
  [Char] -> Mcmc a ()
forall a. [Char] -> Mcmc a ()
mcmcInfoS ([Char] -> Mcmc a ()) -> [Char] -> Mcmc a ()
forall a b. (a -> b) -> a -> b
$ [Char]
"Run chain for " [Char] -> [Char] -> [Char]
forall a. Semigroup a => a -> a -> a
<> Int -> [Char]
forall a. Show a => a -> [Char]
show Int
n [Char] -> [Char] -> [Char]
forall a. Semigroup a => a -> a -> a
<> [Char]
" iterations."
  ByteString -> Mcmc a ()
forall a. ByteString -> Mcmc a ()
mcmcInfoB ByteString
"Initial state."
  Mcmc a ()
forall a. ToJSON a => Mcmc a ()
mcmcMonitorExec

-- Save the status of an MCMC run. See 'saveStatus'.
mcmcSave :: ToJSON a => Mcmc a ()
mcmcSave :: Mcmc a ()
mcmcSave = do
  Status a
s <- StateT (Status a) IO (Status a)
forall (m :: * -> *) s. Monad m => StateT s m s
get
  case Status a -> Maybe Int
forall a. Status a -> Maybe Int
save Status a
s of
    Just Int
n -> do
      ByteString -> Mcmc a ()
forall a. ByteString -> Mcmc a ()
mcmcInfoB (ByteString -> Mcmc a ()) -> ByteString -> Mcmc a ()
forall a b. (a -> b) -> a -> b
$ ByteString
"Save Markov chain with trace of length " ByteString -> ByteString -> ByteString
forall a. Semigroup a => a -> a -> a
<> [Char] -> ByteString
BL.pack (Int -> [Char]
forall a. Show a => a -> [Char]
show Int
n) ByteString -> ByteString -> ByteString
forall a. Semigroup a => a -> a -> a
<> ByteString
"."
      ByteString -> Mcmc a ()
forall a. ByteString -> Mcmc a ()
mcmcInfoB ByteString
"For long traces, or complex objects, this may take a while."
      IO () -> Mcmc a ()
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO () -> Mcmc a ()) -> IO () -> Mcmc a ()
forall a b. (a -> b) -> a -> b
$ [Char] -> Status a -> IO ()
forall a. ToJSON a => [Char] -> Status a -> IO ()
saveStatus (Status a -> [Char]
forall a. Status a -> [Char]
name Status a
s [Char] -> [Char] -> [Char]
forall a. Semigroup a => a -> a -> a
<> [Char]
".mcmc") Status a
s
      ByteString -> Mcmc a ()
forall a. ByteString -> Mcmc a ()
mcmcInfoB ByteString
"Done saving Markov chain."
    Maybe Int
Nothing -> ByteString -> Mcmc a ()
forall a. ByteString -> Mcmc a ()
mcmcInfoB ByteString
"Do not save the Markov chain."

-- | Execute the 'Monitor's of the chain. See 'mExec'.
mcmcMonitorExec :: ToJSON a => Mcmc a ()
mcmcMonitorExec :: Mcmc a ()
mcmcMonitorExec = do
  Status a
s <- StateT (Status a) IO (Status a)
forall (m :: * -> *) s. Monad m => StateT s m s
get
  let i :: Int
i = Status a -> Int
forall a. Status a -> Int
iteration Status a
s
      j :: Int
j = Status a -> Int
forall a. Status a -> Int
iterations Status a
s Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Int -> Maybe Int -> Int
forall a. a -> Maybe a -> a
fromMaybe Int
0 (Status a -> Maybe Int
forall a. Status a -> Maybe Int
burnInIterations Status a
s)
      m :: Monitor a
m = Status a -> Monitor a
forall a. Status a -> Monitor a
monitor Status a
s
      (Int
ss, UTCTime
st) = (Int, UTCTime) -> Maybe (Int, UTCTime) -> (Int, UTCTime)
forall a. a -> Maybe a -> a
fromMaybe ([Char] -> (Int, UTCTime)
forall a. HasCallStack => [Char] -> a
error [Char]
"mcmcMonitorExec: Starting state and time not set.") (Status a -> Maybe (Int, UTCTime)
forall a. Status a -> Maybe (Int, UTCTime)
start Status a
s)
      tr :: Trace a
tr = Status a -> Trace a
forall a. Status a -> Trace a
trace Status a
s
      vb :: Verbosity
vb = Status a -> Verbosity
forall a. Status a -> Verbosity
verbosity Status a
s
  Maybe ByteString
mt <- IO (Maybe ByteString) -> StateT (Status a) IO (Maybe ByteString)
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO (Maybe ByteString) -> StateT (Status a) IO (Maybe ByteString))
-> IO (Maybe ByteString) -> StateT (Status a) IO (Maybe ByteString)
forall a b. (a -> b) -> a -> b
$ Verbosity
-> Int
-> Int
-> UTCTime
-> Trace a
-> Int
-> Monitor a
-> IO (Maybe ByteString)
forall a.
Verbosity
-> Int
-> Int
-> UTCTime
-> Trace a
-> Int
-> Monitor a
-> IO (Maybe ByteString)
mExec Verbosity
vb Int
i Int
ss UTCTime
st Trace a
tr Int
j Monitor a
m
  Maybe ByteString -> (ByteString -> Mcmc a ()) -> Mcmc a ()
forall (t :: * -> *) (m :: * -> *) a b.
(Foldable t, Monad m) =>
t a -> (a -> m b) -> m ()
forM_ Maybe ByteString
mt ByteString -> Mcmc a ()
forall a. ByteString -> Mcmc a ()
mcmcOutB

-- Close the 'Monitor's of the chain. See 'mClose'.
mcmcClose :: ToJSON a => Mcmc a ()
mcmcClose :: Mcmc a ()
mcmcClose = do
  Status a
s <- StateT (Status a) IO (Status a)
forall (m :: * -> *) s. Monad m => StateT s m s
get
  Mcmc a ByteString
forall a. Mcmc a ByteString
mcmcSummarizeCycle Mcmc a ByteString -> (ByteString -> Mcmc a ()) -> Mcmc a ()
forall (m :: * -> *) a b. Monad m => m a -> (a -> m b) -> m b
>>= ByteString -> Mcmc a ()
forall a. ByteString -> Mcmc a ()
mcmcInfoB
  ByteString -> Mcmc a ()
forall a. ByteString -> Mcmc a ()
mcmcInfoB ByteString
"Metropolis-Hastings sampler finished."
  let m :: Monitor a
m = Status a -> Monitor a
forall a. Status a -> Monitor a
monitor Status a
s
  Monitor a
m' <- IO (Monitor a) -> StateT (Status a) IO (Monitor a)
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO (Monitor a) -> StateT (Status a) IO (Monitor a))
-> IO (Monitor a) -> StateT (Status a) IO (Monitor a)
forall a b. (a -> b) -> a -> b
$ Monitor a -> IO (Monitor a)
forall a. Monitor a -> IO (Monitor a)
mClose Monitor a
m
  Status a -> Mcmc a ()
forall (m :: * -> *) s. Monad m => s -> StateT s m ()
put (Status a -> Mcmc a ()) -> Status a -> Mcmc a ()
forall a b. (a -> b) -> a -> b
$ Status a
s {monitor :: Monitor a
monitor = Monitor a
m'}
  Mcmc a ()
forall a. ToJSON a => Mcmc a ()
mcmcSave
  UTCTime
t <- IO UTCTime -> StateT (Status a) IO UTCTime
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO IO UTCTime
getCurrentTime
  let rt :: NominalDiffTime
rt = case Status a -> Maybe (Int, UTCTime)
forall a. Status a -> Maybe (Int, UTCTime)
start Status a
s of
        Maybe (Int, UTCTime)
Nothing -> [Char] -> NominalDiffTime
forall a. HasCallStack => [Char] -> a
error [Char]
"mcmcClose: Start time not set."
        Just (Int
_, UTCTime
st) -> UTCTime
t UTCTime -> UTCTime -> NominalDiffTime
`diffUTCTime` UTCTime
st
  ByteString -> Mcmc a ()
forall a. ByteString -> Mcmc a ()
mcmcInfoB (ByteString -> Mcmc a ()) -> ByteString -> Mcmc a ()
forall a b. (a -> b) -> a -> b
$ ByteString
"Wall clock run time: " ByteString -> ByteString -> ByteString
forall a. Semigroup a => a -> a -> a
<> NominalDiffTime -> ByteString
renderDuration NominalDiffTime
rt ByteString -> ByteString -> ByteString
forall a. Semigroup a => a -> a -> a
<> ByteString
"."
  [Char] -> Mcmc a ()
forall a. [Char] -> Mcmc a ()
mcmcInfoS ([Char] -> Mcmc a ()) -> [Char] -> Mcmc a ()
forall a b. (a -> b) -> a -> b
$ [Char]
"End time: " [Char] -> [Char] -> [Char]
forall a. Semigroup a => a -> a -> a
<> UTCTime -> [Char]
forall t. FormatTime t => t -> [Char]
fTime UTCTime
t
  case Status a -> Maybe Handle
forall a. Status a -> Maybe Handle
logHandle Status a
s of
    Just Handle
h -> IO () -> Mcmc a ()
forall (m :: * -> *) a. MonadIO m => IO a -> m a
liftIO (IO () -> Mcmc a ()) -> IO () -> Mcmc a ()
forall a b. (a -> b) -> a -> b
$ Handle -> IO ()
hClose Handle
h
    Maybe Handle
Nothing -> () -> Mcmc a ()
forall (m :: * -> *) a. Monad m => a -> m a
return ()

-- | Run an MCMC algorithm.
mcmcRun :: ToJSON a => Mcmc a () -> Status a -> IO (Status a)
mcmcRun :: Mcmc a () -> Status a -> IO (Status a)
mcmcRun Mcmc a ()
algorithm = Mcmc a () -> Status a -> IO (Status a)
forall (m :: * -> *) s a. Monad m => StateT s m a -> s -> m s
execStateT (Mcmc a () -> Status a -> IO (Status a))
-> Mcmc a () -> Status a -> IO (Status a)
forall a b. (a -> b) -> a -> b
$ do
  Mcmc a ()
forall a. Mcmc a ()
mcmcInit
  Mcmc a ()
algorithm
  Mcmc a ()
forall a. ToJSON a => Mcmc a ()
mcmcClose