{-# LANGUAGE DeriveDataTypeable #-}
{-# LANGUAGE DeriveGeneric #-}

-- |
--   Module      :  ELynx.Tree.Distribution.BirthDeath
--   Description :  Birth and death distribution
--   Copyright   :  2021 Dominik Schrempf
--   License     :  GPL-3.0-or-later
--
--   Maintainer  :  dominik.schrempf@gmail.com
--   Stability   :  unstable
--   Portability :  portable
--
-- Creation date: Tue Feb 13 13:16:18 2018.
--
-- See Gernhard, T. (2008). The conditioned reconstructed process. Journal of
-- Theoretical Biology, 253(4), 769–778. http://doi.org/10.1016/j.jtbi.2008.04.005.
--
-- Distribution of the values of the point process such that it corresponds to
-- reconstructed trees under the birth and death process.
module ELynx.Tree.Distribution.BirthDeath
  ( BirthDeathDistribution (..),
    cumulative,
    density,
    quantile,
  )
where

import Data.Data
  ( Data,
    Typeable,
  )
import ELynx.Tree.Distribution.Types
import GHC.Generics (Generic)
import qualified Statistics.Distribution as D

-- | Distribution of the values of the point process such that it corresponds to
-- a reconstructed tree of the birth and death process.
data BirthDeathDistribution = BDD
  { -- | Time to origin of the tree.
    BirthDeathDistribution -> Double
bddTOr :: Time,
    -- | Birth rate.
    BirthDeathDistribution -> Double
bddLa :: Rate,
    -- | Death rate.
    BirthDeathDistribution -> Double
bddMu :: Rate
  }
  deriving (BirthDeathDistribution -> BirthDeathDistribution -> Bool
forall a. (a -> a -> Bool) -> (a -> a -> Bool) -> Eq a
/= :: BirthDeathDistribution -> BirthDeathDistribution -> Bool
$c/= :: BirthDeathDistribution -> BirthDeathDistribution -> Bool
== :: BirthDeathDistribution -> BirthDeathDistribution -> Bool
$c== :: BirthDeathDistribution -> BirthDeathDistribution -> Bool
Eq, Typeable, Typeable BirthDeathDistribution
BirthDeathDistribution -> DataType
BirthDeathDistribution -> Constr
(forall b. Data b => b -> b)
-> BirthDeathDistribution -> BirthDeathDistribution
forall a.
Typeable a
-> (forall (c :: * -> *).
    (forall d b. Data d => c (d -> b) -> d -> c b)
    -> (forall g. g -> c g) -> a -> c a)
-> (forall (c :: * -> *).
    (forall b r. Data b => c (b -> r) -> c r)
    -> (forall r. r -> c r) -> Constr -> c a)
-> (a -> Constr)
-> (a -> DataType)
-> (forall (t :: * -> *) (c :: * -> *).
    Typeable t =>
    (forall d. Data d => c (t d)) -> Maybe (c a))
-> (forall (t :: * -> * -> *) (c :: * -> *).
    Typeable t =>
    (forall d e. (Data d, Data e) => c (t d e)) -> Maybe (c a))
-> ((forall b. Data b => b -> b) -> a -> a)
-> (forall r r'.
    (r -> r' -> r) -> r -> (forall d. Data d => d -> r') -> a -> r)
-> (forall r r'.
    (r' -> r -> r) -> r -> (forall d. Data d => d -> r') -> a -> r)
-> (forall u. (forall d. Data d => d -> u) -> a -> [u])
-> (forall u. Int -> (forall d. Data d => d -> u) -> a -> u)
-> (forall (m :: * -> *).
    Monad m =>
    (forall d. Data d => d -> m d) -> a -> m a)
-> (forall (m :: * -> *).
    MonadPlus m =>
    (forall d. Data d => d -> m d) -> a -> m a)
-> (forall (m :: * -> *).
    MonadPlus m =>
    (forall d. Data d => d -> m d) -> a -> m a)
-> Data a
forall u.
Int -> (forall d. Data d => d -> u) -> BirthDeathDistribution -> u
forall u.
(forall d. Data d => d -> u) -> BirthDeathDistribution -> [u]
forall r r'.
(r -> r' -> r)
-> r
-> (forall d. Data d => d -> r')
-> BirthDeathDistribution
-> r
forall r r'.
(r' -> r -> r)
-> r
-> (forall d. Data d => d -> r')
-> BirthDeathDistribution
-> r
forall (m :: * -> *).
Monad m =>
(forall d. Data d => d -> m d)
-> BirthDeathDistribution -> m BirthDeathDistribution
forall (m :: * -> *).
MonadPlus m =>
(forall d. Data d => d -> m d)
-> BirthDeathDistribution -> m BirthDeathDistribution
forall (c :: * -> *).
(forall b r. Data b => c (b -> r) -> c r)
-> (forall r. r -> c r) -> Constr -> c BirthDeathDistribution
forall (c :: * -> *).
(forall d b. Data d => c (d -> b) -> d -> c b)
-> (forall g. g -> c g)
-> BirthDeathDistribution
-> c BirthDeathDistribution
forall (t :: * -> *) (c :: * -> *).
Typeable t =>
(forall d. Data d => c (t d)) -> Maybe (c BirthDeathDistribution)
forall (t :: * -> * -> *) (c :: * -> *).
Typeable t =>
(forall d e. (Data d, Data e) => c (t d e))
-> Maybe (c BirthDeathDistribution)
gmapMo :: forall (m :: * -> *).
MonadPlus m =>
(forall d. Data d => d -> m d)
-> BirthDeathDistribution -> m BirthDeathDistribution
$cgmapMo :: forall (m :: * -> *).
MonadPlus m =>
(forall d. Data d => d -> m d)
-> BirthDeathDistribution -> m BirthDeathDistribution
gmapMp :: forall (m :: * -> *).
MonadPlus m =>
(forall d. Data d => d -> m d)
-> BirthDeathDistribution -> m BirthDeathDistribution
$cgmapMp :: forall (m :: * -> *).
MonadPlus m =>
(forall d. Data d => d -> m d)
-> BirthDeathDistribution -> m BirthDeathDistribution
gmapM :: forall (m :: * -> *).
Monad m =>
(forall d. Data d => d -> m d)
-> BirthDeathDistribution -> m BirthDeathDistribution
$cgmapM :: forall (m :: * -> *).
Monad m =>
(forall d. Data d => d -> m d)
-> BirthDeathDistribution -> m BirthDeathDistribution
gmapQi :: forall u.
Int -> (forall d. Data d => d -> u) -> BirthDeathDistribution -> u
$cgmapQi :: forall u.
Int -> (forall d. Data d => d -> u) -> BirthDeathDistribution -> u
gmapQ :: forall u.
(forall d. Data d => d -> u) -> BirthDeathDistribution -> [u]
$cgmapQ :: forall u.
(forall d. Data d => d -> u) -> BirthDeathDistribution -> [u]
gmapQr :: forall r r'.
(r' -> r -> r)
-> r
-> (forall d. Data d => d -> r')
-> BirthDeathDistribution
-> r
$cgmapQr :: forall r r'.
(r' -> r -> r)
-> r
-> (forall d. Data d => d -> r')
-> BirthDeathDistribution
-> r
gmapQl :: forall r r'.
(r -> r' -> r)
-> r
-> (forall d. Data d => d -> r')
-> BirthDeathDistribution
-> r
$cgmapQl :: forall r r'.
(r -> r' -> r)
-> r
-> (forall d. Data d => d -> r')
-> BirthDeathDistribution
-> r
gmapT :: (forall b. Data b => b -> b)
-> BirthDeathDistribution -> BirthDeathDistribution
$cgmapT :: (forall b. Data b => b -> b)
-> BirthDeathDistribution -> BirthDeathDistribution
dataCast2 :: forall (t :: * -> * -> *) (c :: * -> *).
Typeable t =>
(forall d e. (Data d, Data e) => c (t d e))
-> Maybe (c BirthDeathDistribution)
$cdataCast2 :: forall (t :: * -> * -> *) (c :: * -> *).
Typeable t =>
(forall d e. (Data d, Data e) => c (t d e))
-> Maybe (c BirthDeathDistribution)
dataCast1 :: forall (t :: * -> *) (c :: * -> *).
Typeable t =>
(forall d. Data d => c (t d)) -> Maybe (c BirthDeathDistribution)
$cdataCast1 :: forall (t :: * -> *) (c :: * -> *).
Typeable t =>
(forall d. Data d => c (t d)) -> Maybe (c BirthDeathDistribution)
dataTypeOf :: BirthDeathDistribution -> DataType
$cdataTypeOf :: BirthDeathDistribution -> DataType
toConstr :: BirthDeathDistribution -> Constr
$ctoConstr :: BirthDeathDistribution -> Constr
gunfold :: forall (c :: * -> *).
(forall b r. Data b => c (b -> r) -> c r)
-> (forall r. r -> c r) -> Constr -> c BirthDeathDistribution
$cgunfold :: forall (c :: * -> *).
(forall b r. Data b => c (b -> r) -> c r)
-> (forall r. r -> c r) -> Constr -> c BirthDeathDistribution
gfoldl :: forall (c :: * -> *).
(forall d b. Data d => c (d -> b) -> d -> c b)
-> (forall g. g -> c g)
-> BirthDeathDistribution
-> c BirthDeathDistribution
$cgfoldl :: forall (c :: * -> *).
(forall d b. Data d => c (d -> b) -> d -> c b)
-> (forall g. g -> c g)
-> BirthDeathDistribution
-> c BirthDeathDistribution
Data, forall x. Rep BirthDeathDistribution x -> BirthDeathDistribution
forall x. BirthDeathDistribution -> Rep BirthDeathDistribution x
forall a.
(forall x. a -> Rep a x) -> (forall x. Rep a x -> a) -> Generic a
$cto :: forall x. Rep BirthDeathDistribution x -> BirthDeathDistribution
$cfrom :: forall x. BirthDeathDistribution -> Rep BirthDeathDistribution x
Generic)

instance D.Distribution BirthDeathDistribution where
  cumulative :: BirthDeathDistribution -> Double -> Double
cumulative = BirthDeathDistribution -> Double -> Double
cumulative

-- | Cumulative distribution function Eq. (3).
cumulative :: BirthDeathDistribution -> Time -> Double
cumulative :: BirthDeathDistribution -> Double -> Double
cumulative (BDD Double
t Double
l Double
m) Double
x
  | Double
x forall a. Ord a => a -> a -> Bool
<= Double
0 = Double
0
  | Double
x forall a. Ord a => a -> a -> Bool
> Double
t = Double
1
  | Bool
otherwise = Double
t1 forall a. Num a => a -> a -> a
* Double
t2
  where
    d :: Double
d = Double
l forall a. Num a => a -> a -> a
- Double
m
    t1 :: Double
t1 = (Double
1.0 forall a. Num a => a -> a -> a
- forall a. Floating a => a -> a
exp (-Double
d forall a. Num a => a -> a -> a
* Double
x)) forall a. Fractional a => a -> a -> a
/ (Double
l forall a. Num a => a -> a -> a
- Double
m forall a. Num a => a -> a -> a
* forall a. Floating a => a -> a
exp (-Double
d forall a. Num a => a -> a -> a
* Double
x))
    t2 :: Double
t2 = (Double
l forall a. Num a => a -> a -> a
- Double
m forall a. Num a => a -> a -> a
* forall a. Floating a => a -> a
exp (-Double
d forall a. Num a => a -> a -> a
* Double
t)) forall a. Fractional a => a -> a -> a
/ (Double
1.0 forall a. Num a => a -> a -> a
- forall a. Floating a => a -> a
exp (-Double
d forall a. Num a => a -> a -> a
* Double
t))

instance D.ContDistr BirthDeathDistribution where
  density :: BirthDeathDistribution -> Double -> Double
density = BirthDeathDistribution -> Double -> Double
density
  quantile :: BirthDeathDistribution -> Double -> Double
quantile = BirthDeathDistribution -> Double -> Double
quantile

-- | Density function Eq. (2).
density :: BirthDeathDistribution -> Time -> Double
density :: BirthDeathDistribution -> Double -> Double
density (BDD Double
t Double
l Double
m) Double
x
  | Double
x forall a. Ord a => a -> a -> Bool
< Double
0 = Double
0
  | Double
x forall a. Ord a => a -> a -> Bool
> Double
t = Double
0
  | Bool
otherwise = Double
d forall a. Floating a => a -> a -> a
** Double
2 forall a. Num a => a -> a -> a
* Double
t1 forall a. Num a => a -> a -> a
* Double
t2
  where
    d :: Double
d = Double
l forall a. Num a => a -> a -> a
- Double
m
    t1 :: Double
t1 = forall a. Floating a => a -> a
exp (-Double
d forall a. Num a => a -> a -> a
* Double
x) forall a. Fractional a => a -> a -> a
/ ((Double
l forall a. Num a => a -> a -> a
- Double
m forall a. Num a => a -> a -> a
* forall a. Floating a => a -> a
exp (-Double
d forall a. Num a => a -> a -> a
* Double
x)) forall a. Floating a => a -> a -> a
** Double
2)
    t2 :: Double
t2 = (Double
l forall a. Num a => a -> a -> a
- Double
m forall a. Num a => a -> a -> a
* forall a. Floating a => a -> a
exp (-Double
d forall a. Num a => a -> a -> a
* Double
t)) forall a. Fractional a => a -> a -> a
/ (Double
1.0 forall a. Num a => a -> a -> a
- forall a. Floating a => a -> a
exp (-Double
d forall a. Num a => a -> a -> a
* Double
t))

-- | Inverted cumulative probability distribution 'cumulative'. See also
-- 'D.ContDistr'.
quantile :: BirthDeathDistribution -> Double -> Time
quantile :: BirthDeathDistribution -> Double -> Double
quantile (BDD Double
t Double
l Double
m) Double
p
  | Double
p forall a. Ord a => a -> a -> Bool
>= Double
0 Bool -> Bool -> Bool
&& Double
p forall a. Ord a => a -> a -> Bool
<= Double
1 =
      Double
res
  | Bool
otherwise =
      forall a. HasCallStack => [Char] -> a
error forall a b. (a -> b) -> a -> b
$
        [Char]
"BirthDeath.quantile: p must be in range [0,1] but got "
          forall a. [a] -> [a] -> [a]
++ forall a. Show a => a -> [Char]
show Double
p
          forall a. [a] -> [a] -> [a]
++ [Char]
"."
  where
    d :: Double
d = Double
l forall a. Num a => a -> a -> a
- Double
m
    t2 :: Double
t2 = (Double
l forall a. Num a => a -> a -> a
- Double
m forall a. Num a => a -> a -> a
* forall a. Floating a => a -> a
exp (-Double
d forall a. Num a => a -> a -> a
* Double
t)) forall a. Fractional a => a -> a -> a
/ (Double
1.0 forall a. Num a => a -> a -> a
- forall a. Floating a => a -> a
exp (-Double
d forall a. Num a => a -> a -> a
* Double
t))
    res :: Double
res = (-Double
1.0 forall a. Fractional a => a -> a -> a
/ Double
d) forall a. Num a => a -> a -> a
* forall a. Floating a => a -> a
log ((Double
1.0 forall a. Num a => a -> a -> a
- Double
p forall a. Num a => a -> a -> a
* Double
l forall a. Fractional a => a -> a -> a
/ Double
t2) forall a. Fractional a => a -> a -> a
/ (Double
1.0 forall a. Num a => a -> a -> a
- Double
p forall a. Num a => a -> a -> a
* Double
m forall a. Fractional a => a -> a -> a
/ Double
t2))

instance D.ContGen BirthDeathDistribution where
  genContVar :: forall g (m :: * -> *).
StatefulGen g m =>
BirthDeathDistribution -> g -> m Double
genContVar = forall d g (m :: * -> *).
(ContDistr d, StatefulGen g m) =>
d -> g -> m Double
D.genContinuous