core

Authors: Chase Coleman, Spencer Lyon, Daisuke Oyama, Tom Sargent,
John Stachurski

Filename: core.py

This file contains some useful objects for handling a finite-state discrete-time Markov chain.

Definitions and Some Basic Facts about Markov Chains

Let \(\{X_t\}\) be a Markov chain represented by an \(n \times n\) stochastic matrix \(P\). State \(i\) has access to state \(j\), denoted \(i \to j\), if \(i = j\) or \(P^k[i, j] > 0\) for some \(k = 1, 2, \ldots\); \(i\) and j communicate, denoted \(i \leftrightarrow j\), if \(i \to j\) and \(j \to i\). The binary relation \(\leftrightarrow\) is an equivalent relation. A communication class of the Markov chain \(\{X_t\}\), or of the stochastic matrix \(P\), is an equivalent class of \(\leftrightarrow\). Equivalently, a communication class is a strongly connected component (SCC) in the associated directed graph \(\Gamma(P)\), a directed graph with \(n\) nodes where there is an edge from \(i\) to \(j\) if and only if \(P[i, j] > 0\). The Markov chain, or the stochastic matrix, is irreducible if it admits only one communication class, or equivalently, if \(\Gamma(P)\) is strongly connected.

A state \(i\) is recurrent if \(i \to j\) implies \(j \to i\); it is transient if it is not recurrent. For any \(i, j\) contained in a communication class, \(i\) is recurrent if and only if \(j\) is recurrent. Therefore, recurrence is a property of a communication class. Thus, a communication class is a recurrent class if it contains a recurrent state. Equivalently, a recurrent class is a SCC that corresponds to a sink node in the condensation of the directed graph \(\Gamma(P)\), where the condensation of \(\Gamma(P)\) is a directed graph in which each SCC is replaced with a single node and there is an edge from one SCC \(C\) to another SCC \(C'\) if \(C \neq C'\) and some node in \(C\) has access to some node in \(C'\). A recurrent class is also called a closed communication class. The condensation is acyclic, so that there exists at least one recurrent class.

For example, if the entries of \(P\) are all strictly positive, then the whole state space is a communication class as well as a recurrent class. (More generally, if there is only one communication class, then it is a recurrent class.) As another example, consider the stochastic matrix \(P = [[1, 0], [0,5, 0.5]]\). This has two communication classes, \(\{0\}\) and \(\{1\}\), and \(\{0\}\) is the only recurrent class.

A stationary distribution of the Markov chain \(\{X_t\}\), or of the stochastic matrix \(P\), is a nonnegative vector \(x\) such that \(x' P = x'\) and \(x' \mathbf{1} = 1\), where \(\mathbf{1}\) is the vector of ones. The Markov chain has a unique stationary distribution if and only if it has a unique recurrent class. More generally, each recurrent class has a unique stationary distribution whose support equals that recurrent class. The set of all stationary distributions is given by the convex hull of these unique stationary distributions for the recurrent classes.

A natural number \(d\) is the period of state \(i\) if it is the greatest common divisor of all \(k\)‘s such that \(P^k[i, i] > 0\); equivalently, it is the GCD of the lengths of the cycles in \(\Gamma(P)\) passing through \(i\). For any \(i, j\) contained in a communication class, \(i\) has period \(d\) if and only if \(j\) has period \(d\). The period of an irreducible Markov chain (or of an irreducible stochastic matrix) is the period of any state. We define the period of a general (not necessarily irreducible) Markov chain to be the least common multiple of the periods of its recurrent classes, where the period of a recurrent class is the period of any state in that class. A Markov chain is aperiodic if its period is one. A Markov chain is irreducible and aperiodic if and only if it is uniformly ergodic, i.e., there exists some \(m\) such that \(P^m[i, j] > 0\) for all \(i, j\) (in this case, \(P\) is also called primitive).

Suppose that an irreducible Markov chain has period \(d\). Fix any state, say state \(0\). For each \(m = 0, \ldots, d-1\), let \(S_m\) be the set of states \(i\) such that \(P^{kd+m}[0, i] > 0\) for some \(k\). These sets \(S_0, \ldots, S_{d-1}\) constitute a partition of the state space and are called the cyclic classes. For each \(S_m\) and each \(i \in S_m\), we have \(\sum_{j \in S_{m+1}} P[i, j] = 1\), where \(S_d = S_0\).

class quantecon.markov.core.MarkovChain(P, state_values=None)[source]

Bases: object

Class for a finite-state discrete-time Markov chain. It stores useful information such as the stationary distributions, and communication, recurrent, and cyclic classes, and allows simulation of state transitions.

Parameters:

P : array_like or scipy sparse matrix (float, ndim=2)

The transition matrix. Must be of shape n x n.

state_values : array_like(default=None)

Array_like of length n containing the values associated with the states, which must be homogeneous in type. If None, the values default to integers 0 through n-1.

Notes

In computing stationary distributions, if the input matrix is a sparse matrix, internally it is converted to a dense matrix.

Attributes

P (ndarray or scipy.sparse.csr_matrix (float, ndim=2)) See Parameters
stationary_distributions (array_like(float, ndim=2)) Array containing stationary distributions, one for each recurrent class, as rows.
is_irreducible (bool) Indicate whether the Markov chain is irreducible.
num_communication_classes (int) The number of the communication classes.
communication_classes_indices (list(ndarray(int))) List of numpy arrays containing the indices of the communication classes.
communication_classes (list(ndarray)) List of numpy arrays containing the communication classes, where the states are annotated with their values (if state_values is not None).
num_recurrent_classes (int) The number of the recurrent classes.
recurrent_classes_indices (list(ndarray(int))) List of numpy arrays containing the indices of the recurrent classes.
recurrent_classes (list(ndarray)) List of numpy arrays containing the recurrent classes, where the states are annotated with their values (if state_values is not None).
is_aperiodic (bool) Indicate whether the Markov chain is aperiodic.
period (int) The period of the Markov chain.
cyclic_classes_indices (list(ndarray(int))) List of numpy arrays containing the indices of the cyclic classes. Defined only when the Markov chain is irreducible.
cyclic_classes (list(ndarray)) List of numpy arrays containing the cyclic classes, where the states are annotated with their values (if state_values is not None). Defined only when the Markov chain is irreducible.

Methods

get_index(value) Return the index (or indices) of the given value (or values) in state_values.
simulate(ts_length[, init, num_reps, …]) Simulate time series of state transitions, where the states are annotated with their values (if state_values is not None).
simulate_indices(ts_length[, init, …]) Simulate time series of state transitions, where state indices are returned.
cdfs
cdfs1d
communication_classes
communication_classes_indices
cyclic_classes
cyclic_classes_indices
digraph
get_index(value)[source]

Return the index (or indices) of the given value (or values) in state_values.

Parameters:

value

Value(s) to get the index (indices) for.

Returns:

idx : int or ndarray(int)

Index of value if value is a single state value; array of indices if value is an array_like of state values.

is_aperiodic
is_irreducible
num_communication_classes
num_recurrent_classes
period
recurrent_classes
recurrent_classes_indices
simulate(ts_length, init=None, num_reps=None, random_state=None)[source]

Simulate time series of state transitions, where the states are annotated with their values (if state_values is not None).

Parameters:

ts_length : scalar(int)

Length of each simulation.

init : scalar or array_like, optional(default=None)

Initial state values(s). If None, the initial state is randomly drawn.

num_reps : scalar(int), optional(default=None)

Number of repetitions of simulation.

random_state : int or np.random.RandomState, optional

Random seed (integer) or np.random.RandomState instance to set the initial state of the random number generator for reproducibility. If None, a randomly initialized RandomState is used.

Returns:

X : ndarray(ndim=1 or 2)

Array containing the state values of the sample path(s). See the simulate method for more information.

simulate_indices(ts_length, init=None, num_reps=None, random_state=None)[source]

Simulate time series of state transitions, where state indices are returned.

Parameters:

ts_length : scalar(int)

Length of each simulation.

init : int or array_like(int, ndim=1), optional

Initial state(s). If None, the initial state is randomly drawn.

num_reps : scalar(int), optional(default=None)

Number of repetitions of simulation.

random_state : int or np.random.RandomState, optional

Random seed (integer) or np.random.RandomState instance to set the initial state of the random number generator for reproducibility. If None, a randomly initialized RandomState is used.

Returns:

X : ndarray(ndim=1 or 2)

Array containing the sample path(s), of shape (ts_length,) if init is a scalar (integer) or None and num_reps is None; of shape (k, ts_length) otherwise, where k = len(init) if (init, num_reps) = (array, None), k = num_reps if (init, num_reps) = (int or None, int), and k = len(init)*num_reps if (init, num_reps) = (array, int).

state_values
stationary_distributions
quantecon.markov.core.mc_compute_stationary(P)[source]

Computes stationary distributions of P, one for each recurrent class. Any stationary distribution is written as a convex combination of these distributions.

Returns:

stationary_dists : array_like(float, ndim=2)

Array containing the stationary distributions as its rows.

quantecon.markov.core.mc_sample_path(P, init=0, sample_size=1000, random_state=None)[source]

Generates one sample path from the Markov chain represented by (n x n) transition matrix P on state space S = {{0,…,n-1}}.

Parameters:

P : array_like(float, ndim=2)

A Markov transition matrix.

init : array_like(float ndim=1) or scalar(int), optional(default=0)

If init is an array_like, then it is treated as the initial distribution across states. If init is a scalar, then it treated as the deterministic initial state.

sample_size : scalar(int), optional(default=1000)

The length of the sample path.

random_state : int or np.random.RandomState, optional

Random seed (integer) or np.random.RandomState instance to set the initial state of the random number generator for reproducibility. If None, a randomly initialized RandomState is used.

Returns:

X : array_like(int, ndim=1)

The simulation of states.