{-# LANGUAGE OverloadedStrings #-} {-# LANGUAGE RankNTypes #-} -- TODO: Allow execution of moves in order of appearance in the cycle. -- TODO: Moves and monitors both use lenses and names for what they move and -- monitor. Should a data structure be used combining the lens and the name, so -- that things are cohesive? -- TODO: Moves on simplices: SimplexElementScale (?). -- TODO: Moves on tree branch lengths. -- - Slide a node on the tree. -- - Scale a tree. -- TODO: Moves on tree topologies. -- - NNI -- - Narrow (what is this, see RevBayes) -- - FNPR (dito) -- TODO: Bactrian moves; https://www.ncbi.nlm.nih.gov/pmc/articles/PMC3845170/. -- -- slideBactrian -- -- scaleBactrian -- | -- Module : Mcmc.Move -- Description : Moves and cycles -- Copyright : (c) Dominik Schrempf 2020 -- License : GPL-3.0-or-later -- -- Maintainer : dominik.schrempf@gmail.com -- Stability : unstable -- Portability : portable -- -- Creation date: Wed May 20 13:42:53 2020. module Mcmc.Move ( -- * Move Move (..), MoveSimple (..), Tuner (tParam, tFunc), tuner, tune, autotune, -- * Cycle Cycle (fromCycle), fromList, autotuneC, summarizeCycle, -- * Acceptance Acceptance (..), emptyA, pushA, resetA, acceptanceRatios, ) where import Data.Aeson import Data.Function import Data.List import qualified Data.Map.Strict as M import Data.Map.Strict (Map) import Data.Maybe import qualified Data.Text.Lazy as T import Data.Text.Lazy (Text) import qualified Data.Text.Lazy.Builder as B import qualified Data.Text.Lazy.Builder.Int as B import qualified Data.Text.Lazy.Builder.RealFloat as B import Numeric.Log hiding (sum) import System.Random.MWC -- | A 'Move' is an instruction about how the Markov chain will traverse the -- state space @a@. Essentially, it is a probability density conditioned on the -- current state. -- -- A 'Move' may be tuneable in that it contains information about how to enlarge -- or shrink the step size to tune the acceptance ratio. data Move a = Move { -- | Name (no moves with the same name are allowed in a 'Cycle'). mvName :: String, -- | The weight determines how often a 'Move' is executed per iteration of -- the Markov chain. mvWeight :: Int, -- | Simple move without tuning information. mvSimple :: MoveSimple a, -- | Tuning is disabled if set to 'Nothing'. mvTune :: Maybe (Tuner a) } instance Show (Move a) where show m = show $ mvName m instance Eq (Move a) where m == n = mvName m == mvName n instance Ord (Move a) where compare = compare `on` mvName -- XXX: One could also use a different type for 'mvSample', so that -- 'mvDensity' can be avoided. In detail, -- -- @ -- mvSample :: a -> GenIO -> IO (a, Log Double, Log, Double) -- @ -- -- where the densities describe the probability of going there and back. -- However, we may need more information about the move for other MCMC samplers -- different from Metropolis-Hastings. -- | Simple move without tuning information. -- -- In order to calculate the Metropolis-Hastings ratio, we need to know the -- probability (density) of jumping forth, and the probability (density) of -- jumping back. data MoveSimple a = MoveSimple { -- | Instruction about randomly moving from the current state to a new -- state, given some source of randomness. mvSample :: a -> GenIO -> IO a, -- | The density of going from one state to another. Set to 'Nothing' for -- symmetric moves. mvDensity :: Maybe (a -> a -> Log Double) } -- | Tune the acceptance ratio of a 'Move'; see 'tune', or 'autotune'. data Tuner a = Tuner { tParam :: Double, tFunc :: Double -> MoveSimple a } -- | Create a 'Tuner'. The tuning function accepts a tuning parameter, and -- returns a corresponding 'MoveSimple'. The larger the tuning parameter, the -- larger the 'Move', and vice versa. tuner :: (Double -> MoveSimple a) -> Tuner a tuner = Tuner 1.0 -- | Tune a 'Move'. Return 'Nothing' if 'Move' is not tuneable. If the parameter -- @dt@ is larger than 1.0, the 'Move' is enlarged, if @0