Modelling state-transition dynamics in resting-state brain signals by the hidden Markov and Gaussian mixture models.

Recent studies have proposed that one can summarize brain activity into dynamics among a relatively small number of hidden states and that such an approach is a promising tool for revealing brain function. Hidden Markov models (HMMs) are a prevalent approach to inferring such neural dynamics among discrete brain states. However, the impact of assuming Markovian structure in neural time series data has not been sufficiently examined. Here, to address this situation and examine the performance of the HMM, we compare the model with the Gaussian mixture model (GMM), which is with no temporal regularization and thus a statistically simpler model than the HMM, by applying both models to synthetic time series generated from empirical resting-state functional magnetic resonance imaging (fMRI) data. We compared the GMM and HMM for various sampling frequencies, lengths of recording per participant, numbers of participants and numbers of independent component signals. We find that the HMM attains a better accuracy of estimating the hidden state than the GMM in a majority of cases. However, we also find that the accuracy of the GMM is comparable to that of the HMM under the condition that the sampling frequency is reasonably low (e.g., TR = 2.88 or 3.60 s) or the data are relatively short. These results suggest that the GMM can be a viable alternative to the HMM for investigating hidden-state dynamics under this condition.

expectation–maximization

echo planar imaging

field of view

functional magnetic resonance imaging

Gaussian mixture model

Human Connectome Project

hidden Markov model

independent component analysis

magnetoencephalography

principal component analysis

region of interest

repetition time

INTRODUCTION

Brain dynamics are a product of large‐scale networks realized by interaction of functionally specialized regions in the brain (Deco et al., 2015; Fox et al., 2005; Raichle et al., 2001; Sporns, 2011). Such dynamics have been considered to underpin the integration of information (Tononi et al., 1998), cognitive functions (Bressler & Menon, 2010) and their impairments (i.e., neuropsychiatric disorders) (Menon, 2011). Understanding dynamical coordination of brain regions necessitates data‐analysis methods that reduce the dimension of large‐scale neural data, which are often provided in the form of multivariate time series, without losing much information. Widely used examples include independent component analysis (ICA) (Calhoun & Adalı, 2012) and network analysis (Bassett & Sporns, 2017; Sporns, 2011).

One approach to investigating integrated dynamics of multivariate time‐varying neural signals is to assume a relatively small number of latent states and summarize the multidimensional brain activity data at each time point into one of these states. One can estimate time series of the latent state by, for example, the hidden Markov models (HMMs) (Baker et al., 2014; Brookes et al., 2018; Nielsen et al., 2018; Ryali et al., 2016; Taghia et al., 2017; Vidaurre, 2021; Vidaurre, Abeysuriya, et al., 2018; Vidaurre et al., 2016; Vidaurre et al., 2017; Vidaurre, Hunt, et al., 2018; Warnick et al., 2017), dynamic functional connectivity (Allen et al., 2014; Calhoun et al., 2014; Nielsen et al., 2018) and energy landscape analysis (Ezaki et al., 2017; Watanabe, Hirose, et al., 2014; Watanabe, Masuda, et al., 2014). This strategy allows us to continue to work on the same time domain as the original data and therefore to, for example, compute transition rates between the latent states and interpret state transition events, rather than to reduce the data to static measures (e.g., functional connectivity) or transform the data to the frequency domain implicitly assuming stationarity of the time series. State‐transition dynamics have been reported to be closely related to various functions of the brain, including executive function (Ezaki et al., 2018), decision‐making (Taghia et al., 2018), and to psychiatric conditions such as autism (Watanabe & Rees, 2017) and schizophrenia (Kottaram et al., 2019). For example, Ezaki et al. (2018) analysed resting‐state functional magnetic resonance imaging (fMRI) data obtained from healthy humans using the energy landscape analysis. They showed that the ease of state transitions between synchronized activity patterns of specific regions of interest (ROIs) explained age‐related changes in executive functions. Taghia et al. (2018) applied a Bayesian switching linear dynamical systems model to fMRI data obtained from participants performing an n‐back working memory task. They found a task‐specific hidden state and dynamical switching path of the estimated hidden states. Using the energy landscape analysis, Watanabe and Rees (2017) showed that high functioning autistic adults had atypically stable brain dynamics with lower transition rates among different brain systems and longer dwelling time and that such over‐stability was predictive of both their symptom severity and unique cognitive skills. For resting‐state fMRI data obtained from individuals with schizophrenia, Kottaram et al. (2019) reported aberrant transition dynamics among latent states estimated by an HMM. Relevance of hidden states to task‐related data in a supervised setting has also been reported. Vidaurre, Hunt, et al. (2018) found hidden states representing task‐related brain states using magnetoencephalography (MEG) data during a button press task.

Among these methods, the HMMs have been widely used for studying fMRI (Ryali et al., 2016; Vidaurre, 2021; Vidaurre et al., 2017; Vidaurre, Hunt, et al., 2018; Warnick et al., 2017) and MEG (Baker et al., 2014; Brookes et al., 2018; Taghia et al., 2017; Vidaurre, Abeysuriya, et al., 2018; Vidaurre et al., 2016; Vidaurre, Hunt, et al., 2018) data recorded from the human brain. The HMM is a model comprising a set of probability distributions of the observables each of which corresponds to a latent (i.e., hidden) state and the transition probabilities between the pairs of latent states. By assumption, the state‐transition dynamics of an HMM are Markovian, that is, with no memory effect longer than a single time step. HMMs have been useful in modelling neural dynamics for the following reasons. First, they can be applied to relatively high‐dimensional time series (Baker et al., 2014; Brookes et al., 2018; Vidaurre, Abeysuriya, et al., 2018; Vidaurre, Hunt, et al., 2018). Second, they can detect changes in signals without delay in the form of changes in the latent state, which is not straightforward with dynamic functional connectivity calculated with sliding time windows (Vidaurre, 2021; Vidaurre, Abeysuriya, et al., 2018).

The HMM modelling is based on the implicit assumption that the stochastic rule of state transitions depends on the last state and not on the states in the further past. In fact, the practical impacts of this temporal assumption have not been well studied at least with fMRI data. The potential deviation of the data from the assumed Markovian temporal structure may be detrimental to the HMM modelling.

The same analysis pipeline to infer hidden states of the time series data (Figure 1) can be realized by a simple static mixture model that does not impose any temporal (HMM‐like) regularization. Mixture models have been used in neuroimaging research for detecting activation of brain regions (Everitt & Bullmore, 1999; Hartvig & Jensen, 2000) and clustering them into larger ROIs (Górriz et al., 2009). In Everitt and Bullmore (1999) and Hartvig and Jensen (2000), for example, a mixture model with two probability distributions corresponding to activation and deactivation of each voxel, respectively, was used to classify the state of the voxels. The mixture models have also been used for capturing the state transitions in neuroimaging data. For example, Nielsen et al. (2018) used mixture models for clustering the dynamic functional connectivity patterns. Vidaurre (2021) proposed mixture models for the principal component analysis (PCA) (as well as their HMM variants) to capture dynamic functional connectivity. Because mixture models do not assume any temporal structure, they are not influenced by the sampling frequency and therefore may serve as useful baseline models with which to assess the efficiency of fitting HMMs. If state transitions in the given data are considerably influenced by the current state, an HMM is expected to perform better than a mixture model. In contrast, if state transitions do not depend on the current state, mixture models may outperform HMMs because of their relative simplicity. In general, complex models may overfit to the data, and their model estimation algorithms are often computationally costly and may end up converging to local optima.

FIGURE 1

Overview of the estimation of hidden‐state dynamics using Gaussian mixture models (GMMs) and hidden Markov models (HMMs). (a) A multivariate time series in discrete time such as functional magnetic resonance imaging (fMRI) data. (b) One fits a GMM with two components to the multivariate time series data shown in (a). The case of N = 2 is schematically shown. The estimation of the GMM enables us to associate one of the hidden states (shown in colour) to the data point at each discrete time, x . Using the estimated GMM, one can estimate the time course of the hidden state. (c) One fits an HMM with two hidden states to the same data. In general, how the data points are clustered into two hidden states is different between the GMM and HMM. Using the estimated HMM, one can estimate the time course of the hidden state

On these grounds, we compare the HMM and Gaussian mixture model (GMM) in the present study. We used the GMM because, in HMMs applied to neuroimaging data, the Gaussian distribution is widely used as the probability distribution conditioned on the hidden state (Baker et al., 2014; Brookes et al., 2018; Ryali et al., 2016; Taghia et al., 2017; Vidaurre, Abeysuriya, et al., 2018; Vidaurre et al., 2016; Vidaurre et al., 2017; Vidaurre, Hunt, et al., 2018; Warnick et al., 2017). Then, if we ignore state‐transition dynamics as described by a hidden Markov process, the distribution of signals as estimated by the fitting of an HMM is a GMM.

We compare HMMs and GMMs on synthetic data sets that we generate from human resting‐state fMRI data. Empirical resting‐state fMRI data vary in the time resolution, corresponding to the TR (i.e., repetition time), the length of the recording time (i.e., number of image volumes) and the number of participants. Therefore, we examine how these three factors influence the performance of the HMMs and GMMs in estimating the hidden states of fMRI signals over time. We show that, while the HMM outperforms the GMM when the sampling frequency is high (e.g., TR = 0.72 s), the advantage of the HMM diminishes as the sampling frequency decreases (i.e., larger TR). In addition, when the recording time is short (i.e., a small number of volumes or participants), the GMM more robustly estimates the hidden states than the HMM especially when the dimension of the signal is large. Our results provide quantitative guidance on when HMMs work better than GMMs or vice versa.

MATERIALS AND METHODS

An overview of the analysis pipeline is shown in Figure 1.

Gaussian mixture model

Assume that there are t max observations of N‐dimensional data, x (t = 1, …, t max) (Figure 1a). Our GMM assumes that each observation, x , is generated from one of the two Gaussian distributions (Figure 1b). We set the number of states to two for simplicity. The probability density of the observed data conditioned on the state is given by , where s  ∈ {1, 2} is the hidden (i.e., latent) state, denotes a N‐dimensional multivariate Gaussian distribution, and and are the mean and covariance matrix of the Gaussian distribution under hidden state s , respectively. The marginal probability distribution of x is given by

where π is the probability that hidden state s is taken. One estimates π , and Σ (s=1, 2) by maximizing the log‐likelihood function given by

We used the expectation–maximization (EM) algorithm, which is typically used for maximizing Equation (2) (Dempster et al., 1977; Lindsay, 1995).

Then, the time course of the hidden state, (t = 1, …, t max), given the observations is estimated by where and are the maximum likelihood estimator obtained by the EM algorithm.

We analyse the data using the GMM package in scikit‐learn (Pedregosa et al., 2011). We used the default setting of scikit‐learn for determining the initial conditions for the EM algorithm.

Hidden Markov model

We consider HMMs with Gaussian components (Ephraim & Merhav, 2002) (Figure 1c). The model assumes that each of the N‐dimensional observations x (t = 1, ..., t max) is generated from one of the two Gaussian distributions, as in GMMs, and that s (t = 1, ⋯, t max) obeys first‐order Markovian dynamics given by

To estimate the HMM, we used an EM algorithm known as the Baum–Welch algorithm (Baum et al., 1970). We used the Viterbi algorithm (Viterbi, 1967) to estimate the time course of the hidden state given the observations. We estimated HMMs for our data using a python package hmm‐learn (https://hmmlearn.readthedocs.io/), which was originally developed as part of scikit‐learn (Lindsay, 1995). We used the default setting of hmm‐learn for determining the initial conditions for the EM algorithm.

The EM algorithm does not guarantee the exact optimization due to local minima. Therefore, for both GMM and HMM, we carried out the optimization procedure 10 times and adopted the model that attained the largest likelihood. In the present study, all the estimations successfully converged.

Resting‐state fMRI data

We used ‘S1200 extensively processed rfMRI data’ provided by the Human Connectome Project (HCP, https://www.humanconnectome.org/) (Van Essen et al., 2012). The data collection was approved by the local ethics committees that participated in HCP, and all participants provided written consent. The data set includes the group ICA data for 1003 participants (22–35 years old, 534 females). In this section, we briefly describe the preprocessing procedure to obtain this data set. See Smith et al. (2013), Griffanti et al. (2014) and Salimi‐Khorshidi et al. (2014) for the details of the preprocessing.

The participants completed four sessions of 15‐min echo planar imaging (EPI) sequence on a 3T Siemens Connectome‐Skyra (TR = 0.72 s, TE = 33.1 ms, 72 slices, 2.0 mm isotropic; field of view (FOV) = 208 × 180 mm) and a single T1‐weighted sequence (TR = 2.4 s, TE = 2.14 ms, 0.7 mm isotropic, FOV = 224 × 224 mm). Each session yielded 1200 volumes (i.e., observations) of EPI images. The fMRI data were first minimally preprocessed according to Smith et al. (2013). Then, artefacts were removed using ICA + FIX (Griffanti et al., 2014; Salimi‐Khorshidi et al., 2014) and inter‐participant registration of cerebral cortex using MSMAII (Glasser et al., 2016; Robinson et al., 2014). Then, the group ICA was performed by the MELODIC's incremental group PCA algorithm (Beckmann & Smith, 2004; Hyvärinen, 1999). The group ICA was carried out for dimensions (i.e., number of independent components) N = 15, 25, 50, 100, 200 and 300. We did not use group‐ICA data generated with N=100, 200 and 300 because of high computational cost for estimating the models used in the present study.

We fed the fMRI data to algorithms for estimating GMMs or HMMs after concatenating the observed signals obtained from all the sessions from a single participant into one sequence and then the sequences obtained from all the participants into one sequence. In the concatenated data, the final volume of Sessions 1, 2 and 3 is followed by the first volume of Sessions 2, 3 and 4, respectively, although different sessions are not causally related to each other. In practice, the influence of the concatenation on the estimation of the HMM is considered to be negligible because each session is sufficiently long (i.e., 1200 volumes) (Vidaurre, Abeysuriya, et al., 2018; Vidaurre et al., 2017).

GMM‐based synthetic time series

We prepared synthetic fMRI data with an underlying hidden‐state dynamics for each of N (= 15, 25 and 50) as follows. First, we fit the GMM with two hidden states to the entire N‐dimensional fMRI data that are composed of 1200 volumes × 4 sessions × 1003 participants. As a result, we obtained the mean vector and covariance matrix of each Gaussian distribution, which we denote by and (k= 1, 2) and the probability of appearance of each state, π (see Equation 1). Then, by maximizing the likelihood given by Equation (2) given the estimated model, we estimated the hidden states for the first session (but not for the second to fourth sessions) of the fMRI recording, which we denote by (with t = 1, …, 1200), for each participant p (with p = 1, …, 1003).

Second, using , and (k=1, 2), we generated N‐dimensional (= 15, 25 or 50) synthetic signals, for p = 1, …, 1003 and t = 1, …, 1200, that is,

Then, we subsampled the signals generated by Equation (5) to prepare synthetic signals with typical lengths of a resting‐state fMRI recording, that is, T = 5, 10 and 14.4 min, n pat participants and a sampling frequency of TR = 0.72, 1.44, 2.16, 2.88 and 3.60. Note that T = 14.4 min is the original length of a single session in the HCP data. We used the subsampled data to test the performance of the GMM and HMM.

We generated synthetic fMRI data with the given TR value and the number of participants, n pat, by subsampling as follows. Consider the case T = 10 min. First, for each participant p (with p = 1, …, n pat), the samples (t = 1, 2, …, 833) with the associated hidden state labels provide test data with TR = 0.72 s. Note that, given TR = 0.72 s, time t = 833 corresponds to 10 min starting from the first volume, that is, t = 1. To generate test data with TR = 1.44 s, we use subsampled data , which is of length 416. Similarly, the subsampled series defines synthetic data of length 278 with TR = 2.16 s; the series defines synthetic data of length 209 with TR = 2.88 s; and the series defines synthetic data of length 167 with TR = 3.60 s. We apply the same subsampling method to generate synthetic data with a range of TR when the assumed recording time is different (i.e., T = 5 or 14.4 min). For example, if T = 5 min, which corresponds to t = 416, the data for TR = 1.44 s are given by .

HMM‐based synthetic time series

We also used another set of synthetic time series data, which we call the HMM‐based synthetic time series. Similarly to the generation of the GMM‐based synthetic time series, we fitted the HMM with two hidden states and the corresponding Gaussian distributions to the entire fMRI data composed of 1200 volumes × 4 sessions × 1003 participants. The fitting yielded the mean vector and covariance matrix of the N‐dimensional Gaussian distribution corresponding to each state, which we denoted by and (k = 1, 2), the transition rates between the two hidden states, and the probability distribution of the initial state. Using the estimated HMM and Viterbi algorithm (Viterbi, 1967), we estimated the hidden state labels for the first session of each participant, which we denote by (with p = 1, …, 1003 and t = 1, …, 1200). The rest of the procedure is same as that for the GMM‐based synthetic time series.

RESULTS

Hidden states estimated from the fMRI data

First, we compared between the GMM and HMM the Gaussian distributions and the time courses of the hidden state inferred from the empirical fMRI data. The Gaussian distributions estimated for the GMM, parametrized by and (with k = 1, 2), and those for the HMM, parametrized by and , were similar to each other, but the means estimated by the HMM were farther from 0 than those estimated by the GMM (Figure 2a). Figure 2b compares the estimated time courses of the hidden state for a single participant between the GMM and HMM. Although the hidden state is the same between the GMM and HMM in a majority of times (77% for N = 15; 77% for N = 25; 82% for N = 50), the GMM yielded more frequent state flips (the fraction of times at which the state flipped was 0.17 for N = 15, 0.19 for N = 25 and 0.17 for N = 50) than the HMM (0.05 for N = 15; 0.04 for N = 25; 0.04 for N = 50). Consistent with this observation, the probability distribution of the duration of each hidden state has a higher mass at small values for GMMs than HMMs (Figure 2c). Figure 2c also indicates that the distribution of the duration of each hidden state has a peak for the HMM but not for the GMM.

FIGURE 2

Comparisons between the Gaussian mixture model (GMM) and hidden Markov model (HMM) fitted to the functional magnetic resonance imaging (fMRI) data. (a) Mean vectors ( vs. ) and covariance matrices ( vs. ). A circle represents each entry of the mean vector or covariance matrix in the estimated GMM and HMM. The solid lines represent the diagonal. (b) A sample time course of the estimated hidden state labels, and , for 10 min and a single participant. (c) Distribution of the duration of a hidden state, computed based on the entire sequences of the hidden state, and (1 ≤ p ≤ 1003,1 ≤ t ≤ 1200)

In the next subsections, we compare the performance of the GMM and HMM using the two types of synthetic data (i.e., the GMM‐based synthetic time series and the HMM‐based synthetic time series).

Results for the GMM‐based synthetic time series

We fit the GMM and HMM to GMM‐based synthetic time series with various TR values, three values of the recording time per participant (i.e., T), various numbers of participants pooled (i.e., n pat) and three values of N and compared the accuracy of fitting. Note that the true hidden state at each time is available in this numerical experiment. We defined the accuracy of estimation as the fraction of times at which the estimated state is the same as the true state. We show the accuracy of estimation for T = 10 min in Figure 3. For both GMM and HMM, the accuracy was low for a small number of participants (n pat < 10), increased as n pat increased and saturated before n pat ≈ 20 at > 83%. In general, more participants were required for attaining a given accuracy value when N is larger or the sampling frequency is lower (i.e., the TR is larger). For the HMM but not for the GMM, the accuracy increased as the sampling frequency increased. When the sampling frequency is high (i.e., small TR) and n pat is relatively large, the HMM outperformed the GMM (Figure 3a). However, this advantage of the HMM diminishes when the TR is large, n pat is small or N is large (Figure 3a–e). For example, when N=25 and TR = 2.88 s, the GMM outperformed the HMM for 1 ≤ n pat ≤ 11 (Figure 3e).

FIGURE 3

Accuracy of estimating the hidden state for the Gaussian mixture model (GMM)‐based synthetic time series for T = 10 min. (a) TR = 0.72. (b) TR = 1.44. (c) TR = 2.16. (d) TR = 2.88. (e) TR = 3.60. The shaded regions represent one standard deviation

For T = 5 and 14.4 min, the results were qualitatively the same as those for T = 10 min (Figures S1 and S2). Because T is proportional to the length of the data as n pat is, the accuracy is higher for a larger T value. It should also be noted that, with T = 5 min, the GMM outperformed the HMM in wider parameter regions than with T = 10 min.

Results for the HMM‐based synthetic time series

We applied the same analysis as that in the previous section to HMM‐based synthetic time series. The accuracy of estimating the true hidden state is shown in Figure 4 for the GMM and HMM estimators, various values of TR, T = 10 min, various values of n pat and three values of N. The results were similar to those obtained for the GMM‐based synthetic time series. When TR is large, N is large and n pat is small, the accuracy of estimation for the GMM was comparable to that for the HMM. Otherwise, the HMM outperformed the GMM. With the HMM‐based synthetic time series, the HMM yielded a higher accuracy than the GMM in a wider parameter region than with the GMM‐based synthetic time series.

FIGURE 4

Accuracy of estimating the hidden state for the hidden Markov model (HMM)‐based synthetic time series for T = 10 min. (a) TR = 0.72. (b) TR = 1.44. (c) TR = 2.16. (d) TR = 2.88. (e) TR = 3.60. The shaded regions represent one standard deviation

We also carried out the same analysis with T= 5 and 14.4 min to confirm that the results were qualitatively the same as those with T = 10 min and that the accuracy was larger for a larger T in general (Figures S3 and S4).

DISCUSSION

We compared the performance of estimating the time course of the hidden state between the GMM and HMM. On the two types of synthetic time series that we tested, the HMM inferred the true hidden state more accurately than the GMM in the majority of cases. However, when the length (i.e., n pat × T) and sampling frequency of the data are limited, the accuracy of estimation for the GMM was comparable to or even better than that for the HMM. The performance of the GMM and HMM estimators also depended on the underlying dynamics of the hidden state, that is, the type of synthetic time series.

The accuracy of estimating the hidden state using HMMs depended on the sampling frequency of the data. In practice, HMMs are efficient with a reasonable sampling frequency with which the contribution of the first‐order dynamics dominates (Martinez‐Diaz et al., 2007). In our results, the HMM outperformed the GMM in the case of the highest sampling frequency (i.e., TR = 0.72 s) on both the GMM‐ and HMM‐based synthetic time series. Moreover, for the HMM‐based synthetic time series, the accuracy score for the HMM was larger when TR was smaller (blue lines in Figure 3). These results suggest that the HMM will probably perform better than the GMM if TR is less than 0.72 s. Note that TR = 0.72 s or similar sampling frequencies are becoming common in fMRI experiments. Furthermore, with TR = 2.88 and 3.60 s, which also fall in a range of TR that remains common in human fMRI experiments, the HMM and GMM performed similarly given a sufficient amount of data. This result suggests that, in this situation, we may benefit from using the GMM, which is conceptually, mathematically and algorithmically simpler than the HMM.

When the given data are short and sampled with a low frequency (e.g., TR = 3 s), GMMs are probably more advantageous than HMMs in that a GMM with the same number of components as an HMM is better at avoiding overfitting. Furthermore, GMMs seem to be also advantageous over HMMs when the given data are long or the estimation procedure has to be run many times. Up to our numerical effort, estimation of GMMs was at least a couple of times faster than that of HMMs including when we attempted to estimate models with more than two states. This observation is consistent with a previous study reporting that estimation of the HMM is often computationally challenging when the data is long or the number of participants is large (Vidaurre, Abeysuriya, et al., 2018).

We refrained from optimizing the number of hidden states. This is because there is no established way to do so (Celeux & Durand, 2008; Pohle et al., 2017), although some methods based on, for example, the free energy in the variational Bayes (VB) algorithms (Rezek & Roberts, 2005) and the number of appearances of each state in the estimated hidden‐state time courses (Vidaurre, Hunt, et al., 2018) have practically been used. (Also see Nielsen et al., 2018 for an infinite HMM approach, which finds the appropriate number of hidden states but is not computationally feasible in many practical cases.) In previous studies using fMRI data, the estimated numbers of hidden states are distributed in a wide range, that is, between 5 and 19 (Kottaram et al., 2019; Ryali et al., 2016; Scofield et al., 2019; Stevner et al., 2019; Vidaurre, Abeysuriya, et al., 2018; Vidaurre et al., 2017). In contrast to these studies, we assumed two hidden states for the sake of simplicity. This choice was also motivated by a previous study reporting that the hidden states were robustly agglomerated into two clusters in human fMRI data and that the frequency of the two states was heritable and related to cognitive measures (Vidaurre, Hunt, et al., 2018). We also found in our previous work with energy landscape analysis that transitions among two or three macroscopic states were correlated with participants' behaviour in a bistable visual perception task (Watanabe, Masuda, et al., 2014) and executive function (Ezaki et al., 2018). Therefore, we believe that characterizing brain dynamics by transitions among an a priori determined small number of states, as we have done in the present study, is a useful approach.

The performance of the GMMs and HMMs may depend on various factors, which we have not examined in this study, such as the preprocessing method, scanner and tasks. In particular, the type of task is expected to modulate the hidden‐state dynamics and may inform the choice of the model with state dynamics. We observed that the duration of the hidden state is qualitatively different between the GMM and HMM (Figure 3c). Note that various generalizations of the HMM (i.e., hidden semi‐Markov models) (Yu, 2010) have been proposed to incorporate different types of the distribution of the duration. In the absence of the ground‐truth data for the hidden state dynamics, which is generally the case, we do not have a particular hypothesis regarding the shape of the distribution of the hidden state's duration. Instead, in this study, we focused on the effect of other key parameters on the estimation accuracy of the GMM and HMM. Additional information that accompanies recorded brain signals, such as behavioural switches during a task, may help us to better compare the GMM, HMM and their variants.

CONFLICT OF INTEREST

Authors declare no conflict of interest.

AUTHOR CONTRIBUTIONS

TE and NM designed the research. TE and YH performed numerical experiments. TE, TW and NM discussed the results. TE, TW and NM wrote the paper.

PEER REVIEW

The peer review history for this article is available at https://publons.com/publon/10.1111/ejn.15386.

Figure S1. Accuracy of estimating the hidden state for the GMM‐based synthetic time series for T = 5 min. (a) TR = 0.72. (b) TR = 1.44. (c) TR = 2.16. (d) TR = 2.88. (e) TR = 3.60. The shaded regions represent one standard deviation.

Figure S2. Accuracy of estimating the hidden state for the GMM‐based synthetic time series for T = 14.4 min. (a) TR = 0.72. (b) TR = 1.44. (c) TR = 2.16. (d) TR = 2.88. (e) TR = 3.60. The shaded regions represent one standard deviation.

Figure S3. Accuracy of estimating the hidden state for the HMM‐based synthetic time series for T = 5 min. (a) TR = 0.72. (b) TR = 1.44. (c) TR = 2.16. (d) TR = 2.88. (e) TR = 3.60. The shaded regions represent one standard deviation.

Figure S4. Accuracy of estimating the hidden state for the HMM‐based synthetic time series for T = 14.4 min. (a) TR = 0.72. (b) TR = 1.44. (c) TR = 2.16. (d) TR = 2.88. (e) TR = 3.60. The shaded regions represent one standard deviation.

Click here for additional data file.