- 1CoCo Lab, Department of Psychology, Université de Montréal, Montréal, QC, Canada
- 2Institute for Applied Mathematics Mauro Picone, National Research Council, Roma, Italy
- 3Department of Electrical Engineering, École de Technologie Supérieure, Montréal, QC, Canada
- 4Mathematical Research Center, Université de Montréal, Montréal, QC, Canada
- 5Centre UNIQUE, Union Neurosciences et Intelligence Artificielle - Québec, Montréal, QC, Canada
- 6CUBRIC, School of Psychology, College of Biomedical and Life Sciences, Cardiff University, Cardiff, United Kingdom
- 7Division of Psychological Medicine and Clinical Neurosciences, MRC Centre for Neuropsychiatric Genetics and Genomics, School of Medicine, College of Biomedical and Life Sciences, Cardiff University, Cardiff, United Kingdom
- 8MEG Center, Université de Montréal, Montréal, QC, Canada
Schizophrenia has a complex etiology and symptomatology that is difficult to untangle. After decades of research, important advancements toward a central biomarker are still lacking. One of the missing pieces is a better understanding of how non-linear neural dynamics are altered in this patient population. In this study, the resting-state neuromagnetic signals of schizophrenia patients and healthy controls were analyzed in the framework of criticality. When biological systems like the brain are in a state of criticality, they are thought to be functioning at maximum efficiency (e.g., optimal communication and storage of information) and with maximum adaptability to incoming information. Here, we assessed the self-similarity and multifractality of resting-state brain signals recorded with magnetoencephalography in patients with schizophrenia patients and in matched controls. Schizophrenia patients had similar, although attenuated, patterns of self-similarity and multifractality values. Statistical tests showed that patients had higher values of self-similarity than controls in fronto-temporal regions, indicative of more regularity and memory in the signal. In contrast, patients had less multifractality than controls in the parietal and occipital regions, indicative of less diverse singularities and reduced variability in the signal. In addition, supervised machine-learning, based on logistic regression, successfully discriminated the two groups using measures of self-similarity and multifractality as features. Our results provide new insights into the baseline cognitive functioning of schizophrenia patients by identifying key alterations of criticality properties in their resting-state brain data.
Introduction
The global prevalence of schizophrenia is reported to be close to 21 million individuals (Charlson et al., 2018). The symptoms and poor prognosis of those affected can deeply impact their daily functioning, and weigh on those close to them. Unfortunately, progress in therapeutic development is slow in the field of psychiatry due to the extreme complexity of the brain, the heterogeneity of patients’ symptoms and difficulties in translational research. More knowledge is needed to better understand what alterations occur in the neural activity of patients. Among the missing pieces, further characterization of the resting neural dynamics of schizophrenia, and their relationship to patients’ symptoms, is needed. Alterations in the rhythmic (oscillatory) neural activity of schizophrenia patients have been widely reported in the neuroimaging literature (reviews: Uhlhaas and Singer, 2010; Maran et al., 2016; Alamian et al., 2017). In addition, an emerging body of research has reported changes in the arrhythmic properties of brain dynamics in schizophrenia (Breakspear, 2006; Fernández et al., 2013). A powerful concept that has so far remained under-exploited and poorly understood in neuropsychiatry is criticality.
What Is Criticality?
The dynamics of many complex systems, such as the human brain, appear to reside around the critical point of a phase transition (Beggs and Plenz, 2003; Stam and De Bruin, 2004; Fraiman and Chialvo, 2012; Palva and Palva, 2018). At this point of criticality, these systems are in a wavering state, at the cusp of a new phase, between the states of order and disorder (Beggs and Timme, 2012; Cocchi et al., 2017; Souza França et al., 2018). The brain requires such a balance of regularity (i.e., structure) on the one hand, to maintain coherent behavior, and flexibility (i.e., local variability) on the other hand, to adapt to ongoing changes in the environment (Chialvo, 2004; Beggs and Timme, 2012). Indeed, critical brain dynamics have been shown to be optimal for fast switching between metastable brain states, for maximizing information transfer and information storage within neural networks (Socolar and Kauffman, 2003; Haldeman and Beggs, 2005), and for optimizing phase synchrony (Yang et al., 2012). Importantly, it is within a critical state that neural communication can span the greatest distance and achieve maximal correlational length (Fraiman and Chialvo, 2012). Thus, the brain’s state of criticality is thought to affect the functional properties of oscillations, local synchronization and signal processing (Palva and Palva, 2018). Changes to this state, due to psychiatric illness for instance, can alter certain properties of this balance (e.g., in terms of strength and number of synaptic connections) (Beggs and Timme, 2012). Some of the tuning parameters of criticality appear to be embedded in the balance between neural excitation and inhibition (e.g., through NMDA receptors; Mazzoni et al., 2007; Shew et al., 2009; Hobbs et al., 2010; Poil et al., 2012), in neural network connection strengths, and synaptic plasticity (Rubinov et al., 2011; Beggs and Timme, 2012).
Measures of Criticality
Self-Similarity and Multifractality
Within the framework of criticality, local and large-scale fluctuations arise from excitatory post-synaptic potentials (EPSPs) and modulate brain states by facilitating or suppressing neuronal firing (Palva and Palva, 2018), with long-range spatial spread (He et al., 2010; Zilber, 2014). Systems in this state are characterized by power-law distributions, fractal geometry and fast metastable state transitions (Plenz and Chialvo, 2009; Cocchi et al., 2017; Chialvo, 2018; Palva and Palva, 2018). These features of a critical state are said to be scale-free or scale invariant. Power-law distributions of a given signal can be recognized as a linear slope in the log-log plot of the feature distribution, and they imply that the signal’s statistics and structural characteristics are preserved across spatiotemporal scales—in other words, that the signal has fractal properties (Beggs and Plenz, 2003; Chialvo, 2018). Fractal architectures describe objects that contain identical, or statistically equivalent, repetitive patterns at different magnifying scales (Mandelbrot, 1983, 1985; Van Orden et al., 2012; Fetterhoff et al., 2015).
Scale invariant dynamics of systems at criticality (i.e., power-law distributions and fractal architecture) have often been described using a 1/fβ power law fitted to Fourier-based spectral estimations. On the other hand, self-similarity is a well-accepted model for scale-free dynamics and is richer than the sole measure of β, as it captures fractional Gaussian noise and fractional Brownian motion. Self-similarity can be measured by the Hurst exponent, H. In the brain, H is thought to index how well neural activity is temporally structured (via its autocorrelation). The smoother the signal, the higher the value of H (Zilber, 2014). However, self-similarity alone does not fully account for scale-free dynamics or criticality, since it can only capture additive processes (La Rocca et al., 2018). Combining self-similarity with multifractality improves on this framework to better capture criticality in a system. Multifractality can account for the remaining non-additive, non-Gaussian processes. The multifractality parameter, M, quantifies the diversity of H’s (singularities) and the overarching geometry of spatiotemporal fluctuations (Leonarduzzi et al., 2016; La Rocca et al., 2018). Generally, fractals are evaluated using the topological dimension, D, which describes the complexity and structure of an object by measuring the change in detail based on the change in scale (Di Ieva, 2016). In multifractal analysis, the local regularity of a signal is quantified using the Hölder exponent, D(h) (Jaffard et al., 2016), allowing a more realistic characterization of phenomena that are too complex to be explained solely by Euclidian models. In sum, the brain’s degree of criticality is defined by its scale-free dynamics, which are best quantified by combining measures of self-similarity and multifractality.
Common Measures of Criticality
Numerous metrics have been developed to measure the scale-free properties that define criticality, such as Detrended Fluctuation Analysis (DFA) applied to oscillatory envelopes (Linkenkaer-Hansen et al., 2001; Hardstone et al., 2012) and neuronal avalanche detection (Beggs and Plenz, 2003). Non-linear dynamics, and specifically multifractal analysis, has been used to address questions of self-similarity and multifractality. Multifractal analysis can characterize both the amount of global self-similarity in a system and the amount of local fluctuations (i.e., number of singularities) (Zilber et al., 2012). This approach allows for more in-depth interpretations of the electrophysiological data compared to more conventional analytical approaches. A number of mathematical frameworks have tapped into this, such as the Multifractal Detrended Fluctuation Analysis (MFDFA; Kantelhardt et al., 2002; Ihlen, 2012) and the Wavelet Leaders-based Multifractal Analysis (WLMA; Wendt and Abry, 2007; Serrano and Figliola, 2009). For reviews of scale-free and multifractal analytical approaches (see Lopes and Betrouni, 2009; Zilber, 2014).
Application to Psychiatry
A scoping review of alterations of brain criticality changes in clinical populations was recently discussed in Zimmern (2020). An insightful illustration of reported changes to the state of criticality across multiple neurological and psychiatric disorders, from the perspective of self-similarity, are illustrated in Figure 6 of that article (Zimmern, 2020). The application of criticality models to psychiatry, and in particular to the study of schizophrenia (SZ), is well in line with leading theories for this pathology, which are centered around dysconnectivity and altered information processing and transfer (Weinberger et al., 1992; Friston and Frith, 1995; Fernández et al., 2013). So far, most of the empirical evidence for dysconnectivity theory in SZ has come from functional magnetic resonance imaging studies, which highlight several important alterations in anatomical and functional connectivity that exist in SZ patients, as well as from electroencephalography (EEG) and magnetoencephalography (MEG) connectivity studies (review: Alamian et al., 2017). However, we still lack a complete, in-depth understanding of the brain alterations inherent to this pathology in the temporal domain.
In terms of scale-free analyses in psychiatry, power spectral densities (PSD) of resting-state fMRI scans have shown SZ patients to have reduced complexity and disrupted scale invariant dynamics compared to controls in the precuneus, inferior frontal gyrus and temporal gyrus, and these changes correlated with their symptoms (Lee et al., 2021). Electrophysiological studies have found altered dimensional complexity and increased variability in SZ patients’ signal (Koukkou et al., 1993). A number of studies have applied different versions of multifractal analysis on electrophysiological (Slezin et al., 2007; Racz et al., 2020) or white-matter MRI data in SZ (Takahashi et al., 2009). One of these used the multifractal analysis on resting-state EEG data, and found increased long-range autocorrelation and multifractality in patients compared to controls (Racz et al., 2020).
In addition, two insightful reviews have examined how non-linear methods could improve our understanding of SZ (Breakspear, 2006; Fernández et al., 2013). They highlighted conflicting results among studies reporting on complexity changes in SZ, which they proposed were attributable to participants’ symptomatic state, the method of imaging or medication. Complexity as measured by Lempel–Ziv complexity (LZC) or correlation dimension (D2) was typically found to be increased in SZ in studies that recruited younger, first-episode patients who were drug-naïve and symptomatic, while studies reporting SZ-related reductions in complexity tended to recruit older, chronic, patients who were on medication and hence less symptomatic (Lee et al., 2008; Fernández et al., 2013). Although these measures have been widely applied to neuroscientific data, they each come with caveats that affect their precision or generalizability. Moreover, these reviews highlight the importance of controlling for factors such as age and medication when studying complex pathologies, such as SZ.
Goals of the Study
The brain is functionally optimal when in a state of criticality—in other words, when neural activity can spread equally well at long and short distances in time and space and information is processed and stored efficiently (Shew et al., 2009)—and multifractality analysis is an efficient indicator of criticality. Meanwhile, leading neural theories of SZ emphasize a pathological connectivity among neural signals across both space and time. It follows that multifractal analysis of brain signals in SZ may provide important insights into the nature of the pathological alterations that are associated with the disease and that underlie the severity of its symptoms.
Based on previous research that used self-similarity metrics (e.g., DFA) among the SZ population (Nikulin et al., 2012; Alamian et al., 2020), we expected altered self-similarity and multifractality values compared to healthy controls. Moreover, based on the literature on altered complexity in SZ (e.g., Lee et al., 2008, 2021; Fernández et al., 2013) we hypothesize that our patient group would show reduced multifractality compared to controls. We also predict significant correlations between measures of criticality and patients’ clinical symptom scores. The aim of the present study is to test these hypotheses by examining how criticality is altered in the neural activity of chronic SZ patients. More specifically, we set out to address this question by using a multimodal neuroimaging approach, combining resting-state MEG and structural MRI, and wavelet-based estimations of multifractality and self-similarity.
Materials and Methods
Participants
Participant data collection was conducted at the Cardiff University Brain Research Imaging Centre in Wales, United Kingdom., and the data analyses were conducted at the University of Montreal, QC, Canada. Ethical approval was obtained for the data collection according to the guidelines of the United Kingdom National Health Service ethics board, and the Cardiff University School of Psychology ethics board (EC.12.07.03.3164). Ethical approval was also obtained for these analyses from the research committee of the University of Montréal (CERAS-2018-19-069-D).
Behavioral and neuroimaging data from 25 chronic SZ patients (average age = 44.96 ± 8.55, 8 females) and 25 healthy controls (average age = 44.04 ± 9.20, 8 females) were included in this study. Healthy controls had no history of psychiatric or neurological disorders. The collected demographic information from all participants included: age, gender, depression score on the Beck Depression Inventory—II (BDI-II, Beck et al., 1996), and mania score on the Altman Self-Rating Mania Scale (ASRM, Altman et al., 1997). For the SZ patient group, additional information was collected: scores on the Scale of the Assessment of Positive Symptoms (SAPS) and the Scale of the Assessment of Negative Symptoms (SANS) (Kay et al., 1987), and information on antipsychotic doses standardized using olanzapine equivalents (Gardner et al., 2010). All of these data were anonymized, such that no identifiable information of participants was associated with their data nor with data from subsequent analyses. Patients were overall fairly asymptomatic on the testing day. No statistically significant group differences were observed across these demographic and clinical metrics, except for BDI-II scores, where SZ patients had on average mild depression (14.83 ± 9.11), compared to controls (4.50 ± 4.67). Additional details on participant information (i.e., recruitment procedure, exclusions, inclusions, and sample size calculation) can be found in Alamian et al. (2020).
Magnetoencephalography Experimental Design
The brain imaging data used for this study comes from 5-min of resting-state MEG recorded during an eyes-closed condition, with a 275-channel CTF machine. Reference electrodes were placed on each participant to account for cardiac, ocular, and other potential artifacts (Messaritaki et al., 2017). The MEG signal was initially recorded at a sampling frequency of 1,200 Hz. A 3 Tesla General Electric Signa HDx scanner with an eight-channel receive-only head RF coil was used to acquire MRI data. Each participant had a 5-min weighted 3D T1 anatomical scan (TR/TE/TI = 7.8/3.0/450 ms, flip angle = 20°FOV = 256*192*172 mm, 1 mm isotropic resolution) that was later used for source-reconstruction of the MEG data.
Data Preprocessing and Magnetoencephalography Source Reconstruction
Reference electrodes were placed on each participant above and below the center of the left eye, on the left and right pre-auricular, under the left and right temples and behind the left ear, to account for cardiac, ocular, and other potential artifacts (Messaritaki et al., 2017). The MEG signal was initially recorded at a sampling frequency of 1,200 Hz. NeuroPycon (Meunier et al., 2020), an open-source python toolbox, was used for the preprocessing and source-reconstruction analyses. First, the continuous raw data was down-sampled from 1,200 to 600 Hz, and band-pass filtered between 0.1 and 150 Hz using a finite impulse response filtering (FIR 1, order = 3) and a Hamming window. Next, independent component analysis (ICA) was used to remove artifacts (i.e., blinks, horizontal eye movements, heartbeat) from the MEG signal using MNE-python (Hyvarinen, 1999; Gramfort et al., 2013). ICs related to heart and ocular artifacts were identified based on the correlation with ECG and EoG channels. ICs were visually inspected to check the reliability of the automatic procedure implemented in MNE. On average we removed 1–2 ICs related to cardiac artifacts and 1–2 ICs related to ocular artifacts.
Since it has been reported that the values of the Hurst exponent, H, are unusually low in sensor-space, and tend to increase when moving from sensor to source space (based on simulations and real data: Blythe et al., 2014), source-reconstruction steps were taken to present cortical-level results in multifractal analysis. To generate individual anatomical source-spaces, the anatomical T1-MRI information of each subject was segmented with FreeSurfer (Fischl, 2012). However, given that this process would produce different source-space dimensions for each participant, individual source spaces were morphed and projected onto a standardized space from FreeSurfer (fsaverage) (Greve et al., 2013). The resulting source-space comprised 8,196 nodes on the cortical surface, where dipoles were 5 mm apart. The single layer model boundary element method implemented in MNE-python was used to compute the lead field matrix (Gramfort et al., 2013). Weighted Minimum Norm Estimate (Dale and Sereno, 1993; Hämäläinen and Ilmoniemi, 1994; Hincapié et al., 2016), implemented in the MNE-python package (Hyvarinen, 1999; Gramfort et al., 2013), was used to compute the inverse solution with a Tikhonov regularization parameter of lambda = 1.0 (Hincapié et al., 2016). Dipoles of the source-space were constrained to have an orientation perpendicular to the cortical surface. Thus, for this study, 8,196 time series were extracted at the cortical level.
Characterization of Criticality Through Self-Similarity and Multifractality
Measuring Self-Similarity and Multifractality
The singularity spectrum is a concise way to summarize information about scale-free dynamics. It allows the plotting of the Hölder exponents (h) about local variability in a time series, against the Fractional (Hausdorff) Dimensions, D(h), as can be seen in Figure 1.
Figure 1. Sketch of a singularity spectrum. These sketches illustrate the multifractal scaling function, which depicts a singularity spectrum. Local variability in the signal is represented by Hölder exponents, h, on the x-axis, while the amount of singularities is represented by the Fractal Dimension, D(h), on the y-axis. The apex of the curve reveals the most common h exponent, while the width of the curve reveals the multifractal spectrum. Using log-cumulants from the WLBMF (described in section “Defining Parameters of Log-Cumulants”) to describe the singularity spectrum, C1 informs on the apex, while C2 informs on the width of the function. (A) Shows a monofractal function, where C1 = H, the Hurst exponent, and C2 = 0. (B) Shows a multifractal function, where the concavity shows the distribution of h singularities.
Multifractal analysis builds on measures of self-similarity (e.g., slope of the PSD, DFA) to provide information about local fluctuations (singularities) in time. The multifractality spectrum and the scaling function ζ(q) (in terms of statistical moments q) are related, and can be described using the Legendre transformation:
When a signal is monofractal, this becomes a linear function, where ζ(q) = qH, as it would only have a single singularity (one unique property, Figure 1A). Here, the self-similarity parameter would be equal to H, the Hurst exponent (Wendt and Abry, 2007). When a signal is multifractal, the function ζ(q) has a curvature, as in Figure 1B, which shows the global spectrum of singularities. The Hölder exponent (h) with the largest Fractal dimension, D (apex of the curve), is said to be the most common singularity in the time-series. The width of the curve can be described with the multifractality parameter, M (Wendt and Abry, 2007).
In this study, to meaningfully estimate self-similarity and multifractality, we used the Wavelet p-Leader and Bootstrap based MultiFractal analysis (WLBMF). This approach builds on the Wavelet leaders-based multifractal analysis (WLMA) method that has been thoroughly described elsewhere (Wendt and Abry, 2007; Wendt et al., 2007; Serrano and Figliola, 2009; Ciuciu et al., 2012; Fetterhoff et al., 2015). Briefly, this WLMA method of estimating the singularity spectrum was shown to be efficient in untangling the scaling properties of neuronal signal, and more robust than other algorithms in addressing non-stationarity issues (Wendt, 2008). The curved shape of the scaling function ζ(q) can be written in its polynomial expansion around its maximum to allow the evaluation of Cp, log-cumulants:
The singularity spectrum can be thus derived from the series-expansion of Cp. The first two log-cumulants are the most informative, with C1, the first log-cumulant, reflecting self-similarity [and the location of the maximum of D(h), similar to H]. Its values approximate those of the H, and typically range between 0 and 1, although values above 1 have been observed (Samoradnitsky and Taqqu, 1994). C1 values above 0.5 indicate positive correlation (signal has memory), values below 0.5 indicate negative correlation, and a value of 0.5 indicates lack of correlation (random signal). Meanwhile, C2, the second log-cumulant, reflects multifractality (and the width of the singularity spectrum, like M) (Wendt and Abry, 2007; Wendt et al., 2009; Zilber, 2014; Diallo and Mendy, 2019). Given the concavity of the scaling function, C2 is always negative, and when C2 equals 0, it is said to indicate monofractality. Typically, the few studies that have applied this novel analytical approach have observed values between 0 and −0.02 (Zilber, 2014) or 0 and −0.07 (Ciuciu et al., 2012).
Hölder exponents cannot take on negative values. Thus, most multifractal analyses are constrained to scaling functions that have only positive local regularities, implying that there is a continuous temporal positive correlation in the signal (i.e., locally bound everywhere in the function). However, this is not true of all brain signals, which can present with discontinuities in the signal and can thus take on negative regularities. Thus, p-leaders have been proposed as a way to circumvent this limitation (Jaffard et al., 2016). The p-leader formalism has been proposed as an extension of and improvement on older mathematical frameworks of multifractal analysis (e.g., MFDFA) using wavelet-projections, by allowing the analysis of negative local regularities and by providing more accurate and detailed characterization of singularities in the signal. Different p-leader values change the regularity exponents, where p = infinity corresponds to the original wavelet-leaders analysis, p = 2 brings about similar exponents as observed using DFA. For a deeper understanding of the mathematical details, we refer the reader to Jaffard et al. (2016) and Leonarduzzi et al. (2016).
Defining Parameters of Log-Cumulants
One method to detect criticality in the brain is through the Wavelet p-Leader and Bootstrap based MultiFractal (WLBMF) analysis and, more specifically, through the evaluation of log-cumulants (Wendt and Abry, 2007; Wendt et al., 2007). This MATLAB-implemented technique uses the discrete wavelet domain for the analysis of self-similarity and multifractality in signals. In order to compute C1 and C2 in our study, we first plotted the PSD of each participant group (SZ patients, controls) in log-log space and identified the portion of the PSD function exhibiting a log-linear relationship. In our data, the log-linear portion of the PSD belonged to j1 = 7 and j2 = 10, which correspond to 3.5 and 0.4 Hz, respectively, as deduced by the following equation: Scale = , where Sf represents the sampling frequency, j1 and j2 represent the start and end points of the log-linear portion, respectively, and the scale represents the frequency bin to which it corresponds. This frequency range is similar to those of other researchers who have used the same multifractal analysis (Zilber, 2014). For a step-by-step illustration of the method, we direct the reader to Figure 7.1 in Zilber (2014) for an illustration of these steps. The PSD was calculated at the overall cortical level and also at the ROI level, using the Destrieux Atlas (Destrieux et al., 2010), to ensure that the linear part of the spectrum was comparable across brain regions. For the purposes of this study, we used second order statistics in the evaluation of the log-cumulants (i.e., p-leader of p = 2), which is comparable to long-range temporal correlations computed with DFA (Leonarduzzi et al., 2016). For the ROI-based investigations, the C1 and C2 log-cumulants were first computed for each node (n = 8,196 sources) in cortical source-space, and then averaged across ROIs (n = 148 ROIs based on the Destrieux atlas, Destrieux et al., 2010). Although we calculate group differences across all individual nodes, we chose to also run ROI based analysis to help with the interpretability of the brain regions involved.
Statistics and Machine-Learning Analyses
Conventional Statistics and Correlation Analyses
Group statistical analyses were conducted between SZ patients and matched-controls to test for group-level differences in C1, C2, and demographic and clinical data. This was done at the ROI and source levels. To do so, we used non-parametric statistical tests (two-tailed, unpaired, pseudo t-tests), corrected with maximum statistics using permutations (n = 1000, p < 0.001) (Nichols and Holmes, 2001; Pantazis et al., 2005).
Moreover, Pearson correlations with False Discovery Rate (FDR) correction (Genovese et al., 2002) were used to explore the relationship between cortex-level C1/C2 values and scores on the SANS, SAPS and medication-dosage, in patients. FDR correction (Benjamini-Hochberg) was applied to each p-value (computed for each of the 8,196 nodes) to account for the multiple comparisons in order to achieve a significance threshold of p < 0.05, corrected.
Machine Learning Analyses
MEG signal classification was conducted using a logistic regression model and a stratified 10-fold cross-validation scheme to evaluate the discriminative power of the log-cumulants C1 and C2 in classifying SZ patients and controls. First, at each of the 8,196 nodes, the feature vector (either C1 or C2 values), computed for each participant, was split into 10-folds, while maintaining a balance between the two classes (SZ and controls). Next, the classifier was trained on the data from nine of the 10-folds and tested on the remaining fold (test set). The classification performance was assessed using the decoding accuracy (DA) on the test set (i.e., percentage of correctly classified participants across the total number of participants in the test set). This operation was repeated iteratively until all the folds were used as test sets. The mean DA was used as the classification performance metric. In order to infer the statistical significance of the obtained DAs, permutation tests were applied to derive a statistical threshold as described in Combrisson and Jerbi (2015). This method consists of generating a null-distribution of DAs obtained by running multiple instances of the classification (n = 1,000), and randomly shuffling class labels each time. Maximum statistics were applied in order to control for multiple comparisons across all the nodes (Nichols and Holmes, 2001; Pantazis et al., 2005). Visbrain was used for all the ROI and cortical-level visualizations (Combrisson et al., 2019).
Results
Alterations in Self-Similarity and Multifractality
The group averages of C1 and C2 values for SZ patients and healthy controls can be seen in Figure 2. Across both participant groups, a clear gradient in C1 values was observed, where self-similarity values increase gradually from the frontal lobe to the occipital lobe. Interestingly, a similar gradient, but in the opposite direction, is observed in terms of C2 values in both groups, with C2 values gradually increasing from the occipital lobe to the frontal lobe. Moreover, the magnitude of this gradient appears less pronounced in patients than in controls.
Figure 2. Group averages of C1 and C2 values in SZ patients and controls. Averaged C1 and C2 values were computed for each of the 8,196 nodes, within each group. P-leader p = 2 was used. SZ, schizophrenia.
Conventional unpaired t-tests between the two subject groups did not yield any statistically significant differences in terms of C1 or C2 values (p < 0.05, two-tailed t-test). Figure 3A shows t-values for the direction and magnitude of group differences for C1 and C2 values, where positive (red) t-values indicate brain areas where SZ patients have smaller C1 or C2 values compared to controls, and negative (blue) t-values indicate brain areas where patients have larger C1 or C2 values compared to controls.
Figure 3. Group differences and machine-learning results. (A) Shows t-values from the unpaired t-tests (non-significant), showing (controls—patients). Positive (red) t-values illustrate brain regions where patients show smaller C1/C2 values than controls, while negative (blue) t-values illustrate regions where patients have larger C1/C2 values than controls. (B) Shows unthresholded DA values based on logistic regression, using C1/C2 as a single feature. (C) Shows the same DA values, thresholded at p < 0.05. (D) Shows the DA values corrected for multiple comparisons using maximum statistics (p < 0.05), thresholded at the chance level of 70%. P-leader p = 2 was used. DA, decoding accuracy.
By contrast, when using a machine-learning approach to test for out-of-sample generalization in the same data, we found that C1 and C2 in multiple brain regions led to statistically significant classification of the two subject groups, with up to 77% decoding accuracy (Figure 3D, max statistics correction, p < 0.05). More specifically, using source-space C1 values as a decoding feature led to statistically significant discrimination of SZ and controls in the subcallosal gyrus, middle fontal gyrus and anterior part of the cingulate gyrus, bilaterally. The left superior frontal gyrus, the left inferior frontal gyrus and sulci, and the right orbital, straight and frontomarginal gyri were also significant. The maximum decoding occurred in the left superior frontal gyrus (77%, compared to the chance level of 70%). Meanwhile, using source-space C2 values as a decoding feature led to statistically significant classification of SZ patients and controls in the superior parietal lobule, precuneus and posterior-ventral part of the cingulate gyrus in the right hemisphere. The left post-central gyrus, and superior temporal gyrus and occipital gyrus, bilaterally, were also significant. The maximum decoding accuracy took place in the right temporal gyrus (76%, compared to the chance level of 70%). Figures 3B–D show the unthresholded DA values for C1 and C2, as well as the uncorrected results at p < 0.05, and the corrected classification results at p < 0.05, with multiple comparisons correction using max statistics.
Figure 4 shows the classification results based on C1 and C2 values computed at the ROI-level (p < 0.05, corrected for multiple comparisons). The ROIs involved in the significant discrimination of patients and controls were the left straight gyrus, the triangular part of the inferior frontal gyrus and the medial transverse frontopolar gyrus and sulcus for C1, and the superior occipital gyrus, the right cuneus and the left angular gyrus for C2. To illustrate how the classifier was able to successfully separate SZ patients from healthy controls, individual C1 and C2 values were computed and averaged across all brain sites that had shown significant decoding at the source-level. These values are presented in a scatter plot in Figure 5. The distribution of the individual C1 and C2 values (averaged over all sources with significant decoding accuracy) shows that C1 values are higher in patients than in controls (i.e., a trend toward more self-similarity) and C2 values also shift upwards in patients (i.e., a trend toward less multifractality).
Figure 4. ROI-based classification of SZ and controls using C1 and C2. Machine-learning classification of SZ patients and healthy controls using logistic regression and the features of C1 or C2 at the ROI-level. The ROI analysis was based on the Destrieux Atlas, p < 0.05, corrected for multiple comparisons. DA, decoding accuracy; SZ, schizophrenia.
Figure 5. Scatter plot visualization of individual C1 and C2 values. This figure shows individual C1 and C2 values, averaged across all the nodes that showed statistically significant patient vs. controls decoding (n = 50). This scatter plot illustrates that patients exhibit overall higher self-similarity (higher C1) and less multifractality (higher, less negative, C2).
It is noteworthy that this scatter plot reveals the presence of positive C2 values in the dataset, primarily in patients. Although mathematically ill-defined, the observation of positive C2 is not unprecedented. Positive C2 values in some individuals can be attributed to numerical instabilities (and might be statistically undistinguishable from 0) or to the fact that the data in these participants cannot be modeled using the multifractal formalism. The safest interpretation for the positive C2 values observed in Figure 5 (primarily in patients), is that data in these individuals were neither multifractal (C2 < 0) nor monofractal (C2 = 0). Given that this specific type of multifractal analysis has never been conducted on clinical data before, we explored how the results would change when using a p-leader of p = 4 (as opposed to the p = 2 we have used up to now). This analysis found fewer participants to have positive C2 values compared to p = 2, and generally allowed for a better modeling of multifractality in the resting neuromagnetic signal of participants. Figures of C1/C2 group averages and classification patterns based on p = 4 can be found in Supplementary Material 1. In summary, we observed a similar albeit stronger decoding of patients and controls based on C2 values in p = 4 than p = 2. Interestingly, C1 values were smaller (Supplementary Figure 1), and the strong frontal lobe classification results based on C1 values at p = 2 diminished at p = 4 (Supplementary Figures 2C,D). Taken together, the results of C1 estimation (self-similarity) were more reliable in our data when using a p-leader of p = 2, while C2 estimation (multifractality) provided more robust results with p = 4. Most importantly, the trends in terms of increasing C1 and C2 values in patients compared to controls was present irrespective of the choice of p.
Correlations With Clinical Scores/Information
The investigation of potential correlations between C1/C2 and clinical information resulted in a number of interesting results. Specifically, the correlations between C1 values and patients’ SANS scores (maximum r = 0.78, p < 0.05) in the left inferior frontal gyrus and sulcus (Figure 6A), and between C2 values and patients’ SAPS scores (maximum r = 0.78, p < 0.05) in the circular sulcus of the insula (Figure 6B) were statistically significant. In addition, the relationship between C1 and medication dosage yielded a statistically significant positive correlation (maximum r = 0.79, p < 0.05, after correcting across all nodes). Figures 6C, 7C illustrate that patients with higher medication dosage exhibited higher C1 values. This was especially significant in the superior frontal gyri, the right middle temporal gyrus, left mid-anterior cingulate gyrus and left inferior temporal sulcus (see Figure 6C). The positive correlations in these analyses are shown in the scatter plots in Figures 7A–C. These plots depict the relationship between individually averaged C1 and C2 values (based on the significant nodes), and patients’ symptom severity and medication dosages. To further clarify the C1 × SANS correlational results, a Pearson correlation was conducted between SANS scores and medication dosage, revealing a low-to-moderate correlation coefficient. The r2 of the regression model suggested that this relationship explained 27–40% of the data, meaning that the correlation of C1 × SANS was only partially mediated by medication.
Figure 6. Correlational results between C1 and C2 values and patients’ clinical information. Pearson correlation results between patients’ (A) C1 values and negative symptom scores on the SANS (p < 0.05), (B) C2 values and positive symptom scores on the SAPS (p < 0.05), and (C) C1 values and medication dosages (olanzapine equivalent in mg), p < 0.05.
Figure 7. Scatter plots showing the positive correlations between C1/C2 values and clinical information. These scatter plots depict correlations (p < 0.05) between individually averaged C1 and C2 values and subjects’ clinical information. The averaging of C1 and C2 values was over significant nodes. (A) Shows the correlation between C1 values and patients’ SANS scores, (B) shows the correlation between C2 values and SAPS scores, and (C) shows the correlation between C1 values and patients medication dosage (olanzapine equivalent in mg).
Uncorrected Spearman correlations were also computed and reasonable overlap was observed between the Pearson-based correlations and the Spearman findings. Similar to the Pearson correlations, the Spearman analysis revealed a positive correlation between C1 values and patients’ SANS scores in the left mid-anterior part of the cingulate gyrus and sulcus (r = 0.65, p < 0.0005), in the left precentral gyrus (r = 0.67, p < 0.0005), in the left temporal pole (r = 0.66, p < 0.0005), and in the right middle frontal sulcus (r = 0.71, p < 0.0005). Future studies with larger cohorts would be critical to probe the robustness of these results and should take into account covariates such as age, sex and illness duration.
Discussion
The central goal of this study was to examine and characterize criticality features in the baseline neural dynamics of schizophrenia. To do so, we evaluated the first two log-cumulants of the Wavelet p-Leader and Bootstrap based MultiFractal (WLBMF) analysis on the resting-state neuromagnetic signals of chronic SZ patients and healthy controls. This allowed us to determine the values of C1 (reflective of self-similarity) and C2 (reflective of multifractality) on the linear, scale-free portion of participants’ arrhythmic MEG signal in source-space. In brief, our findings partially supported our initial hypotheses about self-similarity and multifractality changes in SZ, whilst also revealing unexpected alterations in criticality.
Specifically, the findings of this study show that there are clear opposite gradients in the values of C1 and C2, along the rostro-caudal axis. A progression from low to high values of C1 were observed from anterior to posterior poles (i.e., frontal to occipital lobes), while C2 values showed the reverse progression. For both of these metrics, the gradient was less clear in SZ patients than in healthy controls. The t-values of the unpaired t-tests showed that patients had higher C1 values in the fronto-temporal area, and lower C1 values in the parieto-occipital areas compared to controls. In contrast, patients appeared to have higher C2 values in the temporal, parietal, and occipital areas than controls. Conventional t-test statistics failed to reach significance after multiple comparisons correction. However, a machine-learning approach based on logistic regression yielded statistically significant decoding (up to 77%) of patients and controls in a number of brain regions. Indeed, SZ patients and controls were categorized using C1 values in the anterior part of the cingulate gyrus (ACC), the left inferior gyrus, and the mid and superior frontal gyri, among other brain regions. Meanwhile, using C2 as a feature, we were able to statistically significantly classify patients and controls in the right temporal gyrus, precuneus, and occipital gyrus, among other brain regions.
In terms of the first log-cumulant, patients had a range of C1 values of [0.07, 1.44] in significant regions. In controls, this range was of [0.18, 1.16]. Typically, C1 (and thus H) values would be expected to be between 0 and 1 (where 0 < C1 < 0.5 implies negatively autocorrelated signal, C1 = 0.5 implies uncorrelated signal, and 0.5 < C1 < 1 implies positively autocorrelated signal), although values above 1 have been observed within the theory of generalized processes and tempered distributions (Samoradnitsky and Taqqu, 1994). In terms of the second log-cumulant, patients had a range of C2 values of [−0.01, 0.015] in significant brain regions, while controls had a range of [−0.02, 0.011]. These values fall within the same ranges reported by previous researchers (e.g., Zilber, 2014). As a reminder, higher C1 values are indicative of more self-similarity and memory in the signal, while lower (more negative) C2 values are indicative of more complexity in the form of multifractality. From our results, we infer that SZ patients exhibited more self-similar neural dynamics than healthy controls, and thus more regularity in the frontal and temporal brain areas. In addition, patients had fewer singularities (less diverse h) in the parietal and occipital brain regions, compared to healthy controls whose neural signals were more multifractal.
Further investigation of this analysis revealed that a subportion of participants (predominantly patients) had some positive C2 values. Theoretically, only [C2 < 0] (multifractal signals) or [C2 = 0] (monofractal signals) are expected. Observing positive C2 values implies that the multifractal formalism could not properly model the neuromagnetic data recorded in these patients. So, what does this tell us about the success of the classifier in using C2 to distinguish between patients and controls? The simplest explanation is that individuals with more negative C2 (stronger multifractal properties) were identified as healthy, whereas individuals with C2 values closer to zero (monofractal), or even higher than zero (neither multifractal no monfractal), were classified as patients. As a side note, we found that using an alternate p-leader of p = 4 improved C2 values, and the classifier reaffirmed the diminished multifractality characteristics of patients’ resting neuromagnetic signal. Taken together, we observe clear rostro-caudal gradients of ascending self-similarity and multifractality across both participant groups, albeit more clearly in controls. The reduced multifractality and increased self-similarity might reflect a certain rigidity in the temporal dynamics of SZ patients’ neural activity.
Our findings are consistent with recent publications that have characterized complexity in SZ in the same regions in which we observed alteration in the log-cumulants C1 and C2 [i.e., precuneus, inferior frontal gyrus and temporal gyrus, (e.g., Lee et al., 2021)]. Interestingly, a recent resting-state MEG-based study of SZ patients by La Rocca et al. (2018) also found a gradient in C1 values along the longitudinal axis; however, in contrast to our own finding of an ascending anterior-posterior gradient, they instead found an opposite, descending anterior-posterior gradient (La Rocca et al., 2018). In addition, La Rocca et al. (2018) compared how criticality features changed during a perceptual task. They reported that in healthy individuals, global self-similarity decreased, while focal multifractality increased when switching from rest to task. Moreover, the changes in multifractality correlated with brain regions implicated in the task. This finding could suggest that the metric of C2 has a functional role in cognitive processes (La Rocca et al., 2018). Of note, there are some methodological differences between our studies, such as the choice of scale (j1 and j2) for the linear portion of the PSD. Differences could also be due to age differences. Indeed, the authors reported the mean age of their participants to be 22 years old, while our group’s mean age was of 44 years old. In the complexity literature, it has been often reported that the properties of scale-free dynamics change with age (e.g., Fernández et al., 2011; Churchill et al., 2016), and so it is possible that there is a reversal of the self-similarity gradient with age. More work is needed to elucidate this.
Positive correlations were observed between the metrics of self-similarity and multifractality and patients’ clinical information. In particular, we observed an increase in C1 values in patients with increasing severity of scores on the negative symptoms scale (SANS) in the inferior frontal gyrus, as well as with patient’s medication dosage, the latter of which was especially strong (r = 0.79). The left frontal gyrus plays an important role in cognitive functioning (Swick et al., 2008) and language (Klaus and Hartwigsen, 2019). At the structural level, cortical thinning has been observed in the inferior frontal gyrus in SZ patients compared to healthy controls, which correlated with cognitive dysfunction (Kuperberg et al., 2003; Oertel-Knöchel et al., 2013). Correlation between inferior frontal gyrus volume and negative symptoms in SZ patients have been previously observed, but not in their non-affected siblings (Harms et al., 2010). At the functional level, higher cluster coefficients have been observed in the left inferior frontal compared to bipolar patients or controls (Kim et al., 2020), as well as weaker connectivity within the language network (Jeong et al., 2009). In addition to the reported structural alterations in this language processing center, the reduction in the temporal flexibility and enhanced regularity in the signal might explain why patients’ have poorer speech understanding, such as difficulty detecting metaphors, sarcasm or jokes (Rossetti et al., 2018). A correlational trend was also observed between multifractality and patients’ scores on the positive symptom scale (SAPS) in the circular sulcus of the insula. In past studies, negative correlations have been observed between reduced gray matter volume of the insula and SZ patients’ positive symptoms (Wylie and Tregellas, 2010; Cascella et al., 2011). It is interesting to note that self-similarity and multifractality were oppositely (and perhaps complementarily) correlated with symptom severity scores.
Taking into account the correlational findings, it is not surprising that, in our dataset of chronic and medicated patients, antipsychotic medication dosage was related to symptom severity, which itself was related to scale-free neural properties. Psychiatrists typically increase pharmaceutical dosage, gradually and as needed, to help manage symptoms. Sometimes, certain drug combinations that help manage positive symptoms (hallucinations, delusions) can worsen negative symptoms (Schooler, 1994; Goff et al., 1996). Evidence from other studies (Koukkou et al., 1993; Saito et al., 1998; Raghavendra et al., 2009) suggests that drug-naïve and first-episode patients may display a different pattern of criticality, thus the generalizability of our results is limited to other medicated, chronic SZ patients.
Another parallel can be drawn between this study’s results and findings from DFA analyses. The log-cumulants (C1 and C2) derived from WLBMF analysis using a p-leader of p = 2, as was used in the present study, are similar to scaling exponents obtained using DFA (Leonarduzzi et al., 2016), in that they both reflect temporal autocorrelations. In one of our recent publications, we computed DFA exponents on oscillatory envelopes in this same dataset of SZ patients and healthy controls (Alamian et al., 2020). The scale used for the computation of the log-cumulants (j1,2: 0.4–3.5 Hz) overlaps with the delta oscillatory band (0.5–3.5 Hz). Comparing delta DFA exponents and C1 between the studies reveals a good agreement: DFA exponents were reduced in patients compared to controls in the occipital and parietal lobes and increased values in the prefrontal and temporal lobes, similar to the C1 topology. The overlap was remarkably good considering that DFA was computed on band-limited rhythmic brain signal, while the log-cumulants of the singularity spectrum were computed on the arrhythmic raw brain signal. This comparison shows that while DFA is an adequate measure of the self-similarity aspect of criticality, it does not, however, provide any information on the multifractality of a signal, as does the second log-cumulant, C2. In this respect, they capture different properties of the neural signal, and should be treated as such. Several studies have examined the alterations of DFA across a number of psychiatric and neurological disorders. They found that a drop in DFA exponents occurs in SZ as well as in Alzheimer’s and Parkinson’s disease, whereas other conditions, such as depression, insomnia and epilepsy are typically associated with increases in DFA exponents (Zimmern, 2020). These findings reveal that reduced temporal autocorrelations observed in SZ are not disease-specific, but capture alterations that might be common to multiple psychiatric or neurological conditions. This again highlights the need for more elaborate measures of brain criticality, such as through the WLBMF analysis carried out in the present manuscript.
Criticality in the brain likely informs on the spatiotemporal organization and functioning of neural networks at the micro- and macroscopic levels (Hesse and Gross, 2014; Cocchi et al., 2017). While the origins of criticality are still debated, many agree that scale-free neural fluctuations are the signature of a brain in a state of criticality. A right balance of scale invariant properties (self-similarity, multifractality) is thought to be needed to adapt and respond to ever changing environments (Linkenkaer-Hansen et al., 2001; Plenz and Chialvo, 2009; Beggs and Timme, 2012; Palva et al., 2013; Shew and Plenz, 2013). Consequently, we propose that a change in this equilibrium could disrupt optimal brain functioning. When self-similarity is strong in a signal, as in the brain signals of our SZ cohort, the signal’s temporal autocorrelation decays slowly, such that signal memory lasts a long time. While still the subject of debate, it has been proposed that this enhanced temporal persistence (or redundancy) may make the brain less efficient in information processing (Zilber et al., 2013). Lower levels of self-similarity in signals, as in those of our healthy controls, are thought to reflect enhanced neural excitability and more efficient processing (He, 2011, 2014; Zilber et al., 2013). Interpretations of multifractality are still unclear, but it appears that a richer repertoire of singularities (multifractality > monofractality) suggests more variability and flexibility in the neural signal (Beggs and Timme, 2012), and thus in behavior. In our dataset, patients exhibited reduced multifractality in certain areas, thus suggesting a decrease in complexity and flexibility in their resting neuromagnetic signal. The observed alterations in these criticality metrics in SZ could explain the long, sustained nature of patients’ positive symptoms (delusions, hallucinations) and their difficulty in breaking away from them.
Conclusion
The overarching scale invariance of brain activity is thought to be a useful indicator of its organization across both temporal and anatomical scales (Werner, 2007; Zilber, 2014). Indeed, many have suggested that biological systems optimally process, adapt to and communicate information over long neural distances when in a state of criticality. This critical state involves a balance between regularity (structure) and flexibility (variability, local fluctuations). Disruption of this equilibrium may reduce the efficiency with which the system responds to changes in the environment. In this study, we applied WLBMF analysis to resting MEG signals and observed clear deviations in both the self-similarity and multifractality of these signals in chronic SZ patients compared to healthy controls. These changes in the state of criticality of patients lend further support to the theory of dysconnectivity in SZ from the perspective of temporal dynamics, as it characterizes a different way in which information interruption occurs in patients. This study also demonstrated that alterations in neural criticality can be used to accurately differentiate between chronic SZ patients and controls. We expect that these findings will fuel the search for strong biomarkers in SZ, borrowing a new, largely uncharted path.
Data Availability Statement
The datasets presented in this article are not readily available because the dataset for this manuscript is not publicly available due to ethical restrictions. Requests to access the datasets should be directed to corresponding author.
Ethics Statement
The studies involving human participants were reviewed and approved by the United Kingdom National Health Service Ethics Board, Cardiff University School of Psychology, and CERAS of University of Montreal. The patients/participants provided their written informed consent to participate in this study.
Author Contributions
GA wrote the manuscript, designed and performed the analyses, and interpreted the results. TL designed the analyses and edited the manuscript. AP designed the analyses. J-ML provided theoretical background and edited the manuscript. LK and JW designed the data collection protocol and acquired the data. KS designed the data collection protocol and edited the manuscript. KJ conceptualized the study, interpreted the data, and edited the manuscript. All authors contributed to the article and approved the submitted version.
Funding
This research was supported in part by the FRQNT Strategic Clusters Program (grant no. 2020-RS4-265502—Centre UNIQUE—Union Neurosciences and Artificial Intelligence—Quebec). The scanning was supported by CUBRIC and the Schools of Psychology and Medicine at Cardiff University, together with the MRC/EPSRC funded United Kingdom MEG Partnership Grant (grant no. MR/K005464/1). GA was supported by the Fonds de Recherche du Québec en Santé. TL was supported by the postdoctoral fellowship from IVADO. KJ was supported by funding from the Canada Research Chairs Program and a Discovery Grant (grant no. RGPIN-2015-04854) from NSERC (Canada), a New Investigators Award from FRQNT (grant no. 2018-NC-206005), and the IVADO-Apogée fundamental research project grant.
Conflict of Interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Publisher’s Note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
Acknowledgments
We would like to thank Herwig Wendt for helping with the interpretations of our findings, and Jordan O’Byrne for helpful comments on the manuscript.
Supplementary Material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fncir.2022.630621/full#supplementary-material
References
Alamian, G., Hincapié, A.-S., Pascarella, A., Thiery, T., Combrisson, E., Saive, A.-L., et al. (2017). Measuring alterations in oscillatory brain networks in schizophrenia with resting-state MEG: state-of-the-art and methodological challenges. Clin. Neurophysiol. 128, 1719–1736. doi: 10.1016/j.clinph.2017.06.246
Alamian, G., Pascarella, A., Lajnef, T., Knight, L., Walters, J., Singh, K. D., et al. (2020). Patient, interrupted: MEG oscillation dynamics reveal temporal dysconnectivity in schizophrenia. Neuroimage Clin. 28:102485. doi: 10.1016/j.nicl.2020.102485
Altman, E. G., Hedeker, D., Peterson, J. L., and Davis, J. M. (1997). The altman self-rating mania scale. Biol. Psychiatry 42, 948–955. doi: 10.1016/S0006-3223(96)00548-3
Beck, A. T., Steer, R. A., and Brown, G. K. (1996). Manual for the Beck Depression Inventory-II. San Antonio, TX: Psychological Corporation.
Beggs, J. M., and Plenz, D. (2003). Neuronal avalanches in neocortical circuits. J. Neurosci. 23, 11167–11177. doi: 10.1523/JNEUROSCI.23-35-11167.2003
Beggs, J. M., and Timme, N. (2012). Being critical of criticality in the brain. Front. Physiol. 3:163. doi: 10.3389/fphys.2012.00163
Blythe, D. A. J., Haufe, S., Müller, K. R., and Nikulin, V. V. (2014). The effect of linear mixing in the EEG on Hurst exponent estimation. Neuroimage 99, 377–387. doi: 10.1016/j.neuroimage.2014.05.041
Breakspear, M. (2006). The nonlinear theory of schizophrenia. Aust. N. Z. J. Psychiatry 40, 20–35. doi: 10.1111/j.1440-1614.2006.01737.x
Cascella, N. G., Gerner, G. J., Fieldstone, S. C., Sawa, A., and Schretlen, D. J. (2011). The insula-claustrum region and delusions in schizophrenia. Schizophr. Res. 133, 77–81. doi: 10.1016/j.schres.2011.08.004
Charlson, F. J., Ferrari, A. J., Santomauro, D. F., Diminic, S., Stockings, E., Scott, J. G., et al. (2018). Global epidemiology and burden of schizophrenia: findings from the global burden of disease study 2016. Schizophr. Bull. 44, 1195–1203. doi: 10.1093/schbul/sby058
Chialvo, D. R. (2004). Critical brain networks. Phys. A Stat. Mech. Appl. 340, 756–765. doi: 10.1016/j.physa.2004.05.064
Chialvo, D. R. (2018). Life at the edge: complexity and criticality in biological function. Acta Phys. Pol. B 49, 1955–1979. doi: 10.5506/APhysPolB.49.1955
Churchill, N. W., Spring, R., Grady, C., Cimprich, B., Askren, M. K., Reuter-Lorenz, P. A., et al. (2016). The suppression of scale-free fMRI brain dynamics across three different sources of effort: aging, task novelty and task difficulty. Sci. Rep. 6:30895. doi: 10.1038/srep30895
Ciuciu, P., Varoquaux, G., Abry, P., Sadaghiani, S., and Kleinschmidt, A. (2012). Scale-free and multifractal time dynamics of fMRI signals during rest and task. Front. Physiol. 3:186. doi: 10.3389/fphys.2012.00186
Cocchi, L., Gollo, L. L., Zalesky, A., and Breakspear, M. (2017). Criticality in the brain: a synthesis of neurobiology, models and cognition. Prog. Neurobiol. 158, 132–152. doi: 10.1016/j.pneurobio.2017.07.002
Combrisson, E., and Jerbi, K. (2015). Exceeding chance level by chance: the caveat of theoretical chance levels in brain signal classification and statistical assessment of decoding accuracy. J. Neurosci. Methods 250, 126–136. doi: 10.1016/J.JNEUMETH.2015.01.010
Combrisson, E., Vallat, R., O’Reilly, C., Jas, M., Pascarella, A., Saive, A., et al. (2019). Visbrain: a multi-purpose GPU-accelerated open-source suite for multimodal brain data visualization. Front. Neuroinform. 13:14. doi: 10.3389/fninf.2019.00014
Dale, A. M., and Sereno, M. I. (1993). Improved localization of cortical activity by combining EEG and MEG with MRI cortical surfarce reconstruction: a linear approach. J. Cogn. Neurosci. 5, 162–176. doi: 10.1162/jocn.1993.5.2.162
Destrieux, C., Fischl, B., Dale, A., and Halgren, E. (2010). Automatic parcellation of human cortical gyri and sulci using standard anatomical nomenclature. Neuroimage 53, 1–15. doi: 10.1016/j.neuroimage.2010.06.010
Di Ieva, A. (2016). The Fractal Geometry of the Brain (Springer Series in Computational Neuroscience). New York, NY: Springer.
Diallo, O. K., and Mendy, P. (2019). Wavelet leader and multifractal detrended fluctuation analysis of market efficiency: evidence from the WAEMU market index. World J. Appl. Econ. 5, 1–23. doi: 10.22440/wjae.5.1.1
Fernández, A., Gómez, C., Hornero, R., and López-Ibor, J. J. (2013). Complexity and schizophrenia. Prog. Neuro Psychopharmacol. Biol. Psychiatry 45, 267–276. doi: 10.1016/j.pnpbp.2012.03.015
Fernández, A., López-Ibor, M.-I., Turrero, A., Santos, J.-M., Morón, M.-D., Hornero, R., et al. (2011). Lempel–Ziv complexity in schizophrenia: a MEG study. Clin. Neurophysiol. 122, 2227–2235. doi: 10.1016/j.clinph.2011.04.011
Fetterhoff, D., Opris, I., Simpson, S. L., Deadwyler, S. A., Hampson, R. E., and Kraft, R. A. (2015). Multifractal analysis of information processing in hippocampal neural ensembles during working memory under δ9-tetrahydrocannabinol administration. J. Neurosci. Methods 244, 136–153. doi: 10.1016/j.jneumeth.2014.07.013
Fraiman, D., and Chialvo, D. R. (2012). What kind of noise is brain noise: anomalous scaling behavior of the resting brain activity fluctuations. Front. Physiol. 3:307. doi: 10.3389/fphys.2012.00307
Friston, K. J., and Frith, C. D. (1995). Schizophrenia: a disconnection syndrome? Clin. Neurosci. 3, 89–97.
Gardner, D., Murphy, A., O’Donnell, H., Centorrino, F., and Baldessarini, R. (2010). International consensus study of antipsychotic dosing. Am. J. Psychiatry 167, 686–693. doi: 10.1176/appi.ajp.2009.09060802
Genovese, C. R., Lazar, N. A., and Nichols, T. (2002). Thresholding of statistical maps in functional neuroimaging using the false discovery rate. Neuroimage 15, 870–878. doi: 10.1006/NIMG.2001.1037
Goff, D. C., Tsai, G., Manoach, D. S., Flood, J., Darby, D. G., and Coyle, J. T. (1996). D-cycloserine added to clozapine for patients with schizophrenia. Am. J. Psychiatry 153, 1628–1630. doi: 10.1176/ajp.153.12.1628
Gramfort, A., Luessi, M., Larson, E., Engemann, D. A., Strohmeier, D., Brodbeck, C., et al. (2013). MEG and EEG data analysis with MNE-Python. Front. Neurosci. 7:267. doi: 10.3389/fnins.2013.00267
Greve, D. N., Van der Haegen, L., Cai, Q., Stufflebeam, S., Sabuncu, M. R., Fischl, B., et al. (2013). A surface-based analysis of language lateralization and cortical asymmetry. J. Cogn. Neurosci. 25, 1477–1492. doi: 10.1162/jocn_a_00405
Haldeman, C., and Beggs, J. M. (2005). Critical branching captures activity in living neural networks and maximizes the number of metastable states. Phys. Rev. Lett. 94:058101. doi: 10.1103/PhysRevLett.94.058101
Hämäläinen, M. S., and Ilmoniemi, R. J. (1994). Interpreting magnetic fields of the brain: minimum norm estimates. Med. Biol. Eng. Comput. 32, 35–42. doi: 10.1007/BF02512476
Hardstone, R., Poil, S.-S., Schiavone, G., Jansen, R., Nikulin, V. V., Mansvelder, H. D., et al. (2012). Detrended fluctuation analysis: a scale-free view on neuronal oscillations. Front. Physiol. 3:450. doi: 10.3389/fphys.2012.00450
Harms, M. P., Wang, L., Campanella, C., Aldridge, K., Moffitt, A. J., Kuelper, J., et al. (2010). Structural abnormalities in gyri of the prefrontal cortex in individuals with schizophrenia and their unaffected siblings. Br. J. Psychiatry 196, 150–157. doi: 10.1192/bjp.bp.109.067314
He, B. J. (2011). Scale-free properties of the functional magnetic resonance imaging signal during rest and task. J. Neurosci. 31, 13786–13795. doi: 10.1523/JNEUROSCI.2111-11.2011
He, B. J. (2014). Scale-free brain activity: past, present, and future. Trends Cogn. Sci. 18, 480–487. doi: 10.1016/j.tics.2014.04.003
He, B. J., Zempel, J. M., Snyder, A. Z., and Raichle, M. E. (2010). The temporal structures and functional significance of scale-free brain activity. Neuron 66, 353–369. doi: 10.1016/j.neuron.2010.04.020
Hesse, J., and Gross, T. (2014). Self-organized criticality as a fundamental property of neural systems. Front. Syst. Neurosci. 8:166. doi: 10.3389/fnsys.2014.00166
Hincapié, A.-S., Kujala, J., Mattout, J., Daligault, S., Delpuech, C., Mery, D., et al. (2016). MEG connectivity and power detections with minimum norm estimates require different regularization parameters. Comput. Intell. Neurosci. 2016:3979547. doi: 10.1155/2016/3979547
Hobbs, J. P., Smith, J. L., and Beggs, J. M. (2010). Aberrant neuronal avalanches in cortical tissue removed from juvenile epilepsy patients. J. Clin. Neurophysiol. 27, 380–386. doi: 10.1097/WNP.0b013e3181fdf8d3
Hyvarinen, A. (1999). Fast and robust fixed-point algorithms for independent component analysis. IEEE Trans. Neural Netw. 10, 626–634. doi: 10.1109/72.761722
Ihlen, E. A. F. (2012). Introduction to multifractal detrended fluctuation analysis in Matlab. Front. Physiol. 3:141. doi: 10.3389/fphys.2012.00141
Jaffard, S., Melot, C., Leonarduzzi, R., Wendt, H., Abry, P., Roux, S. G., et al. (2016). P-exponent and p-leaders, part I: negative pointwise regularity. Phys. A Stat. Mech. Appl. 448, 300–318. doi: 10.1016/j.physa.2015.12.061
Jeong, B., Wible, C. G., Hashimoto, R. I., and Kubicki, M. (2009). Functional and anatomical connectivity abnormalities in left inferior frontal gyrus in schizophrenia. Hum. Brain Mapp. 30, 4138–4151. doi: 10.1002/hbm.20835
Kantelhardt, J. W., Zschiegner, S. A., Koscielny-Bunde, E., Havlin, S., Bunde, A., and Stanley, H. E. (2002). Multifractal detrended fluctuation analysis of nonstationary time series. Phys. A Stat. Mech. Appl. 316, 87–114. doi: 10.1016/S0378-4371(02)01383-3
Kay, S. R., Fiszbein, A., and Opler, L. A. (1987). The positive and negative syndrome scale (PANSS) for schizophrenia. Schizophr. Bull. 13, 261–276. doi: 10.1093/schbul/13.2.261
Kim, S., Kim, Y.-W., Shim, M., Jin, M. J., Im, C.-H., and Lee, S.-H. (2020). Altered cortical functional networks in patients with schizophrenia and bipolar disorder: a resting-state electroencephalographic study. Front. Psychiatry 11:661. doi: 10.3389/fpsyt.2020.00661
Klaus, J., and Hartwigsen, G. (2019). Dissociating semantic and phonological contributions of the left inferior frontal gyrus to language production. Hum. Brain Mapp. 40, 3279–3287. doi: 10.1002/hbm.24597
Koukkou, M., Lehmann, D., Wackermann, J., Dvorak, I., and Henggeler, B. (1993). Dimensional complexity of EEG brain mechanisms in untreated schizophrenia. Biol. Psychiatry 33, 397–407. doi: 10.1016/0006-3223(93)90167-C
Kuperberg, G. R., Broome, M. R., McGuire, P. K., David, A. S., Eddy, M., Ozawa, F., et al. (2003). Regionally localized thinning of the cerebral cortex in schizophrenia. Arch. Gen. Psychiatry 60, 878–888. doi: 10.1001/archpsyc.60.9.878
La Rocca, D., Zilber, N., Abry, P., van Wassenhove, V., and Ciuciu, P. (2018). Self-similarity and multifractality in human brain activity: a wavelet-based analysis of scale-free brain dynamics✩. bioRxiv [Preprint] 315853. doi: 10.1101/315853
Lee, S. H., Choo, J. S., Im, W. Y., and Chae, J. H. (2008). Nonlinear analysis of electroencephalogram in schizophrenia patients with persistent auditory hallucination. Psychiatry Investig. 5, 115–120. doi: 10.4306/pi.2008.5.2.115
Lee, Y. J., Huang, S. Y., Lin, C. P., Tsai, S. J., and Yang, A. C. (2021). Alteration of power law scaling of spontaneous brain activity in schizophrenia. Schizophr. Res. 238, 10–19. doi: 10.1016/j.schres.2021.08.026
Leonarduzzi, R., Wendt, H., Abry, P., Jaffard, S., Melot, C., Roux, S. G., et al. (2016). P-exponent and p-leaders, part II: multifractal analysis. Relations to detrended fluctuation analysis. Phys. A Stat. Mech. Appl. 448, 319–339. doi: 10.1016/j.physa.2015.12.035
Linkenkaer-Hansen, K., Nikouline, V. V., Palva, J. M., and Ilmoniemi, R. J. (2001). Long-range temporal correlations and scaling behavior in human brain oscillations. J. Neurosci. 21, 1370–1377. doi: 10.1523/JNEUROSCI.21-04-01370.2001
Lopes, R., and Betrouni, N. (2009). Fractal and multifractal analysis: a review. Med. Image Anal. 13, 634–649. doi: 10.1016/j.media.2009.05.003
Mandelbrot, B. B. (1985). Self-affine fractals and fractal dimension. Phys. Script. 32, 257–260. doi: 10.1088/0031-8949/32/4/001
Maran, M., Grent-‘t-Jong, T., and Uhlhaas, P. J. (2016). Electrophysiological insights into connectivity anomalies in schizophrenia: a systematic review. Neuropsychiatr. Electrophysiol. 2:6. doi: 10.1186/s40810-016-0020-5
Mazzoni, A., Broccard, F. D., Garcia-Perez, E., Bonifazi, P., Ruaro, M. E., and Torre, V. (2007). On the dynamics of the spontaneous activity in neuronal networks. PLoS One 2:e439. doi: 10.1371/journal.pone.0000439
Messaritaki, E., Koelewijn, L., Dima, D. C., Williams, G. M., Perry, G., and Singh, K. D. (2017). Assessment and elimination of the effects of head movement on MEG resting-state measures of oscillatory brain activity. Neuroimage 159, 302–324. doi: 10.1016/j.neuroimage.2017.07.038
Meunier, D., Pascarella, A., Altukhov, D., Jas, M., Combrisson, E., Lajnef, T., et al. (2020). NeuroPycon: an open-source Python toolbox for fast multi-modal and reproducible brain connectivity pipelines. Neuroimage 219:117020. doi: 10.1016/j.neuroimage.2020.117020
Nichols, T. E., and Holmes, A. P. (2001). Nonparametric permutation tests for functional neuroimaging: a primer with examples. Hum. Brain Mapp. 15, 1–25. doi: 10.1002/hbm.1058
Nikulin, V. V. Jönsson, E. G., and Brismar, T. (2012). Attenuation of long-range temporal correlations in the amplitude dynamics of alpha and beta neuronal oscillations in patients with schizophrenia. Neuroimage 61, 162–169. doi: 10.1016/j.neuroimage.2012.03.008
Oertel-Knöchel, V., Knöchel, C., Rotarska-Jagiela, A., Reinke, B., Prvulovic, D., Haenschel, C., et al. (2013). Association between psychotic symptoms and cortical thickness reduction across the schizophrenia spectrum. Cereb. Cortex 23, 61–70. doi: 10.1093/cercor/bhr380
Palva, J. M., Zhigalov, A., Hirvonen, J., Korhonen, O., Linkenkaer-Hansen, K., and Palva, S. (2013). Neuronal long-range temporal correlations and avalanche dynamics are correlated with behavioral scaling laws. Proc. Natl. Acad. Sci. U.S.A. 110, 3585–3590. doi: 10.1073/pnas.1216855110
Palva, S., and Palva, J. M. (2018). Roles of brain criticality and multiscale oscillations in temporal predictions for sensorimotor processing. Trends Neurosci. 41, 729–743. doi: 10.1016/j.tins.2018.08.008
Pantazis, D., Nichols, T. E., Baillet, S., and Leahy, R. M. (2005). A comparison of random field theory and permutation methods for the statistical analysis of MEG data. Neuroimage 25, 383–394. doi: 10.1016/j.neuroimage.2004.09.040
Plenz, D., and Chialvo, D. R. (2009). Scaling properties of neuronal avalanches are consistent with critical dynamics. arXiv [Preprint]. arXiv:0912.5369
Poil, S.-S., Hardstone, R., Mansvelder, H. D., and Linkenkaer-Hansen, K. (2012). Critical-state dynamics of avalanches and oscillations jointly emerge from balanced excitation/inhibition in neuronal networks. J. Neurosci. 32, 9817–9823. doi: 10.1523/JNEUROSCI.5990-11.2012
Racz, F. S., Stylianou, O., Mukli, P., and Eke, A. (2020). Multifractal and entropy-based analysis of delta band neural activity reveals altered functional connectivity dynamics in schizophrenia. Front. Syst. Neurosci. 14:49. doi: 10.3389/fnsys.2020.00049
Raghavendra, B. S., Dutt, D. N., Halahalli, H. N., and John, J. P. (2009). Complexity analysis of EEG in patients with schizophrenia using fractal dimension. Physiol. Meas. 30, 795–808. doi: 10.1088/0967-3334/30/8/005
Rossetti, I., Brambilla, P., and Papagno, C. (2018). Metaphor comprehension in schizophrenic patients. Front. Psychol. 9:670. doi: 10.3389/fpsyg.2018.00670
Rubinov, M., Sporns, O., Thivierge, J. P., and Breakspear, M. (2011). Neurobiologically realistic determinants of self-organized criticality in networks of spiking neurons. PLoS Comput. Biol. 7:e1002038. doi: 10.1371/journal.pcbi.1002038
Saito, N., Kuginuki, T., Yagyu, T., Kinoshita, T., Koenig, T., Pascual-Marqui, R. D., et al. (1998). Global, regional, and local measures of complexity of multichannel electroencephalography in acute, neuroleptic-naive, first-break schizophrenics. Biol. Psychiatry 43, 794–802. doi: 10.1016/S0006-3223(97)00547-7
Samoradnitsky, G., and Taqqu, M. S. (1994). Stable Non-Gaussian Random Processes: Stochastic Models with Infinite Variance, 1st Edn. New York, NY: Chapman and Hall/CRC.
Schooler, N. R. (1994). Deficit symptoms in schizophrenia: negative symptoms versus neuroleptic-induced deficits. Acta Psychiatr. Scand. 89, 21–26. doi: 10.1111/j.1600-0447.1994.tb05827.x
Serrano, E., and Figliola, A. (2009). Wavelet leaders: a new method to estimate the multifractal singularity spectra. Phys. A Stat. Mech. Appl. 388, 2793–2805. doi: 10.1016/j.physa.2009.03.043
Shew, W. L., and Plenz, D. (2013). The functional benefits of criticality in the cortex. Neuroscientist 19, 88–100. doi: 10.1177/1073858412445487
Shew, W. L., Yang, H., Petermann, T., Roy, R., and Plenz, D. (2009). Neuronal avalanches imply maximum dynamic range in cortical networks at criticality. J. Neurosci. 29, 15595–15600. doi: 10.1523/JNEUROSCI.3864-09.2009
Slezin, V. B., Korsakova, E. A., Dytjatkovsky, M. A., Schultz, E. A., Arystova, T. A., and Siivola, J. R. (2007). Multifractal analysis as an aid in the diagnostics of mental disorders. Nord. J. Psychiatry 61, 339–342. doi: 10.1080/08039480701643175
Socolar, J. E. S., and Kauffman, S. A. (2003). Scaling in ordered and critical random Boolean networks. Phys. Rev. Lett. 90:068702. doi: 10.1103/PhysRevLett.90.068702
Souza França, L. G., Vivas Miranda, J. G., Leite, M., Sharma, N. K., Walker, M. C., Lemieux, L., et al. (2018). Fractal and multifractal properties of electrographic recordings of human brain activity: toward its use as a signal feature for machine learning in clinical applications. Front. Physiol. 9:1767. doi: 10.3389/fphys.2018.01767
Stam, C. J., and De Bruin, E. A. (2004). Scale-free dynamics of global functional connectivity in the human brain. Hum. Brain Mapp. 22, 97–109. doi: 10.1002/hbm.20016
Swick, D., Ashley, V., and Turken, A. U. (2008). Left inferior frontal gyrus is critical for response inhibition. BMC Neurosci. 9:102. doi: 10.1186/1471-2202-9-102
Takahashi, T., Kosaka, H., Murata, T., Omori, M., Narita, K., Mitsuya, H., et al. (2009). Application of a multifractal analysis to study brain white matter abnormalities of schizophrenia on T2-weighted magnetic resonance imaging. Psychiatry Res. Neuroimaging 171, 177–188. doi: 10.1016/j.pscychresns.2008.03.009
Uhlhaas, P. J., and Singer, W. (2010). Abnormal neural oscillations and synchrony in schizophrenia. Nat. Rev. Neurosci. 11, 100–113. doi: 10.1038/nrn2774
Van Orden, G., Hollis, G., and Wallot, S. (2012). The blue-collar brain. Front. Physiol. 3:207. doi: 10.3389/fphys.2012.00207
Weinberger, D. R., Berman, K. F., Suddath, R., and Fuller Torrey, E. (1992). Evidence of dysfunction of a prefrontal-limbic network in schizophrenia: a magnetic resonance imaging and regional cerebral blood flow study of discordant monozygotic twins. Am. J. Psychiatry 149, 890–897. doi: 10.1176/ajp.149.7.890
Wendt, H. (2008). Contributions of Wavelet Leaders and Bootstrap to Multifractal Analysis: Images, Estimation Performance, Dependence Structure and Vanishing Moments. Confidence Intervals and Hypothesis Tests. Available online at: http://citeseerx.ist.psu.edu/viewdoc/download?doi=10.1.1.372.8718&rep=rep1&type=pdf (accessed September 1, 2020).
Wendt, H., and Abry, P. (2007). Multifractality tests using bootstrapped wavelet leaders. IEEE Trans. Signal Process. 55, 4811–4820. doi: 10.1109/TSP.2007.896269
Wendt, H., Abry, P., and Jaffard, S. (2007). Bootstrap for empirical multifractal analysis. IEEE Signal Process. Mag. 24, 38–48. doi: 10.1109/MSP.2007.4286563
Wendt, H., Abry, P., Jaffard, S., Ji, H., and Shen, Z. (2009). “Wavelet leader multifractal analysis for texture classification,” in Proceedings of the International Conference on Image Processing, ICIP (Washington, DC: IEEE Computer Society), 3829–3832. doi: 10.1109/ICIP.2009.5414273
Werner, G. (2007). Brain dynamics across levels of organization. J. Physiol. 101, 273–279. doi: 10.1016/j.jphysparis.2007.12.001
Wylie, K. P., and Tregellas, J. R. (2010). The role of the insula in schizophrenia. Schizophr. Res. 123, 93–104. doi: 10.1016/j.schres.2010.08.027
Yang, H., Shew, W. L., Roy, R., and Plenz, D. (2012). Maximal variability of phase synchrony in cortical networks with neuronal avalanches. J. Neurosci. 32, 1061–1072. doi: 10.1523/JNEUROSCI.2771-11.2012
Zilber, N. (2014). ERF and Scale-Free Analyses of Source-Reconstructed MEG Brain Signals During a Multisensory Learning Paradigm. Available online at: https://tel.archives-ouvertes.fr/tel-00984990 (accessed July 6, 2020).
Zilber, N., Ciuciu, P., Abry, P., and Van Wassenhove, V. (2012). “Modulation of scale-free properties of brain activity in MEG,” in Proceedings of the International Symposium on Biomedical Imaging, Barcelona, 1531–1534. doi: 10.1109/ISBI.2012.6235864
Zilber, N., Ciuciu, P., Abry, P., and Van Wassenhove, V. (2013). “Learning-induced modulation of scale-free properties of brain activity measured with MEG,” in Proceedings of the 10th IEEE International Symposium on Biomedical Imaging, San Francisco, CA, 998–1001. doi: 10.1109/ISBI.2013.6556645
Keywords: complexity, criticality, multifractal analysis, machine-learning, magnetoencephalography, resting-state, scale-free dynamics
Citation: Alamian G, Lajnef T, Pascarella A, Lina J-M, Knight L, Walters J, Singh KD and Jerbi K (2022) Altered Brain Criticality in Schizophrenia: New Insights From Magnetoencephalography. Front. Neural Circuits 16:630621. doi: 10.3389/fncir.2022.630621
Received: 18 November 2020; Accepted: 03 March 2022;
Published: 28 March 2022.
Edited by:
Axel Sandvig, Norwegian University of Science and Technology, NorwayReviewed by:
Ken-ichi Amemori, Kyoto University, JapanRiccardo Iandolo, Norwegian University of Science and Technology, Norway
Copyright © 2022 Alamian, Lajnef, Pascarella, Lina, Knight, Walters, Singh and Jerbi. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Golnoush Alamian, golnoush.alamian@ubc.ca