class KoopmanReweightedMSM

class deeptime.markov.msm.KoopmanReweightedMSM(transition_matrix: ndarray, stationary_distribution: Optional[ndarray] = None, reversible: Optional[bool] = None, n_eigenvalues: Optional[int] = None, ncv: Optional[int] = None, twostep_count_matrices: Optional[ndarray] = None, oom_components: Optional[ndarray] = None, count_model: Optional[TransitionCountModel] = None, oom_eigenvalues: Optional[ndarray] = None, oom_evaluator: Optional[ndarray] = None, oom_information_state_vector: Optional[ndarray] = None)

This class belongs to a markov state model which was estimated by Koopman reweighting.

Parameters:
  • transition_matrix ((n, n) ndarray) – The transition matrix, see MarkovStateModel.__init__() .

  • stationary_distribution (ndarray(n), optional, default=None) – Stationary distribution if already computed, see MarkovStateModel.__init__() .

  • reversible (bool, optional, default=None) – Whether MSM is reversible, see MarkovStateModel.__init__() .

  • n_eigenvalues (int, optional, default=None) – The number of eigenvalues / eigenvectors to be kept, see MarkovStateModel.__init__() .

  • ncv (int optional, default=None) – Performance parameter for reversible MSMs, see MarkovStateModel.__init__() .

  • count_model (TransitionCountModel, optional, default=None) – In case the model was estimated with a OOMReweightedMSM estimator, this contains a count matrix based on lagged data, i.e., data with the last lag frames removed and histogram information based on the full dtrajs.

  • twostep_count_matrices ((n, n, n) ndarray, optional, default=None) – Two-step count matrices for all states, each twostep_count_matrices[:, i, :] is a count matrix for each \(i=1,\ldots,n\).

  • oom_components ((m, n, m) ndarray, optional, default=None) – Matrix of set-observable operators, where m is the rank based on the rank decision made during estimation.

  • oom_eigenvalues ((m,) ndarray, optional, default=None) – The eigenvalues from OOM.

  • oom_evaluator ((m,) ndarray, optional, default=None) – Evaluator of OOM.

  • oom_information_state_vector ((m,) ndarray, optional, default=None) – Information state vector of OOM.

Attributes

count_model

Returns a transition count model, can be None.

empirical_koopman_model

Yields a CovarianceKoopmanModel based on the count matrix of this Markov state model.

has_count_model

Yields whether this Markov state model has a count model.

is_real

Checks if all eigenvalues as well as eigenvectors/functions are real.

koopman_model

Yields a CovarianceKoopmanModel based on the transition matrix and stationary

lagtime

The lagtime this model was estimated at.

n_eigenvalues

number of eigenvalues to compute.

n_states

Number of active states on which all computations and estimations are done

ncv

Number of Lanczos vectors used when computing the partial eigenvalue decomposition

oom_components

Return OOM components.

oom_eigenvalues

System eigenvalues estimated by OOM.

oom_evaluator

Return OOM evaluator vector.

oom_information_state_vector

Return OOM initial state vector.

oom_rank

Return OOM model rank.

oom_timescales

System timescales estimated by OOM.

reversible

Returns whether the MarkovStateModel is reversible

sparse

Returns whether the MarkovStateModel is sparse

stationary

Whether the MSM is stationary, i.e. whether the initial distribution is the stationary distribution of the hidden transition matrix.

stationary_distribution

The stationary distribution on the MarkovStateModel states

transition_matrix

The transition matrix on the active set.

transition_matrix_tolerance

The tolerance under which a matrix is considered a transition matrix.

twostep_count_matrices

] is a count matrix for each n.

Methods

ck_test(models, n_metastable_sets[, ...])

Validates a model estimated at lag time tau by testing its predictions for longer lag times.

committor_backward(A, B)

Backward committor from set A to set B

committor_forward(A, B)

Forward committor (also known as p_fold or splitting probability) from set A to set B.

compute_state_indices(dtrajs)

Generates a trajectory/time indices for the given list of states.

compute_trajectory_weights(dtrajs)

Uses the MarkovStateModel to assign a probability weight to each trajectory frame.

copy()

Makes a deep copy of this model.

correlation(a[, b, maxtime, k, ncv])

Time-correlation for equilibrium experiment.

eigenvalues([k])

Compute or fetch the transition matrix eigenvalues.

eigenvectors_left([k])

Compute the left transition matrix eigenvectors

eigenvectors_right([k])

Compute the right transition matrix eigenvectors

expectation(a)

Equilibrium expectation value of a given observable.

fingerprint_correlation(a[, b, k, ncv])

Dynamical fingerprint for equilibrium time-correlation experiment.

fingerprint_relaxation(p0, a[, k, ncv])

Dynamical fingerprint for perturbation/relaxation experiment.

get_params([deep])

Get the parameters.

hmm(dtrajs, nhidden[, return_estimator])

Estimates a hidden Markov state model as described in [1].

mfpt(A, B)

Mean first passage times from set A to set B, in units of the input trajectory time step.

pcca(n_metastable_sets)

Runs PCCA+ [2] to compute a metastable decomposition of MarkovStateModel states.

propagate(p0, k)

Propagates the initial distribution p0 k times

reactive_flux(source_states, target_states)

A->B reactive flux from transition path theory (TPT)

relaxation(p0, a[, maxtime, k, ncv])

Simulates a perturbation-relaxation experiment.

score([dtrajs, r, dim])

Scores the MSM using the dtrajs using the variational approach for Markov processes.

set_params(**params)

Set the parameters of this estimator.

simulate(n_steps[, start, stop, dt, seed])

Generates a realization of the Markov Model.

submodel(states)

Restricts this markov state model to a subset of states by taking a submatrix of the transition matrix and re-normalizing it, as well as restricting the stationary distribution and count model if given.

timescales([k])

Relaxation timescales corresponding to the eigenvalues.

to_koopman_model([empirical, epsilon])

Computes the SVD of the symmetrized Koopman operator in the analytical or empirical distribution, returns as Koopman model.

update_stationary_distribution(value)

Explicitly sets the stationary distribution, re-normalizes

update_transition_matrix(value)

Sets the transition matrix and invalidates all cached and derived properties.

ck_test(models, n_metastable_sets, include_lag0=True, err_est=False, progress=None)

Validates a model estimated at lag time tau by testing its predictions for longer lag times. This is known as the Chapman-Kolmogorov test as it is based on the Chapman-Kolmogorov equation. The test is performed on metastable sets of states rather than the micro-states themselves.

Parameters:
  • models (list of MarkovStateModel) – A list of models which were estimated at different and in particular longer lagtimes.

  • n_metastable_sets (int) – Number of metastable sets, estimated via pcca().

  • include_lag0 (bool, optional, default=True) – Whether to include a lagtime 0 in the test, which corresponds to an identity transition matrix for MSMs.

  • err_est (bool, optional, default=False) – Whether to compute errors on observable evaluations of the models in the models parameter.

  • progress (object, optional, default=None) – Optional progress bar. Tested for tqdm.

Returns:

ck_test – A Chapman-Kolmogorov test object that can be used with plot_ck_test.

Return type:

ChapmanKolmogorovTest

committor_backward(A, B)

Backward committor from set A to set B

Parameters:
  • A (int or int array) – set of starting states

  • B (int or int array) – set of target states

committor_forward(A, B)

Forward committor (also known as p_fold or splitting probability) from set A to set B.

Parameters:
  • A (int or int array) – set of starting states

  • B (int or int array) – set of target states

compute_state_indices(dtrajs) List[ndarray]

Generates a trajectory/time indices for the given list of states. If a count model is provided in this MSM and it does not represent the full state space, the discrete trajectories are first mapped to the active state space, inactive states are mapped to -1.

Parameters:

dtrajs (array_like or list of array_like) – Discretized trajectories.

Returns:

state_indices – A list of arrays with trajectory/time indices for the provided discretized trajectories.

Return type:

List of ndarray

compute_trajectory_weights(dtrajs)

Uses the MarkovStateModel to assign a probability weight to each trajectory frame.

This is a powerful function for the calculation of arbitrary observables in the trajectories one has started the analysis with. The stationary probability of the MarkovStateModel will be used to reweigh all states. Returns a list of weight arrays, one for each trajectory, and with a number of elements equal to trajectory frames. Given \(N\) trajectories of lengths \(T_1\) to \(T_N\), this function returns corresponding weights:

\[(w_{1,1}, ..., w_{1,T_1}), (w_{N,1}, ..., w_{N,T_N})\]

that are normalized to one:

\[\sum_{i=1}^N \sum_{t=1}^{T_i} w_{i,t} = 1\]

Suppose you are interested in computing the expectation value of a function \(a(x)\), where \(x\) are your input configurations. Use this function to compute the weights of all input configurations and obtain the estimated expectation by:

\[\langle a \rangle = \sum_{i=1}^N \sum_{t=1}^{T_i} w_{i,t} a(x_{i,t})\]

Or if you are interested in computing the time-lagged correlation between functions \(a(x)\) and \(b(x)\) you could do:

\[\langle a(t) b(t+\tau) \rangle_t = \sum_{i=1}^N \sum_{t=1}^{T_i} w_{i,t} a(x_{i,t}) a(x_{i,t+\tau})\]
Returns:

weights – The normalized trajectory weights. Given \(N\) trajectories of lengths \(T_1\) to \(T_N\), returns the corresponding weights:

\[(w_{1,1}, ..., w_{1,T_1}), (w_{N,1}, ..., w_{N,T_N})\]

Return type:

list of ndarray

copy() Model

Makes a deep copy of this model.

Returns:

A new copy of this model.

Return type:

copy

correlation(a, b=None, maxtime=None, k=None, ncv=None)

Time-correlation for equilibrium experiment.

In order to simulate a time-correlation experiment (e.g. fluorescence correlation spectroscopy [3], dynamical neutron scattering [4], …), first compute the mean values of your experimental observable \(a\) by MarkovStateModel state:

\[a_i = \frac{1}{N_i} \sum_{x_t \in S_i} f(x_t)\]

where \(S_i\) is the set of configurations belonging to MarkovStateModel state \(i\) and \(f()\) is a function that computes the experimental observable of interest for configuration \(x_t\). If a cross-correlation function is wanted, also apply the above computation to a second experimental observable \(b\).

Then the accurate (i.e. without statistical error) autocorrelation function of \(f(x_t)\) given the Markov model is computed by correlation(a), and the accurate cross-correlation function is computed by correlation(a,b). This is done by evaluating the equation

\[\begin{aligned} acf_a(k\tau) &= \mathbf{a}^\top \mathrm{diag}(\boldsymbol{\pi}) \mathbf{P(\tau)}^k \mathbf{a} \\ ccf_{a,b}(k\tau) &= \mathbf{a}^\top \mathrm{diag}(\boldsymbol{\pi}) \mathbf{P(\tau)}^k \mathbf{b} \end{aligned}\]

where \(acf\) stands for autocorrelation function and \(ccf\) stands for cross-correlation function, \(\mathbf{P(\tau)}\) is the transition matrix at lag time \(\tau\), \(\boldsymbol{\pi}\) is the equilibrium distribution of \(\mathbf{P}\), and \(k\) is the time index.

Note that instead of using this method you could generate a long synthetic trajectory from the MarkovStateModel and then estimating the time-correlation of your observable(s) directly from this trajectory. However, there is no reason to do this because the present method does that calculation without any sampling, and only in the limit of an infinitely long synthetic trajectory the two results will agree exactly. The correlation function computed by the present method still has statistical uncertainty from the fact that the underlying MarkovStateModel transition matrix has statistical uncertainty when being estimated from data, but there is no additional (and unnecessary) uncertainty due to synthetic trajectory generation.

Parameters:
  • a ((n,) ndarray) – Observable, represented as vector on state space

  • maxtime (int or float) – Maximum time (in units of the input trajectory time step) until which the correlation function will be evaluated. Internally, the correlation function can only be computed in integer multiples of the Markov model lag time, and therefore the actual last time point will be computed at \(\mathrm{ceil}(\mathrm{maxtime} / \tau)\) By default (None), the maxtime will be set equal to the 5 times the slowest relaxation time of the MarkovStateModel, because after this time the signal is almost constant.

  • b ((n,) ndarray (optional)) – Second observable, for cross-correlations

  • k (int (optional)) – Number of eigenvalues and eigenvectors to use for computation. This option is only relevant for sparse matrices and long times for which an eigenvalue decomposition will be done instead of using the matrix power.

  • ncv (int (optional)) – Only relevant for sparse matrices and large lag times where the relaxation will be computed using an eigenvalue decomposition. The number of Lanczos vectors generated, ncv must be greater than k; it is recommended that ncv > 2*k.

Returns:

  • times (ndarray (N)) – Time points (in units of the input trajectory time step) at which the correlation has been computed

  • correlations (ndarray (N)) – Correlation values at given times

Examples

This example computes the autocorrelation function of a simple observable on a three-state Markov model and plots the result using matplotlib:

>>> import numpy as np
>>> import deeptime.markov.msm as msm
>>>
>>> P = np.array([[0.99, 0.01, 0], [0.01, 0.9, 0.09], [0, 0.1, 0.9]])
>>> a = np.array([0.0, 0.5, 1.0])
>>> M = msm.MarkovStateModel(P)
>>> times, acf = M.correlation(a)
>>>
>>> import matplotlib.pylab as plt 
>>> plt.plot(times,acf)  
eigenvalues(k=None)

Compute or fetch the transition matrix eigenvalues.

Parameters:

k (int) – number of eigenvalues to be returned. By default will return all available eigenvalues

Returns:

ts – transition matrix eigenvalues \(\lambda_i, i = 1, ..., k\)., sorted by descending norm.

Return type:

ndarray(k,)

eigenvectors_left(k=None)

Compute the left transition matrix eigenvectors

Parameters:

k (int) – number of eigenvectors to be returned. By default uses value of n_eigenvalues.

Returns:

L – left eigenvectors in a row matrix. l_ij is the j’th component of the i’th left eigenvector

Return type:

ndarray(k,n)

eigenvectors_right(k=None)

Compute the right transition matrix eigenvectors

Parameters:

k (int) – number of eigenvectors to be computed. By default uses value of n_eigenvalues.

Returns:

R – right eigenvectors in a column matrix. r_ij is the i’th component of the j’th right eigenvector

Return type:

ndarray(n,k)

expectation(a: ndarray)

Equilibrium expectation value of a given observable.

Parameters:

a ((n,) ndarray) – Observable vector on the MarkovStateModel state space

Returns:

val – Equilibrium expectation value fo the given observable

Return type:

float

Notes

The equilibrium expectation value of an observable \(a\) is defined as follows

\[\mathbb{E}_{\mu}[a] = \sum_i \pi_i a_i\]

\(\pi=(\pi_i)\) is the stationary vector of the transition matrix \(P\).

fingerprint_correlation(a, b=None, k=None, ncv=None)

Dynamical fingerprint for equilibrium time-correlation experiment.

Parameters:
  • a ((n,) ndarray) – Observable, represented as vector on MarkovStateModel state space

  • b ((n,) ndarray, optional) – Second observable, for cross-correlations

  • k (int, optional) – Number of eigenvalues and eigenvectors to use for computation. This option is only relevant for sparse matrices and long times for which an eigenvalue decomposition will be done instead of using the matrix power

  • ncv (int, optional) – Only relevant for sparse matrices and large lag times, where the relaxation will be computed using an eigenvalue decomposition. The number of Lanczos vectors generated, ncv must be greater than k; it is recommended that ncv > 2*k

Returns:

  • timescales ((k,) ndarray) – Time-scales (in units of the input trajectory time step) of the transition matrix

  • amplitudes ((k,) ndarray) – Amplitudes for the correlation experiment

  • Spectral densities are commonly used in spectroscopy. Dynamical

  • fingerprints are a useful representation for computational

  • spectroscopy results and have been introduced in [3].

References

fingerprint_relaxation(p0, a, k=None, ncv=None)

Dynamical fingerprint for perturbation/relaxation experiment.

Parameters:
  • p0 ((n,) ndarray) – Initial distribution for a relaxation experiment

  • a ((n,) ndarray) – Observable, represented as vector on state space

  • k (int, optional) – Number of eigenvalues and eigenvectors to use for computation

  • ncv (int, optional) – Only relevant for sparse matrices and large lag times, where the relaxation will be computes using an eigenvalue decomposition. The number of Lanczos vectors generated, ncv must be greater than k; it is recommended that ncv > 2*k

Returns:

  • timescales ((k,) ndarray) – Time-scales (in units of the input trajectory time step) of the transition matrix

  • amplitudes ((k,) ndarray) – Amplitudes for the relaxation experiment

  • Spectral densities are commonly used in spectroscopy. Dynamical fingerprints are a useful representation

  • for computational spectroscopy results and have been introduced in [3].

get_params(deep=False)

Get the parameters.

Returns:

params – Parameter names mapped to their values.

Return type:

mapping of string to any

hmm(dtrajs, nhidden: int, return_estimator=False)

Estimates a hidden Markov state model as described in [1].

Parameters:
  • dtrajs (list of int-array like) – discrete trajectories to use for estimation of the HMM.

  • nhidden (int) – number of hidden (metastable) states

  • return_estimator (boolean, optional) – if False only the Model is returned, if True both the Estimator and the Model is returned.

Returns:

hmm – A hidden markov model.

Return type:

deeptime.markov.hmm.HiddenMarkovModel

mfpt(A, B)

Mean first passage times from set A to set B, in units of the input trajectory time step.

Parameters:
  • A (int or int array) – set of starting states

  • B (int or int array) – set of target states

pcca(n_metastable_sets: int) PCCAModel

Runs PCCA+ [2] to compute a metastable decomposition of MarkovStateModel states.

After calling this method you can access metastable_memberships(), metastable_distributions(), metastable_sets() and metastable_assignments() of the returned object.

Parameters:

n_metastable_sets (int) – Number of metastable sets

Returns:

pcca_obj – An object containing all PCCA+ quantities.

Return type:

PCCAModel

Notes

If you coarse grain with PCCA+, the order of the obtained memberships might not be preserved.

propagate(p0, k: int)

Propagates the initial distribution p0 k times

Computes the product

\[p_k = p_0^T P^k\]

If the lag time of transition matrix \(P\) is \(\tau\), this will provide the probability distribution at time \(k \tau\).

Parameters:
  • p0 (ndarray(n,)) – Initial distribution. Vector of size of the active set.

  • k (int) – Number of time steps

Returns:

pk – Distribution after k steps. Vector of size of the active set.

Return type:

ndarray(n,)

reactive_flux(source_states, target_states) ReactiveFlux

A->B reactive flux from transition path theory (TPT)

The returned object can be used to extract various quantities of the flux, as well as to compute A -> B transition pathways, their weights, and to coarse-grain the flux onto sets of states.

Parameters:
  • source_states (array_like) – List of integer state labels for set A

  • target_states (array_like) – List of integer state labels for set B

Returns:

tptobj – An object containing the reactive A->B flux network and several additional quantities, such as the stationary probability, committors and set definitions.

Return type:

ReactiveFlux

See also

ReactiveFlux

Reactive Flux model

relaxation(p0, a, maxtime=None, k=None, ncv=None)

Simulates a perturbation-relaxation experiment.

In perturbation-relaxation experiments such as temperature-jump, pH-jump, pressure jump or rapid mixing experiments, an ensemble of molecules is initially prepared in an off-equilibrium distribution and the expectation value of some experimental observable is then followed over time as the ensemble relaxes towards equilibrium.

In order to simulate such an experiment, first determine the distribution of states at which the experiment is started, \(p_0\) and compute the mean values of your experimental observable \(a\) by MarkovStateModel state:

\[a_i = \frac{1}{N_i} \sum_{x_t \in S_i} f(x_t)\]

where \(S_i\) is the set of configurations belonging to MarkovStateModel state \(i\) and \(f()\) is a function that computes the experimental observable of interest for configuration \(x_t\).

Then the accurate (i.e. without statistical error) time-dependent expectation value of \(f(x_t)\) given the Markov model is computed by relaxation(p0, a). This is done by evaluating the equation

\[E_a(k\tau) = \mathbf{p_0}^{\top} \mathbf{P(\tau)}^k \mathbf{a}\]

where \(E\) stands for the expectation value that relaxes to its equilibrium value that is identical to expectation(a), \(\mathbf{P(\tau)}\) is the transition matrix at lag time \(\tau\), \(\boldsymbol{\pi}\) is the equilibrium distribution of \(\mathbf{P}\), and \(k\) is the time index.

Note that instead of using this method you could generate many synthetic trajectory from the MarkovStateModel with starting points drawn from the initial distribution and then estimating the time-dependent expectation value by an ensemble average. However, there is no reason to do this because the present method does that calculation without any sampling, and only in the limit of an infinitely many trajectories the two results will agree exactly. The relaxation function computed by the present method still has statistical uncertainty from the fact that the underlying MarkovStateModel transition matrix has statistical uncertainty when being estimated from data, but there is no additional (and unnecessary) uncertainty due to synthetic trajectory generation.

Parameters:
  • p0 ((n,) ndarray) – Initial distribution for a relaxation experiment

  • a ((n,) ndarray) – Observable, represented as vector on state space

  • maxtime (int or float, optional) – Maximum time (in units of the input trajectory time step) until which the correlation function will be evaluated. Internally, the correlation function can only be computed in integer multiples of the Markov model lag time, and therefore the actual last time point will be computed at \(\mathrm{ceil}(\mathrm{maxtime} / \tau)\). By default (None), the maxtime will be set equal to the 5 times the slowest relaxation time of the MarkovStateModel, because after this time the signal is constant.

  • k (int (optional)) – Number of eigenvalues and eigenvectors to use for computation

  • ncv (int (optional)) – Only relevant for sparse matrices and large lag times, where the relaxation will be computes using an eigenvalue decomposition. The number of Lanczos vectors generated, ncv must be greater than k; it is recommended that ncv > 2*k

Returns:

  • times (ndarray (N)) – Time points (in units of the input trajectory time step) at which the relaxation has been computed

  • res (ndarray) – Array of expectation value at given times

score(dtrajs=None, r=2, dim=None)

Scores the MSM using the dtrajs using the variational approach for Markov processes.

Implemented according to [5] and [6].

Currently only implemented using dense matrices - will be slow for large state spaces.

Parameters:
  • dtrajs (list of arrays, optional, default=None) – Test data (discrete trajectories). Note that if the test data contains states which are not represented in this model, they are ignored.

  • r (float or str, optional, default=2) –

    Overwrite scoring method to be used if desired. Can be any float greater or equal 1 or ‘E’ for VAMP-r score or VAMP-E score, respectively.

    Special cases [5] [6]:

    • ’VAMP1’: Sum of singular values of the symmetrized transition matrix. If the MSM is reversible, this is equal to the sum of transition matrix eigenvalues, also called Rayleigh quotient [7].

    • ’VAMP2’: Sum of squared singular values of the symmetrized transition matrix [6]. If the MSM is reversible, this is equal to the kinetic variance [8].

    • ’VAMPE’: Approximation error of the estimated Koopman operator with respect to the true Koopman operator up to an additive constant [6].

  • dim (int or None, optional, default=None) – The maximum number of eigenvalues or singular values used in the score. If set to None, all available eigenvalues will be used.

set_params(**params)

Set the parameters of this estimator.

The method works on simple estimators as well as on nested objects (such as pipelines). The latter have parameters of the form <component>__<parameter> so that it’s possible to update each component of a nested object.

Parameters:

**params (dict) – Estimator parameters.

Returns:

self – Estimator instance.

Return type:

object

simulate(n_steps: int, start: Optional[int] = None, stop: Optional[Union[int, List[int]]] = None, dt: int = 1, seed: int = -1)

Generates a realization of the Markov Model.

Parameters:
  • n_steps (int) – trajectory length in steps of the lag time

  • start (int, optional, default = None) – starting state

  • stop (int or int-array-like, optional, default = None) – stopping hidden set. If given, the trajectory will be stopped before N steps once a hidden state of the stop set is reached

  • dt (int, default=1) – trajectory will be saved every dt time steps. Internally, the dt’th power of P is taken to ensure a more efficient simulation.

  • seed (int, default=-1) – If non-negative, this fixes the seed used for sampling to the provided value so that results are reproducible.

Returns:

The state trajectory with length N/dt

Return type:

(N/dt,) ndarray

Examples

>>> msm = MarkovStateModel(transition_matrix=np.array([[.7, .3], [.3, .7]]))
>>> trajectory = msm.simulate(n_steps=15)
>>> print(trajectory)  
[...]
submodel(states)

Restricts this markov state model to a subset of states by taking a submatrix of the transition matrix and re-normalizing it, as well as restricting the stationary distribution and count model if given.

Parameters:

states (array_like of int) – states to restrict to

Returns:

submodel – A onto the given states restricted MSM.

Return type:

MarkovStateModel

timescales(k=None)

Relaxation timescales corresponding to the eigenvalues.

Parameters:

k (int, optional, default=None) – Number of timescales to be returned. By default, this uses all available eigenvalues except for the first (stationary) eigenvalue.

Returns:

ts – Relaxation timescales in units of the input trajectory time step, defined by \(-\tau / \mathrm{ln} | \lambda_i |, i = 2,...,k+1\).

Return type:

ndarray(m)

to_koopman_model(empirical: bool = True, epsilon: float = 1e-10)

Computes the SVD of the symmetrized Koopman operator in the analytical or empirical distribution, returns as Koopman model.

Parameters:
  • empirical (bool, optional, default=True) – Determines whether the model should refer to the analytical distributions based on the transition matrix or the empirical distributions based on the count matrix.

  • epsilon (float, optional, default=1e-10) – Regularization parameter for computing inverses of covariance matrices.

Returns:

model – The model.

Return type:

CovarianceKoopmanModel

update_stationary_distribution(value: Optional[ndarray])

Explicitly sets the stationary distribution, re-normalizes

update_transition_matrix(value: ndarray)

Sets the transition matrix and invalidates all cached and derived properties.

property count_model: Optional[TransitionCountModel]

Returns a transition count model, can be None. The transition count model statistics about data that was used for transition counting as well as a count matrix.

Returns:

model – The transition count model or None.

Return type:

TransitionCountModel or None

property empirical_koopman_model: CovarianceKoopmanModel

Yields a CovarianceKoopmanModel based on the count matrix of this Markov state model.

Returns:

model – The model.

Return type:

CovarianceKoopmanModel

Notes

If \(C\in\mathbb{R}^{n\times n}\) denotes the count matrix and \(P\) the transition matrix, we define covariance matrices based on transition count statistics

\[\begin{aligned} C_{00} &= \text{diag} \left( \sum_i C_{i1}, \ldots, \sum_i C_{in} \right) \\ C_{11} &= \text{diag} \left( \sum_i C_{1i}, \ldots, \sum_i C_{ni} \right) \\ C_{01} &= C, \end{aligned}\]

and reweight the operator \(P\) to the empirical distribution via \(C_{01\text{, re}} = C_{00}P\). Based on these we define the Koopman operator \(K = C_{00}^{-1/2}C_{01\text{, re}}C_{11}^{-1/2}\).

See also

to_koopman_model

property has_count_model: bool

Yields whether this Markov state model has a count model.

Type:

bool

property is_real

Checks if all eigenvalues as well as eigenvectors/functions are real.

property koopman_model: CovarianceKoopmanModel
Yields a CovarianceKoopmanModel based on the transition matrix and stationary

distribution of this Markov state model.

Returns:

model – The model.

Return type:

CovarianceKoopmanModel

Notes

If \(P\in\mathbb{R}^{n\times n}\) denotes the transition matrix and \(\mu\in\mathbb{R}^n\) denotes the stationary distribution, we define covariance matrices

\[\begin{aligned} C_{00} &= C_{11} = \text{diag} (\mu_1,\ldots,\mu_n)\\ C_{01} &= C_{00}P \end{aligned} \]

and based on these a Koopman operator \(K = C_{00}^{-1/2}C_{01}C_{11}^{-1/2}\).

See also

to_koopman_model

property lagtime: int

The lagtime this model was estimated at. In case no count model was provided, this property defaults to the lagtime set in the constructor or a lagtime of \(\tau = 1\) if it was left None.

Returns:

lagtime – The lagtime.

Return type:

int

property n_eigenvalues: int

number of eigenvalues to compute.

property n_states: int

Number of active states on which all computations and estimations are done

property ncv

Number of Lanczos vectors used when computing the partial eigenvalue decomposition

property oom_components

Return OOM components.

property oom_eigenvalues

System eigenvalues estimated by OOM.

property oom_evaluator

Return OOM evaluator vector.

property oom_information_state_vector

Return OOM initial state vector.

property oom_rank

Return OOM model rank.

property oom_timescales

System timescales estimated by OOM.

property reversible: bool

Returns whether the MarkovStateModel is reversible

property sparse: bool

Returns whether the MarkovStateModel is sparse

property stationary

Whether the MSM is stationary, i.e. whether the initial distribution is the stationary distribution of the hidden transition matrix.

property stationary_distribution

The stationary distribution on the MarkovStateModel states

property transition_matrix

The transition matrix on the active set.

property transition_matrix_tolerance

The tolerance under which a matrix is considered a transition matrix. This means that all elements must be non-negative and the row sums must be 1.

property twostep_count_matrices

] is a count matrix for each n.

Type:

Two-step count matrices for all states. C2t[

Type:

, n,