-- Note: It is not necessary to add another type @b@ to store supplementary -- information about the chain. The information can just be stored in @a@ -- equally well. -- XXX: Status tuned exclusively to the Metropolis-Hastings algorithm. We should -- abstract the algorithm from the chain. Maybe something like: -- -- @ -- data Status a = Status { Chain a; Algorithm a} -- @ -- | -- Module : Mcmc.Status -- Description : What is an MCMC? -- Copyright : (c) Dominik Schrempf 2020 -- License : GPL-3.0-or-later -- -- Maintainer : dominik.schrempf@gmail.com -- Stability : unstable -- Portability : portable -- -- Creation date: Tue May 5 18:01:15 2020. module Mcmc.Status ( Cleaner (..), Status (..), status, cleanWith, saveWith, force, quiet, debug, noData, ) where import Data.Maybe import Data.Time.Clock import Mcmc.Item import Mcmc.Monitor import Mcmc.Proposal import Mcmc.Trace import Mcmc.Verbosity (Verbosity (..)) import Numeric.Log import System.IO import System.Random.MWC hiding (save) import Prelude hiding (cycle) -- | Clean the state periodically. -- -- The prior and the likelihood will be updated after the cleaning process. -- -- For long chains, successive numerical errors can accumulate such that the -- state diverges from honoring specific required constraints. In these cases, a -- 'Cleaner' can be used to ensure that the required constraints of the state -- are honored. For example, the branches of an ultrametric phylogeny may -- diverge slightly after successful many proposals such that the phylogeny is -- not anymore ultrametric. -- -- Please be aware that the Markov chain will not converge to the true posterior -- distribution if the state is changed substantially! Only apply subtle changes -- that are absolutely necessary to preserve the required properties of the -- state such as specific numerical constraints. data Cleaner a = Cleaner { -- | Clean every given number of iterations. clEvery :: Int, -- | Cleaning function. Executed before monitoring the state. clFunction :: a -> a } -- | The 'Status' contains all information to run an MCMC chain. It is -- constructed using the function 'status'. -- -- The polymorphic type @a@ stores the state of the chain. It can also be used -- to store auxiliary information. data Status a = Status { -- MCMC related variables; saved. -- | The name of the MCMC chain; used as file prefix. name :: String, -- | The current 'Item' of the chain combines the current state and the -- current likelihood. item :: Item a, -- | The iteration is the number of completed cycles. iteration :: Int, -- | The 'Trace' of the Markov chain in reverse order, the most recent -- 'Item' is at the head of the list. trace :: Trace a, -- | For each 'Proposal', store the list of accepted (True) and rejected (False) -- proposals; for reasons of efficiency, the list is also stored in reverse -- order. acceptance :: Acceptance (Proposal a), -- | Number of burn in iterations; deactivate burn in with 'Nothing'. burnInIterations :: Maybe Int, -- | Auto tuning period (only during burn in); deactivate auto tuning with -- 'Nothing'. autoTuningPeriod :: Maybe Int, -- | Number of normal iterations excluding burn in. Note that auto tuning -- only happens during burn in. iterations :: Int, -- -- Auxiliary variables; saved. -- | Overwrite output files? Default is 'False', change with 'force'. forceOverwrite :: Bool, -- | Save the chain with trace of given length at the end of the run? -- Default is no save ('Nothing'). Change with 'saveWith'. save :: Maybe Int, -- | Verbosity. verbosity :: Verbosity, -- | The random number generator. generator :: GenIO, -- -- Auxiliary variables; not saved. -- | Starting time and starting iteration of chain; used to calculate -- run time and ETA. start :: Maybe (Int, UTCTime), -- | Handle to log file. logHandle :: Maybe Handle, -- -- Auxiliary functions; not saved. -- | The prior function. The un-normalized posterior is the product of the -- prior and the likelihood. priorF :: a -> Log Double, -- | The likelihood function. The un-normalized posterior is the product of -- the prior and the likelihood. likelihoodF :: a -> Log Double, -- | Clean the state periodically. cleaner :: Maybe (Cleaner a), -- -- Variables related to the algorithm; not saved. -- | A set of 'Proposal's form a 'Cycle'. cycle :: Cycle a, -- | A 'Monitor' observing the chain. monitor :: Monitor a } -- | Initialize the 'Status' of a Markov chain Monte Carlo run. status :: -- | Name of the Markov chain; used as file prefix. String -> -- | The prior function. (a -> Log Double) -> -- | The likelihood function. (a -> Log Double) -> -- | A list of 'Proposal's executed in forward order. The chain will be logged -- after each cycle. Cycle a -> -- | A 'Monitor' observing the chain. Monitor a -> -- | The initial state in the state space @a@. a -> -- | Number of burn in iterations; deactivate burn in with 'Nothing'. Maybe Int -> -- | Auto tuning period (only during burn in); deactivate auto tuning with -- 'Nothing'. Maybe Int -> -- | Number of normal iterations excluding burn in. Note that auto tuning only -- happens during burn in. Int -> -- | A source of randomness. For reproducible runs, make sure to use -- generators with the same, fixed seed. GenIO -> Status a status n p l c m x mB mT nI g | isJust mT && isNothing mB = error "status: Auto tuning period given, but no burn in." | otherwise = Status n i 0 (singletonT i) (emptyA $ ccProposals c) mB mT nI False Nothing Info g Nothing Nothing p l Nothing c m where i = Item x (p x) (l x) -- | Clean the state every given number of generations using the given function. -- See 'Cleaner'. cleanWith :: Cleaner a -> Status a -> Status a cleanWith c s = s {cleaner = Just c} -- | Save the Markov chain with trace of given length. saveWith :: Int -> Status a -> Status a saveWith n s = s {save = Just n} -- | Overwrite existing files; it is not necessary to use 'force', when a chain -- is continued. force :: Status a -> Status a force s = s {forceOverwrite = True} -- | Do not print anything to standard output. Do not create log file. File -- monitors and batch monitors are executed normally. quiet :: Status a -> Status a quiet s = s {verbosity = Quiet} -- | Be verbose. debug :: Status a -> Status a debug s = s {verbosity = Debug} -- | Set the likelihood function to 1.0. Useful for debugging and testing. noData :: Status a -> Status a noData s = s {likelihoodF = const 1.0}