Semi-automatic Extraction of Functional Dynamic Networks Describing Patient's Epileptic Seizures

Abstract

Intracranial electroencephalography (EEG) studies using stereotactic EEG (SEEG) have shown that during seizures, epileptic activity spreads across several anatomical regions from the seizure onset zone toward remote brain areas. A full and objective characterization of this patient-specific time-varying network is crucial for optimal surgical treatment. Functional connectivity (FC) analysis of SEEG signals recorded during seizures enables to describe the statistical relations between all pairs of recorded signals. However, extracting meaningful information from those large datasets is time consuming and requires high expertise. In the present study, we first propose a novel method named Brain-wide Time-varying Network Decomposition (BTND) to characterize the dynamic epileptogenic networks activated during seizures in individual patients recorded with SEEG electrodes. The method provides a number of pathological FC subgraphs with their temporal course of activation. The method can be applied to several seizures of the patient to extract reproducible subgraphs. Second, we compare the activated subgraphs obtained by the BTND method with visual interpretation of SEEG signals recorded in 27 seizures from nine different patients. As a whole, we found that activated subgraphs corresponded to brain regions involved during the course of the seizures and their time course was highly consistent with classical visual interpretation. We believe that the proposed method can complement the visual analysis of SEEG signals recorded during seizures by highlighting and characterizing the most significant parts of epileptic networks with their activation dynamics.

1. Introduction

About 30–40% of epileptic patients are drug resistant (1). For those patients, surgical resection of the epileptogenic brain structures is considered to promote seizure freedom (1). Intracranial EEG using depth EEG recordings or the stereoencephalography (SEEG) method is often required to guide tailored-surgical resection (2, 3).

The primary aim of SEEG is to delineate precisely the epileptogenic regions. However, since the pioneering works in SEEG, it has been shown that seizures cannot be considered as static phenomena with a single focus activation leading to clinical manifestations (4). SEEG recordings of focal seizures typically show that the epileptic activity spreads during seizures across several anatomical regions. It begins at the seizure onset zone and spreads toward remote brain areas. This dynamical pathological process can be described by several brain states characterized by transient and abnormal connectivity profiles within the epileptogenic network (5). A full and objective characterization of this patient-specific dynamic network is crucial for optimal surgical treatment.

Understanding brain network modifications operating at different time scales in the interictal state and during seizures is a very active stream of research. This is in line with computational neuroscience studies modeling neurological diseases as brain network disease (6–11).

In the field of epileptology, the structure of epileptic (or ictal) networks has been explored using connectivity analysis (functional or effective connectivity) and more recently, graph-theory analysis, (5). Functional connectivity (FC) approaches, based on linear or non-linear measures, quantify the statistical relations between any pair of recorded SEEG signals and their evolution across time. In most of the studies, pairwise correlations between remote sites are computed at specific time points [for a review, see (5)]. Those approaches unveil the organizational interactions of different regions of interest in the brains. They have proved very fruitful to investigate the neurophysiological correlations of symptoms during seizures, or to describe several subtypes of focal epilepsy involving specific networks (12–14). Graph theory is the study of graphs, which are mathematical structures used to model networks, and specifically pairwise relations between objects (8, 15). A graph is made up of vertices (also called nodes) that are connected by edges. Graph theory approaches allow the description of both local and global characteristics (16). For epilepsy, the nodes usually represent electrode contacts and the edges represent the FC measures. The resulting graph structure is known to contain relevant fingerprints of the seizure dynamic. Studying network structures with graph theory provides mathematical tools to investigate different subtypes of epilepsy (17). For example, it has been shown that the properties of networks' topology are different for temporal lobe, mesial temporal lobe, and neocortical epilepsy (18). Moreover, some studies suggest that investigating local properties of the network structure through graph theory concepts provides biomarkers for epileptogenic focus localization (19–21). Lastly, dynamic graph theory can also provide step-by-step modeling of the propagation of the seizure in the brain (19, 22).

A major challenge for the study of ictal networks is that seizures are a highly dynamical process with rapid transitions between network states (5, 23). To track network changes during the seizures, the network organization has to be described at short time scales. In the most straightforward approach, FC approaches forming the backbone of the network are estimated at different time steps of the seizure. Therefore, they are prone to create spurious connections and constitute a noisy estimation of the pathologic dynamic network. Moreover, as the brain activity is commonly monitored with more than 100 electrode contacts distributed along the stereotactic rods, the number of FC scales quadratically (100 electrode contacts providing around 5000 FC measures), making the study of FC over time a resources consuming task that requires high expertise.

One advantage of SEEG monitoring is to allow the recording of several seizures from the same patient, with the same measurement points (or nodes), yielding thus as many realizations of identically distributed dynamical networks of FC. Hence, several seizures from the same patient can be investigated that may share some connectivity features but with different dynamics. Methodological tools have to be proposed to extract relevant information from those large datasets. Ideally, such methods should summarize the dynamics of the seizures by providing several epileptogenic networks, with stable network structures, that characterize the different steps of the seizures and that are involved in the production and the propagation of the ictal events. Those network states may be common to the different seizures of the patient, but the timeline activation should follow a pattern that remains specific to each seizure (22).

Along these lines, we propose a novel semi-automatic method to characterize the dynamic epileptogenic network quantitatively and across time using SEEG signals. We propose to perform the joint analysis of all seizures of the same patient, first to reduce the measurement uncertainty in the calculation of FC indices, and second to make robust the identification of a time-functional pattern, systematic in all seizures and characteristic of a patient's pathology. First, the full FC matrix is computed for each time step using a classical FC measure, namely the Phase Locking Value (PLV) (24, 25). Then, the method extracts several pathological subgraphs with their own activation score during seizures. We expect each subgraph to comprise several brain nodes with high connectivity values. The paper is organized as follows: We extend our previous work (26) to seizures with different durations, we call this method the Brain-wide Time-varying Network Decomposition (BTND). On the application side, we validate the clinical use of the method on a larger clinical dataset.

2. Materials and Methods

2.1. Description of the Brain-Wide Time-Varying Network Decomposition Method

For each patient, the dataset is composed of several SEEG recordings for different seizures. Seizures can have different durations. The proposed strategy can be summarized in four steps: (a) We chop each recording into short segments; (b) for each segment, we estimate via FC measures the connectivity for each pair of electrode contacts; (c) we rearrange the FC measures into a list of matrices representing the time evolution of FC for each seizure of a patient; (d) the list of matrices representing the multi-seizures brain-wide time-varying network is decomposed into FC subgraphs characteristic of one patient but common to all his seizures, along with their activation profile specific to each seizure. The main steps of the method are illustrated in Figure 1.

Figure 1

2.2. Representation of the Multi-Seizure Brain-Wide Time-Varying Network

Practically, FC measurements are stored in a three-dimensional structure X_lt{s} corresponding to the FC of index l at time segment t and for the seizure s. The list of matrix X{s} ∈ ℝ^L×T(s) ∀s (1, ..., S) is the mathematical representation of the multi-seizures brain-wide time-varying network. Here, S is the total number of seizures recorded for the same patient. We remind that the number of time steps T(s) can vary for each seizure. Figure 1C illustrates an example of list of matrices X{s}.

2.3. The Optimization Problem Related to the BTND Method

The core of the BTND method is then to seek, for each seizure, for the following decomposition:

where F ∈ ℝ^L×K contains the K FC subgraphs and V{s} ∈ ℝ^T(s)×K are their respective temporal activations corresponding to the specific seizure s. This approach directly entails the requested decomposition since the columns of the matrix F contains the weights of the edges in the sub-graphs, as it can be seen in Figure 1D and matrices V{s} directly correspond to the activation profiles of each seizure depicted in Figure 1D.

However, the solution for the decomposition (1) is not unique, and to favor handily interpretability from the medical viewpoint, we impose several constraints on the components F and V{s}:

• because most of the functional connectivity measures and activation indices are naturally positive values, we impose F and V to be non-negative matrices;

• to limit the complexity of the inferred subgraphs, we restrict the number of non-zero significant FC values, yielding a sparse matrix F. The sparsity constraint is meaningful in the context of epilepsy, where a large number of functional connectivities can be passively implied in the neurological process;

• to promote FC subgraphs that are continuously activated over specific periods, we impose sparsity and compactness on V{s}. These two constraints drastically improve the interpretability of the solution, as they prompt sparse and piecewise continuous activation periods that are close to cluster-like solutions. Formally, they correspond to the fused lasso constraint that reads:

  where the parameter γ and the parameter η compel each activation profile v = V_:k{s}, to be sparse and compact, respectively. As for the third term in expression (2), it prevents incoherent solutions due to scaling indeterminacy (27).

Finally, the BTND boils down to an instance of joint non-negative matrix factorization (28–31) that takes on the following form:

The ζ_s are free parameters to balance the relative importance of each seizure (in our study, we consider all seizures evenly and select the ζ_s parameter according to the energy of each seizure, see the Supplementary Materials for more precision). The sparsity factor γ_s and the compactness factor η_s are chosen to adapt to the particular duration of each seizure s.

Then, the method is associated with three hyperparameters: λ and γ, respectively, controlling the sparsity level of the FC subgraphs and the activation profiles, and η regulating the temporal compactness of the activation profiles.

2.4. Additional Comments on BTND

Let us stress that the optimization problem of Equation (3) is non-convex, and therefore, different initial conditions yield different solutions corresponding to local maxima. As it is common practice for non-convex methods in machine learning [e.g., k-means for data clustering, (32)], we repeatedly solve (3) with a different initialization and retain the one that reached the smaller minimum cost function. In practice, we empirically chose 20 trials. As for the other hyper-parameters, we observed that η = 0.2 produces coherent activation profiles. Also, fixing λ = γ simplifies the procedure without altering the results significantly. Then, λ is tuned so as to identify the 20% most resilient activated FC. To select the most pertinent number of subgraphs, we successively compute the decomposition for K ranging from 3 to 10. Based on an Elbow criterion and on visual inspection, we compare the quality of the resulting temporal activation and subgraphs, identifying thus the best value for K. Finally, we normalize each subgraph such that their connectivity strengths are between 0 and 1. We only retain connections above some threshold (empirically set to 0.2 in our experiments) to eliminate non-significant interactions.

For more details on the practical use of the method, see the Supplementary Materials. We provide the URL¹ for a Github repository with Matlab implementation of the proposed BTND method.

3. Application on a Real Dataset of Epileptic Patients

3.1. Patients

To illustrate the clinical relevance of the BTND method, we applied the method on seizures recorded with intracranial EEG in 9 epileptic patients.

We included 9 adult patients suffering from drug-resistant focal epilepsy, followed in the Department of Functional Neurology and Epileptology at Lyon's University Hospital, who underwent intracranial EEG with SEEG according to the following criteria: (i) at least 1 seizure recorded during long-term monitoring; and (ii) conventional visual analysis of the SEEG signals identified clearly the seizure-onset zone.

Clinical details of all the patient included are listed in Table 1.

Among them, 8 patients presented seizures suggesting a temporal lobe involvement but with clinical features or morphological alterations on brain magnetic resonance imaging (MRI) not typical for medial temporal epilepsy requiring intracranial EEG. For one patient, clinical semiology suggested an involvement of operculo-insular cortex. Three patients underwent surgical resection of the epileptogenic cortex and all had a good surgical outcome (Engel class Ia for all patients with a follow-up duration between 4 and 48 months). For 2 patients, a focal thermolesion using SEEG electrodes was performed, resulting in a dramatic improvement of epilepsy (Engel Ia for both patients with 4 and 5 months of follow-up). For 2 patients, surgery was contraindicated because of the involvement of both temporal lobes during seizures. For 2 patients, surgical resection is planned based on SEEG findings but has not been performed at the time of the present study.

This study, involving human participants, was reviewed and approved by Ethics Committee CPP Lyon Sud EST IV (24/05/2012 N2012-A00516-37). The patients/participants provided their written informed consent to participate in this study.

3.2. SEEG Recordings

Intracerebral multi-contact electrodes (5–15 contacts, diameter 0.8 mm, length 2 mm, and 1.5 mm apart) were implanted according to Talairach's stereotactic method (3). Electrode location was verified with post-implantation MRI. Prolonged extra-operative recordings were performed to capture each patient's habitual seizures.

SEEG data were acquired with a 256 channel video EEG monitoring system, Micromed video EEG acquisition system (SD LTM express, Micromed, Treviso, Italy), using the following parameters: sampling rate 256 Hz, high-pass filter 0.15 Hz, low-pass filter 200 Hz, notch filter at 50 Hz.

The median number of bipolar contacts recorded per patient was 101 (range 64–130). The main cerebral structures targeted by intracranial electrodes for each patient are listed in the Supplementary Table 1. As a whole, for all patients, the medial temporal lobe (anterior hippocampus, posterior hippocampus, amygdala, entorhinal cortex), the temporal lateral neocortex and the insular cortex were always targeted. Depending on the electroclinical findings of each patient, frontal lobe, parietal lobe, and occipital cortex were targeted.

For 9 patients, all available seizures were extracted, forming a dataset composed of a total of 27 seizures. For 6 patients, 3 seizures were analyzed. For 2 patients, 4 seizures were analyzed and for 1 patient, a single seizure was available for analysis. For each seizure, we extracted signals for at least 1 min before the onset of the seizure and the whole course of the seizure.

The duration of the seizure was largely heterogeneous both at the inter-individual and intra-individual level. The median length of the seizures across patients was 96 s (range 18–337).

3.3. SEEG Signal Analysis

SEEG signals are considered in bipolar derivations; hence each signal is referenced to its closest neighbor. A high-pass filter, with a cut-off frequency equal to 20 Hz at -3 dB was applied on SEEG signals in order to highlight the high-frequency activity typical of seizure activity, particularly at seizure onset (33).

For each seizure, we computed the FC matrices during the pre-seizure period and the course of the whole seizure using a classical connectivity measure, the PLV (24, 25). Briefly, the PLV quantifies the synchronization in phase between two signal. The phase of each signal is computed by the Hilbert transform. The PLV is the time average of the relative phase difference. To compute the FC matrices, the SEEG signals were windowed with 4 s sliding windows moving by steps of 1 s. For each 4 s window, the PLV between all pairs of bipolar contacts was computed to produce the overall network at each time step.

Finally, the BTND method is applied to the set of several seizures for each patient, to decompose all seizures in several subgraphs, activating through time.

3.4. Comparison of Network Dynamics Estimated Through Conventional Visual Analysis and the BTND Method

Each seizure was pragmatically segmented in three main periods defined by visual analysis: seizure-onset, seizure propagation, and seizure ending. Visual analysis of the seizures was performed by an expert in clinical SEEG interpretation (JJ). Seizure onset corresponded to the time period with a dramatic change of SEEG signals with either low-voltage fast activity (typically above 20 Hz) or rhythmic spikes in a subset of electrode contacts. Seizure propagation corresponded to an extended time period where ictal SEEG discharge spread to several brain structures either locally or remotely from the seizure onset zone. The recruitment of these regions in the propagation zone can happen either by independent activation of the single areas or by activating multiple areas at the same time. Lastly, a seizure was supposed to have ended when the activation across brain structures was mostly synchronous (typically synchronous spikes) and stable in time and resolved ultimately.

For each time period, we determined the electrode contacts that were involved in the ictal wave with a conventional visual inspection. The electrode contacts were then pooled in several anatomical predefined subregions.

For each patient, the output of the BTND method provided the temporal profile of activation of several common subgraphs of the whole network during each seizure. The list of activated subgraphs at each time period of the seizure was then collected. Each subgraph included several contacts with strong functional connectivity. At each time period, we determined which anatomical subregions were connected based on the activated subgraphs.

Lastly, a qualitative comparison between the set of activated structures determined by visual analysis and the BTND method was performed for each seizure.

3.5. Results for the Real Dataset of Epileptic Patients

The overall functional connectivity organization was extracted in all 27 seizures using the BTND method. This means that the seizures (from 1 to 4) of one patient are processed together according to BTND. Table 2 provides the qualitative comparison between the set of activated structures between visual analysis and the BTND method for patients 1 and 2. Also, Tables 3, 4 show the same qualitative comparison for, respectively, patients 2–4 and patients 5–9.

Using this method, we found, in the case of 6 patients, that 6 distinct functional subgraphs characterized the organization of seizures; 7 subgraphs for 1 patient; 5 subgraphs for another one; and only 4 subgraphs characterized seizures for the last patient. However, for each patient, some subgraphs were more strongly activated before seizure onset or were continuously activated before seizure onset and remained active during the course of the seizures. Those subgraphs were considered as non-specific subgraphs for the ictal events. For 5 patients, 2 subgraphs were non-specific while for 4 patients a single subgraph was non-specific.

At seizure onset, a single subgraph was activated for 1 patient, two subgraphs were activated for 4 patients, three subgraphs were activated for 3 patients, and 5 subgraphs were activated for 1 patient. For 24 seizures in 6 patients, the seizure onset determined by visual analysis overlapped closely with the network disclosed by the BTND method. For those patients, the cortical regions underlying the seizure-onset zone determined through visual analysis were included in the seizure-onset subgraphs. However, the seizure-onset subgraph also included other regions with strong functional connectivity not directly outside of the seizure-onset zone. For one of those 6 patients (Pt 2), the seizure involved either the left or the right medial temporal lobe at seizure onset. The seizure-onset subgraph was different for each seizure, and the lateralization of the activated structures was concordant with visual analysis. For 3 seizures in one patient (Pt 9), the seizure-onset subgraph was discordant from the seizure onset-zone. For this patient, the seizure-onset zone involved either right or left medial temporal lobe depending on the seizure. The seizure-onset subgraph for this patient was wrongly lateralized to the right or the left temporal lobe.

During seizure propagation, there was always a close spatial overlap between the activated subgraphs and the brain regions involved at each part of the seizure in all patients. For the 27 seizures in the 9 patients, the brain regions involved during seizure propagation were included in activated subgraphs. However, the congruence between activated subgraphs and regions disclosed by visual analysis was not perfect: a minority of regions were revealed by the BTND method but was not detected by visual analysis.

During seizure ending, a tight spatial overlap was also observed between activated subgraphs and brain regions determined by visual analysis. For 23 seizures, the brain regions involved at seizure ending were included in activated subgraphs. For 4 seizures in 1 patient, visual analysis disclosed more activated regions than the BTND method (Pt 4).

The detailed results are now presented for two cases (Pt 1 and Pt 2).

CASE 1:

Pt 1 is a 49 years old male patient. Presurgical non-invasive investigations suggested left temporal lobe epilepsy but some radiological features were considered as atypical for mesial temporal lobe epilepsy syndrome and prompted invasive EEG with SEEG. SEEG targeted several regions within left temporal lobe (anterior hippocampus, posterior hippocampus, amygdala, temporal pole, anterior temporal neocortex, posterior temporal lobe), left orbito frontal cortex, and right temporal lobe (right amygdala, right anterior temporal neocortex).

Three seizures were recorded during SEEG. During the 3 seizures, the initial seizure-onset activity developed in left anterior and posterior hippocampus with secondary involvement of the temporal pole, amygdalar nucleus, and left anterior temporal neocortex at the end of the seizures.

The BTND method applied to the three seizures decomposed the connectivity pattern in 6 subgraphs. Figure 2 shows the recording of two seizures (seizure 1 and 2) of the patient 1 for selected electrode contacts. Below each recording, we provide the activation profiles of all subgraphs obtained by the BTND for this specific seizure. On top is represented the main cerebral structures targeted by intracranial electrodes for this patient. Figure 3 shows the 6 FC subgraphs revealed by the BTND.

Figure 2

Figure 3

One subgraph was active before the seizure and during the whole course of the seizures. This subgraph was mostly composed of local connections within temporal lobe (mostly within anterior hippocampus, posterior hippocampus, amygdala). At the seizure onset, during the first seconds of the seizure, there was a reproducible activation of one subgraph, that involved mostly connections between anterior hippocampus and amygdala, posterior hippocampus, and posterior temporal neocortex. A few seconds later, a strong activation of another subgraph was observed that involved mostly connections between medial temporal lobe and temporal pole. During the course of the seizures, there was a consistent activation of three other subgraphs, with a very similar pattern between seizures. Those subgraphs contained mostly connections between anterior lateral temporal neocortex, temporal posterior neocortex, and temporal pole. As a whole, the pattern of activations was very similar for the three seizures and was very consistent with visual analysis of the seizures.

From a clinical point of view, at seizure onset, the patient was asymptomatic and the seizure remained clinically silent for almost 80 s. First clinical symptoms (behavioral arrest and loss of consciousness) occurred more than 1 min after seizure onset during the course propagation. Secondary clinical manifestations included oro alimentary automatisms (seizure 1) and right arm dystonia (seizure 2). Those symptoms occurred while several modules were simultaneously activated.

The patient underwent left anterior lobectomy that resulted in seizure freedom with more than 24 months of follow-up.

CASE 2:

Pt 2 is a 37 years old female patient. Presurgical non-invasive investigations suggested that both temporal lobes could trigger habitual epileptic seizures of the patient. Intracranial EEG using SEEG was thus required to evaluate the intrinsic epileptogenicity of each temporal lobe. Intracranial SEEG electrodes targeted mostly both medial and lateral temporal lobes (left and right anterior hippocampus, left and right temporal pole, right amygdala, left and right anterior temporal neocortex, left and right posterior temporal neocortex, left and right insula), but also left and right orbitofrontal cortex.

Three seizures were recorded during SEEG. For the seizure 1, seizure-onset was characterized by a rapid discharge in the left anterior and posterior hippocampus with secondary spread to left anterior lateral temporal and left orbito frontal cortex fast activity and at the end of the seizure propagation to right temporal lobe (right hippocampus and right amygdala). For two seizures (seizures 2 and 3), the initial seizure-onset activity developed in right hippocampus, right amygdala, and right entorhinal cortex with a secondary ictal spread to right anterior temporal cortex and with a propagation to left temporal lobe at seizure ending.

The BTND method applied to the three seizures decomposed the connectivity pattern in 7 subgraphs. Figure 4 shows the recording of two seizures (seizures 1 and 2) of the Patient 2 for selected electrode contacts. Below each recording, we provide the activation profiles of all subgraphs obtained by the BTND for this specific seizure. On top is represented the main cerebral structures targeted by intracranial electrodes for this patient. Figure 5 shows the 7 FC subgraphs revealed by the BTND.

Figure 4

Figure 5

One subgraph was mostly composed of connections within left medial temporal lobe and left orbitofrontal cortex. Another subgraph was mostly composed of connections within left anterior temporal neocortex. The remaining five subgraphs consisted of regions connecting mostly right medial temporal structures (hippocampus, amygdala, entorhinal cortex) and/or right lateral temporal neocortex. The time course of activation of those subgraphs was closely related to the ictal involvement of both temporal lobes revealed by visual analysis of the seizures: involvement of the left (right) medial temporal lobe was paralleled by an activation of the left (right) subgraphs in a timely fashion, depending on the seizure.

Clinically, ictal semiology for seizure 1 consisted of nausea and olfactory hallucinations with preserved consciousness reported consciously by the patient 30 s after EEG onset followed by behavioral arrest with loss of consciousness more 60 seconds later. Loss of consciousness occurred lately during the course of the seizure when both temporal lobes were involved. At that time, the BTND method showed activation of modules in both temporal lobes. Seizure 2 was a nocturnal seizure with mild clinical semiology, mostly consisting of nocturnal arousal and confusion.

Since SEEG revealed an intrinsic epileptogenicity of both temporal lobes, surgical resection was contraindicated.

4. Discussion

This study investigated a new method named BTND to decompose the multi-seizure brain-wide time-varying network obtained by means of FC. The dataset is visualized as FC sub-graphs characterizing the dynamic of all seizures from the same patient. The FC measure used was the PLV, estimated at different time steps of the seizure and applied on large band signal (20–100 Hz). We compared the obtained decomposition of ictal events from 9 patients who have drug-resistant focal epilepsy (observation of a total of 27 seizures) to the visual interpretation from the clinician. Overall, for every patient, results consisting of spatially localized FC subgraphs with stepwise activation were easily interpretable. For 8 of the 9 patients, the decomposition matched with the clinical observation entailing the BTND method as a relevant tool to visualize the multi-seizure brain-wide time-varying network.

4.1. Investigating Seizure Dynamics With Brain-Wide Time-Varying Network

It is well-established that the brain is a complex network, with a sophisticated structural connectivity architecture and specific anatomical networks shaping sensory and cognitive processes (34, 35). Functional connectivity measures are a statistical way to investigate the interrelations between brain regions, forming a physiological or pathological brain network (36). When applied on static processes, FC network analysis can provide an instantaneous picture of a stable network. In the field of clinical neuroscience, this proved to be useful for studying stable disease traits like the effects of specific lesions (10) or the effects of drugs on the brain (6, 11, 37). Graph theory-based measures can then provide quantitative tools to explore the overall topology of the network (35). Thus, specific metrics such as modularity, clustering coefficients, or efficiency originating from graph theory (16) were proposed as a useful strategy, for example, to classify patients with Alzheimer disease from standard patients (9) or to explain the effect of physical activity on relations between brain regions (38).

However, static network analysis does not capture one fundamental property of epileptic seizures, the dynamic propagation of the ictal wave. Analysis of time-varying network inferred by dynamic FC measures is a recent topic. Newly emerging dynamic measures that quantify how community organization evolved in time have been proposed in recent years (39). In this context, Kerr et al. (40) shows that simple metrics using the first eigenvector of each FC network lead to the separation of ictal, pre-ictal, or non-ictal events of a recording. A similar metric demonstrated in (22) used to describe the seizure as a succession of states. The strategy can be used to decompose the time axes in states, and then to extract the major FC network of each state (41). Despite producing intelligible results, this strategy hampers the identification of FC subgraphs with interconnected temporal activation, and prevents highlighting some complex relations between brain regions. A complementary approach has been proposed in (23): the authors identify first the main subgraphs of the seizures by modularity optimization (42), and the evolution in time of the main subgraphs leads to a decomposition of the seizure in time states. It should be noticed that this approach is not fundamentally different from the study of static graphs since the subgraphs and their time evolution are determined independently. However, in addition to finding a measure characterizing each FC graph, one must know how to analyze the temporal evolution of the proposed scores. It becomes even more complicated when several modalities are used, as in (43), where several seizures are analyzed for different frequency bands over time.

The present study proposes another strategy, consisting of decomposing all modalities characterizing the dataset simultaneously. Thus, the main advantages of the proposed method are that the analysis pipeline provides both the subgraphs and their temporal activations. The output of the process highlights the main components of the connectivity structure and summarizes a large amount of data with an automatic approach. A similar approach has already been employed in the context of epilepsy in (44) to decompose FC matrices from seizures into FC subgraphs with their respective activation profile. Khambhati et al. (44) demonstrate that inferred FC subgraphs during interictal periods can predict brain regions that generate seizures, and that those subgraphs undergo slower and more coordinated fluctuations during the ictal events compared to interictal states. However, this kind of simultaneous decompositions can produce results that are difficult to interpret. Therefore, for an application different from epilepsy, a sparsity constraint applied on activation profile was shown to enforce intelligibility of the results (45), providing FC subgraphs discriminating the brain network's dynamic in neurodevelopment. Moreover, contrary to our study, the decomposition from (44) does not integrate the different seizures of the same patient to obtain more reproducible FC subgraphs. Tools to decompose several modalities, like tensor decomposition (46–48), were already applied in neuroscience (41, 49, 50) and recently in the context of epilepsy (26). In (26), we proposed a specific tensor decomposition with relevant constraints to encourage the inference of interpretable clusters of FC common to several seizures of the same patient. However, this method is only applicable if all seizures from the same patient have similar durations. In the present study, we developed a new method offering the possibility of decomposing several seizures with different durations from the same patient, which is a more realistic situation in a clinical setting.

4.2. Limitations

The BTND method can produce relevant FC, but several limitations have to be addressed. First, the method requires 3 parameters that are directly related to the obtained FC results. We expose a simple procedure to choose each parameter that produces relevant results for each of the 27 seizures. However, the optimal selection of these parameters may be dependent on the FC measure and might be adapted for other clinical applications. Second, for one patient, the BTND method identified a seizure-onset subgraph that was discordant from the seizure onset-zone. For that patient, at seizure-onset, there was a rapid discharge within the seizure-onset zone (without any focal change of synchrony) and a spiking activity within the contralateral hemisphere (accompanied with an increase in synchrony). The increase of FC evaluated by the phase-locking value was wrongly lateralized, emphasizing that the BTND method is dependent on the FC measure used to infer the brain-wide time-varying network. Thus, using another FC measure could be envisaged for this patient. In addition, several studies have already shown that seizure-onset is often marked by a dramatic decrease in seizure connectivity among recorded brain structures. At the same time, synchrony increases progressively during the seizure (51, 52). It might be thus expected that seizure-onset should be characterized by a decrease of activation of subgraphs located within the seizure-onset zone with the BTND method in some patients. It is important to consider that our procedure only highlights functional subgraphs associated with high values of functional connectivity. Thus, it is ideal for showing activations of synchrony in different areas of the brain. However, it is not able to explicitly demonstrate the deactivations that may occur in the brain at the early start of the seizure when SEEG activities in different areas of the brain are suddenly decorrelated. Lastly, direct validation of the method is out of reach since there is no perfect gold standard for the estimation of the connectivity pattern in epileptic patients. We chose to make a correspondence between the propagation patterns disclosed by classical visual interpretation of SEEG signals and the BTND method. Simulation studies generating neural models of epileptic activity with a known connectivity pattern could represent an alternative in future studies.

4.3. Clinical Application of the Method

In the present study, we presented an automatic method to describe the connectivity structures of epileptic seizures using an original algorithm with several constraints optimized for that clinical context.

As a whole, we found that the method produces several subgraphs of connections with their activation time course that parallel the patterns of propagation of the seizures closely. In 24/27 seizures, one or two subgraphs overlapped clearly with the seizure onset zone. This suggests that the method can help to localize the seizure-onset zone if there is an increase in synchrony at seizure-onset revealed by FC. This finding confirms several studies evidencing that focal modulations of synchrony help to localize the seizure-onset zone to be surgically resected (5, 53). Moreover, thanks to the sparsity and temporal coherence constraints, the method effortlessly reveals how FC synchrony propagates to the different regions of the brain. As a whole, the method summarizes a vast amount of complex interactions with high readability. We believe that the method is thus highly valuable in clinical studies focusing on connectivity in epileptic patients. For example, this may help to unravel the neural bases of clinical semiology of seizures, which are often related to ictal dysfunction of widespread brain networks involving cortical or subcortical structures. The BTND method provides an exhaustive view of the structure of functional networks at each period of the seizure enabling a fine-grained correlation ictal semiology and network analysis. Lastly, clinical interpretation of SEEG signals is mostly focused on changes of power in large frequency bands (typically high-frequency above 20 Hz at seizure-onset and low band or large band frequency during propagation) revealed by visual inspection. The BTND method enables the detection of changes of synchrony revealed by functional connectivity measures, which are not necessarily paralleled by changes of power. In that respect, we believe that the method may thus contribute to bringing into clinical practice computational measures complementary to a visual interpretation of seizures.

5. Conclusion

We present here a novel approach to decompose epileptic seizures in several time-varying subgraphs with high and reproducible functional connectivity values. The method extracts the most significant subgraphs and their corresponding time course of activation. We suggest that this represents a first step to simplify the interpretation of large datasets of functional connectivity for clinical practice. We believe that this will enable further studies investigating the clinical relevance of networks identification for epilepsy surgery.

References

• 1.

  RyvlinPCrossJHRheimsS. Epilepsy surgery in children and adults. Lancet Neurol. (2014) 13:1114–26. 10.1016/S1474-4422(14)70156-5

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 2.

  GuenotMIsnardJRyvlinPFischerCOstrowskyKMauguiereFet al. Neurophysiological monitoring for epilepsy surgery: the Talairach SEEG method. Stereotact Funct Neurosurg. (2001) 77:29–32. 10.1159/000064595

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 3.

  IsnardJTaussigDBartolomeiFBourdillonPCatenoixHChassouxFet al. French guidelines on stereoelectroencephalography (SEEG). Neurophysiol Clin. (2018) 48:5–13. 10.1016/j.neucli.2017.11.005

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 4.

  ChauvelPGonzalez-MartinezJBulacioJ. Presurgical intracranial investigations in epilepsy surgery. Handb Clin Neurol. (2019) 161:45–71. 10.1016/B978-0-444-64142-7.00040-0

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 5.

  BartolomeiFLagardeSWendlingFMcGonigalAJirsaVGuyeMet al. Defining epileptogenic networks: contribution of SEEG and signal analysis. Epilepsia. (2017) 58:1131–47. 10.1111/epi.13791

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 6.

  HallettMde HaanWDecoGDenglerRDi IorioRGalleaCet al. Human brain connectivity: Clinical applications for clinical neurophysiology. Clin Neurophysiol. (2020) 131:1621–51. 10.1016/j.clinph.2020.03.031

  • CrossRef
  • Google Scholar

• 7.

  van MierloPHöllerYFockeNKVulliemozS. Network perspectives on epilepsy using EEG/MEG source connectivity. Front Neurol. (2019) 10:721. 10.3389/fneur.2019.00721

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 8.

  van den HeuvelMPSpornsO. A cross-disorder connectome landscape of brain dysconnectivity. Nat Rev Neurosci. (2019) 20:435–46. 10.1038/s41583-019-0177-6

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 9.

  RubinovMKnockSAStamCJMicheloyannisSHarrisAWWilliamsLMet al. Small-world properties of nonlinear brain activity in schizophrenia. Hum Brain Mapp. (2009) 30:403–16. 10.1002/hbm.20517

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 10.

  SponheimSRMcGuireKAKangSSDavenportNDAviyenteSBernatEMet al. Evidence of disrupted functional connectivity in the brain after combat-related blast injury. Neuroimage. (2011) 54:S21–9. 10.1016/j.neuroimage.2010.09.007

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 11.

  AlonsoJFMananasMARomeroSHoyerDRibaJBarbanojMJ. Drug effect on EEG connectivity assessed by linear and nonlinear couplings. Hum Brain Mapp. (2010) 31:487–97. 10.1002/hbm.20881

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 12.

  BartolomeiFTrebuchonAGavaretMRegisJWendlingFChauvelP. Acute alteration of emotional behaviour in epileptic seizures is related to transient desynchrony in emotion-regulation networks. Clin. Neurophysiol. (2005) 116:2473–9. 10.1016/j.clinph.2005.05.013

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 13.

  BartolomeiFGavaretMHewettRValtonLAubertSRegisJet al. Neural networks underlying parietal lobe seizures: a quantified study from intracerebral recordings. Epilepsy Res. (2011) 93:164–76. 10.1016/j.eplepsyres.2010.12.005

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 14.

  BartolomeiFBarbeauEJNguyenTMcGonigalARegisJChauvelPet al. Rhinal-hippocampal interactions during déjà vu. Clin. Neurophysiol. (2012) 123:489–95. 10.1016/j.clinph.2011.08.012

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 15.

  StaceyWKramerMGunnarsdottirKGonzalez-MartinezJZaghloulKInatiSet al. Emerging roles of network analysis for epilepsy. Epilepsy Res. (2020) 159:106255. 10.1016/j.eplepsyres.2019.106255

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 16.

  BullmoreESpornsO. Complex brain networks: graph theoretical analysis of structural and functional systems. Nat. Rev. Neurosci. (2009) 10:186–98. 10.1038/nrn2575

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 17.

  HongSJBernhardtBCGillRSBernasconiNBernasconiA. The spectrum of structural and functional network alterations in malformations of cortical development. Brain. (2017) 140:2133–43. 10.1093/brain/awx145

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 18.

  BartolomeiFBettusGStamCJGuyeM. Interictal network properties in mesial temporal lobe epilepsy: a graph theoretical study from intracerebral recordings. Clin. Neurophysiol. (2013) 124:2345–53. 10.1016/j.clinph.2013.06.003

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 19.

  Van MierloPCarretteEHallezHRaedtRMeursAVandenbergheSet al. Ictal-onset localization through connectivity analysis of intracranial EEG signals in patients with refractory epilepsy. Epilepsia. (2013) 54:1409–18. 10.1111/epi.12206

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 20.

  WangMyWangJZhouJGuanYgZhaiFLiuCqet al. Identification of the epileptogenic zone of temporal lobe epilepsy from stereo-electroencephalography signals: a phase transfer entropy and graph theory approach. Neuroimage Clin. (2017) 16:184–95. 10.1016/j.nicl.2017.07.022

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 21.

  WilkeCVan DrongelenWKohrmanMHeB. Neocortical seizure foci localization by means of a directed transfer function method. Epilepsia. (2010) 51:564–72. 10.1111/j.1528-1167.2009.02329.x

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 22.

  BurnsSPSantanielloSYaffeRBJounyCCCroneNEBergeyGKet al. Network dynamics of the brain and influence of the epileptic seizure onset zone. Proc Natl Acad Sci USA. (2014) 111:E5321–30. 10.1073/pnas.1401752111

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 23.

  KhambhatiANDavisKAOommenBSChenSHLucasTHLittBet al. Dynamic network drivers of seizure generation, propagation and termination in human neocortical epilepsy. PLoS Comput Biol. (2015) 11: e1004608. 10.1371/journal.pcbi.1004608

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 24.

  LachauxJPRodriguezEMartinerieJVarelaFJ. Measuring phase synchrony in brain signals. Hum Brain Mapp. (1999) 8:194–208. 10.1002/(sici)1097-0193(1999)8:4<194::aid-hbm4>3.0.co;2-c

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 25.

  AydoreSPantazisDLeahyRM. A note on the phase locking value and its properties. Neuroimage. (2013) 74:231–44. 10.1016/j.neuroimage.2013.02.008

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 26.

  FrusqueGMJungJBorgnatPGoncalvesP. Multiplex network inference with sparse tensor decomposition for functional connectivity. IEEE Trans Signal Inform Process Netw. (2020) 6:316–28. 10.1109/TSIPN.2020.2984853

  • CrossRef
  • Google Scholar

• 27.

  MairalJBachFPonceJSapiroG. Online learning for matrix factorization and sparse coding. J Mach Learn Res. (2010) 11:19–60. 10.1145/1553374.1553463

  • CrossRef
  • Google Scholar

• 28.

  WangHQZhengCHZhaoXM. j NMFMA: a joint non-negative matrix factorization meta-analysis of transcriptomics data. Bioinformatics. (2015) 31:572–80. 10.1093/bioinformatics/btu679

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 29.

  LahotiPGarimellaKGionisA. Joint non-negative matrix factorization for learning ideological leaning on twitter. In: Proceedings of the Eleventh ACM International Conference on Web Search and Data Mining. Marina Del Rey, CA (2018). p. 351–9. 10.1145/3159652.3159669

  • CrossRef
  • Google Scholar

• 30.

  KimHChooJKimJReddyCKParkH. Simultaneous discovery of common and discriminative topics via joint nonnegative matrix factorization. In: Proceedings of the 21th ACM SIGKDD International Conference on Knowledge Discovery and Data Mining. Sydney, NSW (2015). p. 567–76. 10.1145/2783258.2783338

  • CrossRef
  • Google Scholar

• 31.

  JiangXHuXXuW. Microbiome data representation by joint nonnegative matrix factorization with Laplacian regularization. IEEE/ACM Trans Comput Biol Bioinform. (2015) 14:353–9. 10.1109/TCBB.2015.2440261

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 32.

  JainAK. Data clustering: 50 years beyond K-means. Pattern Recogn Lett. (2010) 31:651–66. 10.1016/j.patrec.2009.09.011

  • CrossRef
  • Google Scholar

• 33.

  BartolomeiFChauvelPWendlingF. Epileptogenicity of brain structures in human temporal lobe epilepsy: a quantified study from intracerebral EEG. Brain. (2008) 131(Pt 7):1818–30. 10.1093/brain/awn111

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 34.

  BassettDSSpornsO. Network neuroscience. Nat Neurosci. (2017) 20:353–64. Number:3. 10.1038/nn.4502

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 35.

  StamCJ. Modern network science of neurological disorders. Nat Rev Neurosci. (2014) 15:683–95. 10.1038/nrn3801

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 36.

  HeBAstolfiLValdes-SosaPAMarinazzoDPalvaSBenarCGet al. Electrophysiological brain connectivity: theory and implementation. IEEE Trans Biomed Eng. (2019) 66:2115–37. 10.1109/TBME.2019.2913928

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 37.

  DouwLvan DellenEGouwAAGriffaAde HaanWvan den HeuvelMet al. The road ahead in clinical network neuroscience. Netw Neurosci. (2019) 3:969–93. 10.1162/netn_a_00103

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 38.

  KamijoKTakedaYHillmanCH. The relation of physical activity to functional connectivity between brain regions. Clin Neurophysiol. (2011) 122:81–9. 10.1016/j.clinph.2010.06.007

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 39.

  KhambhatiANSizemoreAEBetzelRFBassettDS. Modeling and interpreting mesoscale network dynamics. Neuroimage. (2018) 180:337–49. 10.1016/j.neuroimage.2017.06.029

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 40.

  KerrMSBurnsSPGaleJGonzalez-MartinezJBulacioJSarmaSV. Multivariate analysis of SEEG signals during seizure. In: 2011 Annual International Conference of the IEEE Engineering in Medicine and Biology Society. Boston, MA: IEEE (2011). p. 8279–82. 10.1109/IEMBS.2011.6092041

  • CrossRef
  • Google Scholar

• 41.

  MahyariAGZoltowskiDMBernatEMAviyenteS. A tensor decomposition-based approach for detecting dynamic network states from EEG. IEEE Trans Biomed Eng. (2017) 64:225–37. 10.1109/TBME.2016.2553960

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 42.

  NewmanM. Networks. Oxford University Press (2018). 10.1093/oso/9780198805090.001.0001

  • CrossRef
  • Google Scholar

• 43.

  Le Van QuyenMSossJNavarroVRobertsonRChavezMBaulacMet al. Preictal state identification by synchronization changes in long-term intracranial EEG recordings. Clin Neurophysiol. (2005) 116:559–68. 10.1016/j.clinph.2004.10.014

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 44.

  KhambhatiANBassettDSOommenBSChenSHLucasTHDavisKAet al. Recurring functional interactions predict network architecture of interictal and ictal states in neocortical epilepsy. eNeuro. (2017) 4:ENEURO.0091-16.2017. 10.1523/ENEURO.0091-16.2017

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 45.

  ChaiLRKhambhatiANCiricRMooreTMGurRCGurREet al. Evolution of brain network dynamics in neurodevelopment. Netw Neurosci. (2017) 1:14–30. 10.1162/NETN_a_00001

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 46.

  CichockiA. Tensor decompositions: a new concept in brain data analysis?arXiv preprint arXiv:13050395. (2013).

  • Google Scholar

• 47.

  KoldaTBaderB. Tensor decompositions and applications. SIAM Rev. (2009) 51:455–500. 10.1137/07070111X

  • CrossRef
  • Google Scholar

• 48.

  ComonP. Tensors: a brief introduction. IEEE Signal Process Mag. (2014) 31:44–53. 10.1109/MSP.2014.2298533

  • CrossRef
  • Google Scholar

• 49.

  LeonardiNVan De VilleD. Identifying network correlates of brain states using tensor decompositions of whole-brain dynamic functional connectivity. In: 2013 International Workshop on Pattern Recognition in Neuroimaging (PRNI). IEEE (2013). p. 74–7. 10.1109/PRNI.2013.28

  • CrossRef
  • Google Scholar

• 50.

  OzdemirABernatEMAviyenteS. Recursive tensor subspace tracking for dynamic brain network analysis. IEEE Trans Signal Inform Process Netw. (2017) 3:669–82. 10.1109/TSIPN.2017.2668146

  • CrossRef
  • Google Scholar

• 51.

  WendlingFBartolomeiFBellangerJJBourienJChauvelP. Epileptic fast intracerebral EEG activity: evidence for spatial decorrelation at seizure onset. Brain. (2003) 126:1449–59. 10.1093/brain/awg144

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 52.

  JiruskaPde CurtisMJefferysJGRSchevonCASchiffSJSchindlerK. Synchronization and desynchronization in epilepsy: controversies and hypotheses. J Physiol. (2013) 591:787–97. Number:4. 10.1113/jphysiol.2012.239590

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 53.

  van MierloPPapadopoulouMCarretteEBoonPVandenbergheSVonckKet al. Functional brain connectivity from EEG in epilepsy: seizure prediction and epileptogenic focus localization. Prog Neurobiol. (2014) 121:19–35. 10.1016/j.pneurobio.2014.06.004

  • Pubmed Abstract
  • CrossRef
  • Google Scholar