1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
//! Simulation of (sub-)stochastic processes.
//!
//! # Goal
//!
//! Serve as an extension of the [rand crate](https://crates.io/crates/rand) for sub-stochastic processes.
//! 
//! # Examples
//!
//! ## Discrete time
//! 
//! Construction of a random walk in the integers.
//! ```
//! # #![allow(unused_mut)]
//! # use markovian::prelude::*;
//! # use rand::prelude::*;
//! let init_state: i32 = 0;
//! let transition = |state: &i32| raw_dist![(0.5, state + 1), (0.5, state - 1)];
//! let rng = thread_rng();
//! let mut mc = markovian::MarkovChain::new(init_state, transition, rng);
//! ```  
//! 
//! ## Branching process
//!
//! Construction using density p(0) = 0.3, p(1) = 0.4, p(2) = 0.3. 
//! ```
//! # #![allow(unused_mut)]
//! # use markovian::prelude::*;
//! # use rand::prelude::*;
//! let init_state: u32 = 1;
//! let base_distribution = raw_dist![(0.3, 0), (0.4, 1), (0.3, 2)];
//! let rng = thread_rng();
//! let mut branching_process = markovian::BranchingProcess::new(init_state, base_distribution, rng);
//! ``` 
//! ## Continuous time
//! 
//! Construction of a random walk in the integers, with expponential time for each transition.
//! ```
//! # #![allow(unused_mut)]
//! # use rand::prelude::*;
//! # use rand_distr::{Exp, Uniform};
//! # use markovian::prelude::*;
//! let init_state: i32 = 0;
//! struct MyTransition;
//! impl markovian::Transition<i32, (f64, i32)> for MyTransition {
//!     fn sample_from<R: ?Sized>(&self, state: &i32, rng: &mut R) -> (f64, i32)
//!     where
//!         R: Rng
//!     {
//!         let time = Exp::new(2.0).unwrap().sample(rng);
//!         let step = Uniform::from(0..=1).sample(rng) * 2 - 1;
//!         (time, state + step)
//!     }
//! }
//! let transition = MyTransition;
//! let rng = thread_rng();
//! let mut mc = markovian::TimedMarkovChain::new(init_state, transition, rng);
//! ``` 
//!
//! # Remarks
//!
//! All methods are `inline`, by design.
//! 
//! Non-trivial ways to use the crate are described below, including time dependence, 
//! continuous space and non-markovian processes.
//!
//! ## Time dependence
//!
//! Include the time as part of the state of the process.
//!
//! ### Examples
//! 
//! A random walk on the integers that tends to move more to the right as time goes by.  
//! ```
//! # #![allow(unused_mut)]
//! # use rand::prelude::*;
//! # use rand_distr::{Exp, Uniform};
//! # use markovian::prelude::*;
//! let init_state: (usize, i32) = (0, 0);
//! let transition = |(time, state): &(usize, i32)| raw_dist![
//!     (0.6 - 1.0 / (time + 2) as f64, (time + 1, state + 1)),
//!     (0.4 + 1.0 / (time + 2) as f64, (time + 1, state - 1))
//! ];
//! let rng = thread_rng();
//! let mut mc = markovian::MarkovChain::new(init_state, &transition, rng);
//! 
//! // Take a sample of 10 elements 
//! mc.take(10).map(|(_, state)| state).collect::<Vec<i32>>();
//! ```
//! 
//! ## Continuous space
//!
//! Randomize the transition: return a random element together with a probability one
//!
//! ### Examples
//! 
//! A random walk on the real line with variable step size.  
//! ```
//! # use rand_distr::Exp;
//! # use rand::prelude::*;
//! # use markovian::prelude::*;
//! let init_state: f64 = 0.0;
//! struct MyTransition;
//! impl markovian::Transition<f64, f64> for MyTransition {
//!     fn sample_from<R: ?Sized>(&self, state: &f64, rng: &mut R) -> f64
//!     where
//!         R: Rng
//!     {
//!         let step = Exp::new(2.0).unwrap().sample(rng);
//!         state + step
//!     }
//! }
//! let transition = MyTransition;
//! let rng = thread_rng();
//! let mut mc = markovian::MarkovChain::new(init_state, transition, rng);
//! mc.next();
//!  
//! // current_state is positive 
//! assert!(mc.state().unwrap() > &0.0);
//! ```
//! 
//! ## Non markovian
//!
//! Include history in the state. For example, instead of `i32`, use `Vec<i32>`. 
//!
//! ### Examples
//! 
//! A random walk on the integers that is atracted to zero in a non markovian
//! fashion. 
//! ```
//! # use rand::prelude::*;
//! # use markovian::prelude::*;
//! let init_state: Vec<i32> = vec![0];
//! let transition = |state: &Vec<i32>| {
//!     // New possible states
//!     let mut right = state.clone();
//!     right.push(state[state.len() - 1] + 1);
//!     let mut left = state.clone();
//!     left.push(state[state.len() - 1] - 1);
//! 
//!     // Some non markovian transtion
//!     let path_stadistic: i32 = state.iter().sum();
//!     if path_stadistic.is_positive() {
//!         raw_dist![
//!             (1.0 / (path_stadistic.abs() + 1) as f64, right), 
//!             (1.0 - 1.0 / (path_stadistic.abs() + 1) as f64, left)
//!         ]
//!     } else {
//!         raw_dist![
//!             (1.0 - 1.0 / (path_stadistic.abs() + 1) as f64, right), 
//!             (1.0 / (path_stadistic.abs() + 1) as f64, left)
//!         ]
//!     }
//! };
//! let rng = thread_rng();
//! let mut mc = markovian::MarkovChain::new(init_state, transition, rng);
//!  
//! // state has history
//! mc.next();
//! assert_eq!(mc.state().unwrap().len(), 2);
//! ```
//! 
pub use self::branching_process::BranchingProcess;
pub use self::continuous_finite_markov_chain::ContFiniteMarkovChain;
pub use self::finite_markov_chain::FiniteMarkovChain;
pub use self::markov_chain::MarkovChain;
pub use self::timed_markov_chain::TimedMarkovChain;
pub use self::traits::{State, StateIterator, Transition};

mod branching_process;
mod continuous_finite_markov_chain;
mod finite_markov_chain;
mod markov_chain;
mod timed_markov_chain;
mod traits;
mod macros;

/// Ease interoperability with rand_distr crate.
pub mod distributions;
/// Errors of this crate.
pub mod errors;


/// Ease of use of this crate in general.
pub mod prelude {
    pub use crate::traits::*;
    pub use crate::{raw_dist};
    pub use crate::distributions::Raw;
}


/// Testing random variables.
#[cfg(test)]
pub mod tests {
    // Notes on testing
    //
    // Testing random number distributions correctly is hard. The following
    // testing is desired:
    //
    // - Construction: test initialisation with a few valid parameter sets.
    // - Erroneous usage: test that incorrect usage generates an error.
    // - Vector: test that usage with fixed inputs (including RNG) generates a
    //   fixed output sequence on all platforms.
    // - Correctness at fixed points (optional): using a specific mock RNG,
    //   check that specific values are sampled (e.g. end-points and median of
    //   distribution).
    // - Correctness of PDF (extra): generate a histogram of samples within a
    //   certain range, and check this approximates the PDF. These tests are
    //   expected to be expensive, and should be behind a feature-gate.
    //
    // TODO: Vector and correctness tests are largely absent so far.
    // NOTE: Some distributions have tests checking only that samples can be
    // generated. This is redundant with vector and correctness tests.

    /// Construct a deterministic RNG with the given seed
    pub fn rng(seed: u64) -> impl rand::RngCore {
        // For tests, we want a statistically good, fast, reproducible RNG.
        // PCG32 will do fine, and will be easy to embed if we ever need to.
        const INC: u64 = 11634580027462260723;
        rand_pcg::Pcg32::new(seed, INC)
    }
}