Traceability and dynamical resistance of precursor of extreme events

Extreme events occur in a variety of natural, technical, and societal systems and often have catastrophic consequences. Their low-probability, high-impact nature has recently triggered research into improving our understanding of generating mechanisms, providing early warnings as well as developing control strategies. For the latter to be effective, knowledge about dynamical resistance of a system prior to an extreme event is of utmost importance. Here we introduce a novel time-series-based and non-perturbative approach to efficiently monitor dynamical resistance and apply it to high-resolution observations of brain activities from 43 subjects with uncontrollable epileptic seizures. We gain surprising insights into pre-seizure dynamical resistance of brains that also provide important clues for success or failure of measures for seizure prevention. The novel resistance monitoring perspective advances our understanding of precursor dynamics in complex spatio-temporal systems with potential applications in refining control strategies.

We develop an approach that allows a time-resolved monitoring of dynamical resistance from multivariant observations of a system's interaction dynamics only. Estimating resistance usually requires knowledge about a system's different dynamical regimes and their response to some exogenous or endogenous perturbation 36,37 , an information which may not be accessible for the aforementioned systems. Here, we assume that induced and/or spontaneous changes in the system's underlying coupling structure-as assessed by a time-series-based analysis of interactions-suffice to properly identify the system's different dynamical regimes. An important estimate for dynamical resistance is the minimum "distance" between these operationally-defined regimes.
At the example of precursors of epileptic seizures in the human brain 32 , we demonstrate the potentials of our approach which eventually could be used to help maintain or move the epileptic brain and other vulnerable systems towards more desirable and sustainable dynamics, track thresholds of potential concern and help with evaluations on how such systems are being managed.

Results
A time-series-based approach to tracking dynamical resistance. Our approach (Fig. 1) is based on probing with high temporal resolution the dynamical coupling structure between interacting, even non-stationary and non-linear elements (or subsystems) of a system taking into account the subsystems' individual signals only. We assess the dynamical coupling structure with an analysis of the strength of interactions between all pairs of subsystems (Methods), which results in a temporal sequence of interaction matrices ρ.
To illustrate the possible existence of various dynamical regimes, we define similarity between two matrices ρ(t k ) and ρ(t l ) at times t k and t l as ξ ρ ρ where … denotes the Euclidean norm 38 . We emphasize that the similarity matrix ξ -estimated for all times t k and t l -contains pertinent information about the system's dynamics, and that recurrent patterns in the similarity matrix indicate dynamical regimes 39 . Now, in order to identify these regimes and their number we employ a time-resolved hierarchical clustering analysis of coupling structures in an abstract space spanned by all pairwise interactions (Methods). The minimum "distance" between the resulting clusters, that represent different dynamical regimes, is a central element of our approach and we take it as an estimate for dynamical resistance R: the larger this distance between regimes the higher is the capacity of a system to absorb disturbances and to reorganize while undergoing dynamical changes on the verge of an extreme event so as to still retain essentially the same functionality.
Dynamical resistance of the human epileptic brain prior to seizures. We illustrate this new approach with an important example from clinical medicine. Epilepsy is a very common neurological disorder affecting approximately 65 million people worldwide 40 . It can not be controlled sufficiently with anti-epileptic drugs in approximately 30% of people with epilepsy 41 . Central to the burden of uncontrollable epilepsy for the person with epilepsy is the seemingly unpredictability of seizures. A novel approach to control previously uncontrollable seizures in people with epilepsy would consist in identifying seizure precursors combined with delivering a counteracting influence (e.g., neurostimulation, local cooling, local drug perfusion, or behavioral intervention) [42][43][44] to prevent the generation of the extreme event 45,46 .
It is now well established that seizures in many people with epilepsy are preceded by a measurable changes in brain dynamics, which constitute a precursor of sufficient duration 32 . Lacking knowledge about the resistance of brain dynamics during such a precursor, it is unclear whether perturbing the brain with a counteracting influence would indeed prevent seizure generation. We address this problem and track dynamical resistance of brains of 43 subjects with epilepsy over days during which individual brains transit into and out of seizures (Methods).
We study the individual brains' dynamical coupling structures (Methods, Fig. 1) and derive a similarity matrix for each subject. Zooms into examples of such matrices from two subjects are shown in Fig. 2. This representation comprises changes of the coupling structure over 24 h in a single figure and allows to compare the similarity of coupling structures at different times. We find various recurrent patterns in each subject's similarity matrix that suggests the presence of different dynamical regimes. Changes in the brains' dynamical coupling structures occur on various time scales, and similarity to coupling structures at previous times appears to depend on the time of day. Prior to seizures, we observe time intervals of rather constant similarity that are rarely interrupted by changes in the coupling structure (Fig. 2, left panel) or intervals of comparably low similarity but with frequent intermissions (Fig. 2, right panel). We find similar changes in brains' coupling structures preceding the other seizures. Nevertheless, we also observe similar patterning during time intervals far off seizure-related activities.
Having identified different dynamical regimes for each subject (Methods), we find that dynamical resistance fluctuates over time to a lesser or greater extent (Fig. 3, panel a; Supplementary Fig. S1), however, with some peculiarities. First, in the hours prior to the vast majority of seizures (≈70%) we find dynamical resistance of brains to increase (Fig. 3, panel b). At first glance, this observation is quite surprising, since intuitively we would expect a diminished resistance of the epileptic brain in order to facilitate the generation of the extreme event seizure. Having said that, we here investigate seizures from people with drug-resistant epilepsies. We thus speculate that an increased dynamical resistance of the epileptic brain prior to seizures accounts for the reduced effectiveness of anticonvulsants.
Second, we observe some partly periodic temporal structure in the data (Fig. 3, panel c), which indicates dynamical resistance of brains to be influenced by ultradian, circadian, and probably even infradian rhythms (Fig. 3, panel d). When comparing resistance values from recordings in the daytime (6 am to 10 pm) with those at night-time (10 pm to 6 am; since no sleep-scoring was available for the patients investigated here, we cannot evaluate the influence of different sleep stages), we find the latter to be significantly increased (Fig. 4).
Although an influence of biological rhythms can to some extent be expected [47][48][49][50] , it is conceivable that the resulting waxing and waning of brain's dynamical resistance but not a distinct seizure-generating mechanism leads to the observed pre-seizure changes. We therefore investigate whether fluctuations in brain's dynamical resistance prior to seizures differ from the ones seen in between seizures. We borrow an approach commonly Tracking dynamical resistance of a system when transiting into and out of an extreme event. (a) The system's dynamical coupling structure is probed from a sliding-window analysis of the strength of interactions between all pairs (i, j) of subsystems (here: phase-synchronization analysis of intracranial electroencephalographic signals from various brain regions). (b) Temporal sequence of symmetric interaction matrices ρ. (c) Similar coupling structures with small Euclidean separation in an abstract space spanned by all pairwise interactions ρ ij are subsumed (hierarchical clustering analysis) to a finite number  of the system's accessible dynamical regimes C n (grey spheres). Black dots inside spheres indicate similar coupling structures at different, not necessarily consecutive times during a preset time window (here 1 h). Dynamical resistance R is defined as the minimum distance (dotted lines) between all accessible dynamical regimes. (d) Exemplary temporal evolution of dynamical resistance R when transiting into and out of an extreme event (here: epileptic seizure; onset marked with a bolt on x-axis). The inter-seizure period is colored black, the pre-seizure period (assumed duration 4 h) is colored red. The seizure and the post-seizure period (assumed duration 1.5 h) are colored orange. R int denote resistance values from the inter-seizure period and R pre those from the pre-seizure period.
Scientific RepoRts | (2019) 9:1744 | https://doi.org/10.1038/s41598-018-38372-y employed in seizure-prediction studies 22 that is based on a surrogate-based comparison of the distributions of resistance values from the presumed pre-seizure periods and from inter-seizure periods (Methods and Fig. 5 panel a). We find that in more than two thirds of cases (68.8% of seizures; = 6  accessible dynamical regimes), pre-seizure fluctuations in dynamical resistance clearly differ from those seen during inter-seizure periods ( Fig. 5 panel b). We perform similar investigations on the level of individual brains and find pre-seizure fluctuations of dynamical resistance from all but one subjects to clearly differ from inter-seizure fluctuations ( Supplementary  Fig. S4). This makes our approach comparable to previous, high-performing 51 precursor-indexing approaches.
We hypothesize that such abnormal pre-seizure fluctuations in dynamical resistance not only reflect the emergence of a seizure-permissive state but also indicate success/failure of potential seizure prevention strategies. We also hypothesize that tracking brain's dynamical resistance over time provides important clues for targeted and personalized seizure prevention techniques, particularly about how and when to counteract. Both of these hypotheses are now testable by applying our approach to monitor dynamical resistance.

Discussion
We have introduced a time-series-based approach to efficiently monitor dynamical resistance of complex systems prone to extreme events. Our approach differs from approaches that aim at quantifying stability. These typically require finding fixed points of a system for every possible value of control parameters. This ansatz, which is suitable for low-dimensional systems only 36 , is not applicable to complex systems that consist of a large number of coupled elements interacting with each other in a non-linear way on multiple scales. More recent developments 37 require precise knowledge of the governing equations of motion. In contrast, our novel approach exploits-in a data-driven way-the system's underlying dynamical coupling structure to define accessible dynamical regimes and derives an estimate for resistance from the minimum distance between regimes.
We have demonstrated the high suitability of our approach by a time-resolved tracking of dynamical resistance of individual brains, that transit into and out of extreme events epileptic seizures. The analysis goes substantially beyond state-of-the-art investigations of epileptic brain dynamics and also provides, for the first time, an objective measure for brain resistance, which is of paramount importance to optimize existing and to develop novel therapeutic interventions for controlling seizures.
Like any data-driven approach our method is limited by several assumptions: lacking an appropriate model for testing our approach's general suitability, we introduced the method using a simple and robust measure for phase synchronization 52 to quantify couplings, but the approach can to some extent also be implemented with other linear and non-linear quantifiers, for example using information-theoretic measures 53 . Phase-coherence analysis between all pairs of interacting subsystems assumes the dynamical coupling structure to mainly reflect phase-synchronization phenomena. As regards our brain analysis and despite there being strong evidence for phase synchronization to underlie various physiologic and pathophysiologic brain functions [54][55][56] , other forms of synchronization 57 are certainly relevant for brain dynamics and discounting them can confound the assessment of resistance of brain dynamics. Likewise, there are other confounders that arise from the data recording 35 .
For the second major step of the method, namely the identification of dynamical regimes, we employed an hierarchical, centroid-based clustering algorithm. Our choice was based on the well-accepted notion that no clustering algorithm suitable to all problems and universally applicable to arbitrary datasets exists, and in view of our data analysis, we deliberately chose an easy-to-handle and widely-used algorithm. Since many of the existing algorithms-including the one employed here-fail to identify a large number of dynamical regimes from high-dimensional data due to the curse of dimensionality 58 , we restricted our downstream analyses to a lower number of accessible regimes. Given these limitations it is important to interpret dynamical resistance only relative to the variables that were taken into account. A further complication are potentially predictive or other pertinent changes in dynamical resistance that are not captured with the temporal resolution considered here. Our hourly time resolution reflects a balance between a robust probing of the brain's dynamical coupling structure based on non-stationary non-linear signals and an efficient centroidal Voronoi tessellation of the data space, thereby taking into account the duration of a potential precursor and the low-probability nature of epileptic seizures. Investigating dynamical resistance of other vulnerable systems may require an adaptation of the temporal resolution together with an allowance of a possible interplay of different time scales or time-varying time delays of couplings. More detailed research questions that aim at tracking dynamical resistance in prospective studies can take into account improvements towards real-time computations as well as adaptations of statistical approaches to judge sensitivity and specificity of critical changes in dynamical resistance. Our novel resistance monitoring perspective provides a complementary approach to improve both detection and understanding of precursor dynamics in complex spatio-temporal systems, with potential applications in developing timely and targeted prevention techniques to facilitate better preparedness and to diminish loss of lives and damage to health, ecosystems, infrastructure, and property. We envisage advancements in various disciplines that are concerned with extreme events, ranging from the earth and climate sciences via the life and social sciences to finance and energy sciences.

Materials and Methods
Data. We analyse multiday, multichannel intracranial electroencephalographic (iEEG) signals from 43 subjects with epilepsy (25 female, 18 male, with ages between 9-65, average 34 years) that underwent presurgical evaluation for drug-resistant focal epilepsies. The data were part of previous studies 45,59 , and all subjects had signed informed consent that their clinical data might be used and published for research purposes. The study protocol had previously been approved by the ethics committee of the University of Bonn, and methods were carried out in accordance with the approved guidelines.
Signals were recorded from chronically implanted intrahippocampal depth electrodes and subdural grid-and/ or strip-electrodes. Decisions regarding placement of electrodes were purely clinically driven and were made independently of this study. Signals were band-pass-filtered between 1-45 Hz, sampled at 200 Hz (sampling interval 5 ms) using a 16 bit analog-to-digital converter, and referenced against the average of signals of two electrode contacts outside the presumed focal region. Reference contacts were chosen individually for each subject. The recordings with, on average, 56 electrode contacts lasted between 18-258 h (average 107 h) during which 112 clinical seizures (3 seizures/subject, range 1-7) occurred. The time of seizure onset was visually identified on the iEEG as the time of earliest clear change from the subject's baseline or normal background iEEG that eventually led to an electrographic seizure. We neglect subclinical seizures in our analyses.
Probing the dynamical coupling structure. To quantify the strength of coupling between two sub-systems X and Y, we consider their output signals {x} and {y} each of length T. Since synchronization plays an important role in brain function and dysfunction, here we employ an analysis approach that is based on the concept of phase synchronization. We note that other analysis approaches like e.g. based on information theory provide similar results ( Supplementary Figs S5 to S11).
We calculate in a time-resolved manner (non-overlapping windows of 20.48 s duration; T = 4096 data points) a phase-synchronization index (mean phase coherence) 52,60 x y between phase time series {Φ x } and {Φ y } that we derive by Hilbert-transforming iEEG signals. As an estimator for the strength of interactions, ρ xy takes on values between 0 and 1, indicating either complete asynchrony (zero coupling) or complete synchrony (strongest coupling) between sub-systems X and Y. We perform these steps of analysis for every possible combination of pairs of sampled brain regions (Fig. 1). Identifying dynamical regimes and estimating dynamical resistance. We investigate changes in the similarity of coupling structures by means of a time-resolved clustering analysis (time base of 1 h; k-means algorithm 61 ). With this technique, similar coupling structures from not necessarily consecutive points in time but with small Euclidean separation are subsumed to clusters in an abstract space spanned by all pairwise interactions. We identify these different clusters with different dynamical regimes (Fig. 1). The minimum distance d between regimes is an estimate for the dynamical resistance R of a system (Fig. 6): the larger this distance the higher is the system's capacity to absorb disturbances, given that a larger distance between dynamical regimes requires a stronger perturbation to induce changes. We initialize the clustering algorithm by assuming all coupling structures to belong to the same cluster. Next we divide this initial cluster into two clusters by randomly choosing two initial centroids of these sub-clusters. We then assign coupling structures to either the first or the second cluster according to a minimization of the variance of Euclidean separations of cluster constituents. We use the resulting centroids of these clusters as new initial centroids and repeat the process in an iterative manner 62 until the current centroid equals the previous one. We minimize a possible bias due to the choice of the initial centroids by only considering those centroids for which we achieve a minimum variance of Euclidean separations of cluster constituents after reiterating the aforementioned steps of analysis for randomly chosen initial centroids.
The exact number of a system's dynamical regimes is a priori unknown. We therefore repeat the clustering analysis for each of the two clusters separately in a divisive hierarchical manner 63 until the Euclidean separation between all resulting clusters falls below numerically resolvable values. This approach allows us to identify up to = 10  accessible dynamical regimes when applied to a temporal sequence of coupling matrices built from electroencephalographic signals and to create-for each analysis time frame-a dendrogram-like pattern of dynamical regimes 38 that is sorted by the Euclidean separation between cluster centroids. The minimum Euclidean separation between cluster centroids is taken as the minimum distance d between dynamical regimes (Fig. 6). Testing for differences between dynamical resistances prior to and in between seizures. For our retrospective analyses, we use a statistical approach that is commonly employed in seizure-prediction studies 22 . In doing so, we assume existence of a critical transitional period (pre-seizure period) that is distinguishable from other, non-critical transitional periods that can be observed in between seizures (inter-seizure period). Given that seizure generation is likely to take place over minutes to hours 15 , we here assume the pre-seizure period to last for 4 h. For our analyses, we discard data from the 90 min interval after the onset of a seizure to minimize the risk of biasing our findings with effects related to brain dynamics during and after the seizure 64 . We then inspect separability  of populations of resistance estimates from the inter-seizure periods and from the pre-seizure periods, merging data from all subjects. To this end, we employ the Kolmogorov-Smirnov test to evaluate the difference between the two resistance populations, as it makes no assumption on the underlying distributions. We take the Kolmogorov-Smirnov statistic as an estimate for separability  and evaluate whether it deviates from randomness by testing against the null hypothesis of the non-existence of a critical transitional period. To do so, we employ the concept of seizure time surrogates (STS) 65 that also allows us to account for possible confounding variables such as seizure clustering, daily rhythms, and changes in anticonvulsants. We derive STS from a random permutation of the original inter-seizure intervals and the time interval from the first seizure back to an arbitrarily defined starting point (Fig. 7 left). We reject the null hypothesis if separability for original seizure times exceeds the maximum one obtained with 19 STS (p < 0.05). This was the case for data derived with a number of accessible dynamical regimes  ∈ {2, 6, 7, 8} (Fig. 7 right). We present our findings obtained with  = 6; for the other numbers of accessible dynamical regimes, we obtained comparable findings.  . Geometrical illustration of dynamical resistance. Consider a system whose dynamical regimes C i and C j are separated in an abstract space by the minimum distance d (C i , C j ). (a) If d is large, small endogenous and/ or exogenous perturbations of the system that lead to small modifications of the system dynamics (grey-shaded areas) do not suffice to induce changes (e.g., a switch from regime C i to C j ) or a merging of different regimes: the system's dynamics has a high resistance and can even absorb larger perturbations in the sense that dynamical regimes remain distinguishable. (b) With a small distance between dynamical regimes even tiny perturbations can induce changes or a merging of regimes; the system's dynamics has a low resistance and can not absorb perturbations in the sense that dynamical regimes remain distinguishable.