glhmm.glhmm

Gaussian Linear Hidden Markov Model @author: Diego Vidaurre 2023

class glhmm.glhmm.glhmm(K=10, covtype='shareddiag', model_mean='state', model_beta='state', dirichlet_diag=10, connectivity=None, Pstructure=None, Pistructure=None)

Bases: object

Gaussian Linear Hidden Markov Model class to decode stimulus from data.

Attributes:

Kint, default=10

number of states in the model.

covtypestr, {‘shareddiag’, ‘diag’,’sharedfull’,’full’}, default ‘shareddiag’

Type of covariance matrix. Choose ‘shareddiag’ to have one diagonal covariance matrix for all states, or ‘diag’ to have a diagonal full covariance matrix for each state, or ‘sharedfull’ to have a shared full covariance matrix for all states, or ‘full’ to have a full covariance matrix for each state.

model_meanstr, {‘state’, ‘shared’, ‘no’}, default ‘state’

Model for the mean. If ‘state’, the mean will be modelled state-dependent. If ‘shared’, the mean will be modelled globally (shared between all states). If ‘no’ the mean of the timeseries will not be used to drive the states.

model_betastr, {‘state’, ‘shared’, ‘no’}, default ‘state’

Model for the beta. If ‘state’, the regression coefficients will be modelled state-dependent. If ‘shared’, the regression coefficients will be modelled globally (shared between all states). If ‘no’ the regression coefficients will not be used to drive the states.

dirichlet_diagfloat, default=10

The value of the diagonal of the Dirichlet distribution for the transition probabilities. The higher the value, the more persistent the states will be. Note that this value is relative; the prior competes with the data, so if the timeseries is very long, the dirichlet_diag may have little effect unless it is set to a very large value.

connectivityarray_like of shape (n_states, n_states), optional

Matrix of binary values defining the connectivity of the states. This parameter can only be used with a diagonal covariance matrix (i.e., covtype=’diag’).

Pstructurearray_like, optional

Binary matrix defining the allowed transitions between states. The default is a (n_states, n_states) matrix of all ones, allowing all possible transitions between states.

Pistructurearray_like, optional

Binary vector defining the allowed initial states. The default is a (n_states,) vector of all ones, allowing all states to be used as initial states.

Notes:

This class requires the following modules: numpy, math, scipy, sys, warnings, copy, and time.

decode(X, Y, indices=None, files=None, viterbi=False, set=None, serial=False, gpu_acceleration=0, gpuChunks=1)

Calculates state time courses for all the data using either parallel or sequential processing.

Parameters:

Xarray-like of shape (n_samples, n_parcels)

The timeseries of set of variables 1.

Yarray-like of shape (n_samples, n_parcels)

The timeseries of set of variables 2.

indicesarray-like of shape (n_sessions, 2), optional, default=None

The start and end indices of each trial/session in the input data.

fileslist of str, optional, default=None

List of filenames corresponding to the indices.

viterbibool, optional, default=False

Whether or not the Viterbi algorithm should be used.

setint, optional, default=None

Index of the sessions set to decode.

Returns:

If viterbi=True:

vpatharray of shape (n_samples,)

The most likely state sequence.

If viterbi=False:

Gammaarray of shape (n_samples, n_states)

The state probability timeseries.

Xiarray of shape (n_samples - n_sessions, n_states, n_states)

The joint probabilities of past and future states conditioned on data.

scalearray-like of shape (n_samples,)

The scaling factors from the inference, used to compute the free energy. In normal use, we would do

Gamma,Xi,_ = hmm.decode(X,Y,indices)

Raises:

Exception

If the model has not been trained. If both ‘files’ and ‘Y’ arguments are provided.

dual_estimate(X, Y, indices=None, Gamma=None, Xi=None, for_kernel=False)

Dual estimation of HMM parameters.

Parameters:

Xarray-like of shape (n_samples, n_variables_1)

The timeseries of set of variables 1.

Yarray-like of shape (n_samples, n_variables_2)

The timeseries of set of variables 2.

indicesarray-like of shape (n_sessions, 2), optional

The start and end indices of each trial/session in the input data. If None, a single segment spanning the entire sequence is used.

Gammaarray-like of shape (n_samples, n_states), optional

The state probabilities. If None, it is computed from the input observations.

Xiarray-like of shape (n_samples - n_sessions, n_states, n_states), optional

The joint probabilities of past and future states conditioned on data. If None, it is computed from the input observations.

for_kernelbool, optional

Whether purpose of dual estimation is kernel (gradient) computation, or not If True, function will also return Gamma and Xi (default False)

Returns:

hmm_dualobject

A copy of the HMM object with updated dynamics and observation distributions.

get_P()

Returns transition probability matrix

Returns:

P: ndarray of shape (K,K), where K is the number of states

Raises:

Exception

If the model has not yet been trained.

get_Pi()

Returns initial probabilities

Returns:

Pi: ndarray of shape (K,), where K is the number of states

Raises:

Exception

If the model has not yet been trained.

get_active_K()

Returns the number of active states

Returns:

K_activeint

Number of active states.

get_beta(k=0)

Returns the regression coefficients (beta) for the specified state.

Parameters:

kint, optional, default=0

The index of the state for which to retrieve the beta value.

Returns:

beta: ndarray of shape (n_variables_1 x n_variables_2)

The regression coefficients of each variable in X on each variable in Y for the specified state.

Raises:

Exception

If the model has not yet been trained. If the model has no beta.

get_betas()

Returns the regression coefficients (beta) for all states.

Returns:

betas: ndarray of shape (n_variables_1 x n_variables_2 x n_states)

The regression coefficients of each variable in X on each variable in Y for all states.

Raises:

Exception

If the model has not yet been trained. If the model has no beta.

get_covariance_matrix(k=0)

Returns the covariance matrix for the specified state.

Parameters:

kint, optional

The index of the state. Default=0.

Returns:

array of shape (n_parcels, n_parcels)

The covariance matrix for the specified state.

Raises:

Exception

If the model has not been trained.

get_fe(X, Y, Gamma, Xi, scale=None, indices=None, todo=None, non_informative_prior_P=False)

Computes the Free Energy of an HMM depending on observation model.

Parameters:

Xarray of shape (n_samples, n_parcels)

The timeseries of set of variables 1.

Yarray of shape (n_samples, n_parcels)

The timeseries of set of variables 2.

Gammaarray of shape (n_samples, n_states), default=None

The state timeseries probabilities.

Xiarray-like of shape (n_samples - n_sessions, n_states, n_states)

The joint probabilities of past and future states conditioned on data.

scalearray-like of shape (n_samples,), default=None

The scaling factors used to compute the free energy of the dataset. If None, scaling is automatically computed.

indicesarray-like of shape (n_sessions, 2), optional, default=None

The start and end indices of each trial/session in the input data.

todo: bool of shape (n_terms,) or None, default=None

Whether or not each of the 5 elements (see fe_terms) should be computed. Only for internal use.

non_informative_prior_P: array-like of shape (n_states, n_states), optional, default=False

Prior of transition probability matrix Only for internal use.

Returns:

fe_termsarray of shape (n_terms,)

The variational free energy, separated into different terms: - element 1: Gamma Entropy - element 2: Data negative log-likelihood - element 3: Gamma negative log-likelihood - element 4: KL divergence for initial and transition probabilities - element 5: KL divergence for the state parameters

Raises:

Exception

If the model has not been trained.

Notes:

This function computes the variational free energy using a specific algorithm. For more information on the algorithm, see [^1].

References:

[^1] Smith, J. et al. “A variational approach to Bayesian learning of switching dynamics in dynamical systems.” Journal of Machine Learning Research, vol. 18, no. 4, 2017.

get_inverse_covariance_matrix(k=0)

Returns the inverse covariance matrix for the specified state.

Parameters:

kint, optional

The index of the state. Default=0.

Returns:

array of shape (n_parcels, n_parcels)

The inverse covariance matrix for the specified state.

Raises:

Exception

If the model has not been trained.

get_mean(k=0)

Returns the mean for the specified state.

Parameters:

kint, optional, default=0

The index of the state for which to retrieve the mean.

Returns:

mean: ndarray of shape (n_variables_2,)

The mean value of each variable in Y for the specified state.

Raises:

Exception

If the model has not yet been trained. If the model has no mean.

get_means()

Returns the means for all states.

Returns:

means: ndarray of shape (n_variables_2, n_states)

The mean value of each variable in Y for all states.

Raises:

Exception

If the model has not yet been trained. If the model has no mean.

get_r2(X, Y, Gamma, indices=None)

Computes the explained variance per session/trial and per column of Y

Parameters:

Xarray of shape (n_samples, n_variables_1)

The timeseries of set of variables 1.

Yarray of shape (n_samples, n_variables_2)

The timeseries of set of variables 2.

Gammaarray of shape (n_samples, n_states), default=None

The state timeseries probabilities.

indicesarray-like of shape (n_sessions, 2), optional, default=None

The start and end indices of each trial/session in the input data.

Returns:

r2array of shape (n_sessions, n_variables_2)

The R-squared (proportion of the variance explained) for each session and each variable in Y.

Raises:

Exception

If the model has not been trained, or if it does not have neither mean or beta

Notes:

This function does not take the covariance matrix into account

initialize(p, q)

Initialize the parameters of the HMM with initial random values. IMPORTANT: This should not be run before training. This is only useful for sampling data, and should be combined with subsequent calls to set_beta, set_mean, set_covariance_matrix, set_P and set_Pi.

Parameters:

p : number of channels in X (not used if only Y is modelled) q : number of channels in Y

loglikelihood(X, Y)

Computes the likelihood of the model per state and time point given the data X and Y.

Parameters:

Xarray-like of shape (n_samples, n_parcels)

The timeseries of set of variables 1.

Yarray-like of shape (n_samples, n_parcels)

The timeseries of set of variables 2.

Returns:

Larray of shape (n_samples, n_states)

The likelihood of the model per state and time point given the data X and Y.

Raises:

Exception

If the model has not been trained.

sample(size, X=None, Gamma=None)

Generates Gamma and Y for timeseries of lengths specified in variable size.

Parameters:

sizearray of shape (n_sessions,) or (n_sessions, 2)

If size is 1-dimensional, each element represents the length of a session. If size is 2-dimensional, each row of size represents the start and end indices of a session in a timeseries.

Xarray of shape (n_samples, n_parcels), default=None

The timeseries of set of variables 1.

Gammaarray of shape (n_samples, n_states), default=None

The state probability timeseries.

Returns:

If X is not None:

Xarray of shape (n_samples, n_parcels)

The timeseries of set of variables 1.

Y: array of shape (n_samples,n_parcels)

The timeseries of set of variables 2.

Gammaarray of shape (n_samples, n_states)

The state probability timeseries.

sample_Gamma(size)

Generates Gamma, for timeseries of lengths specified in variable size.

Parameters:

sizearray

Array of shape (n_sessions,) or (n_sessions, 2). If size is 1-dimensional, each element represents the length of a session. If size is 2-dimensional, each row of size represents the start and end indices of a session in a timeseries.

Returns:

Gammaarray of shape (n_samples, n_states)

The state probability timeseries.

set_P(P)

Set transition probability matrix. Useful to create synthetic data for simulations.

Parameters:

P: ndarray of shape (K,K), where K is the number of states

set_Pi(Pi)

Returns initial probabilities. Useful to create synthetic data for simulations.

Parameters:

Pi: ndarray of shape (K,), where K is the number of states

set_beta(beta, k=0)

Sets the regression coefficients (beta) to specific values. Useful to create synthetic data for simulations.

Parameters:

beta: ndarray of shape (n_variables_1 x n_variables_2)

The regression coefficients of each variable in X on each variable in Y for the specified state k.

kint, optional, default=0

The index of the state for which to retrieve the beta value.

set_covariance_matrix(shape, self, k=0)

Sets the covariance matrix to specific values. Useful to create synthetic data for simulations.

Parameters:

ratendarray of shape (n_variables_2 x n_variables_2),

The rate parameter of the covariance

shape : int, the shape parameter of the covariance k : int, optional

The index of the state. Default=0.

set_mean(mean, k=0)

Sets the mean to specific values. Useful to create synthetic data for simulations.

Parameters:

mean: ndarray of shape (n_variables_2,)

The mean value of each variable in Y for the specified state.

kint, optional, default=0

The index of the state for which to retrieve the beta value.

train(X=None, Y=None, indices=None, files=None, Gamma=None, Xi=None, scale=None, options=None)

Train the GLHMM on input data X and Y, which most general formulation is Y = mu_k + X beta_k + noise where noise is Gaussian with mean zero and standard deviation Sigma_k

It supports both standard and stochastic variational learning; for the latter, data must be supplied in files format

Parameters:

Xarray-like of shape (n_samples, n_variables_1)

The timeseries of set of variables 1.

Yarray-like of shape (n_samples, n_variables_2)

The timeseries of set of variables 2.

indicesarray-like of shape (n_sessions, 2), optional

The start and end indices of each trial/session in the input data. If None, one big segment with no cuts is assumed.

filesstr or list of str, optional

The filename(s) containing the data to load. If not None, X, Y, and indices are ignored.

Gammaarray-like of shape (n_samples, n_states), optional

The initial values of the state probabilities.

Xiarray-like of shape (n_samples - n_sessions, n_states, n_states), optional

The joint probabilities of past and future states conditioned on data.

scalearray-like of shape (n_samples,), optional

The scaling factors used to compute the free energy of the dataset. If None, scaling is automatically computed.

optionsdict, optional

A dictionary with options to control the training process.

Returns:

Gammaarray-like of shape (n_samples, n_states)

The state probabilities. To avoid unnecessary use of memory, Gamma is only returned if learning is non-stochastic; otherwise it is returned as an empty numpy array. To get Gamma after stochastic learning, use the decode method.

Xiarray-like of shape (n_samples - n_sessions, n_states, n_states)

The joint probabilities of past and future states conditioned on data. To avoid unnecessary use of memory, Xi is only returned if learning is non-stochastic; otherwise it is returned as an empty numpy array. To get Xi after stochastic learning, use the decode method.

fearray-like

The free energy computed at each iteration of the training process.

Raises:

Exception

If files and Y are both provided or if neither are provided. If X is not provided and the hyperparameter ‘model_beta’ is True.

If ‘files’ is not provided and stochastic learning is called upon