Quantifying the relationship between the brain’s functional activity patterns and its structural backbone is crucial when relating the severity of brain pathology to disability in multiple sclerosis (MS). Network control theory (NCT) characterizes the brain’s energetic landscape using the structural connectome and patterns of brain activity over time. We applied NCT to investigate brain-state dynamics and energy landscapes in controls and people with MS (pwMS). We also computed entropy of brain activity and investigated its association with the dynamic landscape’s transition energy and lesion volume. Brain states were identified by clustering regional brain activity vectors, and NCT was applied to compute the energy required to transition between these brain states. We found that entropy was negatively correlated with lesion volume and transition energy, and that larger transition energies were associated with pwMS with disability. This work supports the notion that shifts in the pattern of brain activity in pwMS without disability results in decreased transition energies compared to controls, but, as this shift evolves over the disease, transition energies increase beyond controls and disability occurs. Our results provide the first evidence in pwMS that larger lesion volumes result in greater transition energy between brain states and decreased entropy of brain activity.

We investigated the brain-state dynamic and energy landscapes in healthy individuals and people with multiple sclerosis (pwMS). We also investigated the entropy of brain activity and its association with transition energy between brain states and lesion volume. We clustered regional brain activity time series to identify the brain states. Then, we applied network control theory using structural connectivity network to identify the minimum required energy to transition between brain states. We observed the pwMS without disability showed decreased transition energy, while pwMS with evidence of disability showed increased transition energy compared to healthy individuals. Lower entropy of brain activity was associated with greater lesion load and larger transition energy. This study provides a possible mechanism of how MS-related damage to the brain’s structural backbone can impact brain dynamics, entropy, and energetics.

Multiple sclerosis (MS) is a chronic disease characterized by neuroinflammation and, eventually, neurodegeneration in the central nervous system (CNS). The size and location of the lesions in the CNS are very heterogeneous (Barkhof, 2002), resulting in different patterns of structural damage among people with MS (pwMS). Structural damage to the brain may result in changes to its functional activity patterns (Tommasin et al., 2018); however, how different patterns of structural damage can modify dynamics of functional activity have not been fully characterized in pwMS.

Functional MRI (fMRI) measures functional brain activity patterns over time, information which can be used to quantify brain activity dynamics or to estimate patterns of regional co-activation, that is, functional connectivity (FC). Both increased or decreased fMRI-based activity/connectivity has been observed in pwMS, where increases are often interpreted as a compensatory mechanism serving to limit clinical disability (Buyukturkoglu et al., 2020; Filippi et al., 2013; Sbardella et al., 2015). FMRI in pwMS has largely been used to investigate the association between the brain’s static and/or dynamic FC and motor/cognitive impairment. (Buyukturkoglu et al., 2020; Tozlu et al., 2021a; Tozlu et al., 2021b). Most FC studies in pwMS have been static in nature (Filippi et al., 2013; Sbardella et al., 2015), ignoring shorter scale changes in FC that have been shown to occur (Calhoun et al., 2014). These shorter scale changes in FC can be captured via dynamic FC (dFC) analysis that involves clustering FC matrices computed via sliding windows of the fMRI time series data (Allen et al., 2014). This approach was previously used to differentiate between the pwMS and healthy controls (HC) (Leonardi et al., 2013; Rocca et al., 2019; Tozlu et al., 2021a), to compare the dynamics between cognitively impaired versus preserved pwMS (d’Ambrosio et al., 2019), and to investigate the relationship between alterations in dFC dynamics and cognition in MS (Fuchs et al., 2022; van Geest et al., 2018a; van Geest et al., 2018b). However, calculating dFCs with a sliding-window approach is not straightforward since the estimation of the window length and shift size of the sliding windows is required; moreover, correlations are still used to describe coactivations of pairs of regions as in static FC. As an alternative approach, clustering can be directly applied to the time series of regional activity to define distinct brain activity states. This approach was successfully used in HC to identify brain activity states that occur while performing a working memory task, quantify the effects of psychedelics, and identify changing dynamics of activity patterns in stroke patients (Cornblath et al., 2020; Olafson et al., 2022; Singleton et al., 2022). In addition to the identification of the brain’s dynamic states, the same studies used network control theory (NCT) approach (Gu et al., 2015; Stiso et al., 2019; E. Tang et al., 2017) to identify the minimum energy required to transition between these dynamic brain states. However, no study to date has applied brain activity clustering or NCT in a population of pwMS, let alone investigate differences in brain dynamics and energetics across disability subgroups in MS.

In terms of capturing brain activity dynamics, entropy is also an important metric with which to quantify the amount of regularity/unpredictability in brain activity time series (Pincus, 1991). Higher entropy of brain activity has been associated with a larger repertoire of available brain states (Carhart-Harris et al., 2014), and decreased NCT-based transition energy under psychedelics compared to placebo was found to be associated with greater increases in entropy of brain states (Singleton et al., 2022). Studies in disease/disorders have shown (a) greater entropy in controls compared to people with attention-deficit/hyperactivity disorder (ADHD) as well as a correlation between lower entropy and more severe ADHD symptom scores (Sokunbi et al., 2013) and (b) decreased entropy in regions containing a stroke lesion and their contralesional hemisphere homologues (Saenger et al., 2018). One study in relapsing-remitting MS (RRMS) showed entropy was associated with an individual’s disability level, where more disability was related to increased entropy (Sokunbi et al., 2013; Zhou et al., 2016). However, no study to date has investigated the link between the entropy of brain activity and (a) lesion volume or (b) NCT-based transition energy between brain states in MS.

In this study, we first aimed to identify recurrent dynamic brain states in both HC and pwMS by unsupervised clustering of regional resting-state fMRI activity. Note that, although related (see Olafson et al., 2022), this approach differs from dFC (sliding-window FC) in that instead of clustering dFC matrices we are clustering the brain activity vectors from each TR to identify commonly occurring patterns of brain activation (not FC). We calculated metrics quantifying the dynamics of these brain activity states, including transition probability and fractional occupancy, and used NCT and individuals’ structural connectomes (SC) to compute the energy required to transition between pairs of brain states. We hypothesized that greater energy is required for brain-state transitions in pwMS with disability compared to both controls and pwMS with no evidence of disability. Finally, we investigated the association between overall transition energy and the information contained in the brain activity signal, as captured by entropy. We calculated the entropy of each region’s activity over the entire fMRI scan and summarized the global entropy by taking the average of the regional measures. Our hypothesis was that decreased global entropy is associated with both increased transition energy and lesion load. This paper is the first to quantify MS-related shifts in the brain’s energetic landscape and to link entropy/energetic demand of brain activity to lesion burden in pwMS.

Subjects

One hundred pwMS (age: 45.5 [36.7, 56.0], 66% females) with a diagnosis of clinically isolated syndrome (CIS)/MS (CIS = 7, RRMS = 88, primary/secondary progressive primary progressive MS [PPMS] and secondary progressive MS [SPMS] = 5) and 19 HC (age: 45 [35, 49], 55% females) were included in our study. All MS subgroups were included in the analyses as we wanted to investigate these imaging biomarkers across the spectrum of disability severities, which are often related to disease subtype. All subjects who had preprocessed imaging data available included in our study. MRIs and demographic data were collected (age, sex, and race for both HC and pwMS, clinical phenotype and disease duration for pwMS). The expanded disability status scale (EDSS) was used to quantify disability in pwMS, where an EDSS of 2 was used as a threshold to categorize the pwMS into no disability (EDSS < 2) or evidence of disability (EDSS ≥ 2) groups. An EDSS of greater than or equal to 2 is the cutoff for clinical definitions of “disability,” whereas scores below that are not considered to be disabled and scores above it correspond to minimal, moderate, or severe disability. All studies were approved by an ethical standards committee on human experimentation and written informed consent was obtained from all patients. Participants were excluded if they had contraindications to MRI, and controls were further excluded if they had ever been diagnosed with or were currently on medication for a neurological or psychological disorder.

Image Acquisition, Processing, and Structural Connectome Extraction

MRI data were acquired on a 3T Siemens Skyra scanner (Siemens, Erlangen, Germany) with a 20-channel head-neck coil. Anatomical MRI (T1/T2/T2-FLAIR, 1-mm3 isovoxel), resting-state fMRI (6 min, TR = 2.3 s, 3.75 × 3.75 × 4 mm voxels) and diffusion MRI (55 directions HARDI, b = 800, 1.8 × 1.8 × 2.5 mm voxels) acquisitions were performed. Multi-echo 2D GRE field maps were collected for use with both fMRI and diffusion MRI (0.75 × 0.75 × 2 mm voxels, TE1 = 6.69 ms, ΔTE = 4.06 ms, number of TEs = 6).

White matter (WM) and gray matter (GM) were segmented and GM further parcellated into 86 FreeSurfer-based regions (68 cortical and 18 subcortex/cerebellum) (Fischl & Dale, 2000) and a fine-grained atlas cc400 (321 cortical and 71 subcortex/cerebellum) using FreeSurfer (Craddock et al., 2012). The white and gray matter surfaces were manually checked and hand edited with control points and reconstruction editing if necessary. As described elsewhere (Kuceyeski et al., 2016), fMRI preprocessing included simultaneous nuisance regression and removal of WM and cerebrospinal fluid (CSF) effects (Hallquist et al., 2013), followed by high-pass filtering (≥0.008) using the CONN v18b toolbox (Whitfield-Gabrieli & Nieto-Castanon, 2012) and SPM12 in Matlab. Nuisance regressors included 24 motion parameters (6 rotation and translation, temporal derivatives, and squared version of each) and the top 5 eigenvectors from eroded masks of both WM and cerebro-spinal fluid. The mean fMRI signal over all voxels in each of the 86 regions for each TR was calculated to obtain a regional time series matrix of BOLD activity.

Diffusion MRI was interpolated to isotropic 1.8-mm voxels, and then corrected for eddy current, motion, and EPI distortion with the eddy command from FSL 5.0.11 (Andersson & Sotiropoulos, 2016) using the outlier detection and replacement option (Andersson et al., 2016). MRtrix3Tissue (https://3tissue.github.io), a fork of Mrtrix3 (Tournier et al., 2019), was used to estimate a voxel-wise single-shell, three-tissue constrained spherical deconvolution model (SS3T-CSD) and then compute whole-brain tractography for each subject. We performed deterministic (sd-stream) tractography with MRtrix3 (Tournier et al., 2010; Tournier et al., 2012) with uniform seeding at each white matter voxel, after which the SIFT2 global filtering algorithm (Smith et al., 2015) was applied to account for bias that exists in greedy, locally optimal streamline reconstruction algorithms. The SC matrix was constructed by taking the sum of the SIFT2 weights of streamlines connecting pairs of regions and then dividing by the sum of the two regions’ volumes.

Lesion masks were produced for each pwMS in a semiautomated process. The WM hyperintensity lesion masks were created by running the T2 FLAIR images through the lesion segmentation tool (LST) and were further hand edited as needed. T2 FLAIR-based lesion masks were transformed to the individual’s T1 native space using the inverse of the T1 to GRE transform and nearest neighbor interpolation. Individual T1 images were then normalized to MNI space using FSL’s linear (FLIRT) and nonlinear (FNIRT) transformation tools (https://www.fmrib.ox.ac.uk/fsl/index.html); transformations with nearest neighbor interpolation were then applied to transform the native anatomical space T2FLAIR lesion masks to MNI space. The transformations were concatenated (T2FLAIR to T1 to MNI) to minimize interpolation effects. Lesions were manually inspected after the transformation to MNI space to verify the accuracy of coregistration and lesion volume (in mm3) calculated.

Brain-State Identification

K-means clustering was used to identify commonly recurring brain states (see Figure 1 for the workflow of the study) over all subjects’ normalized regional brain activity time series. We used both the elbow criterion and the gain in explained variance to identify the optimal number of clusters, which we identified to be six (see Supporting Information Figure S2). Once the optimal number of clusters was identified, clustering was run 10 times. Within each of the 10 runs, clustering was restarted 50 times with different random initializations to identify the final clusters having the best separation, that is, those that minimized the Pearson correlation between pairs of cluster centroids. Adjusted mutual information (AMI) (Nguyen, 2021) was computed to assess the stability of the clustering results over the 10 repetitions; AMI was greater than 0.89 for all pairs of repetitions indicating robust cluster assignment (see Supporting Information Figure S1).

Figure 1.

The workflow of the study. (A) Regional BOLD time series for all subjects (19 HC and 99 pwMS) were collected and k-means was applied to identify commonly occurring brain states. This approach allows us to assign a state to each TR in each individual’s scan. Metrics such as fractional occupancy and transition probability for each state were then calculated. (B) Network control theory was applied to compute the minimum energy required to transition between pairs of brain states or to remain in the same state. Logistic regression with ridge regularization was used to predict class (HC, pwMS with no disability, pwMS with disability) using state transition/persistence energy. (C) Entropy of regional BOLD time series were also computed, compared across the three groups and global entropy was associated with global transition energy and global disease burden (log lesion volume).

Figure 1.

The workflow of the study. (A) Regional BOLD time series for all subjects (19 HC and 99 pwMS) were collected and k-means was applied to identify commonly occurring brain states. This approach allows us to assign a state to each TR in each individual’s scan. Metrics such as fractional occupancy and transition probability for each state were then calculated. (B) Network control theory was applied to compute the minimum energy required to transition between pairs of brain states or to remain in the same state. Logistic regression with ridge regularization was used to predict class (HC, pwMS with no disability, pwMS with disability) using state transition/persistence energy. (C) Entropy of regional BOLD time series were also computed, compared across the three groups and global entropy was associated with global transition energy and global disease burden (log lesion volume).

Close modal

Finally, we replicated the cluster analysis using a higher dimensional atlas (cc400), where optimal k was identified as 5 (see Supporting Information Figure S7). The results obtained with the cc400 atlas were interpreted in the main text and the figures associated with those results were given in the supplementary materials.

Each of the 86 regions was assigned to one of nine networks, including the seven Yeo functional networks, a subcortical and a cerebellar network (Yeo et al., 2011). The network-level contributions of the six centroids were characterized by calculating the cosine similarity of each state’s centroid (positive and negative parts were analyzed separately) and nine binarized vectors representing each of the networks. Because the mean signal from each scan’s regional time series was removed during high-pass filtering, positive values reflect signal intensity above the mean (high-amplitude) and negative values reflect signal intensity below the mean (low-amplitude). The brain states were named by identifying the network with the largest magnitude activation (either high+ or low- amplitude). The following dynamic brain-state metrics were calculated: (a) transition probability between pairs of states and the persistence probability of remaining in the same state and (b) fractional occupancy, defined as the number of TRs assigned to each cluster out of the total number of TRs (see Table 1).

Table 1.

The MRI-based metrics used in the study

MetricDefinitionInterpretation
Transition probability The probability that a brain switches to a different state or remains in the same state (i.e., persistence probability) at the next time point A higher transition probability between two states means this transition is more likely to occur than others. 
Fractional occupancy Number of TRs assigned to each state out of the total number of TRs Higher fractional occupancy in a specific state means more time is spent in this state compared to another states. 
Transition energy The minimum energy required to transition between pairs of brain states or to remain in the same state over some time Higher transition energy means the transition requires more energy to complete. Transition energy is highly negatively correlated with transition probability. 
Entropy A measure indicating the complexity and predictability of a dynamic system Higher entropy means increased complexity and decreased predictability of a system. 
MetricDefinitionInterpretation
Transition probability The probability that a brain switches to a different state or remains in the same state (i.e., persistence probability) at the next time point A higher transition probability between two states means this transition is more likely to occur than others. 
Fractional occupancy Number of TRs assigned to each state out of the total number of TRs Higher fractional occupancy in a specific state means more time is spent in this state compared to another states. 
Transition energy The minimum energy required to transition between pairs of brain states or to remain in the same state over some time Higher transition energy means the transition requires more energy to complete. Transition energy is highly negatively correlated with transition probability. 
Entropy A measure indicating the complexity and predictability of a dynamic system Higher entropy means increased complexity and decreased predictability of a system. 

Transition Energy Calculations

As previously described (Cornblath et al., 2020; Singleton et al., 2022), NCT can be used to understand how the structural (white matter) connectome constrains dynamic brain-state changes. Specifically, it can be used to compute the minimum energy (Em) required to transition between pairs of brain states or to remain in the same state (see Table 1). To begin, we employed a time-invariant model:
s˙t=Ast+But
where A is the structural connectivity matrix with a dimension of 86 × 86, s(t) is a vector of length 86 containing the regional activity at time t, and u(t) is the external input to the system. In our case, B is an identity matrix with a dimension of 86 × 86 that indicates each region has uniform control over the system.
When computing the transition energy (Em) from the initial activity s0 to final activity sf over some time T, we use an invertible controllability Gramian W to control the structural connectivity network A from the set of network nodes B. The invertible controllability Gramian W is defined as
W=0TeATτBBtreAtrTτdτ
Optimal T was identified by grid searching between 0.001 and 10 and identifying the value that maximized the magnitude of the Spearman correlation between the entries in the transition probability matrix and the entries in the transition energy matrix. The correlation between transition probability and energy is expected to be negative since the brain prefers trajectories through state space requiring minimal input energy given structural constraints, therefore the brain tends to change less between states if the minimum required energy is greater. Here, we identified an optimal T of 1, which yielded a correlation of −0.83 (p value < 0.001) (see Supporting Information Figure S5), while the optimal T was 0.5 for the cc400 atlas with a correlation coefficient of −0.91 (see Supporting Information Figure S12).
After computing the invertible controllability Gramian W and optimizing T, the transition energy (Em) is then calculated as the quadratic product between the inverted controllability Gramian W and the difference between eATs0, which represents the activation at time T assuming brain activity spreads in a network diffusion manner over the SC from the initial state s0, and the final state sT:
Em=eATs0sTtrW1eATs0sT.
Since there were six brain states, the transition energy matrix is of size 6 × 6 where off-diagonal entries contain the energy required to transition between each pair of unique states s0 and sT, where s0sT, and diagonal entries contain persistence energy required to stay in the same state (s0 = sT).

Entropy of Brain Activity

To investigate how disease burden (lesion load) and disability level in pwMS were related to the amount of information contained in the brain activity signal, regional entropy was calculated directly on the BOLD time series. Denoting the regional BOLD time series in a region by s(t) = (s(1), …, s(n)) where n is the number of the BOLD time series acquisition. We define two template vectors of length m such as sm,i = (s(i), s(i + 1), …, s(i + m − 1)) and sm,j = (s(j), s(j + 1), …, s(j + m − 1)). Entropy was calculated as the negative logarithm of the probability that if two template vectors of length m have distance less than r (i.e., d[sm,i, sm,j] < r) then two template vectors of length m + 1 also have distance less than r (i.e., d[sm+1,i, sm+1,j] < r) such as
lnFmrGmr
where Fm(r) is the number of template vector pairs having d[sm+1,i, sm+1,j] < r, Gm(r) is the number of template vector pairs having d[sm,i, sm,j] < r, and where d was the Euclidean distance (Richman & Moorman, 2000) (see Table 1). In our study, m = 3, as suggested in two previous studies that computed entropy from resting-state fMRI data in controls (Z. Wang et al., 2014) and in pwMS (Zhou et al., 2016). The parameter r was defined as 0.2 sd(s) where sd(s) is the standard deviation of the time series (s), as recommended by previous work (Tomčala, 2020). An individual’s global entropy was calculated as the average of their regional entropies.

Mass Univariate Analysis

First, demographics and clinical variables were tested for differences between the three groups (HC, pwMS who had no disability, and pwMS who had disability) using Chi-squared test for qualitative variables and permutation test (David, 2008) for quantitative variables. Second, fractional occupancy of and transition probabilities between brain states and regional entropy were compared between groups via permutation test. P values from permutation test were computed as the number of permutations that had a difference in means bigger than the original difference. Group differences were considered significant when p < 0.05 after Benjamini–Hochberg (BH) correction (Benjamini & Hochberg, 1995) for multiple comparisons. Finally, the association between average regional (global) entropy and average (global) transition energy was computed using Spearman’s rank correlation, while the association between global entropy and log of lesion volume (transformed to ensure normality) was computed using Pearson’s correlation. All statistical analyses and graphs were performed using R (https://www.r-project.org) version 3.4.4 and Matlab (https://www.mathworks.com/) version R2020a.

Classification Analysis

Logistic regression with ridge regularization was used to classify (a) HC and pwMS without disability, (b) HC versus pwMS with disability, and (c) pwMS without versus with disability, using demographics/clinical information (age for HC vs. pwMS classification and age and disease duration only for the classification of pwMS subgroups) and vectorized transition energy matrices.

The models were trained with outer and inner loops of k-fold cross-validation (k = 5) to optimize the hyperparameter (λ) and test model performance. The folds for both inner and outer loops were stratified to ensure that each fold contained the same proportion of subjects in the two classes as the original dataset. The inner loop (repeated over five different partitions of the training dataset only) optimized the set of hyperparameters that maximized the validation set area under curve (AUC). A model was then fitted using the entire training dataset and the optimal hyperparameters, which was then assessed on the hold-out test set from the outer loop. The outer loop was repeated for 100 different random partitions of the data. The median AUC (over all 5 folds × 100 iterations = 500 test sets) was calculated to assess the performance of the models.

When the data contains class imbalance, models tend to favor the majority class. Due to the class imbalance in our data (66 vs. 33 pwMS with no disability vs. evidence of disability), the oversampling approach synthetic majority oversampling technique (SMOTE) (Chawla et al., 2002) was used to obtain a balanced training dataset during cross-validation. SMOTE compensates for imbalanced classes by creating synthetic examples using nearest neighbor information and has been shown to be among the most robust and accurate methods with which to control for imbalanced data (Santos et al., 2018).

The interpretation of the parameter coefficients from linear models is difficult when there are colinearities in the data (Haufe et al., 2014). To mitigate this problem, we applied the Haufe transform to the model coefficients estimated in each outer loop. Specifically, solving logistic regression with ridge regularization yields model coefficient estimates
wridge=Σx+λI1μx+μx
where Σx is the covariance matrix of the input variables x in the training dataset, λ is the hyperparameter of the ridge classifier optimized in the inner loop for the training dataset, I is the identity matrix with the same dimension as Σx, μx+, and μx are the average of the input variables for the positive and negative classes in the training dataset. The Haufe transform for these model coefficients is thus
ΣxΣx+λI1μx+μxvaryˆ
where yˆ is the output of the model for the training dataset. The Haufe-transformed model coefficients were averaged over all 5 folds × 100 iterations = 500 training sets, and the result was used to interpret the relationship between TE and the three groups of individuals.

Data/Code Availability Statement

The de-identified data that support the findings of this study are available upon reasonable request from the corresponding author. The codes that were used to generate the results and figures are publicly available. Please see (a) https://github.com/cerent/MS-NCT/ for the codes that were used for the classification analysis and figures and (b) the paper of Cornblath et al. (2020) for the codes that were used for the clustering and network control theory approach.

Patient Characteristics

Table 2 shows the subjects’ demographic and clinical information including sex, age, disease duration, EDSS, MS phenotype, and lesion volume. Unsurprisingly, pwMS with disability were older than both HCs and pwMS without disability (corrected p < 0.05) and had a trend toward shorter disease duration compared to pwMS without disability (corrected p = 0.161). Sex was not different across any of the groups. Phenotype and disability groups were not independent (corrected p < 0.05), where CIS individuals were included in the without disability group and PPMS individuals were included in the with disability group. The MS subgroups did not have a significant difference in lesion volume (corrected p value > 0.05).

Table 2.

Subject demographics and clinical information

VariableHC (n = 19)pwMS without disability (n = 66)pwMS with disability (n = 33)P value (pwMS without disability vs. pwMS with disability)
Age 45 [35.55, 49.50] 40.50 [35, 50] p value vs. HC = 0.88 56 [46, 58] p value vs. HC = 0.005* 0.0001* 
Female (%) 11 (55%) 46 (70%) p value vs. HC = 0.41 20 (61%) p value vs. HC = 0.96 0.49 
Disease duration   10 [7,14.75] 13 [9, 17] 0.161 
EDSS   0 [0, 1] 2 [2, 3] <2.2e−16* 
Phenotype   7 CIS, 59 RRMS 28 RRMS, 5 PPMS <0.001* 
Lesion volume (mm3  1,904 [726, 4111] 2,482 [453, 7788] 0.0948 
VariableHC (n = 19)pwMS without disability (n = 66)pwMS with disability (n = 33)P value (pwMS without disability vs. pwMS with disability)
Age 45 [35.55, 49.50] 40.50 [35, 50] p value vs. HC = 0.88 56 [46, 58] p value vs. HC = 0.005* 0.0001* 
Female (%) 11 (55%) 46 (70%) p value vs. HC = 0.41 20 (61%) p value vs. HC = 0.96 0.49 
Disease duration   10 [7,14.75] 13 [9, 17] 0.161 
EDSS   0 [0, 1] 2 [2, 3] <2.2e−16* 
Phenotype   7 CIS, 59 RRMS 28 RRMS, 5 PPMS <0.001* 
Lesion volume (mm3  1,904 [726, 4111] 2,482 [453, 7788] 0.0948 

Note. Values are presented as median [1st, 3rd quantile] for the continuous variables and as number (percent) for sex. The HC versus pwMS without and with disability as well as the two groups of pwMS were tested for differences; p values shown are corrected for multiple comparisons; * indicates significance.

Brain States and Comparison of Dynamics

Figure 2 shows the six brain states, consisting of three pairs of anticorrelated states with high/low-amplitude activity in the visual network (VIS+/−), ventral attention network (VAN+/−), and somatomotor network (SOM+/−). The transition probability between dynamic states was not significantly different between HC and any MS subgroup or between MS subgroups (corrected p value > 0.05) (see Supporting Information Figure S3). The pwMS subgroups had a trend toward shorter fractional occupancy in the SOM+ and longer fractional occupancy in the VAN+/− and SOM− compared to HC. The pwMS with evidence of disability had a trend toward greater fractional occupancy in the VIS+/−, VAN+, and SOM+ compared to those without disability. However, the differences were not significant for all comparisons (corrected p value > 0.05) (see Supporting Information Figure S4).

Figure 2.

Brain states are visualized using glass brain and radial plots. The radial plots show the mean of the positive and negative brain-state values over the Yeo functional network; DAN = dorsal attention; VAN = ventral attention; LIM = limbic; FP = fronto-parietal; DMN = default-mode network; SUB = subcortex; CER = cerebellum; VIS = visual; SOM = somatomotor. States were named based on the network having the maximum magnitude in the radial plot.

Figure 2.

Brain states are visualized using glass brain and radial plots. The radial plots show the mean of the positive and negative brain-state values over the Yeo functional network; DAN = dorsal attention; VAN = ventral attention; LIM = limbic; FP = fronto-parietal; DMN = default-mode network; SUB = subcortex; CER = cerebellum; VIS = visual; SOM = somatomotor. States were named based on the network having the maximum magnitude in the radial plot.

Close modal

Entropy of Brain Activity

Figure 3 shows the t-statistics comparing regional entropy between HC and both MS subgroups and between MS subgroups. There was no significant difference in regional entropy between HC and either MS subgroup (corrected p value > 0.05). However, entropy in the right insula, which is in the VAN network, was significantly decreased in the pwMS with disability compared to those without disability (corrected p value = 0.003). Figure 3 also shows the scatter plot of global entropy and (a) global transition energy (average across all entries in the transition energy matrix) and (b) log of lesion volume. Lower entropy was significantly correlated with greater transition energy (r = −0.17, p value = 0.04) and greater lesion volume (r = −0.20, p value = 0.03). The number of state transitions and overall transition energy are inversely related (higher transition energy demand = fewer transitions), so, unsurprisingly, lower entropy was also significantly correlated with fewer overall state transitions (r = 0.30, p = 0.0003) (see Supporting Information Figure S6).

Figure 3.

Group difference in entropy and correlation between entropy versus transition energy and lesion load. (A) Group differences in regional entropy between HC versus pwMS without disability (top panel), HC versus pwMS with evidence of disability (middle panel), and between two MS subgroups (bottom panel). The positive t statistics presented on the figure indicate greater average entropy in HC compared to MS subgroups as well as greater average entropy in the pwMS with no disability compared to those with evidence of disability. The scatterplot of global entropy versus (A) global transition energy and (B) log of lesion volume (in mm3). The entropy and energy were negatively correlated with average transition energy (r = −0.17, p value = 0.04) and lesion volume (r = −0.20, p value = 0.03). Green points represent HC (only present in the global energy plot), blue dots represent pwMS without disability, and red dots represent the pwMs with evidence of disability.

Figure 3.

Group difference in entropy and correlation between entropy versus transition energy and lesion load. (A) Group differences in regional entropy between HC versus pwMS without disability (top panel), HC versus pwMS with evidence of disability (middle panel), and between two MS subgroups (bottom panel). The positive t statistics presented on the figure indicate greater average entropy in HC compared to MS subgroups as well as greater average entropy in the pwMS with no disability compared to those with evidence of disability. The scatterplot of global entropy versus (A) global transition energy and (B) log of lesion volume (in mm3). The entropy and energy were negatively correlated with average transition energy (r = −0.17, p value = 0.04) and lesion volume (r = −0.20, p value = 0.03). Green points represent HC (only present in the global energy plot), blue dots represent pwMS without disability, and red dots represent the pwMs with evidence of disability.

Close modal

Group Classification Using State Transition Energies

Figure 4 shows the Haufe-transformed logistic regression model coefficients for classifying HC versus pwMS without disability, HC versus pwMS with disability, and pwMS without versus pwMS with disability; the models had an average AUC of 0.63, 0.62, and 0.62, respectively. Overall, greater transition energy was associated with pwMS with disability compared to both HC and pwMS without disability, while lower transition energy was associated with pwMS without disability compared to HC. Greater transition energies between VAN+ and SOM− were associated with pwMS with disability compared to HC, while smaller transition energies between SOM− and VAN− as well as between VIS− and VAN+ was associated with pwMS without disability compared to HC. Greater transition energy out of VAN+ was most associated with disability in pwMS compared to those without disability.

Figure 4.

Assessing the relationship between transition energy and disability. Haufe-transformed coefficients for the logistic regression with ridge regularization models classifying HC versus pwMS without disability, HC versus pwMS with disability, and pwMS without disability versus pwMS with disability. VIS = visual; VAN = ventral attention network; SOM = somatomotor.

Figure 4.

Assessing the relationship between transition energy and disability. Haufe-transformed coefficients for the logistic regression with ridge regularization models classifying HC versus pwMS without disability, HC versus pwMS with disability, and pwMS without disability versus pwMS with disability. VIS = visual; VAN = ventral attention network; SOM = somatomotor.

Close modal

Robustness of Results

The analysis was replicated using the cc400 atlas; five optimal states were found: VIS+/−, SOM+/−, and VAN− (see Supporting Information Figure S8), which were overlapping with 4/5 states found with the 86 region atlas. Like the original results, there was no significant difference in fractional occupancy between groups (corrected p value > 0.05) (see Supporting Information Figure S9). There was a significant negative correlation between global transition energy and global entropy (r = −0.16, p value = 0.03), while the correlation between global entropy and lesion load was not significant. Results for the classification model were similar to those obtained with the 86-region atlas include the following: (a) greater overall transition energy was associated with pwMS who had evidence of disability compared to HC and (b) greater transition energies between SOM− and FP+ were associated with pwMS who had disability compared to HC.

In this study, we identified commonly occurring brain states in a group of HC and pwMS individuals, investigated the energy required to transition between them using the NCT approach, and, furthermore, associated these energetic measurements with lesion load and brain activity entropy. Our main findings are (a) there were no group differences (HC vs. pwMS) in the dynamics of the six recurring brain states, made up of high- and low-amplitude activity in visual, ventral attention, and somatomotor networks; (b) lower entropy of global brain activity was associated with greater lesion load and larger global transition energy; and (c) pwMS without disability had decreased transition energies compared to controls, while pwMS with disability had increased transition energy compared to both controls and pwMS without disability. We hypothesized that this latter finding may indicate a shift in the patterns of brain activity in the pwMS without disability that decreases transition energies, but as shift evolves over the disease course, transition energies increase and disability occurs.

A recent study showed HCs had greater brain activity entropy than people with ADHD, and, furthermore, that lower entropy in people with ADHD was related to worse symptom scores (Sokunbi et al., 2013). Additionally, decreased entropy was found in the activity of regions containing a stroke lesion and their contralesional hemisphere homologues (Saenger et al., 2018). Our findings are consistent with these previous findings in that increased lesion volume, and likely a weaker structural backbone that may contribute to increased state transition energy, relates to a decrease in the entropy of brain activity. One study found decreased entropy in the right precentral and left parahippocampus and increased entropy in the left superior temporal, right postcentral, and right transverse temporal in RRMS compared to HC compared to controls, which is largely consistent with our regional results (Zhou et al., 2016). However, the same study showed no significant correlation between overall lesion volume and the averaged entropy computed for each voxel. The difference in our global findings could be due to several factors, including differences in the pwMS cohort characteristics, MRI acquisition and processing differences, and differences in the calculation of entropy (voxel-wise as opposed to region-wise entropy in our study). Finally, greater reductions in global transition energy from the placebo to psychedelic state (in controls) were shown to be associated with larger increases in brain-state entropy (Singleton et al., 2022), which is consistent with our findings of negative correlation between global entropy and transition energy across HC and pwMS.

The brain dynamics have previously been studied by clustering regional fMRI time series directly or sliding-window dynamic functional connectivity in HC (Cohen, 2018; Cornblath et al., 2020; Zalesky et al., 2014) and in patient populations such as stroke (Bonkhoff et al., 2021; L. Wang et al., 2010), schizophrenia (Braun et al., 2021; Damaraju et al., 2014; Rashid et al., 2016), and MS (d’Ambrosio et al., 2019; Rocca et al., 2019; Tozlu et al., 2021a). While many studies investigated the dynamics of sliding-window functional connectivity networks in MS (d’Ambrosio et al., 2019; Rocca et al., 2019; Tozlu et al., 2021a), this approach has some limitations as the window length and shift size need to be selected in a somewhat ad hoc manner. In this study, we use an alternative approach to identify the dynamic brain states in pwMS, for the first time, by clustering the regional time series activity directly. This approach has an advantage compared to other techniques as it is easy to implement, no ad hoc parameter estimation is required, and it uses each of the TRs as an individual data point which is especially helpful with short scans where dFC/FC estimates may be unreliable. NCT and the associated transition energy analysis has been applied previously in HC (Cornblath et al., 2020; Cui et al., 2020; Deng et al., 2022; Gu et al., 2022), in pharmacological states, and in neurological/neuropsychiatric disorders such as schizophrenia and stroke (Braun et al., 2021; Olafson et al., 2022; B. Tang et al., 2022). Using the NCT approach for the pwMS allowed us to jointly capture the effects of MS pathology on SC and co-occurring shifts in brain-state dynamics.

It has been shown in early MS there may be an upregulation of functional activity that initially acts to limit disability, followed by exhaustion of this functional compensation mechanism as disability increases (Tommasin et al., 2018). Previous studies showed these adaptive changes might be most prominent in motor and cognition-related networks, thus limiting clinical measures of disability (Lee et al., 2000) as well as cognitive impairment (Audoin et al., 2005) in MS. Moreover, the extent of the compensatory mechanism was found to be negatively correlated with the degree of tissue damage (Audoin et al., 2005; Rocca et al., 2002), such that more damage was related to decreased compensation. Our study’s finding of overall lower energetic demand in pwMS without disability compared to controls could be reflective of this previously identified compensation mechanism, and our finding of higher energetic demand in pwMS with disability compared to controls could be reflective of the previously identified decrease in compensation mechanism with increased damage (and/or increased damage to SCs).

The SOM and VAN states appeared to be the most important in the three classification models, with largest decreases in transition energy between SOM− and VAN− for pwMS without disability compared to controls and largest decreases in transition energies between SOM− and VAN+ for pwMS with disability compared to controls, and into/out of VAN+ for pwMS without disability versus those with disability. The biggest contributions to the HC versus pwMS with disability classification task were transitions between SOM− and VAN+. Our recent study showed that greater functional connectivity in the SOM and VAN states was associated with pwMS who had no disability compared to those with disability (Tozlu et al., 2021a). In our current study, the SOM− state also had the highest amplitude default mode network (DMN) activity. The VAN and DMN are known to be highly anticorrelated since the DMN is task-negative and the VAN is a task-positive network. Decreased anticorrelation between DMN and attention networks was related to worse behavioral performance in HC (Clare Kelly et al., 2008); moreover, increased functional connectivity between VAN and DMN was observed in pwMS who had cognitive impairment compared to those without cognitive impairment and HC (Huiskamp et al., 2021). Even though we did not have cognitive scores for our cohort, previous studies have identified a strong association between disability and cognitive impairment in MS (Lynch et al., 2005). Our results suggest that dysfunction of the interaction between the VAN and DMN states may result in greater energetic demand between these states in pwMS with evidence of disability compared to the other groups. Finally, our results also showed that greater transition energy between SOM+ and VIS− was associated with pwMS who had disability. The most frequently observed disturbances are the motor and visual impairments in pwMS, which are reflected in worse (higher) EDSS (Fuhr et al., 2001; Jasse et al., 2013). Our study shows, for the first time, greater energetic demand to transition between the motor and visual states compared to HC is associated with disability in MS.

Limitations

One of the limitations of the study was in the quality of the MRI collected, as the fMRI acquisition time was relatively short and had a longer TR (6 min, TR = 2.3 s) and the dMRI acquisition had only a single b-value of 800 and only 55 directions. A longer fMRI acquisition may result in finding other, less prevalent, dynamic states. However, Cornblath et al. (2020) applied the same clustering method to 15-min fMRI scans and found five optimal states that all had an appearance rate of approximately 2 per minute; thus, we believe our relatively shorter scan still represents the underlying state space well. Another limitation was the relatively low number of HCs compared to pwMS, which could result in an imbalance of the dynamic states. We used a permutation test when comparing HC and pwMS to minimize the effect of this group size imbalance. Larger lesion volume in some subjects may negatively impact coregistration to MNI space, which is where the FCs were calculated. In our cohort, lesions were generally on the small side, but we did visually inspect the MNI coregistration for quality and adjust if needed. Another limitation of the study was that most of the MS patients used in our study were of the RRMS phenotype, therefore the results may mostly represent disease processes within this subtype. A future study that includes a similar number of all MS phenotypes such as CIS, RRMS, PPMS, and SPMS is needed to capture the variability due to different level of disability in MS. Lastly, this study was cross-sectional; future studies will investigate how brain dynamics change over time and relate to physical and cognitive decline.

Conclusion

This is the first study using NCT to investigate how dynamic brain states and their associated transition energies change in MS and, furthermore, to associate brain activity entropy with lesion load and NCT-derived brain-state transition energies. As hypothesized, we found that global brain activity entropy decreased with increasing lesion load and increasing transition energies. Greater overall transition energies were associated with pwMS who had evidence of disability compared to both HC and pwMS without disability, while decreased overall transition energies were associated with pwMS without disability compared to the other two groups. Transition energies in ventral attention and somatomotor networks largely drove the group classifications. This study, through NCT, sheds light on a possible mechanism of how MS-related damage to the brain’s structural backbone can impact brain activity dynamics, entropy and energetics.

Supporting information for this article is available at https://doi.org/10.1162/netn_a_00292.

Ceren Tozlu: Data curation; Investigation; Visualization; Writing – original draft. Sophie Card: Validation. Keith Jamison: Data curation; Writing – review & editing. Susan Gauthier: Conceptualization; Data curation; Supervision; Writing – review & editing. Amy Kuceyeski: Conceptualization; Supervision; Writing – review & editing.

Amy Kuceyeski, NIH, Award ID: NS104634-01. Amy Kuceyeski, NIH, Award ID: NS102646-01A1. Susan Gauthier, Weill Cornell Medicine, Award ID: TR000456-06. Ceren Tozlu, National Multiple Sclerosis Society (https://dx.doi.org/10.13039/100000890), Award ID: FG-2008-36976.

Functional connectivity:

Temporal coactivation of the brain’s neural activity between different brain regions of interest.

Network control theory:

An approach to identify the brain’s energetic landscape using the structural connectome and patterns of brain activity over time.

Brain state:

The recurring patterns of the brain’s dynamic neural activity.

Entropy:

A measure indicating the complexity and predictability of a dynamic system.

Transition energy:

The minimum energy required to transition between pairs of brain states or to remain in the same state over some time.

Transition probability:

The probability that a brain switches to a different state or remains in the same state (i.e., persistence probability) at the next time point.

Fractional occupancy:

The number of TRs (repetition time) assigned to each cluster out of the total number of TRs.

EDSS:

Expanded Disability Status Scale that measures the disability level in people with multiple sclerosis. EDSS ranges between 0 and 10, where 0 means that the subject does not have a disability and 10 means death due to multiple sclerosis.

Structural connectivity:

Strength of the anatomical connections between different brain regions of interest.

Allen
,
E. A.
,
Damaraju
,
E.
,
Plis
,
S. M.
,
Erhardt
,
E. B.
,
Eichele
,
T.
, &
Calhoun
,
V. D.
(
2014
).
Tracking whole-brain connectivity dynamics in the resting state
.
Cerebral Cortex
,
24
(
3
),
663
676
. ,
[PubMed]
Andersson
,
J. L. R.
,
Graham
,
M. S.
,
Zsoldos
,
E.
, &
Sotiropoulos
,
S. N.
(
2016
).
Incorporating outlier detection and replacement into a non-parametric framework for movement and distortion correction of diffusion MR images
.
NeuroImage
,
141
,
556
572
. ,
[PubMed]
Andersson
,
J. L. R.
, &
Sotiropoulos
,
S. N.
(
2016
).
An integrated approach to correction for off-resonance effects and subject movement in diffusion MR imaging
.
NeuroImage
,
125
,
1063
1078
. ,
[PubMed]
Audoin
,
B.
,
Van Au Duong
,
M.
,
Ranjeva
,
J. P.
,
Ibarrola
,
D.
,
Malikova
,
I.
,
Confort-Gouny
,
S.
,
Soulier
,
E.
,
Viout
,
P.
,
Ali-Chérif
,
A.
,
Pelletier
,
J.
, &
Cozzone
,
P. J.
(
2005
).
Magnetic resonance study of the influence of tissue damage and cortical reorganization on PASAT performance at the earliest stage of multiple sclerosis
.
Human Brain Mapping
,
24
(
3
),
216
228
. ,
[PubMed]
Barkhof
,
F.
(
2002
).
The clinico-radiological paradox in multiple sclerosis revisited
.
Current Opinion in Neurology
,
15
(
3
),
239
245
. ,
[PubMed]
Benjamini
,
Y.
, &
Hochberg
,
Y.
(
1995
).
Controlling the false discovery rate: A practical and powerful approach to multiple testing
.
Journal of the Royal Statistical Society: Series B (Methodological)
,
57
(
1
),
289
300
.
Bonkhoff
,
A. K.
,
Schirmer
,
M. D.
,
Bretzner
,
M.
,
Etherton
,
M.
,
Donahue
,
K.
,
Tuozzo
,
C.
,
Nardin
,
M.
,
Giese
,
A.-K.
,
Wu
,
O.
,
Calhoun
,
V. D.
,
Grefkes
,
C.
, &
Rost
,
N. S.
(
2021
).
Abnormal dynamic functional connectivity is linked to recovery after acute ischemic stroke
.
Human Brain Mapping
,
42
(
7
),
2278
2291
. ,
[PubMed]
Braun
,
U.
,
Harneit
,
A.
,
Pergola
,
G.
,
Menara
,
T.
,
Schäfer
,
A.
,
Betzel
,
R. F.
,
Zang
,
Z.
,
Schweiger
,
J. I.
,
Zhang
,
X.
,
Schwarz
,
K.
,
Chen
,
J.
,
Blasi
,
G.
,
Bertolino
,
A.
,
Durstewitz
,
D.
,
Pasqualetti
,
F.
,
Schwarz
,
E.
,
Meyer-Lindenberg
,
A.
,
Bassett
,
D. S.
, &
Tost
,
H.
(
2021
).
Brain network dynamics during working memory are modulated by dopamine and diminished in schizophrenia
.
Nature Communications
,
12
(
1
),
3478
. ,
[PubMed]
Buyukturkoglu
,
K.
,
Zeng
,
D.
,
Bharadwaj
,
S.
,
Tozlu
,
C.
,
Mormina
,
E.
,
Igwe
,
K. C.
,
Lee
,
S.
,
Habeck
,
C.
,
Brickman
,
A. M.
,
Riley
,
C. S.
,
De Jager
,
P. L.
,
Sumowski
,
J. F.
, &
Leavitt
,
V. M.
(
2020
).
Classifying multiple sclerosis patients on the basis of SDMT performance using machine learning
.
Multiple Sclerosis Journal
,
27
(
1
),
107
116
. ,
[PubMed]
Calhoun
,
V. D.
,
Miller
,
R.
,
Pearlson
,
G.
, &
Adali
,
T.
(
2014
).
The chronnectome: Time-varying connectivity networks as the next frontier in fMRI data discovery
.
Neuron
,
84
(
2
),
262
274
. ,
[PubMed]
Carhart-Harris
,
R. L.
,
Leech
,
R.
,
Hellyer
,
P. J.
,
Shanahan
,
M.
,
Feilding
,
A.
,
Tagliazucchi
,
E.
,
Chialvo
,
D. R.
, &
Nutt
,
D.
(
2014
).
The entropic brain: A theory of conscious states informed by neuroimaging research with psychedelic drugs
.
Frontiers in Human Neuroscience
,
8
,
20
. ,
[PubMed]
Chawla
,
N. V.
,
Bowyer
,
K. W.
,
Hall
,
L. O.
, &
Kegelmeyer
,
W. P.
(
2002
).
SMOTE: Synthetic minority over-sampling technique
.
Journal of Artificial Intelligence Research
,
16
,
321
357
.
Clare Kelly
,
A. M.
,
Uddin
,
L. Q.
,
Biswal
,
B. B.
,
Castellanos
,
F. X.
, &
Milham
,
M. P.
(
2008
).
Competition between functional brain networks mediates behavioral variability
.
NeuroImage
,
39
(
1
),
527
537
. ,
[PubMed]
Cohen
,
J. R.
(
2018
).
The behavioral and cognitive relevance of time-varying, dynamic changes in functional connectivity
.
NeuroImage
,
180
(
Pt B
),
515
525
. ,
[PubMed]
Cornblath
,
E. J.
,
Ashourvan
,
A.
,
Kim
,
J. Z.
,
Betzel
,
R. F.
,
Ciric
,
R.
,
Adebimpe
,
A.
,
Baum
,
G. L.
,
He
,
X.
,
Ruparel
,
K.
,
Moore
,
T. M.
,
Gur
,
R. C.
,
Gur
,
R. E.
,
Shinohara
,
R. T.
,
Roalf
,
D. R.
,
Satterthwaite
,
T. D.
, &
Bassett
,
D. S.
(
2020
).
Temporal sequences of brain activity at rest are constrained by white matter structure and modulated by cognitive demands
.
Communications Biology
,
3
(
1
),
261
. ,
[PubMed]
Craddock
,
R. C.
,
James
,
G. A.
,
Holtzheimer
,
P. E.
,
Hu
,
X. P.
, &
Mayberg
,
H. S.
(
2012
).
A whole brain fMRI atlas generated via spatially constrained spectral clustering
.
Human Brain Mapping
,
33
(
8
),
1914
1928
. ,
[PubMed]
Cui
,
Z.
,
Stiso
,
J.
,
Baum
,
G. L.
,
Kim
,
J. Z.
,
Roalf
,
D. R.
,
Betzel
,
R. F.
,
Gu
,
S.
,
Lu
,
Z.
,
Xia
,
C. H.
,
He
,
X.
,
Ciric
,
R.
,
Oathes
,
D. J.
,
Moore
,
T. M.
,
Shinohara
,
R. T.
,
Ruparel
,
K.
,
Davatzikos
,
C.
,
Pasqualetti
,
F.
,
Gur
,
R. E.
,
Gur
,
R. C.
, …
Satterthwaite
,
T. D.
(
2020
).
Optimization of energy state transition trajectory supports the development of executive function during youth
.
ELife
,
9
,
e53060
. ,
[PubMed]
Damaraju
,
E.
,
Allen
,
E. A.
,
Belger
,
A.
,
Ford
,
J. M.
,
McEwen
,
S.
,
Mathalon
,
D. H.
,
Mueller
,
B. A.
,
Pearlson
,
G. D.
,
Potkin
,
S. G.
,
Preda
,
A.
,
Turner
,
J. A.
,
Vaidya
,
J. G.
,
van Erp
,
T. G.
, &
Calhoun
,
V. D.
(
2014
).
Dynamic functional connectivity analysis reveals transient states of dysconnectivity in schizophrenia
.
NeuroImage: Clinical
,
5
,
298
308
. ,
[PubMed]
d’Ambrosio
,
A.
,
Valsasina
,
P.
,
Gallo
,
A.
,
De Stefano
,
N.
,
Pareto
,
D.
,
Barkhof
,
F.
,
Ciccarelli
,
O.
,
Enzinger
,
C.
,
Tedeschi
,
G.
,
Stromillo
,
M. L.
,
Arévalo
,
M. J.
,
Hulst
,
H. E.
,
Muhlert
,
N.
,
Koini
,
M.
,
Filippi
,
M.
, &
Rocca
,
M. A.
(
2019
).
Reduced dynamics of functional connectivity and cognitive impairment in multiple sclerosis
.
Multiple Sclerosis Journal
,
26
(
4
),
476
488
. ,
[PubMed]
David
,
H. A.
(
2008
).
The beginnings of randomization tests
.
American Statistician
,
62
(
1
),
70
72
.
Deng
,
S.
,
Li
,
J.
,
Thomas Yeo
,
B. T.
, &
Gu
,
S.
(
2022
).
Control theory illustrates the energy efficiency in the dynamic reconfiguration of functional connectivity
.
Communications Biology
,
5
(
1
),
295
. ,
[PubMed]
Filippi
,
M.
,
Agosta
,
F.
,
Spinelli
,
E. G.
, &
Rocca
,
M. A.
(
2013
).
Imaging resting state brain function in multiple sclerosis
.
Journal of Neurology
,
260
(
7
),
1709
1713
. ,
[PubMed]
Fischl
,
B.
, &
Dale
,
A. M.
(
2000
).
Measuring the thickness of the human cerebral cortex from magnetic resonance images
.
Proceedings of the National Academy of Sciences of the United States of America
,
97
(
20
),
11050
11055
. ,
[PubMed]
Fuchs
,
T. A.
,
Schoonheim
,
M. M.
,
Broeders
,
T. A. A.
,
Hulst
,
H. E.
,
Weinstock-Guttman
,
B.
,
Jakimovski
,
D.
,
Silver
,
J.
,
Zivadinov
,
R.
,
Geurts
,
J. J. G.
,
Dwyer
,
M. G.
, &
Benedict
,
R. H. B.
(
2022
).
Functional network dynamics and decreased conscientiousness in multiple sclerosis
.
Journal of Neurology
,
269
(
5
),
2696
2706
. ,
[PubMed]
Fuhr
,
P.
,
Borggrefe-Chappuis
,
A.
,
Schindler
,
C.
, &
Kappos
,
L.
(
2001
).
Visual and motor evoked potentials in the course of multiple sclerosis
.
Brain
,
124
(
11
),
2162
2168
. ,
[PubMed]
Gu
,
S.
,
Fotiadis
,
P.
,
Parkes
,
L.
,
Xia
,
C. H.
,
Gur
,
R. C.
,
Gur
,
R. E.
,
Roalf
,
D. R.
,
Satterthwaite
,
T. D.
, &
Bassett
,
D. S.
(
2022
).
Network controllability mediates the relationship between rigid structure and flexible dynamics
.
Network Neuroscience
,
6
(
1
),
275
297
. ,
[PubMed]
Gu
,
S.
,
Pasqualetti
,
F.
,
Cieslak
,
M.
,
Telesford
,
Q. K.
,
Yu
,
A. B.
,
Kahn
,
A. E.
,
Medaglia
,
J. D.
,
Vettel
,
J. M.
,
Miller
,
M. B.
,
Grafton
,
S. T.
, &
Bassett
,
D. S.
(
2015
).
Controllability of structural brain networks
.
Nature Communications
,
6
(
1
),
8414
. ,
[PubMed]
Hallquist
,
M. N.
,
Hwang
,
K.
, &
Luna
,
B.
(
2013
).
The nuisance of nuisance regression: Spectral misspecification in a common approach to resting-state fMRI preprocessing reintroduces noise and obscures functional connectivity
.
NeuroImage
,
82
,
208
225
. ,
[PubMed]
Haufe
,
S.
,
Meinecke
,
F.
,
Görgen
,
K.
,
Dähne
,
S.
,
Haynes
,
J. D.
,
Blankertz
,
B.
, &
Bießmann
,
F.
(
2014
).
On the interpretation of weight vectors of linear models in multivariate neuroimaging
.
NeuroImage
,
87
,
96
110
. ,
[PubMed]
Huiskamp
,
M.
,
Eijlers
,
A. J. C.
,
Broeders
,
T. A. A.
,
Pasteuning
,
J.
,
Dekker
,
I.
,
Uitdehaag
,
B. M. J.
,
Barkhof
,
F.
,
Wink
,
A.-M.
,
Geurts
,
J. J. G.
,
Hulst
,
H. E.
, &
Schoonheim
,
M. M.
(
2021
).
Longitudinal network changes and conversion to cognitive impairment in multiple sclerosis
.
Neurology
,
97
(
8
),
e794
e802
. ,
[PubMed]
Jasse
,
L.
,
Vukusic
,
S.
,
Durand-Dubief
,
F.
,
Vartin
,
C.
,
Piras
,
C.
,
Bernard
,
M.
,
Pélisson
,
D.
,
Confavreux
,
C.
,
Vighetto
,
A.
, &
Tilikete
,
C.
(
2013
).
Persistent visual impairment in multiple sclerosis: Prevalence, mechanisms and resulting disability
.
Multiple Sclerosis
,
19
(
12
),
1618
1626
. ,
[PubMed]
Kuceyeski
,
A.
,
Shah
,
S.
,
Dyke
,
J. P.
,
Bickel
,
S.
,
Abdelnour
,
F.
,
Schiff
,
N. D.
,
Voss
,
H. U.
, &
Raj
,
A.
(
2016
).
The application of a mathematical model linking structural and functional connectomes in severe brain injury
.
NeuroImage: Clinical
,
11
,
635
647
. ,
[PubMed]
Lee
,
M.
,
Reddy
,
H.
,
Johansen-Berg
,
H.
,
Pendlebury
,
S.
,
Jenkinson
,
M.
,
Smith
,
S.
,
Palace
,
J.
, &
Matthews
,
P. M.
(
2000
).
The motor cortex shows adaptive functional changes to brain injury from multiple sclerosis
.
Annals of Neurology
,
47
(
5
),
606
613
. ,
[PubMed]
Leonardi
,
N.
,
Richiardi
,
J.
,
Gschwind
,
M.
,
Simioni
,
S.
,
Annoni
,
J.-M.
,
Schluep
,
M.
,
Vuilleumier
,
P.
, &
Van De Ville
,
D.
(
2013
).
Principal components of functional connectivity: A new approach to study dynamic brain connectivity during rest
.
NeuroImage
,
83
,
937
950
. ,
[PubMed]
Lynch
,
S. G.
,
Parmenter
,
B. A.
, &
Denney
,
D. R.
(
2005
).
The association between cognitive impairment and physical disability in multiple sclerosis
.
Multiple Sclerosis
,
11
(
4
),
469
476
. ,
[PubMed]
Nguyen
,
X. V.
(
2021
).
The Adjusted Mutual Information
.
MATLAB Central File Exchange
. https://www.mathworks.com/matlabcentral/fileexchange/33144-the-adjusted-mutual-information
Olafson
,
E.
,
Russello
,
G.
,
Jamison
,
K. W.
,
Liu
,
H.
,
Wang
,
D.
,
Bruss
,
J. E.
,
Boes
,
A. D.
, &
Kuceyeski
,
A.
(
2022
).
Frontoparietal network activation is associated with motor recovery in ischemic stroke patients
.
Communications Biology
,
5
(
1
),
993
. ,
[PubMed]
Pincus
,
S. M.
(
1991
).
Approximate entropy as a measure of system complexity
.
Proceedings of the National Academy of Sciences
,
88
(
6
),
2297
2301
. ,
[PubMed]
Rashid
,
B.
,
Arbabshirani
,
M. R.
,
Damaraju
,
E.
,
Cetin
,
M. S.
,
Miller
,
R.
,
Pearlson
,
G. D.
, &
Calhoun
,
V. D.
(
2016
).
Classification of schizophrenia and bipolar patients using static and dynamic resting-state fMRI brain connectivity
.
NeuroImage
,
134
,
645
657
. ,
[PubMed]
Richman
,
J. S.
, &
Moorman
,
J. R.
(
2000
).
Physiological time-series analysis using approximate and sample entropy
.
American Journal of Physiology: Heart and Circulatory Physiology
,
278
(
6
),
H2039
H2049
. ,
[PubMed]
Rocca
,
M. A.
,
Falini
,
A.
,
Colombo
,
B.
,
Scotti
,
G.
,
Comi
,
G.
, &
Filippi
,
M.
(
2002
).
Adaptive functional changes in the cerebral cortex of patients with nondisabling multiple sclerosis correlate with the extent of brain structural damage
.
Annals of Neurology
,
51
(
3
),
330
339
. ,
[PubMed]
Rocca
,
M. A.
,
Hidalgo de La Cruz
,
M.
,
Valsasina
,
P.
,
Mesaros
,
S.
,
Martinovic
,
V.
,
Ivanovic
,
J.
,
Drulovic
,
J.
, &
Filippi
,
M.
(
2019
).
Two-year dynamic functional network connectivity in clinically isolated syndrome
.
Multiple Sclerosis Journal
,
26
(
6
),
645
658
. ,
[PubMed]
Saenger
,
V. M.
,
Ponce-Alvarez
,
A.
,
Adhikari
,
M.
,
Hagmann
,
P.
,
Deco
,
G.
, &
Corbetta
,
M.
(
2018
).
Linking entropy at rest with the underlying structural connectivity in the healthy and lesioned brain
.
Cerebral Cortex
,
28
(
8
),
2948
2958
. ,
[PubMed]
Santos
,
M. S.
,
Soares
,
J. P.
,
Abreu
,
P. H.
,
Araujo
,
H.
, &
Santos
,
J.
(
2018
).
Cross-validation for imbalanced datasets: Avoiding overoptimistic and overfitting approaches [Research Frontier]
.
IEEE Computational Intelligence Magazine
,
13
(
4
),
59
76
.
Sbardella
,
E.
,
Petsas
,
N.
,
Tona
,
F.
, &
Pantano
,
P.
(
2015
).
Resting-state fMRI in MS: General concepts and brief overview of its application
.
BioMed Research International
,
2015
,
212693
. ,
[PubMed]
Singleton
,
S. P.
,
Luppi
,
A. I.
,
Carhart-Harris
,
R. L.
,
Cruzat
,
J.
,
Roseman
,
L.
,
Nutt
,
D. J.
,
Deco
,
G.
,
Kringelbach
,
M. L.
,
Stamatakis
,
E. A.
, &
Kuceyeski
,
A.
(
2022
).
Receptor-informed network control theory links LSD and psilocybin to a flattening of the brain’s control energy landscape
.
Nature Communications
,
13
(
1
),
5812
. ,
[PubMed]
Smith
,
R. E.
,
Tournier
,
J. D.
,
Calamante
,
F.
, &
Connelly
,
A.
(
2015
).
SIFT2: Enabling dense quantitative assessment of brain white matter connectivity using streamlines tractography
.
NeuroImage
,
119
,
338
351
. ,
[PubMed]
Sokunbi
,
M. O.
,
Fung
,
W.
,
Sawlani
,
V.
,
Choppin
,
S.
,
Linden
,
D. E. J.
, &
Thome
,
J.
(
2013
).
Resting state fMRI entropy probes complexity of brain activity in adults with ADHD
.
Psychiatry Research: Neuroimaging
,
214
(
3
),
341
348
. ,
[PubMed]
Stiso
,
J.
,
Khambhati
,
A. N.
,
Menara
,
T.
,
Kahn
,
A. E.
,
Stein
,
J. M.
,
Das
,
S. R.
,
Gorniak
,
R.
,
Tracy
,
J.
,
Litt
,
B.
,
Davis
,
K. A.
,
Pasqualetti
,
F.
,
Lucas
,
T. H.
, &
Bassett
,
D. S.
(
2019
).
White matter network architecture guides direct electrical stimulation through optimal state transitions
.
Cell Reports
,
28
(
10
),
2554
2566
. ,
[PubMed]
Tang
,
B.
,
Zhang
,
W.
,
Deng
,
S.
,
Liu
,
J.
,
Hu
,
N.
,
Gong
,
Q.
,
Gu
,
S.
, &
Lui
,
S.
(
2022
).
Age-associated network controllability changes in first episode drug-naïve schizophrenia
.
BMC Psychiatry
,
22
(
1
),
26
. ,
[PubMed]
Tang
,
E.
,
Giusti
,
C.
,
Baum
,
G. L.
,
Gu
,
S.
,
Pollock
,
E.
,
Kahn
,
A. E.
,
Roalf
,
D. R.
,
Moore
,
T. M.
,
Ruparel
,
K.
,
Gur
,
R. C.
,
Gur
,
R. E.
,
Satterthwaite
,
T. D.
, &
Bassett
,
D. S.
(
2017
).
Developmental increases in white matter network controllability support a growing diversity of brain dynamics
.
Nature Communications
,
8
(
1
),
1252
. ,
[PubMed]
Tomčala
,
J.
(
2020
).
New fast ApEn and SampEn entropy algorithms implementation and their application to supercomputer power consumption
.
Entropy
,
22
(
8
),
863
. ,
[PubMed]
Tommasin
,
S.
,
de Giglio
,
L.
,
Ruggieri
,
S.
,
Petsas
,
N.
,
Giannì
,
C.
,
Pozzilli
,
C.
, &
Pantano
,
P.
(
2018
).
Relation between functional connectivity and disability in multiple sclerosis: A non-linear model
.
Journal of Neurology
,
265
(
12
),
2881
2892
. ,
[PubMed]
Tournier
,
J.-D.
,
Smith
,
R.
,
Raffelt
,
D.
,
Tabbara
,
R.
,
Dhollander
,
T.
,
Pietsch
,
M.
,
Christiaens
,
D.
,
Jeurissen
,
B.
,
Yeh
,
C. H.
, &
Connelly
,
A.
(
2019
).
MRtrix3: A fast, flexible and open software framework for medical image processing and visualisation
.
NeuroImage
,
202
,
116137
. ,
[PubMed]
Tournier
,
J.-D.
,
Calamante
,
F.
, &
Connelly
,
A.
(
2010
).
Improved probabilistic streamlines tractography by 2nd order integration over fibre orientation distributions
. In
Proceedings of the International Society for Magnetic Resonance in Medicine
,
Abstract 1670
.
Tournier
,
J.-D.
,
Calamante
,
F.
, &
Connelly
,
A.
(
2012
).
MRtrix: Diffusion tractography in crossing fiber regions
.
International Journal of Imaging Systems and Technology
,
22
(
1
),
53
66
. .
Tozlu
,
C.
,
Jamison
,
K.
,
Gauthier
,
S. A.
, &
Kuceyeski
,
A.
(
2021a
).
Dynamic functional connectivity better predicts disability than structural and static functional connectivity in people with multiple sclerosis
.
Frontiers in Neuroscience
,
15
,
763966
. ,
[PubMed]
Tozlu
,
C.
,
Jamison
,
K.
,
Gu
,
Z.
,
Gauthier
,
S. A.
, &
Kuceyeski
,
A.
(
2021b
).
Estimated connectivity networks outperform observed connectivity networks when classifying people with multiple sclerosis into disability groups
.
NeuroImage: Clinical
,
32
,
102827
. ,
[PubMed]
van Geest
,
Q.
,
Douw
,
L.
,
van ‘t Klooster
,
S.
,
Leurs
,
C. E.
,
Genova
,
H. M.
,
Wylie
,
G. R.
,
Steenwijk
,
M. D.
,
Killestein
,
J.
,
Geurts
,
J. J. G.
, &
Hulst
,
H. E.
(
2018a
).
Information processing speed in multiple sclerosis: Relevance of default mode network dynamics
.
NeuroImage: Clinical
,
19
,
507
515
. ,
[PubMed]
van Geest
,
Q.
,
Hulst
,
H. E.
,
Meijer
,
K. A.
,
Hoyng
,
L.
,
Geurts
,
J. J. G.
, &
Douw
,
L.
(
2018b
).
The importance of hippocampal dynamic connectivity in explaining memory function in multiple sclerosis
.
Brain and Behavior
,
8
(
5
),
e00954
. ,
[PubMed]
Wang
,
L.
,
Yu
,
C.
,
Chen
,
H.
,
Qin
,
W.
,
He
,
Y.
,
Fan
,
F.
,
Zhang
,
Y.
,
Wang
,
M.
,
Li
,
K.
,
Zang
,
Y.
,
Woodward
,
T. S.
, &
Zhu
,
C.
(
2010
).
Dynamic functional reorganization of the motor execution network after stroke
.
Brain
,
133
(
Pt 4
),
1224
1238
. ,
[PubMed]
Wang
,
Z.
,
Li
,
Y.
,
Childress
,
A. R.
, &
Detre
,
J. A.
(
2014
).
Brain entropy mapping using fMRI
.
PLoS One
,
9
(
3
),
e89948
. ,
[PubMed]
Whitfield-Gabrieli
,
S.
, &
Nieto-Castanon
,
A.
(
2012
).
Conn: A functional connectivity toolbox for correlated and anticorrelated brain networks
.
Brain Connectivity
,
2
(
3
),
125
141
. ,
[PubMed]
Yeo
,
B. T. T.
,
Krienen
,
F. M.
,
Sepulcre
,
J.
,
Sabuncu
,
M. R.
,
Lashkari
,
D.
,
Hollinshead
,
M.
,
Roffman
,
J. L.
,
Smoller
,
J. W.
,
Zöllei
,
L.
,
Polimeni
,
J. R.
,
Fischl
,
B.
,
Liu
,
H.
, &
Buckner
,
R. L.
(
2011
).
The organization of the human cerebral cortex estimated by intrinsic functional connectivity
.
Journal of Neurophysiology
,
106
(
3
),
1125
1165
. ,
[PubMed]
Zalesky
,
A.
,
Fornito
,
A.
,
Cocchi
,
L.
,
Gollo
,
L. L.
, &
Breakspear
,
M.
(
2014
).
Time-resolved resting-state brain networks
.
Proceedings of the National Academy of Sciences of the United States of America
,
111
(
28
),
10341
10346
. ,
[PubMed]
Zhou
,
F.
,
Zhuang
,
Y.
,
Gong
,
H.
,
Zhan
,
J.
,
Grossman
,
M.
, &
Wang
,
Z.
(
2016
).
Resting state brain entropy alterations in relapsing remitting multiple sclerosis
.
PLoS One
,
11
(
1
),
e0146080
. ,
[PubMed]

Author notes

Competing Interests: The authors have declared that no competing interests exist.

Handling Editor: Shella Keilholz

This is an open-access article distributed under the terms of the Creative Commons Attribution 4.0 International License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited. For a full description of the license, please visit https://creativecommons.org/licenses/by/4.0/legalcode.

Supplementary data