High frequency oscillations and high frequency functional network characteristics in the intraoperative electrocorticogram in epilepsy.

OBJECTIVE: High frequency oscillations (HFOs; > 80 Hz), especially fast ripples (FRs, 250-500 Hz), are novel biomarkers for epileptogenic tissue. The pathophysiology suggests enhanced functional connectivity within FR generating tissue. Our aim was to determine the relation between brain areas showing FRs and 'baseline' functional connectivity within EEG networks, especially in the high frequency bands.
METHODS: We marked FRs, ripples (80-250 Hz) and spikes in the electrocorticogram of 14 patients with refractory temporal lobe epilepsy. We assessed 'baseline' functional connectivity in epochs free of epileptiform events within these recordings, using the phase lag index. We computed the Eigenvector Centrality (EC) per channel in the FR and gamma band network. We compared EC between channels that did or did not show events at other moments in time.
RESULTS: FR-band EC was higher in channels with than without spikes. Gamma-band EC was lower in channels with ripples and FRs.
CONCLUSIONS: We confirmed previous findings of functional isolation in the gamma-band and found a first proof of functional integration in the FR-band network of channels covering presumed epileptogenic tissue. SIGNIFICANCE: 'Baseline' high-frequency network parameters might help intra-operative recognition of epileptogenic tissue without the need for waiting for events. These findings can increase our understanding of the 'architecture' of epileptogenic networks and help unravel the pathophysiology of HFOs.

Introduction

Precise identification of the epileptogenic tissue can improve the success rate of epilepsy surgery. High frequency oscillations (HFOs) are promising electrographic biomarkers that identify the epileptogenic cortex more accurately than the currently used interictal spikes (Bragin et al., 1999, Worrell et al., 2004, Jirsch et al., 2006, Urrestarazu et al., 2006, Schindler et al., 2007, Jacobs et al., 2008, Zijlmans et al., 2012b, van't Klooster et al., 2015). HFOs are brief, sinusoid waveforms in the EEG above 80 Hz, and are typically subdivided into ripples (80–250 Hz) and fast ripples (FRs, 250–500 Hz). Retrospective studies in surgical candidates have shown that complete removal of tissue generating HFOs, especially FRs, correlates better with post-operative seizure freedom than removal of tissue generating interictal spikes (Wu et al., 2010, van Klink et al., 2014, van't Klooster et al., 2015) or removal of the seizure onset zone (SOZ) (Ochi et al., 2007, Jacobs et al., 2008, Jacobs et al., 2010). Although HFOs hold great potential as a biomarker for epileptogenic tissue, their clinical implementation is hampered because the current visual analysis requires expertise and time because of their sparse occurrence. In addition, HFOs have also been recorded in healthy neuronal tissue during cognitive processing (Baker et al., 2003, Blanco et al., 2011, Matsumoto et al., 2013, Kucewicz et al., 2014). Distinguishing these physiological from pathological HFOs is crucial for clinical practice, but challenging.

The neuronal mechanisms underlying physiological and pathological HFOs are thought to involve different degrees of interneuron-mediated inhibition. Physiological HFOs are believed to reflect synchronous interneuron-mediated inhibitory postsynaptic potentials (IPSPs) on healthy principal cells, which coordinate the sparse firing of the latter (Ylinen et al., 1995, Le Van Quyen et al., 2006, Jefferys et al., 2012, Schomburg et al., 2012, Alvarado-rojas et al., 2014, Menendez et al., 2015). In contrast, pathological ripples are thought to reflect brief bursts of population spikes arising from groups of abnormally synchronized principal cells, due to reduced interneuron-mediated IPSPs (Bragin et al., 2011, Alvarado-rojas et al., 2014, Menendez et al., 2015). Pathological FRs are proposed to result from in concert, but slightly out-of-phase firing of independent clusters of these pathologically interconnected neurons (Foffani et al., 2007, Ibarz et al., 2010, Jefferys et al., 2012, Fink et al., 2015, Menendez et al., 2015).

The presumed pathophysiology of epileptic HFOs suggests a disruption of normal, micro-scale, (cortical) neuronal wiring resulting in epileptogenic brain tissue. Complex network theory, based on graph analysis, is a promising method to study connections in the brain (Bullmore and Sporns, 2009, Bullmore and Sporns, 2012, Rubinov and Sporns, 2010, Stam, 2014). Functional network analysis decomposes the brain into a collection of nodes, e.g. EEG channels, and edges, e.g. synchronization between the signals measured in different channels. The organization of the resulting meso/macro-scale functional brain networks can be characterized using network measures (Bullmore and Sporns, 2009, Bullmore and Sporns, 2012, Rubinov and Sporns, 2010, Stam, 2014). Global measures describe the overall network structure, whereas local measures aim to identify regions that play a central role in the network because they are well-connected and/or important for communication between brain areas, i.e. hubs. Intracranial EEG studies investigating local, meso-scale, brain network topology in ictal recordings of patients with focal epilepsy have shown an increased synchronization and hub status in the epileptogenic zone (Ponten, 2007, Kramer et al., 2008, Schindler et al., 2008, Wilke et al., 2011, Varotto et al., 2012, Mierlo et al., 2013). Interictal network studies have reported higher centrality measures in the epileptogenic zone (Wilke et al., 2011, Varotto et al., 2012) and lower centrality measures in the SOZ (Warren et al., 2010, Van Diessen et al., 2013b). Nissen et al., 2016 found higher average centrality measures in virtual MEG electrodes in the affected compared to the non-affected hemisphere, but a lower node centrality the closer to the center of brain areas generating epileptiform events. The first studies comparing the presence and localization of HFOs with meso-scale hub regions reported that channels displaying HFOs are functionally isolated in the gamma frequency band, while results in the frequency bands below 30 Hz were contradictive (Ibrahim et al., 2013a, Van Diessen et al., 2013b, Nissen et al., 2016). Based on the presumed enhanced functional connectivity needed to generate FRs, we would expect local functional integration of channels covering these areas, at least in the FR-band. We therefore aimed to determine the relation between locations showing HFOs and hubs in the ‘baseline’ high frequency functional networks constructed of event-free subdural intraoperative electrocorticography (ioECoG) of patients with refractory temporal lobe epilepsy (TLE).

Methods

Patients

A database, consisting of retrospectively collected data from 54 patients with focal refractory epilepsy, who underwent resective epilepsy surgery with ioECoG recorded at 2048 Hz between 2008 and 2012 at the UMC Utrecht, was used. The ioECoG was recorded for tailoring the epileptogenic area based on interictal spikes. HFO information was not available during surgical decision-making. The database was collected conform the guidelines of the institutional ethical committee of the UMC Utrecht, the Netherlands.

Important constraints to allow comparison of network characteristics are (i) a homogeneous patient population (van Diessen et al., 2014b) and (ii) equal network size and density, i.e. the number of nodes and edges (van Wijk et al., 2010, Boersma et al., 2012). We therefore included (1) patients with epilepsy from the same lobe of origin (TLE) who had (2) pre-resection ioECoG recordings with both anterior and posterior temporal sampling of (3) the same 23 electrode positions (Fig. 1). This last number was determined after a first screening of the data and rejecting electrodes that often showed continuous artifacts (the last two strip, and one of the corner electrodes).

Fig. 1

Illustration of the recording situations and included channels.

The pre-resection ECoG recording consisted of two situations: (i) Situation A, where the grid was placed as far towards the anterior temporal pole as possible, covering the anterolateral temporal cortex, and (ii) Situation B, where the grid was flipped 180 degrees to record activity at the posterolateral temporal cortex. In both situations, the 1 × 6 electrode strip was placed subtemporally at 4–6 cm from the temporal pole to allow sampling of mesiotemporal structures. The last two strip electrodes and one of the corner electrodes (red) were excluded from the analysis, because they often showed continuous artifacts and network size is important for comparison of network characteristics. The analyses were performed using bipolar montages rendering 18 bipolar channel pairs (blue) per situation. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

Intra-operative ECoG recordings

IoECoG was recorded using a 4 × 5 electrode grid and 1 × 6 electrode strip (Ad-Tech, Racine, WI, USA) with a 64-channel EEG system (MicroMed, Veneto, Italy). Sampling frequency was 2048 Hz with an anti-aliasing filter at 538 Hz. Platinum electrodes embedded in silicone with 4.2 mm2 contact surface and 1 cm inter-electrode distances were used. The recordings were referenced to an extracranial reference electrode.

The pre-resection recording consisted of two situations: (i) Situation A, where the grid was placed as far towards the anterior temporal pole as possible, covering the anterolateral temporal cortex, and (ii) Situation B, where the grid was flipped 180 degrees to record activity at the posterolateral temporal cortex (Fig. 1). In both situations, the 1 × 6 electrode strip was placed subtemporally at 4–6 cm from the temporal pole to allow sampling of mesiotemporal structures, which are often involved in seizure generation in TLE.

Propofol anesthesia was temporarily interrupted during ioECoG recording to obtain a continuous background pattern, and enable detection and interpretation of interictal spikes for surgical decision-making. Patients did not wake up.

Identification of spikes and HFOs

A one-minute epoch was selected from each ioECoG recording situation (Fig. 2.1). The epochs were selected near the end of the recording to minimize propofol effects (Zijlmans et al., 2012a). Data were analyzed in a bipolar montage to diminish background noise. The bipolar montage of pairs of neighboring electrodes along the length of the electrode grid or strip reduced the number of channel pairs for analysis to 18 per situation.

Fig. 2

Schematic overview of the methodological steps taken.

1. Selecting one minute of ioECoG near the end of the recording to minimize propofol effects.

2. Visual marking of spikes, ripples and fast ripples with corresponding screen settings.

3. Selecting four epileptic event- and artifact free epochs of 2 s for functional network analysis.

4. Computing the PLI between each pair of bipolar electrode channels in the FR-, ripple-, gamma- and theta frequency band

5. Constructing graphs and calculating the strength and eigenvector centrality of each node for different frequency bands and comparing the spatial relation with epileptic event rates.

Epileptiform events (i.e. interictal spikes, ripples and FRs) were visually marked in Stellate Harmonie Reviewer [MvtK, NvK or WZ] (Fig. 2.2). Spikes were marked using conventional EEG filter settings of 0.5–70 Hz, an amplification of 75–100 μV/mm and a time scale of 10 s/page. A spike was defined as a sharp transient of maximum 80 ms duration, standing out at least twice above the baseline amplitude (Alarcon, 1996). Sharp waves co-occurring with spikes on adjacent bipolar electrode channels were also marked as spikes. For HFO marking, the screen was split vertically to simultaneously visualize ripples and FRs. Different filter and amplitude settings were used for ripples (80 Hz high-pass filter, 5 μV/mm amplitude) and FRs (250 Hz high-pass filter, 1 μV/mm amplitude). Data were displayed at maximum temporal resolution of 0.4 s/page. Events consisting of at least four consecutive oscillations, clearly standing out from the baseline of the channel, were marked as ripples or FRs (Jacobs et al., 2008, Zelmann et al., 2009). Marking of spikes was done blinded for HFOs and vice versa. An experienced clinical neurophysiologist [CF or FL] checked the marked spikes; an experienced HFO researcher [MZ] checked the marked ripples and FRs.

Functional connectivity and network analysis

Within each one-minute of ioECoG in which spikes and HFOs were marked, we selected four epochs of 2 s (4096 samples) for functional network analysis. We selected epochs without spikes, ripples and FRs, propofol-induced burst-suppression patterns or artifacts (Fig. 2.3). We started epoch selection at the end of the one-minute recording and used the first four available epochs to reduce propofol effect and avoid selection bias. It was not possible to select longer event-free epochs. In 13 of the fourteen patients we could select five event-free epochs. We decided to use the first four, because a leave-one-out-analysis did not show significant differences in mean PLI per channel over any combination of four epochs. Epochs were selected by one [WZ] and checked by a second observer [MZ]. Selected epochs were exported as EDF and converted to ASCII using MATLAB (The MathWorks, Natick, MA, U.S.A). The ASCII files were used for further analysis in BrainWave version 0.9.133 (http://home.kpn.nl/stam7883/brainwave.html).

‘Baseline’ functional connectivity was computed for each event-free epoch using the Phase Lag Index (PLI). The PLI determines connectivity by calculating the consistency of the nonzero phase lag between two time-series. By disregarding a phase difference of zero modulus π, it eliminates volume conduction as a confounding factor (Stam et al., 2007). In addition, it has the advantage that it is not sensitive to differences in amplitude of the signal (Muthukumaraswamy and Singh, 2011). The PLI ranges between 0, no phase coupling or coupling with a phase difference centered on 0 modulus π, and 1, maximal phase coupling. In formula: PLI = |⟨sign[∆ φ(t]⟩| (Stam et al., 2007). The mean PLI per channel, i.e. the average of the connectivity of that channel with all other channels, over four epochs was used for further analysis. The PLI was calculated for the same bipolar montage as used for epileptiform event marking, resulting in an 18-by-18 adjacency matrix for each situation (Fig. 2.4). We considered the antero- (Situation A) and posterolateral (Situation B) temporal cortex recordings as separate entities, because they were not measured simultaneously. We calculated the PLI for the FR (250–500 Hz), ripple (80–250 Hz), gamma (25–49 Hz) and theta (4–8 Hz) frequency bands. The choice for the gamma and theta frequency bands was based on previous literature (Van Diessen et al., 2013a, Chiang and Haneef, 2014).

PLI matrices were converted into graphs to characterize local network topology (Fig. 2.5). A graph is a simple topographical representation of a network that consists of nodes and edges. In our case, the bipolar channels represent the nodes and the PLI between them the edges. We constructed undirected weighted networks, which are described by the graph G = (N,W), where N is the set of 18 bipolar channels per recording situation and W = {wij} is the NxN symmetric weight matrix with wii = 0 and wij the value of the PLI between node i and j (Rubinov and Sporns, 2010).

We focused on the hub-measures of each node in the network to investigate how brain areas showing epileptiform events are related to alterations in local, meso-scale, network topology. Several hub-measures exist to characterize node centrality. Methodological choices influence the results of network analysis, and make studies difficult to compare (Wang et al., 2014, Bastos and Schoffelen, 2015, Kida et al., 2016). We decided to characterize the local topology of the constructed networks by calculating the Strength (weighted degree) and Eigenvector Centrality (EC) of each network node. We chose these measures because previous work by our group used these measures and similar methodology in depth recordings (Van Diessen et al., 2013b), so it would be possible to compare results. We restricted our analysis to these measures to limit multiple testing.

The Strength is the sum of the weights of the edges connected to a node (Rubinov and Sporns, 2010). A high Strength means that the node has strong and/or many connections to other nodes in the network. The EC is a more advanced measure. In addition to the number and weights of connections to a node, it takes the connectedness of those other nodes into account. The EC is based on the eigenvector ecw corresponding to the greatest eigenvalue λ of the weight matrix W (Newman, 2007). The eigenvectors of the weighted matrix indicate the nodes that show many connections with most other nodes in the network. The eigenvector corresponding to the greatest eigenvalue has only positive entries and assigns relative scores to all nodes in the network based on the concept that connections to high-scoring nodes contribute more to the score of the node in question than equal connections to low-scoring nodes. In formula: if λ is the greatest eigenvalue and ecw the corresponding eigenvector, then W ∗ ecw = λ ∗ ecw or similar is the proportionality factor so that eciw is proportional to the sum of the connectivity scores of all nodes connected to it.

Resection and postsurgical outcome

The position of electrodes on the cortical surface in relation to the resected tissue was determined from photographs taken during surgery. Bipolar channels were classified into resected or non-resected (van Klink et al., 2014). Postsurgical outcome according to the Engel classification was determined from the most recent follow-up (Engel, 1993).

Statistical analyses

We calculated the percentage of spikes, ripples and FRs per channel (Ch) per situation (Sit) per patient to normalize the number of observed events:

Mean Strength and EC per channel per patient differed significantly across the analyzed frequency bands (three-factor nested (n)ANOVAs, Appendix Table A.1). We therefore examined the relation between spikes, ripples or FRs (the independent variable) and Strength or EC (the dependent variable) separately for each frequency band. The cumulative power in the FR and gamma frequency band of the channel with the highest and lowest FR and gamma band Strength of each patient did not differ significantly (paired t-tests, tFR = 0.69, pFR = 0.49; tgamma = 0.92, pgamma = 0.37).

Table A.1

Three-factor nested ANOVAs to assess differences in mean strength and eigenvector centrality per channel per patient across frequency bands.

Source	Strength					EC
	SS	d.f.	MS	F	p-Value	SS	d.f.	MS	F	p-Value
FB	11.01	3	3.67104	11,238.92	< 0.0001	0.01	3	0.00482	3.68	0.01
Pt	0.05	13	0.00421	12.66	< 0.0001	0.01	13	0.00077	0.43	0.96
Ch (Pt)	0.16	490	0.00033	1.02	0.40	0.87	490	0.00177	1.35	< 0.0001
Error	0.49	1509	0.00033			1.98	1509	0.00131
Total	11.72	2015				2.87	2015

Abbreviations: FB = frequency band. Pt = patient. Ch = channel. EC = eigenvector centrality. SS = sum of squares, d.f. = degrees of freedom. MS = mean square. F = f-statistic.

Strength and EC differed significantly across frequency bands.

We assessed differences in mean percentages of events, Strength and EC between resected and non-resected channels, and differences in mean Strength and EC between channels that did or did not show events using paired t-tests.

We explored the relation between the number of events, and Strength and EC using linear mixed models. We considered our data as hierarchically structured with electrode channels, and thus data points, nested within individual subjects. We therefore included patient as a random intercept in our starting model. We performed log-likelihood tests to determine whether expansion of this model with a random coefficient for patients and a random factor for recording situation, i.e. the antero- or posterolateral temporal recording, yielded a better fit of the PLI data in the lowest (theta) and highest (FR) frequency band analyzed (Appendix Table A.2). The model with the best fit was used in all subsequent analyses of spikes, ripples or FRs and Strength or EC per frequency band. We excluded patients without spikes, ripples or FRs from the analyses for the corresponding event type. All statistical analyses were performed using MATLAB (The MathWorks, Natick, MA, U.S.A). A p-value below 0.05 was considered significant.

Table A.2

Log-likelihood tests to compare goodness of fit of models in the fast ripple and theta frequency bands.

Left column: formula of the model. ‘1 + %Event’ is the fixed factor (intersect and slope). We started our model with this fixed factor and a random intersect for each patient (1 | Pt). We determined whether expansion of the model with other random factors led to a better fit of the model using log-likelihood tests. FR = fast ripple. R = ripple. Sp = spikes. Pt = patient. DF: degrees of freedom, i.e. number of factors. AIC: Akaike's Information Criterion. BIC: Bayesian Information Criterion. The lower the AIC and BIC, the better the fit of the model. LLV: LogLikelihood-value: the higher this value, the better the fit. LR-stat: Likelihood ratio statistic: expresses how many times more likely the data are under this model than the previous one. ΔDF: change in degrees of freedom relative to previous, more simple model. LR-stat and ΔDF are used to calculate the p-value. In italics: expansion of the previous model with a random intersect for %-Events | Pt did not lead to a better fit. In bold: model with the best fit based on AIC, BIC and LR-stat.

Results

Fourteen out of the 54 patients with refractory epilepsy who underwent resective epilepsy surgery with ioECoG sampled at 2048 Hz between 2008 and 2012 in the UMC Utrecht met the inclusion criteria (8 females; 5 left sided TLE; Table 1). We excluded 24 patients with extra-temporal lobe epilepsy, 1 TLE patient who had only one recording situation, and 15 TLE patients with < 23 electrodes, or a different electrode configuration. Mean age at surgery of the 14 included patients was 24 (range 3–61) years. Mean duration of epilepsy at the time of surgery was 8 (range 2–20) years. Thirteen patients were seizure free for 3 or more years after surgery, ten without anti-epileptic drugs. One patient had recurrent seizures.

Table 1

Patient characteristics.

Patient	Gender	Age surgery (yrs)	Disease duration (yrs)	Side	Pathology	# Spikes				# Ripples				# FRs				Outcome (Engel)
						SitA		SitB		SitA		SitB		SitA		SitB
						Gr	Rst	Gr	Rst	Gr	Rst	Gr	Rst	Gr	Rst	Gr	Rst
1	Female	46	10	Right	Iatrogenic cortical lesion	1	7	26	33	3	2	8	4	0	0	0	3	1A
2	Female	11	3	Left	Ganglioglioma	3	0	28	0	0	0	0	0	0	0	0	0	1A
3	Male	61	12	Left	No diagnosis	0	12	1	0	0	1	0	0	0	2	0	0	1A
4	Female	40	8	Right	Cavernoma & MTS	47	61	23	73	17	6	2	0	1	0	1	0	1A
5	Female	24	2	Right	DNET	23	20	0	8	0	1	0	0	0	0	0	0	4B
6	Female	37	20	Right	Cavernoma	0	0	0	0	0	0	0	0	0	0	0	0	1A
7	Male	14	11	Right	MTS	81	61	22	146	10	23	12	40	3	3	0	10	1A
8	Male	13	12	Right	MTS	383	48	73	31	130	0	10	0	11	0	4	0	1A
9	Female	20	5	Right	PXA	0	0	0	0	0	0	1	0	0	0	0	0	1A
10	Male	35	12	Left	Glioma	0	20	0	17	0	5	4	5	0	1	0	2	1A
11	Female	3	2	Left	DIG	0	0	0	0	28	0	18	1	0	0	0	0	1A
12	Male	6	5	Right	DIG	417	46	0	63	50	8	15	2	4	0	0	0	1A
13	Female	12	10	Right	No diagnosis	40	74	14	108	15	23	5	40	0	10	2	20	1A
14	Male	12	2	Left	DIG	7	6	0	4	2	1	1	1	0	0	0	0	1A

Abbreviations: # = absolute number of events, FR = fast ripples, SitA = recording situation A, SitB = recording situation B, Gr = grid, Rst = subtemporal reed, MTS = mesiotemporal sclerosis, DNET = dysembrioplastic neuroepithelial tumor, PXA = pleomorphic xanthoastrocytoma, DIG = desmoplastic infantile ganglioglioma.

Epileptiform events

We found spikes in eleven (mean number per patient: 184, range: 13–535), ripples in twelve (mean number: 41, range: 1–140) and FRs in eight patients (mean number: 10, range 2–32). One patient (Pt 6) showed no epileptiform events. Numbers of events per patient per situation are given in Table 1.

Mean percentages of spikes and FRs per channel were higher in resected than in non-resected channels (Res%Spikes: 7.6%, nonRes%Spikes: 2.0%, p = 0.01; Res%FR: 5.1%, nonRes%FR: 0.7%, p = 0.02). There was no significant difference in the mean percentage of ripples per channel between resected and non-resected channels (Res%Ripples: 5.4%, nonRes%Ripples: 3.3%, p = 0.28).

Network measures

Table 2 shows the mean Strength and EC of channels that did or did not show events and of channels that were or were not covering the resected area, for the FR- and gamma frequency band. We found a higher FR-band EC in channels with, compared to those without spikes (p = 0.04). Gamma-band Strength and EC were lower in channels with compared to without ripples and FRs (p < 0.05 for both). Sub-analysis in patients with FR showed that channels with both spikes and FRs tended to have higher respectively lower centrality values in the FR and gamma frequency band than channels with spikes only (Strength: FR: 0.043 vs 0.042, p = 0.58 and gamma: 0.118 vs 0.125 p = 0.02; EC: FR: 0.232 vs 0.229, p = 0.48 and gamma: 0.223 vs 0.233, p = 0.03).

Table 2

Mean strength and eigenvector centrality of channels with/without events or that were/were not resected in the FR and gamma frequency band.

Band	Event type (N)	Strength				EC
		With events/resected		No events/not resected	p-Value	With events/resected		No events/not resected	p-Value
FR-band250–500 Hz	FRs (8)	0.043	>	0.042	0.86	0.232	>	0.227	0.78
	Ripples (12)	0.050	>	0.044	0.18	0.230	>	0.223	0.18
	Spikes (11)	0.043	>	0.041	0.07	0.231	>	0.224	0.04
	Resected (14)	0.047	>	0.042	0.04	0.236	>	0.217	0.04
γ-band25–49 Hz	FRs (8)	0.118	<	0.124	0.02	0.223	<	0.232	0.04
	Ripples (12)	0.124	<	0.127	0.03	0.225	<	0.233	0.007
	Spikes (11)	0.123	<	0.124	0.16	0.232	<	0.232	0.10
	Resected (14)	0.124	<	0.125	0.50	0.231	<	0.231	0.86

Abbreviations: EC = eigenvector centrality. FR-band = fast ripple band = 250–500 Hz. γ-band = gamma-band = 25–49 Hz. (N) = total number of patients in that subgroup. >: value higher in channels with events/that were resected. <: value lower in channels with events/that were resected. Bold: significant result, p < 0.05.

There is a (trend towards) higher FR-band Strength and EC in channels with spikes and channels that were resected. There is a lower gamma-band strength and EC in channels with events and channels that were resected.

There were no significant differences in theta- and ripple-band Strength and EC between channels with and without events and between channels that were or were not covering the resected area (Appendix Table A.3).

Table A.3

Mean strength and eigenvector centrality of channels with/without events or channels that were or were not resected in the ripple and theta frequency band.

Band	Event type (N)	Strength				EC
		With events/resected		No events/not resected	p-Value	With events/resected		No events/not resected	p-Value
R-band80–250 Hz	FR(8)	0.061	<	0.063	0.47	0.225	<	0.231	0.17
	R(12)	0.067	>	0.063	0.32	0.233	>	0.230	0.25
	Sp (11)	0.064	>	0.061	0.08	0.231	>	0.230	0.41
	Res(14)	0.064	>	0.062	0.09	0.232	>	0.230	0.45
θ-band4–8 Hz	FR(8)	0.235	>	0.235	0.88	0.232	=	0.232	0.89
	R(12)	0.233	<	0.233	0.52	0.232	=	0.232	0.30
	Sp (11)	0.234	>	0.234	0.27	0.232	=	0.232	0.46
	Res (14)	0.235	>	0.233	0.36	0.233	>	0.231	0.15

Abbreviations: EC = eigenvector centrality. R-band = 80–250 Hz. θ -band = theta-band = 4–8 Hz. FR = fast ripples. R = ripples. Sp = spikes. Res = resected. (N) = total number of patients in that subgroup. >: value higher in channels with events or channels that were resected. <: value lower in channels with events or channels that were resected.

Relation between the number of epileptiform events and network measures

The model that best fitted the data was described by the formula (Appendix Table A.2):where Network Measure is the dependent variable, either the Strength or the EC. The fixed intercept α represents the average dependent variable start value. The fixed coefficient β1 represents the average relation between the dependent variable and the %Events (spikes, ripples or FRs). aPatient, bPatient%Events and bPatientSituation represent the random variation in intercept and slope for each patient per recording situation. εPatient represents the residuals.

Using this model, we found a positive association between the percentage of spikes per channel and the FR-band Strength and EC. We found negative associations between the percentages of all event-types per channel and gamma-band Strength and EC (Fig. 3, Appendix Table A.4). There were no significant associations between percentages of events and the Strength and EC in the theta- and ripple frequency band (Appendix Table A.5).

Fig. 3

Relation between the number of spikes (left), ripples (middle) and FRs (right) and the eigenvector centrality (EC) in the FR (upper) and gamma (lower) frequency band.

Colored lines are the patient specific regression lines. The bold black line represents the fixed effect of the model. The fixed coefficient β represents the average relation between the dependent variable and the %Events (spikes, ripples or FRs). We show a negative association between the number of all types of epileptic events and node centrality in the gamma band and the trend towards the opposite, a positive association, between the number of epileptic events and node centrality in the fast ripple frequency band. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

Table A.4

Associations between the number of epileptiform events and centrality network measures per channel in the fast ripple and gamma frequency band.

Band	Type of event	Strength		EC
		β (CI95%)	p-Value	β (CI95%)	p-Value
FR-band250–500 Hz	FR	0.001 (− 0.009 – 0.011)	0.89	0.003 (− 0.049 – 0.056)	0.90
	R	0.006 (− 0.006 – 0.019)	0.34	0.024 (− 0.029 – 0.076)	0.38
	Sp	0.029 (0.001 – 0.058)	0.05	0.111 (0.007 – 0.214)	0.04
γ-band25–49 Hz	FR	− 0.014 (− 0.022 – ‐0.005)	0.002	− 0.024 (− 0.039 – ‐0.009)	0.002
	R	− 0.016 (− 0.025 – ‐0.007)	0.001	− 0.027 (− 0.045 – ‐0.009)	0.002
	Sp	− 0.014 (− 0.023 – 0.004)	0.004	− 0.026 (− 0.043 – ‐0.009)	0.002

Abbreviations: FR-band = fast ripple band = 250–500 Hz, γ-band = gamma band = 25–49 Hz. FR = fast ripples. R = ripples. Sp = spikes. EC = eigenvector centrality. β = regression coefficient of the fixed-effects of the model. Bold: significant result, p < 0.05.

There is a positive association between the %-Spikes and the FR-band strength and EC. There is a negative association between the % of all event types and gamma-band strength and EC.

Table A.5

Associations between the number of epileptiform events and centrality network measures per channel in the ripple and theta frequency band.

Band	Event type	Strength		EC
		β (CI95%)	p-Value	β (CI95%)	p-Value
R-band80–250 Hz	FR	− 0.003 (− 0.008 – 0.003)	0.38	− 0.009 (− 0.029 – 0.009)	0.33
	R	0.002 (− 0.004 – 0.009)	0.52	0.008 (− 0.014 – 0.030)	0.45
	Sp	0.004 (− 0.003 – 0.011)	0.21	0.014 (− 0.008 – 0.036)	0.21
θ-band4–8 Hz	FR	0.002 (− 0.013 – 0.017)	0.78	0.001 (− 0.013 – 0.014)	0.91
	R	0.003 (− 0.013 – 0.019)	0.73	0.001 (− 0.014 – 0.016)	0.88
	S	0.004 (− 0.011 – 0.019)	0.64	0.005 (− 0.009 – 0.019)	0.46

Abbreviations: R-band = ripple band = 80–250 Hz. θ-band = theta-band = 4–8 Hz. FR = fast ripples. R = ripples. Sp = spikes. EC = eigenvector centrality. β = regression coefficient of the fixed-effects of the model.

We projected both variables on the schematic electrode representation to visualize the spatial relation between epileptiform events and centrality network measures found using linear mixed models. Fig. 4 shows three patient examples.

Fig. 4

Schematic overview of the relation between epileptiform events and node centrality in the FR and gamma frequency band.

The nodes represent the bipolar electrodes channels. Lines represent the functional connections between them. Note that we thresholded the networks and only show the connections that were higher than the mean strength of that patient in that recording situation. The color of the node indicates the presence/absence of events. Channels that showed FRs also showed ripples and spikes. Channels that showed ripples also showed spikes. The diameter of the node gives an indication of node centrality based on the eigenvector centrality (EC). What can be seen is A) that channels showing FRs in patient 13 are strongly connected in the FR-band functional network, and relatively isolated (still many, but less strong connections than the rest of the grid) in the γ-band functional network. B) The ioECoG of patient 6 did not show epileptic events, however the FR-band functional network clearly shows a hub node. This patient had a mesiotemporal cavernoma resected, and was seizure free afterwards. C) The epileptogenic events and the functional networks were not localized in patient 5. This patient underwent an anterior temporal lobectomy with amygdalohippocampectomy but did not become seizure free. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

Discussion & conclusion

We studied the relation between brain areas showing epileptiform events and hubs in the meso-scale high frequency functional networks in patients with TLE. Our findings confirm the previously reported functional isolation in the gamma-band network of channels showing epileptiform events (Warren et al., 2010, Ibrahim et al., 2013a, Van Diessen et al., 2013b, Nissen et al., 2016). We found a suggestion of local enhanced connectivity of these channels in the FR-band functional network. This could fit within current thoughts on micro-scale pathophysiology that FRs result from synchronously but slightly out of phase firing clusters of neurons that became hyper-excitable due to pathological reduced interneuron-mediated inhibition (Foffani et al., 2007, Ibarz et al., 2010, Bragin et al., 2011, Fink et al., 2015, Menendez et al., 2015). In that line of thought, the hubs in the meso-scale FR network might cover these pathological microscale (sub)network that generate FRs. Importantly, our results have been extracted from epochs that did not show any spikes, ripples or FRs, which implies that there may be high frequency information hidden in the ‘baseline’ signal that can be used to point out channels that will show FR at other moments in time (Bettus et al., 2011).

We found the FR-band functional integration primarily in channels showing spikes and only non-significantly in channels with FRs. This may be due to sample size, i.e. only eight patients showed FRs in the one-minute selection of ioECoG used for epileptiform event marking, and/or scale-issues, i.e. the low number of channels within those patients showing FRs in general implies that the spatial sampling of clinical grids is not sufficient to reveal the underlying micro-scale FR-generating subnetworks. FR rates vary with different pathologies (Ferrari-Marinho et al., 2015) and disease activity (Zijlmans et al., 2009), but are low in general (van Klink et al., 2014). We corrected for the inter-patient variability in event-rate by normalizing the number of events scored per channel by the total number of events per patient. We found a positive association between the percentage of spikes and FR-band node centrality measures. In other words: the higher the percentage of spikes, the more connected that node was in the FR-band functional network. In addition, FR-band Strength and EC was higher in channels with both spikes and FR, compared to channels with spikes only. FRs often occur in channels with high spike rates (Jacobs et al., 2008, van't Klooster et al., 2015). Electrodes displaying spikes and FR have been found to be more closely related to the epileptogenic tissue than electrodes with spikes without FR (van't Klooster et al., 2015). The relation between channels with high numbers of spikes and node centrality measures may thus be similar to the relation between presence of FR and centrality measures, which seems also the case in our data. Moreover, node centrality in the FR network might be a surrogate marker for the presence of FR in that channel at other moments in time. More patients, more channels with FRs per patient (for example using high-density grids), and combined ECoG-depth recordings with micro-wires are needed to properly test these hypotheses.

Both HFOs and hub nodes have been studied in relation to the resected area. Retrospective studies on ECoG data of epilepsy patients after surgery have linked post-surgical seizure-freedom to complete removal of the brain areas generating HFOs, especially FRs (Ochi et al., 2007, Jacobs et al., 2010, Wu et al., 2010, Cho et al., 2014, van Klink et al., 2014, van't Klooster et al., 2015). Retrospective network studies on ECoG data of surgical candidates have shown that inclusion of the nodes with the highest (inter)ictal centrality measures in the gamma-band network in the resection, correlated with post-surgical seizure-freedom (Ortega et al., 2008, Wilke et al., 2011, Varotto et al., 2012, Mierlo et al., 2013). This suggests that resection of brain areas generating HFOs, which we confirmed to be functionally isolated in the gamma-band network (Ibrahim et al., 2013a, Van Diessen et al., 2013b, Nissen et al., 2016), and resection of the most central nodes in the gamma-band functional network represent two approaches to achieve seizure freedom. Gamma-band activity is generated by networks of predominantly inhibitory interneurons (Fries et al., 2007, Cardin et al., 2009, Sohal, 2012). The reduced interneuron-mediated inhibition in brain tissue producing pathological HFOs may underlie the gamma band functional isolation of brain areas showing interictal HFOs at other instances. The closely surrounding healthy tissue may seem relatively well connected in the gamma frequency band. This is supported by the increase in gamma-band node centrality with increasing distance from the center of brain areas generating spikes and HFOs (Nissen et al., 2016). Both removal of the epileptogenic, FR-generating tissue and removal of the closely surrounding tissue through which a seizure initially spreads may, in this framework, lead to seizure-freedom.

The epileptic brain is a dynamic system. Early network studies on intracranial EEG data of patients with focal epilepsy can be divided into ictal and interictal studies. Studies assessing network changes accompanying the ictal state have consistently reported an increased synchronization and a more regular network topology, that is, an increased local connectivity in combination with disrupted long-distance connections, around seizure onset (Ponten, 2007, Kramer et al., 2008, Schindler et al., 2008, Wilke et al., 2011, Varotto et al., 2012, Mierlo et al., 2013). This is followed by a period of waxing and waning of the centrality in electrodes in the epileptogenic zone, which still remains higher than in non-epileptogenic areas, during the seizure (Ponten, 2007, Kramer et al., 2008, Schindler et al., 2008, Wilke et al., 2011, Varotto et al., 2012). At the end of the seizure there is again a slight increase in synchronization of the epileptogenic zone observed (Kramer et al., 2008, Schindler et al., 2008). Interictal network studies have reported increased centrality measures in the epileptogenic zone (Ortega et al., 2008, Wilke et al., 2011, Varotto et al., 2012) and decreased centrality measures in the SOZ (Warren et al., 2010, Van Diessen et al., 2013b). A recent study by (Burns et al., 2014) looked at the dynamics of brain connectivity in continuous ECoG and showed that the SOZ is functionally isolated interictally, but becomes more connected during seizure progression. A possible explanation is that the gamma-band functional isolation of the SOZ in the interictal state is the result of a compensatory strategy of the brain to prevent spreading of epileptic activity, which breaks down at the start of a seizure (Kramer et al., 2008, Ibrahim et al., 2013b, Khambhati et al., 2015). Seizures may then spread to the rest of the epileptogenic network through the pathological hub in the epileptogenic zone, close to the SOZ. Alternatively, the meso-scale gamma-band functional isolation and possible FR-band functional integration may represent a surrogate biomarker for the underlying truly epileptogenic micro-scale network that generates HFOs, from which activity spreads (through different frequency bands) to a wider brain area during seizures (Korzeniewska et al., 2014, Fuertinger et al., 2016). The confirmed functional gamma-band isolation of channels showing HFOs, in combination with the hint towards a positive association between the number of spikes and FR-band EC, and the increased FR-band Strength and EC in channels with spikes and FRs compared to channels with spikes only, support this theory, but the FR-band results are limited by the sample size and scale issues. The increased FR-band strength and EC in resected compared to non-resected channels is in congruence with this theory. Yet, the latter results are biased due to the relatively standard surgical procedure of anterior temporal lobectomy with or without amygdalohippocampectomy in patients with TLE. Future network studies in epilepsy should include the FR frequency band when studying the spatio-temporal dynamics of (epileptogenic) brain activity (Fuertinger et al., 2016).

Selecting patients based on epilepsy location and number and position of electrodes, which are known to influence calculation of network measures (van Wijk et al., 2010, Zalesky et al., 2010a, Zalesky et al., 2010b, Boersma et al., 2012, van Diessen et al., 2014a, van Diessen et al., 2014b), led to a relatively small but from a methodological network perspective homogeneous patient population. We did not correct for network density, but dealt with the network comparison problem by performing within-patient analyses. We included patients of different ages, with different epilepsy duration and underlying pathologies in our study. These clinical variables are known to affect network topology. Longer epilepsy duration is associated with a decreasing clustering coefficient and a more random network topology (Van Dellen et al., 2014). This has been suggested to be due to underlying pathology: for example brain tumors (Bartolomei et al., 2006) or hippocampal sclerosis (Dyhrfjeld-Johnsen et al., 2007) can result in increased randomization of the network. This randomization of the network over time may be a compensatory response to the regularized epilepsy network (Chiang et al., 2014). Decreases in clustering coefficient have also been linked to cognitive decline (Vaessen et al., 2012), which is often observed in patients with epilepsy. We omitted the influence of age by performing within-patient analyses. We minimized but could not exclude the influence of underlying pathology and disease activity by normalizing the event rates within patients. Yet, the robustness of the finding of functional gamma-band isolation of areas generating epileptic events across studies, and in patients with different underlying etiologies suggests an overcoupling epileptogenic mechanism.

To allow comparison with previous work by our group, we chose to compute functional connectivity using the PLI and characterize node centrality calculating the Strength and EC of each network node. The number of connectivity and network measures is rising exponentially, all with different advantages and disadvantages. Methodological choices are known to influence the results of network analysis, but to what extent is largely unknown (Wang et al., 2014, Bastos and Schoffelen, 2015, Kida et al., 2016). This makes it difficult to compare studies. We feel a future study should compare different connectivity and network measures to each other on the same real patient data with a clear gold standard and a control group (similar to what Wang et al., 2014 did on simulation data) to increase understanding. Eventually, international consensus would be welcome on which measures should be used in which circumstances, especially when network analysis will be implemented clinically.

Another potential limitation of our study is that network construction was derived from bipolar montages. Bipolar recordings consist of the measurement of the potential difference between pairs of neighboring electrodes. It removes the contribution of the reference electrode and approximates the local gradient of the electric potential in the direction between the electrodes. This montage was used to limit background noise for visual HFO analysis and subsequently to enable comparison between HFO and network analysis. A bipolar montage does however introduce an unknown directive component, and inequality in electrode usage: most electrodes are sampled twice, while outermost contacts are sampled once. This could influence the connectivity measured. Two studies investigating the effect of reference choice (montage) and volume conduction on different connectivity measures found however, that the PLI was the least affected measure (Stam et al., 2007, Christodoulakis et al., 2013).

Lastly, we double-checked the signals on epileptiform events and visually selected the event-free epochs, yet, although this is less likely in ECoG than in EEG/MEG recordings, low amplitude event-like activity may have been missed. One could even imagine that a continuous scale exists from ‘background’ activity that can be seen in (epileptic) networks towards event-like activity to true HFOs or spikes that fit with the definition of epileptiform activity.

To conclude, if the combination of meso-scale gamma-band functional isolation and FR-band functional integration indeed represents a surrogate biomarker for the underlying truly epileptogenic, micro-scale network that generates HFOs, this may have clinical implications and improve the success of epilepsy surgery. Importantly, the fact that we could extract these results from ‘silent’ epochs that did not show any spikes, ripples or FRs implies that there may be high frequency information hidden in the ‘baseline’ signal that can be used as a surrogate marker to point out channels that will show FRs at other moments in time, in a less time-consuming manner than the current visual analysis which depends on the occurrence of events. In addition, it may add to our understanding of the ‘architecture’ of epileptogenic networks and help unravel the underlying pathophysiology. More research is needed to confirm and extend our preliminary finding of functional integration in the FR frequency band functional network of channels suspected to cover epileptiform event generating brain tissue.

Funding

WZ is supported by the UMC Utrecht Alexandre Suerman Stipendium 2015. MvtK is supported by the Dutch Epilepsy Foundation grant number 2012–04. NvK is supported by the Dutch Brain Foundation grant number 2013–139 and the Dutch Epilepsy Foundation grant number 2015–09. WMO was supported by the Netherlands Organization for Scientific Research NWO-VENI grant number 016.168.038, and the Dutch Brain Foundation grant number F2014(1)-06. MZ is supported by the Rudolf Magnus Institute Talent fellowship, the ZonMW-VENI grant number 91615149.

Conflicts of interest

None of the authors have any conflicts of interest to disclose.