pyemma.thermo.estimate_multi_temperature

pyemma.thermo.estimate_multi_temperature(energy_trajs, temp_trajs, dtrajs, energy_unit='kcal/mol', temp_unit='K', reference_temperature=None, maxiter=10000, maxerr=1e-15, save_convergence_info=0, estimator='wham', lag=1, dt_traj='1 step', init=None, init_maxiter=10000, init_maxerr=1e-08, **kwargs)

This function acts as a wrapper for tram(), dtram(), and wham() and handles the calculation of bias energies (bias) and thermodynamic state trajectories (ttrajs) when the data comes from multi-temperature simulations.

Parameters:
  • energy_trajs (list of N arrays, each of shape (T_i,)) – List of arrays, each having T_i rows, one for each time step, containing the potential energies time series in units of kT, kcal/mol or kJ/mol.
  • temp_trajs (list of N int arrays, each of shape (T_i,)) – List of arrays, each having T_i rows, one for each time step, containing the heat bath temperature time series (at which temperatures the frames were created) in units of K or C. Alternatively, these trajectories may contain kT values instead of temperatures.
  • dtrajs (list of N int arrays, each of shape (T_i,)) – The integers are indexes in 0,...,n-1 enumerating the n discrete states or the bins the trajectory is in at any time.
  • energy_unit (str, optional, default='kcal/mol') – The physical unit used for energies. Current options: kcal/mol, kJ/mol, kT.
  • temp_unit (str, optional, default='K') – The physical unit used for the temperature. Current options: K, C, kT
  • reference_temperature (float or None, optional, default=None) – Reference temperature against which the bias energies are computed. If not given, the lowest temperature or kT value is used. If given, this parameter must have the same unit as the temp_trajs.
  • maxiter (int, optional, default=10000) – The maximum number of self-consistent iterations before the estimator exits unsuccessfully.
  • maxerr (float, optional, default=1E-15) – Convergence criterion based on the maximal free energy change in a self-consistent iteration step.
  • save_convergence_info (int, optional, default=0) – Every save_convergence_info iteration steps, store the actual increment and the actual loglikelihood; 0 means no storage.
  • estimator (str, optional, default='wham') –

    Specify one of the available estimators

    ‘wham’: use WHAM
    ‘mbar’: use MBAR
    ‘dtram’: use the discrete version of TRAM
    ‘tram’: use TRAM
  • lag (int or list of int, optional, default=1) – Integer lag time at which transitions are counted. Providing a list of lag times will trigger one estimation per lag time.
  • dt_traj (str, optional, default='1 step') –

    Description of the physical time corresponding to the lag. May be used by analysis algorithms such as plotting tools to pretty-print the axes. By default ‘1 step’, i.e. there is no physical time unit. Specify by a number, whitespace and unit. Permitted units are (* is an arbitrary string):

    ‘fs’, ‘femtosecond*’
    ‘ps’, ‘picosecond*’
    ‘ns’, ‘nanosecond*’
    ‘us’, ‘microsecond*’
    ‘ms’, ‘millisecond*’
    ‘s’, ‘second*’
  • init (str, optional, default=None) –

    Use a specific initialization for the self-consistent iteration:

    None: use a hard-coded guess for free energies and Lagrangian multipliers
    ‘wham’: perform a short WHAM estimate to initialize the free energies (only with dtram)
    ‘mbar’: perform a short MBAR estimate to initialize the free energies (only with tram)
  • init_maxiter (int, optional, default=10000) – The maximum number of self-consistent iterations during the initialization.
  • init_maxerr (float, optional, default=1.0E-8) – Convergence criterion for the initialization.
  • **kwargs (dict, optional) – You can use this to pass estimator-specific named parameters to the chosen estimator, which are not already coverd by estimate_multi_temperature().
Returns:

estimator_obj – The requested estimator/model object, i.e., WHAM, MBAR, DTRAM or TRAM. If multiple lag times are given, a list of objects is returned (one MEMM per lag time).
Return type:

MEMM or MultiThermModel or list thereof

Example

We look at 1D simulations at two different kT values 1.0 and 2.0, already clustered data, and we use TRAM for the estimation:

>>> from pyemma.thermo import estimate_multi_temperature as estimate_mt
>>> import numpy as np
>>> energy_trajs = [np.array([1.6, 1.4, 1.0, 1.0, 1.2, 1.0, 1.0]), np.array([0.8, 0.7, 0.5, 0.6, 0.7, 0.8, 0.7])]
>>> temp_trajs = [np.array([1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0]), np.array([2.0, 2.0, 2.0, 2.0, 2.0, 2.0, 2.0])]
>>> dtrajs = [np.array([0, 1, 2, 2, 2, 2, 2]), np.array([0, 1, 2, 2, 1, 0, 1])]
>>> tram = estimate_mt(energy_trajs, temp_trajs, dtrajs, energy_unit='kT', temp_unit='kT', estimator='tram', lag=1)
>>> tram.f 
array([ 2.90...,  1.72...,  0.26...])

Note that alhough we only used one temperature per trajectory, estimate_multi_temperature() can handle temperature changes as well.