PyGT.stats

Calculate first passage statistics between macrostates

Tools to calculate the first passage time distribution and phenomenological rate constants between endpoint macrostates \(\mathcal{A}\) and \(\mathcal{B}\).

Note

Install the pathos package to parallelize MFPT computations.

PyGT.stats.compute_passage_stats(A_sel, B_sel, pi, K, dopdf=True, rt=None)[source]

Compute the A->B and B->A first passage time distribution, first moment, and second moment using eigendecomposition of a CTMC rate matrix.

Parameters:
  • A_sel ((N,) array-like) – boolean array that selects out the A nodes
  • B_sel ((N,) array-like) – boolean array that selects out the B nodes
  • pi ((N,) array-like) – stationary distribution
  • K ((N, N) array-like) – CTMC rate matrix
  • dopdf (bool, optional) – Do we calculate full fpt distribution or just the moments. Defaults=True.
  • rt (array, optional) – Vector of times to evaluate first passage time distribution in multiples of \(\left<t\right>\) for A->B and B->A. If None, defaults to a logscale array from \(0.001\left<t\right>\) to \(1000\left<t\right>\) in 400 steps, i.e. np.logspace(-3,3,400). Only relevant if dopdf=True
Returns:

  • tau ((4,) array-like) – First and second moments of first passage time distribution for A->B and B->A [\(\mathcal{T}_{\mathcal{B}\mathcal{A}}\), \(\mathcal{V}_{\mathcal{B}\mathcal{A}}\), \(\mathcal{T}_{\mathcal{A}\mathcal{B}}\), \(\mathcal{V}_{\mathcal{A}\mathcal{B}}\)]
  • pt (( len(rt),4) array-like) – time and first passage time distribution p(t) for A->B and B->A
PyGT.stats.compute_escape_stats(B_sel, pi, K, tau_escape=None, dopdf=True, rt=None)[source]

Compute escape time distribution and first and second moment from the basin specified by B_sel using eigendecomposition.

Parameters:
  • B_sel ((N,) array-like) – boolean array that selects out the nodes in the active basin
  • pi ((N,) array-like) – stationary distribution for CTMC
  • K ((N, N) array-like) – CTMC rate matrix
  • tau_escape (float) – mean time to escape from B. Used to calculate the escape time distribution in multiple of tau_escape (p(t/tau_escape). If None, uses the first moment in network defined by K.
  • dopdf (bool) – whether to calculate full escape time distribution, defaults to True
  • rt (array, optional) – Vector of times to evaluate first passage time distribution in multiples of \(\left<t\right>\) for A->B and B->A. If None, defaults to a logscale array from \(0.001\left<t\right>\) to \(1000\left<t\right>\) in 400 steps, i.e. np.logspace(-3,3,400). Only relevant if dopdf=True
Returns:

  • tau ((2,) array-like) – First and second moments of escape time distribution, [\(\left<t\right>_{\mathcal{B}}\), \(\left<t^2 \right>_{\mathcal{B}}\)]
  • pt ((2, len(rt)) array-like) – time and escape time distribution \(p(t)\left<t\right>\)
PyGT.stats.compute_rates(A_sel, B_sel, B, tau, pi, initA=None, initB=None, MFPTonly=True, fullGT=False, pool_size=None, block=1, screen=False, **kwargs)[source]

In a total state space partitioned into three sets A,B and I, calculate various approximate A<->B transition rates [Wales09] and the MFPT from a rate matrix K. K can be the matrix of an original Markov chain, or a partially graph-transformed Markov chain. [Swinburne20a]

Rate definitions (all assume A,B in local equilibrium) see [Swinburne20a] for details

kSS : assumes I in steady state, A,B local equilibrium

kNSS : Relaxes counts for non-steady state in I. Tends to exact rate if reactant metastable

k^F : Boltzmann weighted sum of inverse of exact MFPT from each reactant state, i.e. \(\sum_j\pi_j/\mathcal{T}_j\)

k^* : Inverse of exact MFPT from reactant, i.e. \(1/\sum_j\pi_j\mathcal{T}_j\)

compute_rates() differs from compute_passage_stats() in that this function removes all intervening states using GT before computing FPT stats and rates on the fully reduced network with state space \((\mathcal{A} \cup \mathcal{B})\). This implementation also does not rely on a full eigendecomposition of the non-absorbing matrix; it instead performs a matrix inversion, or if fullGT is specified, all nodes in the set \((\mathcal{A} \cup b)^\mathsf{c}\) are removed using GT for each \(b \in \mathcal{B}\) so that the MFPT is given by an average:

\[\begin{equation} \mathcal{T}_{\mathcal{A}\mathcal{B}} = \frac{1}{\sum_{b \in \mathcal{B}} p_b(0)} \sum_{b \in \mathcal{B}} \frac{p_b(0) \tau^\prime_b}{1-P^\prime_{bb}} \end{equation}\]

If the MFPT is less than \(10^{20}\), fullGT should not be needed since the inversion of the non-absorbing matrix should be numerically stable. However, a condition number check is performed regardless, which forces a full GT if the MFPT problem is considered numerically unstable.

Parameters:
  • A_sel (array-like (N,)) – selects the N_A nodes in the \(\mathcal{A}\) set.
  • B_sel (array-like (N,)) – selects the N_B nodes in the \(\mathcal{B}\) set.
  • B (sparse or dense matrix (N,N)) – branching probability matrix for CTMC
  • tau (array-like (N,)) – vector of waiting times from each node.
  • pi (array-like (N,)) – stationary distribution of CTMC
  • initA (array-like (N_A,), optional) – normalized initial occupation probabilities in \(\mathcal{A}\) set. Default= local Boltzmann distribution
  • initB (array-like (N_B,), optional) – normalized initial occupation probabilities in \(\mathcal{B}\) set. Default= local Boltzmann distribution
  • MFPTonly (bool) – If True, only MFPTs are calculated (rate calculations ignored).
  • fullGT (bool) – If True, all source nodes are isolated with GT to obtain the average MFPT.
  • pool_size (int) – Number of cores over which to parallelize fullGT computation.
Returns:

results – dictionary of results, with keys ‘MFPTAB’, ‘kSSAB’, ‘kNSSAB’, ‘kQSDAB’, ‘k*AB’, ‘kFAB’ for A<-B and ‘MFPTBA’, ‘kSSBA’, ‘kNSSBA’, ‘kQSDBA’, ‘k*BA’, ‘kFBA’ for B<-A such that results[‘MFPTAB’] = mean first passage time for A<-B
Return type:

dictionary