Altered spreading of fast aperiodic brain waves relates to disease duration in Amyotrophic Lateral Sclerosis

(cid:1) Avalanche transition matrices are used to investigate the spatio-temporal spreading of neuronal avalanches in ALS; (cid:1) ALS patients present higher values of the nodal strength in both cortical and sub-cortical brain areas, as compared to controls; (cid:1) The nodal strength directly correlates with disease duration, corroborating the clinical relevance of brain dynamics.


Introduction
Amyotrophic Lateral Sclerosis (ALS) is a progressive neurodegenerative disease characterized by the gradual degeneration of both upper (UMN) and lower motor neurons (LMN) (Rusconi et al., 2017;Salameh et al., 2015).However, recent evidences suggest that ALS is a multisystem disorder encompassing the whole brain, likely caused by a combination of complex pathogenic mechanisms (Hardiman et al., 2017;Prado et al., 2018).
Functional Connectivity (FC), defined as the pair-wise statistical dependencies between regional signals, captures the large-scale interaction (averaged over time) between brain areas (Friston, 2011).The presence of FC anomalies in ALS patients has been confirmed by several neuroimaging studies utilizing both magnetic resonance imaging (MRI) and electro/magnetoencephalography (EEG/MEG).Specifically, MRI studies showed alterations in the topological organization of brain networks (Basaia et al., 2020;Zhang et al., 2019), with both increases (Castelnovo et al., 2020;Chen et al., 2021) and decreases (Barry et al., 2021;Trojsi et al., 2021) in FC in both motor and extra-motor areas.By using MEG signals on a large cohort of ALS patients, Sorrentino and colleagues showed that the disease progression shifts the functional brain network toward hyper-connectedness, which means that the activity of each brain region over time is strongly constrained by the rest of the brain (Sorrentino et al., 2018).
However, FC is commonly employed to analyze stationary processes, representative of an aggregate account of the activities over a time interval.This might be a suboptimal setup given the complex, multimodal nature of the human dynamics.In fact, largescale brain dynamics shows oscillatory activities nested with aperiodic, scale-free burst of activities (i.e., neuronal avalanches) (Haldeman and Beggs, 2005;Tagliazucchi et al., 2012) that might be the manifestation of a system operating near criticality (Cocchi et al., 2017).
Accordingly, the healthy brain operates near a regime that allows for optimally flexible dynamics, as captured by the size of the functional repertoire (i.e., the number of distinct and unique avalanche patterns), that is a measure of brain flexibility (Cipriano et al., 2023;Romano et al., 2023;Rucco et al., 2020;Sorrentino et al., 2021a).The size of the functional repertoire is a key predictor of symptoms severity and clinical staging (Polverino et al., 2022).Consistently, in a resting-state fMRI study, Chen and colleagues demonstrated a significant correlation between temporal dynamic indices and disease severity in ALS patient, demonstrating again the relevance of ''healthy" dynamics to the design of biomarkers of disease progression (Chen et al., 2021).
Based on these assumptions, we have recently moved on to characterize the spatio-temporal spreading of neuronal avalanches in health and disease.In fact, it is known that avalanches do not spread randomly but, rather, preferentially along the whitematter bundles (i.e.structural connectome) (Duma et al., 2023;Sorrentino et al., 2021b).As a result, brain regions are recruited differently by the avalanches (Tewarie et al., 2022).In particular, structurally more central regions may be recruited more often than peripheral ones, and the probability of any two brain regions being sequentially activated is proportional to the structural connectivity strengths (Sorrentino et al., 2021b).
In this study, we hypothesized that ALS would be reflected in a rearrangement of neuronal avalanches propagation, leading to an overload of perturbations hinging across regions that are structurally more connected.To test our hypothesis, we characterized the spreading of neuronal avalanches.To this end, neuronal avalanches, operationally defined as events starting when large fluctuations of activity are present in at least one brain region and ending when all the regions return to their baseline, were identified from the source-reconstructed MEG signals of thirtysix ALS patients and forty-two healthy controls.Moreover, an avalanche pattern was defined as the set of brain regions recruited in each individual avalanche.In turn, the totality of the unique avalanche patterns (i.e., discarding repetitions) occurring over time was defined as the functional repertoire, and the number of such patterns defines its size, which was utilized to quantify flexibility.Then, we computed the avalanche transition matrices (ATMs), defined as the probability of two brain regions being consecutively recruited in an avalanche.From the ATMs, we calculated the network parameter nodal strength, which provides key information about the brain regions that are mostly involved in the process of spreading avalanches.Finally, we used Pearson's correlation coefficient to test the hypothesis of a relationship between the spatio-temporal dynamics of neuronal avalanches and the clinical conditions of ALS patients.

Participants
ALS patients were recruited from the ALS Center of the First Division of Neurology of the University of Campania ''Luigi Vanvitelli" (Naples, Italy).Patients were right-handed and native Italian speakers diagnosed with ALS according to the revised El-Escorial criteria (Brooks et al., 2000), and all of them underwent Edinburgh Cognitive and Behavioural ALS Screen (ECAS) (Siciliano et al., 2017) to assess global cognitive functioning.Healthy controls were also included.
Inclusion criteria were: (1) no use of drugs that could interfere with MEG signals; (2) no other major systemic, psychiatric or neurological diseases; and (3) no focal or diffuse brain damage at routine MRI.The study protocol was approved by the Local Ethics Committee (University of Campania ''Luigi Vanvitelli") with protocol number 591/2018, and all participants provided written informed consent in accordance with the Declaration of Helsinki.

MRI acquisition
MRI images of ALS patients and healthy controls were acquired on a 3 T scanner equipped with an 8-channel parallel head coil (General Electric Healthcare, Milwaukee, WI), as previously described (Polverino et al., 2020).MR scans were acquired either after the MEG recording or about one month before to minimize the risk of noise in MEG signals.In particular, three dimensional T1-weighted images (Gradient-echo sequence Inversion Recovery prepared Fast Spoiled Gradient Recalled-echo, time repetition = 6 988 ms, TI = 1100 ms, TE = 3.9 ms, flip angle = 10, voxel size = 1 x 1 x 1.2 mm3) were acquired.Four patients and eleven controls did not complete the MRI because of the difficulty in lying down or because refusing the MRI scan.In these cases, we used a standard MRI template.

MEG data acquisition and pre-processing
MEG system is composed of 163-magnetometers placed in a magnetically shielded room (AtB Biomag UG, Ulm, Germany) (Liparoti et al., 2021).Data acquisition, pre-processing, and source reconstruction were detailed in Sorrentino and colleagues (Sorrentino et al., 2018).
Briefly, before each acquisition, four reference positions (nasion, right and left preauricular, and apex) were digitalized, and electrocardiographic and electrooculographic signals were co-recorded (Gross et al., 2013).The brain activity of each subject was recorded with eyes closed in resting-state condition, in two segments of 3.5 min each.Data were acquired with a sampling frequency of 1024 Hz, applying a 4th-order Butterworth IIR band-pass filter to remove components below 0.5 and above 48.0Hz (Rucco et al., 2019).The filter was implemented offline using Matlab 2019a and the Fieldtrip toolbox 2014 (Oostenveld et al., 2011).
An expert rater evaluated brain signals by visual inspection.Noisy channels and noisy time-segments were removed.Principal component analysis (PCA) and independent component analysis (ICA) were used to orthogonalize the sensors on the base of the references and to remove physiological artefacts, respectively.

Source reconstruction
Source reconstruction was performed using the Fieldtrip toolbox (Oostenveld et al., 2011).MEG data were co-registered with the MRI and we used a beamformer approach (Van Veen et al., 1997) to reconstruct time series related to the centroids of 116 regions of interest (ROIs), derived from the Automated Anatomical Labeling (AAL) atlas (Gong et al., 2009;Tzourio-Mazoyer et al., 2002).We excluded the cerebellum from the analysis because of its low reliability (Brookes et al., 2011) and considered only the first 90 ROIs from the AAL.

Neuronal avalanches and functional repertoire
In order to quantify spatio-temporal fluctuations of brain activity, we estimated neuronal avalanches, defined as events starting when at least one brain region deviates from its baseline activity and ending when all areas are back to normal activity (Shriki et al., 2013).Each of the 90 source-reconstructed signals were ztransformed and thresholded according to a cut-off of 2.5 standard deviations (Sorrentino et al., 2021a).We repeated the analyses setting the threshold to z > |2.25| and z > |2.75| to confirm that results do not depend on the choice of a given threshold.
To select a suitable bin length, we computed the branching ratio r (Haldeman and Beggs, 2005).For each subject, for each time bin, and for each avalanche, the geometrically averaged ratio of the number of events (activations) between the subsequent time bin and the current time bin was calculated as where r i is the branching parameter of the i-th avalanche in the dataset, N bin is the total amount of bins in the i-th avalanche, n events (j) is the total number of events active in the j-th bin, and N aval is the total number of avalanches in each participant's recording.
The bin length equal to three samples yielded a critical process with r = 1, suggesting that neuronal avalanches occur in a dynamical regime near a phase transition.However, we also repeated our analysis varying time bins from 1 to 5. To compare the dynamics of brain activity among the subjects, we took into consideration the same duration for each participant time series (122.79 s).Segments of equal duration were randomly selected from the whole recording.

Avalanche transition matrices
We built an avalanche transition matrix (ATM) in which its elements (i, j) represented the probability that region j was active at time t + d, given that region i was active at time t.Hence, for each avalanche an ATM was built, and then these were averaged element-wise and finally symmetrized to obtain a subjectspecific ATM.In the subject-specific ATM, each ij-th entry is computed from those avalanche-specific ATMs when the ij-th edge was recruited, and not from the others.Hence, the quantity of data that contributes to the estimate is edge-specific.
In order to evaluate the spatial features of the avalanche spreading, for each subject-specific ATM we estimated the nodal strength, a topological parameter obtained from the sum of the ij-edges (with ij) of each node, calculated for each node of the ATM.Hence, the nodal strength is a nodal topological parameter representing the sum of the probabilities that a brain region (node) interacts with all the others during the neuronal avalanche spreading, and thus indicates its ''strength" in the process.A comprehensive graphical representation of the data analysis overview is shown in Fig. 1.

Statistical analysis
We performed the t-test to compare age and educational level between ALS patients and healthy controls, while the chi-square test was used for gender comparison.
The topological comparison of the brain dynamics was performed through a permutation test (Nichols and Holmes, 2002).Specifically, for each brain region, the nodal strength values were randomly shuffled between patients and controls, and absolute difference of the means of the new random groups was evaluated.After repeating this procedure 10,000 times, we obtained a null distribution of random absolute differences.At this point, the actual absolute difference between the mean nodal strength values of patients and control was compared to the null distribution to obtain the statistical significance.The results of the permutation testing were corrected by False Discovery Rate (FDR) correction (Benjamini and Hochberg, 1995).The comparison procedure was repeated for each bin and threshold, and we only took into consideration the brain regions that showed a corrected p-value < 0.05 across all conditions.In order to test the homogeneity of our cohort of patients with respect to the different phenotypes (i.e., predominant lower motor neuron (LMN) and predominant upper motor neuron (UMN)), we compared the two pathological groups for demographic parameters (i.e., age, gender and education), clinical data (i.e., disease duration, ECAS, King's/MiToS scales, ALSFRS-R), and for the network parameter nodal strength.No statistically significant differences were found according to the clinical phenotype, with the exception of the age (p-value < 0.001).
Finally, we used Pearson's correlation coefficient to evaluate a possible relationship between the ATMs nodal strength and the clinical conditions of ALS patients (ECAS, ALSFRS-R, clinical staging and disease duration).The results were corrected across both topological and clinical parameters using the FDR.The corrected significance level was set at p-value < 0.05, and all statistical analyses were performed using custom scripts written in MatLab 2022b.

Clinical and demographic features of ALS patients
Thirty-six ALS patients (27 males and 9 females; mean age 63. 61 ± 12.85; mean education 10.67 ± 4.67) and forty-two age-and education-matched healthy controls (28 males and 14 females; mean age 63.57 ± 10.35; mean education 12.26 ± 4.26) were enrolled in the study.
We used the total Amyotrophic Lateral Sclerosis Functional Rating Scale-Revised (ALSFRS-R) (Cedarbaum et al., 1999) to quantify symptoms severity, and patients were classified according to both the King's (Balendra et al., 2014)  Furthermore, according to the Strong criteria for cognitive and behavioral assessment in ALS (Strong et al., 2017), 4 patients had cognitive impairment (ALSci), 9 patients had behavioral impairment (ALSbi) and 4 patients both conditions (ALSbci).
All clinical information about the cohort, such as ALSFRS-R score, King's and MiToS disease staging systems, ALS phenotype, site of onset, disease progression rate, and neurocognitive assessment are reported in Table 1.The light blue squares indicate that the region i was above threshold three times during the avalanche.Region j was active only in two cases out of the three taken into consideration (green arrows).Thus, the probability that region j was active after the activation of region i was 2/3.The avalanche transition matrix (ATM) of an individual is created by averaging across avalanche-specific transition matrices.(For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Avalanche spreading and nodal strength
In order to evaluate the involvement of the brain areas in the neuronal avalanche propagation along the brain network, for each subject-specific ATM we estimated the nodal strength, obtained from the sum of the ij-edges (with ij) of each node, calculated for each node of the ATM.

Correlations with clinical parameters
Pearson's correlation coefficient was used to determine whether there was a correlation between the nodal strength and clinical parameters.Specifically, we focused on disease duration, clinical staging systems (King's and MiToS), ALSFRS-R and ECAS.
We found a significant positive correlation between the nodal strength of the right middle frontal gyrus and disease duration (r = 0.51, pFDR = 0.021; Fig. 3a).We obtained similar results with the right precentral gyrus (r = 0.44, pFDR = 0.028; Fig. 3b) and the right pallidum (r = 0.48, pFDR = 0.021; Fig. 3c).No significant correlation was found between the nodal strength and the other clinical features of the cohort.

Discussion
In the present work, we used source-reconstructed MEG data from a cohort of thirty-six ALS patients and forty-two healthy controls to test the hypothesis that ALS alters the optimal large-scale brain dynamics of neuronal avalanches.We built the ATMs, a recently developed approach to estimate the spatio-temporal propagation of avalanches along the structural brain network.An ATM provides spatial information on both the propagation of neuronal bursts and their regional localization.Specifically, we used the nodal strength as a topological parameter resulting from the sum of the probabilities that a given brain region (i.e., a node) interacts with all the others during a neuronal avalanche.The nodal strength has been calculated for all brain regions of interest, providing information about the areas that are most relevant to the spread of neuronal avalanches.
We started from the substantial observations that ALS is characterized by alterations in the structural and functional brain connectivity, both at the global and at the nodal level (Basaia et al., 2020;Renga, 2022;Trojsi et al., 2018;Zhang et al., 2019), in line with hypothesized spread-out involvement of the brain.In a large ALS cohort study, we have previously shown that functional brain networks are more integrated and vulnerable in patients, and this hyper-connectivity progresses with the clinical evolution of the disease (Sorrentino et al., 2018).Furthermore, ALS patients show less flexible brain dynamics, and this flexibility reduction predicts disease severity (Chen et al., 2021;Polverino et al., 2022).
Based on the ATMs, we found significantly higher values of the nodal strength in the right middle frontal gyrus, in the right precentral gyrus and in the right pallidum in ALS patients compared to healthy controls.These values correlated directly with disease duration.
It is well known that these brain regions are implicated in pathophysiological processes of the disease, involving not only motor and extra-motor cortical areas but also sub-cortical regions such basal ganglia (Hardiman et al., 2017;Verde et al., 2017).On the one hand, the middle frontal gyrus and the precentral gyrus are part of the motor cortex and are identified as crucial areas in the ALS neuropathological process (Agosta et al., 2007;Goldstein and Abrahams, 2013), characterized by progressive degeneration of both gray and white matter as the disease progresses (Bede and Hardiman, 2018).On the other hand, the pallidum, a sub-cortical area, is known to be closely structurally linked to the cortex, regulating voluntary movements (Gillies et al., 2017).Furthermore, it has been long known that pathophysiological processes of ALS are accompanied by degeneration of the pallidum (Smith, 1960).Sub-cortical regions have a key role in recruiting cortical areas, supporting coordinated activity across the brain and confirming the wide-spread involvement of brain areas in ALS (Ahmed et al., 2016;Goldberg et al., 2004;Polverino et al., 2022).
The increased nodal strength in the right middle frontal gyrus, in the precentral gyrus and in the right pallidum suggests that these regions are more frequently recruited in ALS individuals than in controls and attributes to these hubs a crucial role in facilitating the integration of functionally specialized and anatomically disparate neuronal systems in ALS-related neurodegeneration.It is interesting to note that our results show a remodeling of the frontal-subcortical circuitry with a marked right predominance.This is in line with the previous fMRI, functional near infrared spectroscopy (fNIRS) and EEG studies that investigated the brain connectivity changes and the neuronal network remodeling in ALS pathology (Renga, 2022).Borgheai and colleagues showed a shift towards a more centralized frontal functional network organization in ALS patients, describing the existence of hubs with high degree in the right prefrontal cortex (Borgheai et al., 2020).An increased functional connectivity was reported in several regions of the right hemisphere such as the precentral and superior frontal gyri, the fronto-parietal network, the angular gyrus and the inferior cerebellum (Agosta et al., 2013;Shen et al., 2018;Zhou et al., 2013).It cannot be excluded that the right-sided hyper- connectivity could represent a compensatory mechanism to the left-predominant ALS-related neurodegeneration (Devine et al., 2015;Qiu et al., 2019).In fact, these right frontal hubs are specular to the most commonly described regions to be involved in ALS pathology (Li et al., 2017;Qiu et al., 2019;Trojsi et al., 2018).It could be hypothesized that the right hemisphere supplies the damaged specular left-sided areas in the earlier phases of the disease to maintain a proper brain functioning.In this perspective, the positive correlation between right-sided frontal involvement and disease duration could be interpreted as a compensation that reaches this goal.However, it cannot be excluded that this positive correlation may be the consequence of an intrinsic slower progression (and so a minor aggressiveness) of the disease that gives the time to the brain for the compensation.
At present, we cannot state which one is or whether there are other possible explanations.Future studies, including neuropathology or genetic information, will be useful for better understanding the nature of these remodeling.What is relevant is that, in both cases, the increased nodal strength in some cortical and subcortical areas of the right hemisphere seems to be an interesting biomarker for disease duration, opening the path to an increasingly frequent use of the network analysis in the clinical practice and, eventually, in evaluating the therapeutic response to the new developing treatments.
In fact, the concepts of network-wise degeneration (Ahmed et al., 2016), circuit-specific vulnerability (Warren et al., 2013), and spreading of disease along structural connectivity patterns (Verstraete et al., 2014) are now widely recognized and supported by neuroimaging studies, but also by neuropsychological assessments and post-mortem data (Bede, 2017;Bede et al., 2018).
In this work, we provide a deeper characterization of neuronal avalanches propagation, describing their spatio-temporal trajectories and identifying the brain regions most likely to be involved in the process.This makes it possible to localize the bursts of activation and recognize the brain areas that take part in the pathogenic mechanisms of ALS.Our results are consistent with previous evidence on constrained brain dynamics in multiple neurodegenerative diseases (Cipriano et al., 2023;Polverino et al., 2022;Romano et al., 2023;Rucco et al., 2020;Sorrentino et al., 2021a).
It is noteworthy that the nodal strength of the involved regions correlates directly with disease duration, suggesting a rearrangement of the spatio-temporal dynamics of the neuronal avalanches, in the sense of an increased involvement of the most structurally connected regions as ALS progresses (Stam, 2014).
Although the alteration in the spreading of neuronal avalanches in ALS patients shows an interesting relationship with disease duration, there are no significant correlations with other clinical parameters, mainly with the clinical scales (i.e., King's and MiToS scales, ALSFRS-R, ECAS).This probably depends on multiple aspects.Both the King's and the MiToS scales are ordinal, noncontinuous staging systems in which different stages do not directly reflect the involvement of specific brain regions.Similarly, the ALSFRS-R is composed of 12 items investigating specific functions, many of which are consequences of bulbar or spinal involvement (not directly assessable with a MEG system), and can be influenced by various comorbidities.Concerning the ECAS, the final score is the result of the sum of different cognitive and behavioral assessments, each of which derives from a complex interaction of multiple, and often topologically distant, brain regions.Therefore, it is quite unlikely that few frontal hubs, especially in a not so wide sample, are directly related to complex and high-order cortical functions, such as language, executive functions, memory and visuospatial ability.
In addition, the entire picture became more complex in the context of a compensatory hypothesis (as mentioned above).A possible positive or negative correlation could depend on the efficacy of the compensation as well as on the intrinsic difficulty of the contralateral regions to properly compensate for a specific brain hypo-functioning.
Overall, our findings provide a deep and fine characterization of the spatio-temporal distribution of neuronal avalanches, offering a focus on the brain regions most involved in the process and making it possible to specifically recognize the areas that mainly take part in the pathogenic mechanisms of ALS.In fact, the network parameter nodal strength results from the sum of the probabilities that a given brain region interacts with all the others during a neuronal avalanche, and this is true for all brain regions of interest, providing a ''real-time" scenario of what is happening in the brain in a certain period of time.In conclusion, our framework corroborates the clinical relevance of aperiodic, fast large-scale brain activity as a biomarker of microscopic changes induced by neurophysiological processes.

Conclusions
In the current paper, we applied the novel framework based on the ATMs to investigate the spatio-temporal dynamics of neuronal avalanches in ALS, showing that: 1) ALS patients present higher centrality in both cortical and sub-cortical brain areas, as compared to healthy controls; 2) the higher centrality directly correlates with disease duration.

Ethical standards
The study protocol was approved by the Local Ethics Committee (University of Campania ''Luigi Vanvitelli") with protocol number 591/2018, and all participants provided written informed consent in accordance with the Declaration of Helsinki.

Fig. 1 .
Fig. 1.Pipeline overview.(a) Signal processing: a. Neuronal activity is recorded by magnetoencephalography (MEG); b.MEG and MRI signals are co-registered; c.Beamforming is obtained; d.Source-reconstructed time series.(b) Avalanche signal analysis: The blue lines represent the z-score activity of brain regions (ROIs), the rectangles represent the time frames (msec) in which neuronal avalanches occur, the red dots denote regions with above-threshold activity.(c) Flexibility estimation: In the green boxes, the brain-plot and the set of unique avalanche patterns (i.e. the size of the functional repertoire) are illustrated.(d) Avalanche transition matrix:The light blue squares indicate that the region i was above threshold three times during the avalanche.Region j was active only in two cases out of the three taken into consideration (green arrows).Thus, the probability that region j was active after the activation of region i was 2/3.The avalanche transition matrix (ATM) of an individual is created by averaging across avalanche-specific transition matrices.(For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Fig. 2 .
Fig. 2. Nodal strength comparison.Violin plots represent the comparison between healthy controls (HC, light blue) and patients (ALS, red) for the network parameter nodal strength, obtained from the avalanche transition matrices.ALS patients show higher values of nodal strength (a) in the right middle frontal gyrus, (b) in the right precentral gyrus and (c) in the right pallidum, compared to HC. Significance p-value: *p < 0.05.(For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Table 1
Demographic, clinical features and neurocognitive assessment of the cohort.