Time-dependent degree-degree correlations in epileptic brain networks: from assortative to dissortative mixing

We investigate the long-term evolution of degree-degree correlations (assortativity) in functional brain networks from epilepsy patients. Functional networks are derived from continuous multi-day, multi-channel electroencephalographic data, which capture a wide range of physiological and pathophysiological activities. In contrast to previous studies which all reported functional brain networks to be assortative on average, even in case of various neurological and neurodegenerative disorders, we observe large fluctuations in time-resolved degree-degree correlations ranging from assortative to dissortative mixing. Moreover, in some patients these fluctuations exhibit some periodic temporal structure which can be attributed, to a large extent, to daily rhythms. Relevant aspects of the epileptic process, particularly possible pre-seizure alterations, contribute marginally to the observed long-term fluctuations. Our findings suggest that physiological and pathophysiological activity may modify functional brain networks in a different and process-specific way. We evaluate factors that possibly influence the long-term evolution of degree-degree correlations.

Improving our understanding of the emergence of epileptogenesis and ictogenesis from large-scale epileptic brain networks calls for approaches that take into account the interplay between the dynamic properties of network constituents (i.e., nodes and links) and the network topology. When investigating epileptic networks, nodes are usually assumed to represent distinct brain regions and links represent interactions between them, and these nodes and links constitute a functional brain network. Assessing the long-term dynamics of individual brain regions or nodes is mostly based on scalp or invasively recorded electroencephalograms, and links are derived from quantifying the statistical interdependence between signals (also referred to as functional connectivity) captured from electrodes overlying or within different brain regions (see Rubinov and Sporns, 2010;Lehnertz et al., 2014 for details). In the majority of network studies on epilepsy, global properties of epileptic networks during seizures have been characterized with measures such as clustering coefficient, average shortest path length, or synchronizability (see Lehnertz et al., 2014 for an overview). More recent studies explored the relevance of local network properties-such as the importance of individual nodes (Koschützki et al., 2005;Rubinov and Sporns, 2010)-for the dynamics of seizures (Kramer et al., 2008;Wilke et al., 2011;Varotto et al., 2012;Burns et al., 2014;Geier et al., 2015;Zubler et al., 2015). Findings achieved so far for these seizure networks are quite intriguing, given the similarity of their topological evolution across different types of epilepsies, seizures, medication, age, gender, and other clinical features, which might point to a common biophysical mechanism underlying ictogenesis. This similarity, however, is contrasted by strong intra-and interindividual fluctuations of local and global statistical network properties seen for the temporal evolution of epileptic brain networks over periods of days (Kuhnert et al., 2010;Kramer et al., 2011;Geier et al., 2013).
In Bialonski and Lehnertz (2013), functional brain networks before, during, and after one-hundred epileptic seizures with different anatomical onset locations were shown to exhibit assortative mixing patterns. Assortativity increased during seizures, reached a maximum prior to the end of seizures, and decreased already prior to seizure end. Interestingly, the observed concave-like temporal evolution of assortativity over the seizure period resembled closely the ones seen for the global correlation structure (Schindler et al., 2007), for the largest eigenvalue of adjacency matrices based on coherence and cross power (Müller et al., 2011), for clustering coefficient and average shortest path length (Schindler et al., 2008;Bialonski et al., 2011), and-if inverted-for synchronizability (Schindler et al., 2008). Taken together, these findings not only provide important clues on how the topology of functional brain networks changes during seizures but also on how seizures stop. In order to further improve our understanding of network mechanisms underlying seizure generation paralleling generation and maintenance of normal, physiological brain dynamics during the seizure-free interval, we here follow Kuhnert et al. (2010); Kramer et al. (2011);Geier et al. (2013) and investigate the long-term behavior of degree-degree correlations in evolving epileptic networks.

Patient Data
We analyzed long-term, multichannel invasive EEG (iEEG) recordings from seven patients who suffered from pharmacoresistant focal epilepsy with neocortical and/or hippocampal origin (cf. Table 1). After presurgical evaluation, all patients underwent resective surgery that led to complete seizure control. The patients had signed informed consent that their clinical data might be used and published for research purposes, and the study was approved by the ethics committee of the University of Bonn.
iEEG data were recorded continuously from chronically implanted intrahippocampal depth and/or subdural grid and strip electrodes (all manufactured by AD-TECH, WI, USA) using a Stellate Harmonie recording system (Stellate, Montreal, Canada; amplifiers constructed by Schwarzer GmbH, München, Germany). Decisions regarding electrode placement were purely clinically driven and were made independently of this study. The total number of electrode contacts ranged from N rs = 44 to  Kuhnert et al. (2010Kuhnert et al. ( , 2013, we constructed functional networks from iEEG data by associating network nodes with electrode contacts and inferred network links by estimating interdependencies between iEEG time series from pairs (n, m) (n, m ∈ {1, . . . , N rs }) of brain regions, regardless of their anatomical connectivity. For this purpose we used an established method for studying timevariant changes in phase synchronization, namely the mean phase coherence (Mormann et al., 2000): which is the temporal average of the differences of the instantaneous phases of iEEG time series from nodes n and m, and N denotes the number of data points. By definition, R nm is confined to the interval [0,1] where R nm = 1 indicates fully phase-synchronized systems. We used the analytic signal approach (Gabor, 1946;Panter, 1965) and derived the instantaneous phases of an iEEG time series using the Hilbert transform. An important property of this approach is that the instantaneous frequency relates to the predominant frequency in the Fourier spectrum (Boashash, 1992;Frei et al., 2010), which may be subject to fluctuations in the iEEG time series. In such a case, the instantaneous frequency varies rhythmically around the predominant frequency resulting in spurious estimates of the instantaneous phase, which can, however, be reduced by taking the temporal average (cf. Equation 1). From an electrophysiological point of view, we consider it more reasonable to look adaptively (via the Hilbert transform) at synchronization between predominant rhythms in the iEEG than to look at synchronization in some a priori fixed frequency bands for which there is no power in the time series (cf. Bruns, 2004;Osterhage et al., 2007;Frei et al., 2010).

Deriving a Temporal Sequence of Functional Networks
For further analyses, we split the offline bandpass-filtered (1-45 Hz) iEEG time series into consecutive non-overlapping windows of 20.48 s duration each (corresponding to N = 4096 data points) and estimated, for each window, the elements R mn of the phase synchronization matrix R. From this matrix, we constructed binary networks using a thresholding approach and set the non-diagonal elements of the adjacency matrix A to A mn = 1 if the corresponding entry R mn of R exceeded a threshold , and to A mn = 0 otherwise (A mm = 0 ∀m). For each matrix R, we chose such that the resulting network possessed a predefined link density (Anderson et al., 1999) ρ =k (N rs − 1) , with the number of nodes N rs and the mean degreē k = N −1 rs N rs n=1 k n (k n : = N rs m=1 A nm denotes the degree of node n). To do so, we performed a rank ordering of the entries of the lower triangular part of R (excluding the main diagonal) and chose as the (ρM rs − 1)-largest entry (with the number of all possible links M rs = (N 2 rs − N rs )/2). We set ρ = 0.1 to define links. With the aforementioned steps of analysis, we derived a temporal sequence of functional brain networks during inter-ictal, peri-ictal, and ictal periods spanning several days for each patient.

Assortativity
To assess degree-degree correlations in each functional network, we employed the assortativity coefficient a (Newman, 2002(Newman, , 2003. It quantifies whether links of the network tend to connect nodes of similar degrees with each other (in which case the network is called assortative; a > 0) or whether links preferentially connect high-degree nodes with low-degree nodes (dissortative network; a < 0). This tendency can be quantified by assessing the correlation between the degrees of nodes at both ends of links. To simplify the implementation of the assortativity coefficient, we can reformulate it in terms of the degrees of nodes (see Appendix B in Bialonski and Lehnertz, 2013; see also Lehnertz et al., 2014) which then reads where F u = N rs n = 1 k u n . By definition, a is confined to the interval [−1, 1]. Positive (negative) values of a indicate an assortative (dissortative) network, while a value of zero is indicative of a degree-degree uncorrelated network.
The finite size (i.e., the finite number of nodes) of networks can induce degree-degree correlations even for network ensembles that are degree-degree uncorrelated by construction (Barrat et al., 2008). In order to account for such finite-size effects and to distinguish between these and genuine degree-degree correlations, we considered ensembles of Erdős-Rényi random networks. Erdős-Rényi network ensembles are degree-degree uncorrelated and well investigated in the literature on random graphs (Erdős and Rényi, 1959, 1960, 1961Newman, 2002). For each patient, we created an ensemble of 20 Erdős-Rényi networks, having the same finite number of nodes N rs and the same link density ρ as the functional networks. From these random network ensembles, we determined the mean value as well as the standard deviation of the assortativity coefficient. We consider values of the assortativity coefficient of functional networks that deviate more than one standard deviation from the mean value obtained for Erdős-Rényi networks to indicate genuine degree-degree correlations which cannot be explained by finite-size effects. Results of previous studies pointed toward a coevolution of the assortativity coefficient (Bialonski and Lehnertz, 2013) and the clustering coefficient (Schindler et al., 2008) during seizures, which may indicate that assortative mixing comes along with or may be due to the transitivity of network links, a hypothesis which has been discussed for social networks (Newman and Park, 2003). In order to investigate whether such a coevolution can also be observed on longer time scales encompassing interictal and peri-ictal periods, we determined-in addition to the assortativity coefficient-the clustering coefficient for each functional network.

Clustering Coefficient
The clustering coefficient C (Albert and Barabási, 2002;Boccaletti et al., 2006) characterizes the transitivity in networks. We here define the local clustering coefficient C n of node n as taking care of nodes with no or only one connection (i.e., k n ∈ {0, 1}). By averaging C n over all nodes, we obtained the clustering coefficient C. C n and C are confined to the interval [0, 1] by definition.

Results
In Figure 1, we show time courses of the assortativity coefficient a derived from the temporal sequences of functional brain networks from two patients. We observe exclusively assortative mixing for patient 4, whereas for patient 1 repeated switches between phases of assortative and dissortative mixing can be observed (see Table 2 for the temporal means and their standard deviations of the assortativity coefficient for each patient). Of all patients, patient 1 was the only case for which we could observe long periods (up to several hours) of dissortative mixing. Shorter periods (up to several minutes) of dissortative mixing were observed in two other patients (patient 5 and 6). Since the number of recording sites N rs (i.e., the size of a network) varied across patients, we checked whether this might have led to the large fluctuations seen in time-resolved degree-degree correlations but could not observe a clearcut relationship (the Pearson correlation coefficient between N rs and mean assortativity coefficients amounted to 0.12). The evolutions of the assortativity coefficients shown in Figure 1 exhibit large fluctuations over time and appear to be partly periodic. The power spectral density estimates point to strong contributions from processes acting on timescales of some tens of hours but only small contributions from processes acting on timescales less than 30 min. For four patients, we observe a strong component at about 24 h with less pronounced (about a factor of 8 and more) contributions at the subharmonics at about 12 and 8 h. This may point toward an influence of daily rhythms, while components of even longer time scales might be related to alterations of antiepileptic medication during the presurgical evaluation. Not only do we observe large intra-and interindividual fluctuations for the assortativity coefficient a, but FIGURE 2 | Top: Exemplary frequency distributions of the assortativity coefficient derived from data recorded during day (solid) and night times (dashed) for two patients (left and middle) and for the pooled data from all patients (right). Bottom: Relative changes of assortativity (ā night −ā day )/ā day during day and night times for each patient.ā day andā night denote median values of the respective distributions.
Frontiers in Human Neuroscience | www.frontiersin.org even the extent of the intraindividual fluctuations is very different (see the standard deviations of a in Table 2). Interestingly, the time courses shown in Figure 1 indicate that less assortative (or even dissortative) mixing can be observed preferentially during night times. In order to investigate whether this observation extends beyond exemplary data, we split the data recorded during night times (ranging from 22 p.m. to 6 a.m.) and during day times (ranging from 6 a.m. to 22 p.m.). In the upper part of Figure 2, we show distributions of the assortativity coefficient for patients 1 and 4 as well as for the pooled data from all patients. Given the large interindividual variability, differences between distributions are diverse (see lower part of Figure 2), but a preferentially less assortative mixing during night times can be observed for the group of patients investigated here.
Having identified daily rhythms as a potential factor that may strongly influence the long-term evolution of degree-degree correlations in functional epileptic networks, we next evaluated whether relevant aspects of the epileptic process, particularly possible pre-seizure alterations, contribute to the observed longterm fluctuations. Following Mormann et al. (2003) (2015), we assumed that a pre-ictal phase of 4 h duration exists  and compared the distributions of values of the assortativity coefficient from the pre-ictal periods with those from inter-ictal periods. The latter distribution included all data that were recorded at least 4 h prior to and 30 min after a seizure.
In the upper part of Figure 3, we show the distribution of the assortativity coefficient for data from the pre-ictal and inter-ictal period from patient 1 as well as for the pooled data from all patients. Interestingly, the pre-ictal phase appears to be characterized by a slightly (about 10%) decreased assortative mixing, and in only one patient, we could observe a pre-ictal increase of degree-degree correlations (cf. lower part of Figure 3). Although these findings may help to further improve the understanding on how and which network reconfigurations promote seizure generation, more sophisticated analysis techniques (Andrzejak et al., 2009) applied to a larger dataset would be needed in order to statistically judge the observed pre-ictal changes.
Eventually, we studied whether temporal changes in assortative mixing are correlated with temporal changes in network transitivity (as quantified by the clustering coefficient). The temporal evolutions of the clustering coefficient for exemplary patients are shown in Figure 4. Similar to the case of the assortativity coefficient, the time courses of the clustering coefficient show large fluctuations and periodic structures. These periodicities act on timescales of about 24 h, with less pronounced contributions from subharmonics at about 12 and 8 h, which confirms previous results (Kuhnert et al., 2010). The intra-and interindividual fluctuations of the clustering coefficient (cf. Table 2) are much less pronounced than the fluctuation of the assortativity coefficient. We only observed weak correlations between assortativity and clustering coefficient (as quantified by the Pearson correlation coefficient ̺; FIGURE 3 | Top: Frequency distributions of the assortativity coefficient derived from data recorded during pre-ictal (solid) and inter-ictal periods (dashed) for one patients (left) and for the pooled data from all patients (right). Bottom: Relative changes of assortativity (ā pre −ā inter )/ā inter during pre-ictal and inter-ictal periods for each patient (patients 4 and 6 had no seizures during the recording period).ā inter andā pre denote median values of the respective distributions.
cf. Figure 5) in five patients (Pearson correlation coefficient ranged from −0.20 to 0.28), while in two patients, the correlation vanished (|̺| ≪ 0.1). These results indicate that for the functional brain networks investigated here, the assortativity coefficient provides information about the long-term evolution of the functional brain networks, which is complementary to the information provided by the clustering coefficient.

Discussion
We investigated the long-term variability of degree-degree correlations (assortativity) of functional brain networks constructed from invasive EEG recordings from seven patients suffering from focal epilepsies with neocortical and/or hippocampal origin. We observed large fluctuations in timeresolved degree-degree correlations, ranging from dissortative to assortative functional brain networks. For all patients, we observed the temporal evolution of the assortativity coefficient to possess periodic structures. Power spectral densities of the temporal evolutions were dominated by contributions at a time scale of some tens of hours, which may point toward a potential influence of daily rhythms, and/or to changes of the antiepileptic medication. Potential processes related to much smaller time  scales (<30 min) contributed much less to the overall power, a finding also observed in a previous study (Kuhnert et al., 2010) for the clustering coefficient and the average shortest path length. In the majority of patients (five of seven patients), functional brain networks tended to show higher degree-degree correlations during day times than during night times. However, we also observed the opposite behavior for two patients. As there was no sleep scoring available, we cannot relate the data to specific sleep stages. Still, the overall trend indicates that assortativity seems to be reduced during night times.
Only small differences in degree-degree correlations could be found between functional brain networks during inter-ictal periods and those during periods of an assumed pre-ictal state.
While we observed a pre-ictal decrease of the assortativity in four patients, we also observed a pre-ictal increase in assortativity for one patient (no seizures were recorded for two patients). It remains to be shown whether the observed tendency of functional brain networks to show slightly less assortative mixing pre-ictally can be regarded as predictive of the extreme event epileptic seizure.
In the light of results reported in previous studies (Ponten et al., 2007;Schindler et al., 2008;Kramer et al., 2010Kramer et al., , 2011Kuhnert et al., 2010;Bialonski et al., 2011;Bialonski and Lehnertz, 2013), our findings might point toward the following: We speculate that daily rhythms may be reflected in periodic reorganizations of functional brain networks.
FIGURE 6 | Sketch of how functional brain networks may explore the space of accessible network topologies (here parametrized by the clustering coefficient C, average shortest path length L, and assortativity coefficient a) in a process-dependent way: Daily rhythms may be reflected in periodic reorganizations of functional brain networks. Superimposed on that may be a reorganization in functional network connectivity reflecting pathophysiological activity.
Superimposed on that may be a reorganization in functional network connectivity reflecting pathophysiological activity. The space of accessible network topologies, however, may be explored in a different and process-dependent way (cf. Figure 6).
Changes due to daily rhythms seem to be reflected in functional brain networks, which have a larger transitivity and a larger average shortest path length (Kuhnert et al., 2010) but show less assortative mixing (as observed here) during night times than during day times. In contrast, as reported in previous studies, changes related to ictal activity seem to be associated with functional brain networks which are characterized by a larger transitivity, larger average shortest path length (Ponten et al., 2007;Schindler et al., 2008;Kramer et al., 2010;Bialonski et al., 2011), and more assortative mixing (Bialonski and Lehnertz, 2013) during seizures than before or after seizures. These latter changes suggest that functional brain networks might be more segregated during than before or after seizures (Kramer et al., 2010;Bialonski and Lehnertz, 2013): Nodes form different groups (Kramer et al., 2010) in which they share a similar node degree (leading to an increased assortativity, Bialonski and Lehnertz, 2013) and are tightly connected to each other (leading to an increased transitivity, Ponten et al., 2007;Schindler et al., 2008;Kramer et al., 2010;Bialonski et al., 2011). Between these groups, only sparse connections exist (leading to an increased average shortest path length, Ponten et al., 2007;Schindler et al., 2008;Kramer et al., 2010;Bialonski et al., 2011) which, in turn, may weaken the synchronizability of a network as a whole as suggested in previous studies (Schindler et al., 2008;Bialonski and Lehnertz, 2013). Changes in functional brain networks related to daily rhythms, however, appear to be expressed in a different way: During day times, nodes of the networks may be well connected (leading to a decreased average shortest path lengths, Kuhnert et al., 2010) via a few hub-like structures (i.e., nodes with high degrees; leading to an increased assortativity as observed here) which integrate information from all over the network. Individual nodes of lower degrees tend to be interconnected only sparsely (leading to a decreased clustering coefficient, Kuhnert et al., 2010). During night times, the hublike structures disappear or are less pronounced (leading to a decreased assortativity as observed here), giving way to a segregation of groups (associated with an increased average shortest path length, Kuhnert et al., 2010) in which nodes are mutually well interconnected (leading to an increased clustering coefficient, Kuhnert et al., 2010).
We note that these differences between changes of network topologies related to physiological and pathophysiological activities only become apparent when considering degree-degree correlations, in addition to the previously investigated network measures (cf. Figure 6). This underlines the usefulness of characterizing functional brain networks using different network measures, including assortativity. It remains to be shown whether the differences in network changes according to physiological and pathophysiological activity can also be observed in a larger sample of patients (which will allow for a better statistics).
We consider future research as promising which aims at a more detailed characterization of time evolving functional brain networks. For instance, if these networks indeed possess a more pronounced hub-like structure during day than during night times, we speculate that this will likely be reflected in respective changes in the distributions of betweenness centralities. Hubs that appear during day times and which integrate information from most of the other nodes may be identifiable by the largest betweenness centralities that occur in a network.
Another important aspect that is not yet fully explored is the influence of changes of the antiepileptic medication. Altered drug levels appear to affect functional brain networks only locally (Lehnertz and Elger, 1997;Haneef et al., 2015), but it is not clear if and how they might lead to local and/or global network reorganizations. Eventually, it is in general a challenging and non-trivial problem how to relate local to global network properties. Thus, knowing whether and how exactly the reported periodic changes and reorganizations of functional brain networks can be related to features of the node and/or link dynamics may help to gain deeper insights into the complicated dynamics underlying epileptic networks.

Funding
SB and KL acknowledge support by the Volkswagen Foundation (Grants No. 85390 and No. 85392).