Spectral entropy indicates electrophysiological and hemodynamic changes in drug-resistant epilepsy - A multimodal MREG study.

OBJECTIVE: Epilepsy causes measurable irregularity over a range of brain signal frequencies, as well as autonomic nervous system functions that modulate heart and respiratory rate variability. Imaging dynamic neuronal signals utilizing simultaneously acquired ultra-fast 10 Hz magnetic resonance encephalography (MREG), direct current electroencephalography (DC-EEG), and near-infrared spectroscopy (NIRS) can provide a more comprehensive picture of human brain function. Spectral entropy (SE) is a nonlinear method to summarize signal power irregularity over measured frequencies. SE was used as a joint measure to study whether spectral signal irregularity over a range of brain signal frequencies based on synchronous multimodal brain signals could provide new insights in the neural underpinnings of epileptiform activity.
METHODS: Ten patients with focal drug-resistant epilepsy (DRE) and ten healthy controls (HC) were scanned with 10 Hz MREG sequence in combination with EEG, NIRS (measuring oxygenated, deoxygenated, and total hemoglobin: HbO, Hb, and HbT, respectively), and cardiorespiratory signals. After pre-processing, voxelwise SEMREG was estimated from MREG data. Different neurophysiological and physiological subfrequency band signals were further estimated from MREG, DC-EEG, and NIRS: fullband (0-5 Hz, FB), near FB (0.08-5 Hz, NFB), brain pulsations in very-low (0.009-0.08 Hz, VLFP), respiratory (0.12-0.4 Hz, RFP), and cardiac (0.7-1.6 Hz, CFP) frequency bands. Global dynamic fluctuations in MREG and NIRS were analyzed in windows of 2 min with 50% overlap.
RESULTS: Right thalamus, cingulate gyrus, inferior frontal gyrus, and frontal pole showed significantly higher SEMREG in DRE patients compared to HC. In DRE patients, SE of cortical Hb was significantly reduced in FB (p = .045), NFB (p = .017), and CFP (p = .038), while both HbO and HbT were significantly reduced in RFP (p = .038, p = .045, respectively). Dynamic SE of HbT was reduced in DRE patients in RFP during minutes 2 to 6. Fitting to the frontal MREG and NIRS results, DRE patients showed a significant increase in SEEEG in FB in fronto-central and parieto-occipital regions, in VLFP in parieto-central region, accompanied with a significant decrease in RFP in frontal pole and parietal and occipital (O2, Oz) regions.
CONCLUSION: This is the first study to show altered spectral entropy from synchronous MREG, EEG, and NIRS in DRE patients. Higher SEMREG in DRE patients in anterior cingulate gyrus together with SEEEG and SENIRS results in 0.12-0.4 Hz can be linked to altered parasympathetic function and respiratory pulsations in the brain. Higher SEMREG in thalamus in DRE patients is connected to disturbances in anatomical and functional connections in epilepsy. Findings suggest that spectral irregularity of both electrophysiological and hemodynamic signals are altered in specific way depending on the physiological frequency range.

Introduction

Epilepsy is a brain disease characterized by excess neuronal activity leading to abnormal or synchronized bursts of activity called seizures in the brain (Fisher et al., 2014b). In electroencephalographic (EEG) measurements, epileptic seizures show fast activity and high intensity. Epileptiform activity can be seen in most of the frequencies within the physiological EEG range. Focal seizures predominantly occur in theta (4–7 Hz), delta (0–3 Hz), or alpha (8–12 Hz) frequency bands, but also beta (13–30 Hz), gamma (30–100 Hz) and ripple (100–250 Hz) are related to epileptiform activity (Buzsáki et al., 1992; Engel Jr. et al., 2009; Fisher et al., 2014a). Also the very-low frequencies (VLF) of EEG can be used to detect epileptic seizure foci (Caspers et al., 1987; Ikeda et al., 1999; O'leary and Goldrings, 1964; Vanhatalo et al., 2003a). Besides EEG, heart rate variability (HRV) (up to 0.40 Hz) and respiration (0.15–0.40 Hz) are also known to be altered in epilepsy (Ansakorpi et al., 2011; Lotufo et al., 2012; O'Regan and Brown, 2005). Taken together, alterations in epilepsy in electrophysiological activity can be present in various frequencies of several types of signals.

A measure that sums up frequency changes over a wide range could be helpful in detecting abnormal inter-ictal brain activity in the absence of visually identifiable, diagnostic EEG discharges. Spectral entropy (SE) is a measure of signal irregularity, which sums the normalized signal spectral power (Shannon, 1948). Considering that most physiological signals are nonlinear, entropy as a nonlinear method is ideal to study neural signals (Faust and Bairy, 2012). Entropy in time domain, such as approximate entropy (Srinivasan et al., 2007), sample entropy (Song et al., 2012) and permutation entropy (Bruzzo et al., 2008; Li et al., 2014a) has been used to study epilepsy. SE is a frequency domain entropy measure, and it has been calculated in EEG measures in epilepsy (Aarabi et al., 2009; Mirzaei et al., 2010) and to monitor anesthesia, values 40 to 65 indicating adequate depth of anesthesia (Chhabra et al., 2016; Vakkuri et al., 2005; Viertiö-Oja et al., 2004). In addition to anesthesia, entropy has also been used to calculate vigilance estimates from EEG (Bruzzo et al., 2008; Shi et al., 2013). Bruzzo et al. found in same study that permutation entropy worked as a measure in seizure prediction and estimator of vigilance states in patients with focal drug-resistant epilepsy. In epilepsy, vigilance is known to be altered (Bruzzo et al., 2008; Danielson et al., 2011; Poirel, 1987), and in general the state can change during the scan, which can be seen for instance in global fMRI (Liu et al., 2018; Wong et al., 2016). SE thus offers a comparable, scale free metric that has already been shown to reflect epileptic pathology in intracranial (0.5–100 Hz) (Aarabi et al., 2009) and scalp EEG recordings (0–60 Hz) (Mirzaei et al., 2010).

Although it is most common to calculate entropy from EEG signals, it has been adapted in other physiological signals as well. ECG changes in epilepsy have been calculated by multiscale entropy, which is able to measure complexity of physiological data over different temporal scales (Liu et al., 2017). Also complexity of blood oxygen level dependent (BOLD) signal of functional magnetic resonance imaging (fMRI), have been studied (McDonough and Nashiro, 2014; Smith et al., 2014; Sokunbi et al., 2015), recently also in epilepsy (Gupta et al., 2017; Tavares et al., 2017). Gupta et al. have measured wavelet entropy in Rolandic epilepsy (Gupta et al., 2017) and Tavares et al. have used complexity analysis in epileptogenic focus localization (Tavares et al., 2017) based on BOLD signal. To the best of our knowledge, there are no reports on NIRS signal entropy in epilepsy. Simultaneous multimodal acquisition of all these parameters might enable a more comprehensive view of the pathophysiology of inter-ictal activity.

The most popular simultaneous multimodal technique in epilepsy research is combined EEG-fMRI (Bagshaw et al., 2004; Berman et al., 2010; Gotman, 2008; Ives et al., 1993; Laufs et al., 2006) and several studies have shown correlation between EEG and fMRI (Chang et al., 2016; Goldman et al., 2002; Huang-Hellinger et al., 1995; Ives et al., 1993; Moosmann et al., 2003; Wong et al., 2016). Also, EEG-fMRI studies have shown that large-scale networks play crucial role in both focal and generalized seizures (Centeno and Carmichael, 2014). Thus, several types of epilepsy share common functional abnormalities and labeling different types is not unambiguous. Also, despite epileptic seizures may result from wide area of brain, seizure activity also typically spreads over specific anatomic pathways (Laufs et al., 2011; Löscher and Ebert, 1996; Steriade, 2005). In addition to traditional resting state fMRI (Biswal et al., 1995), recently developed ultra-fast fMRI imaging sequence, 10 Hz magnetic resonance encephalography (MREG), increases the sampled signal frequency range 10-fold, i.e. from 0.5 to 5 Hz (Assländer et al., 2013). MREG further improves the detection of single epileptic events (Jacobs et al., 2014; Korhonen et al., 2014; Lee et al., 2013) by more accurate temporal separation of neuronal events and increased sampling rate that prevents cardiorespiratory aliasing (Hennig et al., 2007; Posse et al., 2012). Thus, MREG improves the spectral accuracy of multimodal neuroimaging in terms of characterizing the frequency ranges of epileptiform activity (Jacobs et al., 2014).

Over recent decades, evidence that the neurovascular coupling is altered in epilepsy has accumulated (Gotman et al., 2004; Gotman et al., 2006; Hawco et al., 2007; Mäkiranta et al., 2005). The commonly detected hemodynamic delay may change to precedence of hemodynamic changes compared to electrophysiological activity. Such reversed couplings have been reported in both EEG/fMRI and ECoG/NIRS studies in human studies and animal models of epileptiform activity (Gotman et al., 2004; Hawco et al., 2007; Jacobs et al., 2009; Osharina et al., 2010; Osharina et al., 2017; Pouliot et al., 2012). Near-infrared spectroscopy (NIRS) provides a direct, non-invasive measure of cortical oxyhemoglobin (HbO) and deoxyhemoglobin (Hb) oxygen delivery levels (McCormick et al., 1991; Mehagnoul-Schipper et al., 2002). Thus, NIRS has been found useful in investigation of epilepsy in both focus localization and monitoring oxygenation of responses during epileptic activity (Adelson et al., 1999; Nguyen et al., 2013; Obrig, 2014; Osharina et al., 2010; Watanabe et al., 2002). In this study we utilized NIRS to further verify spectral irregularity in epilepsy.

We used SE to investigate the irregularity of brain signal frequencies in drug-resistant epilepsy (DRE) with multimodal neuroimaging data combining MREG, NIRS, and DC-EEG. Our aim was to investigate dynamic spatiotemporal features of multimodal SE data searching for a marker to distinguish DRE patients and healthy controls (HC). Our null hypothesis was that the signal characteristics in the measured MREG, EEG and NIRS as measured by SE will be identical with DRE patients and HC. In the results we show multimodal evidence that support the rejection of this hypothesis. The findings suggest that spectral irregularity of both electrophysiological and hemodynamic signals are altered in specific ways depending on the physiological frequency range.

Methods and materials

Subjects

The study was approved by The Regional Ethics Committee of the Northern Ostrobothnia Hospital District. A written informed consent was obtained from all of the participants according to the Declaration of Helsinki. Study subjects were selected from the outpatient clinic of neurology department with the inclusion criteria of having a diagnosis of focal drug-resistant epilepsy. They were candidates for invasive treatment options, such as epilepsy surgery (although the patients were selected in this study with the criteria of not having an obvious structural cause of epilepsy in the MRI) and stimulation therapy. Patients either suffered from focal to impaired awareness or focal to bilateral tonic-clonic seizures. Previous clinical structural MRI was considered normal in seven patients, mild hippocampal sclerosis or amyglada pathology was suspected in two patients and mild focal cortical dysplasia was suspected in one patient. Previous interictal EEG-recording were either normal or showed occasional focal epileptiform discharges in eight patients and in two patients the seizure onset in temporal lobe was documented. EEG data of this study was measured and scanned for the presence of interictal discharges by experienced clinical neurophysiologist (S.R.), from Unit of Clinical Neurophysiology, and no clear discharges were found any of the patients. Patients were studied as a population that suffered complex, intractable epilepsy and compared to healthy subjects.

HC (no continuous medication, no alcohol consumption 12 h prior to scanning) were recruited from the email recruitment lists of the local University. Ten DRE patients (median ± standard deviation (STD): 34.3 ± 10.8 years, male/female: 3/7) (Table 1) and ten HC (median ± STD: 24 ± 8.1 years, male/female: 8/2) underwent multimodal imaging of the brain including MREG sequence, NIRS and EEG.

Table 1

Demographics of patients with drug-resistant epilepsy (DRE).

Patient	Gender	Age (years)	Duration of epilepsy (years)	Type of seizures	Etiology of epilepsy	AEDs		Number of scans
1	M	32	8	focal to impaired awareness	S?	CBZ 800 mg, PGB 300 mg		2
2	F	36	10	focal to impaired awareness	S	LTG 550 mg, TPR 200 mg		4
3	F	26	19	focal to bilateral tonic-clonic	U	LCM 200 mg, LTG 175 mg,		2
						VPA 1500 mg, CLB 10 mg
4	F	52	4	focal to impaired awareness	U	LCM 150 mg, ZNS 100 mg		4
5	M	36	6	focal to bilateral tonic-clonic	U	LEV 2500 mg, OXCB 1200 mg,		4
						LCM 300 mg
6	M	53	13	focal to bilateral tonic-clonic	U	LTG 100 mg, LCM 200 mg,		3
						TPR 125 mg
7	F	22	12	focal to bilateral tonic-clonic	U	VPA 900 mg, LTG 250 mg,		4
						CLB 10 mg
8	F	25	11	focal to impaired awareness	U	OXCB 1200 mg, ZNS 400 mg		4
9	F	35	12	focal to impaired awareness	U	TPR 100 mg, LCM 500 mg,		4
						CLB 30 mg
10	F	26	7	focal to impaired awareness	U	LEV 2000 mg, LCM 150 mg,		4
						LTG 400 mg

AEDs, anti-epileptic drugs; S, structural; U, unknown; CBZ, Carbamazebine; PGB, Pregabalin; LTG, Lamotrigine; TPR, Topimarate; LCM, Lacosamide; VPA, Valproic acid; CLB, Clonazepam; ZNS, Zonisamide; LEV, Levetiracetam; OXCB, Oxcarbamazebine.

Data collection

All subjects were scanned in Oulu (Finland) using a Siemens MAGNETOM Skyra 3 T scanner (Siemens Healthcare GmbH, Germany) with a 32-channel head coil. Scanning was carried out with MREG sequence obtained from Freiburg University via collaboration with Jürgen Hennig's group (Assländer et al., 2013; Lee et al., 2013; Zahneisen et al., 2012). MREG is a three dimensional (3D) spiral, single-shot sequence that undersamples 3D k-space trajectory. It allows approximately 20–25 times faster scanning than the conventional fMRI sequence by sampling the brain at the frequency of 10 Hz. 10-min resting-state scan yielding 5822 brain volumes with multimodal set up (at least seven different modalities) was performed as before (Korhonen et al., 2014), with instructions to relax and focus on a cross on the screen. Our multimodal neuroimaging set-up enables simultaneous scanning of MREG, EEG, NIRS, continuous non-invasive blood pressure (NIBP) and anesthesia monitoring as a synchronous system. Four 10-min scans were carried out in the majority of the DRE patients to observe potential seizure activity (specific information in Table 1). Four of HC were scanned twice and the rest only once. The MREG and EEG data from the first scan was chosen for analysis in 6/10 of DRE patients and HC, and second for the rest of the study subjects. The second scan was used, because MREG data of first scan was not comparable in 4/10 HC. To exclude differences between the first and the second scan, the second scan was randomly chosen in four DRE patients. In NIRS, the best available data was chosen into further analysis, due to the technical problems during the patient scans. Exact information can be found from supplementary table 1.

The following parameters were used for the 3D whole brain MREG sequence: repetition time (TR) 100 ms, echo time (TE) 36 ms, field of view (FOV) 192 mm3, 3 mm cubic voxel, and flip angle (FA) 25°. The data were reconstructed using L2-Tikhonov regularization with lambda 0.1, with the latter regularization parameter determined by the L-curve method with a MATLAB recon-tool from sequence developers (Hugger et al., 2011). 3D MPRAGE T1 anatomical images (TR 1900 ms, TE 2.49 ms, TI 900 ms, FA 9°, FOV 240, slice thickness 0.9 mm) were scanned to register the MREG data onto 4 mm Montreal Neurologic Institute space (MNI).

EEG was recorded with 32-channel, MR-compatible BrainAmp system (Brain Products, Gilching, Germany) according to the international 10–20 system. Ag/AgCl electrodes were placed to the head and skin potential was reduced by stick abrasion technique to allow DC-measurement (Hiltunen et al., 2014; Vanhatalo et al., 2003b). BrainAmp SyncBox was used to ensure that the EEG amplifier and the scanner are in synchrony via a transistor–transistor logic pulse. Electrode impedances were < 5 kΩ, sampling rate 5 kHz, band pass from DC to 250 Hz and ECG was recorded. Signal quality was tested outside the scanner room by recording 30 s eyes open and eyes closed.

Fluctuations in cerebral blood flow and volume, (more specifically concentration changes of HbO and Hb) were measured in synchrony with fMRI using an MR-compatible NIRS device (Sorvoja et al., 2010). The device allows to change the combination of wavelengths (Korhonen et al., 2014; Myllylä et al., 2014), and in this study we used wavelengths of 660 and 830 nm. These wavelengths were selected because they are located on both sides of the isosbestic point of the absorption spectrum of HbO and Hb, that is, at 800 nm, where the extinction coefficient of oxygenated and deoxygenated hemoglobin is the same (Boas et al., 2001; Jöbsis, 1977; Korhonen et al., 2014). Furthermore, we have shown that this wavelength combination is sensitive to changes in hemoglobin concentrations (Myllylä et al., 2014). The sampling rate of NIRS data acquisition was 10 kHz and the sensor was placed on forehead.

Besides MREG, EEG, and NIRS, continuous noninvasive blood pressure (Myllylä et al., 2017; Myllylä et al., 2011), respiration (respiration belt) and fingertip photoplethysmography (PPG) (Siemens MAGNETOM Skyra 3 T scanner devices) and PPG, respCO2, ECG, and cuff-based blood pressure were registered using 3 T MRI-compatible anesthesia monitor (GE Datex-Ohmeda™; Aestiva/5 MRI) as described previously (Korhonen et al., 2014) (Fig. 1).

Fig. 1

Scanning, preprocessing and analysis with multimodal setup. A) MREG, EEG and NIRS methods separated. B) Voxelwise MREG signal after preprocessing and FIX. Temporal signals (global MREG, EEG, 32 channels, NIRS (Hbo, Hb, HbT) used in analysis, 10-min measures of all and 2-min sliding window measures of global MREG and NIRS. C) PSD calculated from all signals. D) SE calculated by Shannon's equation separately from all signals on chosen sub-bands (0.009–0.08 Hz, 0.12–0.4 Hz, 0.7–1.6 Hz, 0.08–5 Hz, 0–5 Hz). As a result SE value for all voxels in MREG and for global MREG signal, for EEG (all 32 leads, visualized in map), and for NIRS (Hbo, Hb, HbT). 0 < SE value <1. In SEMREG: bright yellow = highest SE, dark red = lowest SE. In SEEEG: purple = highest SE, bright blue = lowest SE. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

Preprocessing of MREG data

MREG data were preprocessed with the FSL pipeline (FSL 5.09) (Jenkinson et al., 2012). First, the data were high-pass filtered with a cut-off frequency of 0.008 Hz. To minimize T1-relaxation effects, 180 time points were removed from the beginning. Head motions were corrected with FSL 5.08 MCFLIRT software (Jenkinson et al., 2002). After MCFLIRT, brain extraction was carried out with optimization of the deforming smooth surface model, as implemented in FSL 5.08 BET software (Smith, 2002) using threshold parameters f = 0.3 and g = 0. For 3D MPRAGE volumes, parameters f = 0.25 and g = 0.22 with neck clean-up and bias field correction options were used. Spatial smoothing was done with fslmaths 5 mm FWHM Gaussian kernel. Motion left despite correction was controlled by SE calculation and observing movement (mm) from MCFLIRT data. Possible effect of respiration was considered by calculation of SE from respiration belt data (Skyra). FSL-tool MELODIC was driven to set fMRI data into time-courses and spatial maps using single-session, probabilistic independent component analysis (ICA) (Beckmann and Smith, 2004). Three-dimensional MPRAGE images were used to register the MREG data into 4-mm MNI space prior to ICA as a standard procedure in MELODIC. Additionally, preprocessing step FIX (FMRIB's ICA-based X-noiseifier) was driven for MREG data to separate artifacts from neural activity signals in resting-state fMRI data (Griffanti et al., 2014; Salimi-Khorshidi et al., 2014).

Before using FIX in the study data (DRE and HC), separate data set was needed. ICA MELODIC was driven similarly for the separate data set (10 earlier scanned healthy controls) as for the study data set. After ICA the independent components were visually labeled to noise or signal by hand using Melodic Result Viewer (MELVIEW). Noise components were further sub-labeled to cardiac, respiratory, white matter, motion or sinus sagittalis/transversus when relevant. Independent components were classified as resting state networks when several criteria were satisfied (Griffanti et al., 2014; Zerbi et al., 2015). FIX was performed for MREG data of DRE patients and HC.

The global MREG signal was computed in FSL by fslmeants command that counts average time series (intensities).

Preprocessing of EEG data

DC-EEG was preprocessed off-line using Brain Vision Analyzer (version 2.0, Brain Products). Gradient artifacts and ballistocardiographic artifacts were removed using standard average artifact subtraction method (Allen et al., 1998; Allen et al., 2000). DC-EEG data were downsampled to 10 Hz using MATLAB for comparable analysis with MREG and NIRS. All 32 EEG channels were analyzed individually.

Preprocessing of NIRS data

For quality control, NIRS data was visually checked in relation to motion artifacts, baseline shifts and unwanted noise. Sensitivity of the NIRS is regularly tested and validated with optical head phantoms, designed for the purpose (Korhonen et al., 2014). Raw NIRS time series were converted into time series of Hb, HbO, and total hemoglobin (HbT) using Homer2 software (Huppert et al., 2009) and downsampled to 10 Hz using MATLAB. The method for converting raw NIRS data to obtain temporal changes in HbO and Hb concentration is based on modified Beer-Lambert law (Boas et al., 2001; Strangman et al., 2003).

SE analysis

Before analysis, a linear trend was removed in MATLAB (vR2016b; The Math Work, Natick, MA) from all of the signals.

SE analysis for MREG, EEG and NIRS was performed in MATLAB. SE is the measure of irregularity: high entropy tends to be found from broader spectrum with a lot of frequencies, when in low entropy there are only a few spikes where energy is concentrated (Inouye et al., 1991; Urigüen et al., 2017). First we calculated power spectral density (PSD) through FFT analysis taking the signal's length into account, and then, normalized PSD by (Quinquis, 2008). SE was calculated using Shannon's entropy:where pi is probability distribution (normalized PSD), f1 denotes the starting frequency given in number of bins and f2 is ending frequency as bins (5 Hz corresponds to 4096 in 10 min length and 1024 in 2 min length).”

Like described earlier, epilepsy affects in various frequency bands in the signals having different origins. In addition the conventional frequencies (>0.5 Hz) in EEG, interictal and ictal activity have been seen in VLF (Caspers et al., 1987; Ikeda et al., 1999; O'leary and Goldrings, 1964; Vanhatalo et al., 2003a), and also respiratory and cardiac changes seen in ECG (e.g heart rate variability) may be distracted by epilepsy (Ansakorpi et al., 2011; Lotufo et al., 2012; O'Regan and Brown, 2005). In addition we have also shown that MREG signal in DRE has increased variance in very-low and respiratory frequency, while the cardiac frequency has less variance (Kananen et al., 2018). Thus we wanted to further investigate this multimodally with separate sub-bands for the analysis. SE values were computed for whole 10-min global MREG, EEG and NIRS signals (5822 time points) in 0–5 Hz (fullband, FB), 0.08–5 Hz (near-fullband, NFB), brain pulsations in 0.009–0.08 Hz (very-low frequency, VLFP), 0.12–0.4 Hz (respiratory frequency, RFP), and 0.7–1.6 Hz (cardiac frequency, CFP) by choosing the sub-band from the power spectrum. For a dynamic analysis of the global MREG signal and NIRS (1164 time points) for all chosen bandwidths, SE was calculated for 2-min windows with 50% overlap to see possible effect of vigilance level.

SE was calculated for all sub-bands of 10-min MREG for each voxel individually, which allows whole-brain localization of SE values. Three-dimensional standard sequence was used to register the SE data into 2-mm MNI space for comparable analysis. Data outside the brain was removed by multiplying the data with the 2-mm MNI standard brain mask. Moreover, SE values of EEG channels in separate sub-bands were visualized spatially by inhouse Matlab script based on 10–20 system EEG.

Statistical analyses

Two-tailed unpaired samples t-test, randomize analysis with 5000 permutations was performed using FSL to produce a statistical map of differences between DRE patients and HC as groups (Winkler et al., 2014). Threshold-free cluster enhancement (TFCE) correction in both directions (HC > DRE, HC < DRE) separately (Smith and Nichols, 2009) and family-wise error rate control (FWE-corrected) (Holmes et al., 1996) were implemented in the analysis. For the two-dimensional signals, two-sample, unpaired Wilcoxon test (MATLAB ranksum) was used to compare SE-values of DRE patients and HC.

Results

Motion control

Mean and SE of MREG relative movement data (obtained from MCFLIRT) were compared between groups. No significant differences (SE of relative movement 0.907 ± 0.015 vs. 0.916 ± 0.018, p = 0.2 and relative movement mean 0.053 ± 0.011 vs. 0.051 ± 0.010, p = 0.9) between DRE patients vs. HC were found. Additionally, SE of relative movement signal in 2-min time windows (50% overlap) was calculated, but showed no significant group differences.

Respiration control

No significant differences were found between DRE patients and HC in FB (0.526 ± 0.066, n = 8 vs. 0.554 ± 0.063, n = 6, p = 0.49) nor in respiratory frequency (0.721 ± 0.089, n = 8 vs. 0.752 ± 0.092, n = 6, p = 0.57).

Spatial voxel-level SE analysis of MREG signal

In FIX-corrected SEMREG map at FB, DRE patients had significantly higher (p < 0.05) SE located in right thalamus, cingulate gyrus (mostly anterior), inferior frontal gyrus and frontal pole (Fig. 2B). No difference was found without FIX correction (Fig. 2A) (Table 2).

Fig. 2

A) Average spectral entropy (SE) maps of FB MREG and FB MREG including FIX for DRE patients (n = 10) and HC (n = 10). Red illustrates the lowest and yellow the highest SEMREG values (cut off 0.5). B) SEMREG of FB MREG FIX was significantly higher in DRE patients than in HC in right thalamus, mid cingulate gyrus, inferior frontal gyrus, and frontal pole.

Table 2

Group comparison of spatial SEMREG between DRE patients and HC.

DRE > HC	Volume (μl)	Mean t-score	SEM	Max int.	MNI coordinates
					X	Y	Z
Right Thalamus	608	3.9548	0.0312	4.7676	14	−14	2
Cingulate Gyrus, anterior division	496	3.8828	0.0467	5.1316	4	−2	38
Inferior Frontal Gyrus, pars triangularis	352	4.0919	0.041	4.6421	44	22	12
Frontal pole, Left Cerebral Cortex	136	4.2881	0.0375	4.6113	−4	62	−6

FB SEMREG was significantly (p < 0.05) higher in DRE (n = 10) than HC (n = 10) in several clusters. Volume, cluster size as voxels, Mean t-score, standard error of mean (SEM), maximum intensity (Max int.) and Max (X,Y;Z) as MNI coordinates are listed.

SEEEG

SE of DC-EEG was visually more variable in DRE patients than in HC. Significantly higher values of SE were found in DRE patients in DC-EEG signal (channels FC2 and Poz) and VLF band (CP2) (Fig. 3A). These were located around the vertex on the cranium. Significantly lower SE was found in RFP in DRE patients compared to HC (channels Fp1, Fp2, P3, O2 and Oz, p < 0.05), values located in electrodes surrounding lower parts of the cranium (Fig. 3B). Furthermore, in DRE patients compared to HC, RFP changes in SE tended to be reduced in all electrodes at the lower edge of the cranium. Spatially the DC-EEG SE values in DRE patients were more irregular between neighboring electrodes (individually and as a group), whereas results of HC followed a more controlled pattern; lowest values were found in the center vertex of the brain, and highest values in the periphery.

Fig. 3

Group-differences between DRE patients and HC in different SE measures, presented for 0–5 Hz, FB. Orange and yellow indicate a significant group-difference at p < 0.05(*), and p < 0.01 (**), respectively. A) Static SEEEG for separate channels (mean ± STD) and spatial representation of EEG channels with indicators for increases in DRE patients. NIRS sensor indicated. B) SEMREG map indicating an increase in frontal pole (where also NIRS sensor placed) and cingulate gyrus in DRE patients. C) SENIRS group differences in HbO, Hb (significantly lower in DRE patients), and HbT.

SENIRS

SENIRS was significantly lower (Table 3) in DRE patients compared to HC in FB (p = 0.045), in NFB (p = 0.017) and CFP (p = 0.038) in Hb and in RFP in HbO (p = 0.038) and HbT (p = 0.045).

Table 3

Group comparison of 10-min SENIRS between DRE patients and HC, significant results.

Band		SENIRS, Average ± Standard deviation		p-value
		DRE patients	HC	* < 0.05
		n = 9	n = 10
0–5 Hz	Hb	0.40 ± 0.06	0.46 ± 0.10	0.045*
0.08–5 Hz	Hb	0.44 ± 0.13	0.58 ± 0.14	0.017*
0.12–0.4 Hz	HbO	0.85 ± 0.05	0.88 ± 0.02	0.038*
	HbT	0.85 ± 0.05	0.88 ± 0.03	0.045*
0.7–1.6 Hz	Hb	0.39 ± 0.20	0.54 ± 0.24	0.038*

HbO, oxygenated; Hb, deoxygenated; HbT, total hemoglobin.

Dynamic sliding-window analysis of SENIRS (Table 4) revealed that DRE patients tended to have lower SE values compared to HC (Fig. 3B). In FB, SENIRS (HbO, Hb, and HbT) in Hb was lower in DRE patients compared to HC in Hb in minutes 5 to 7 (p = 0.045), in 0.08–5 Hz bandwidth in Hb in minutes 0 to 2 (p = 0.045), 1 to 3 (p = 0.038), 2 to 4 (p = 0.017), 3 to 5 (p = 0.031) and 4 to 6 (p = 0.014) and HbT 0 to 2 (p = 0.045). SENIRS was significantly lower in DRE patients also in VLFP in Hb in minutes 8 to 10 (p = 0.045), in RFP in HbO in minutes 4 to 6 (p = 0.009) and 8 to 10 (p = 0.038), Hb in minutes 4 to 6 (p = 0.021) and 6 to 8 (p = 0.031) and HbT 2 to 4 (p = 0.038), 3 to 5 (p = 0.011), 4 to 6 (p = 0.0046) and 8 to 10 (p = 0.038) (Fig. 4).

Table 4

List of significant differences between DRE patients and HC in dynamic 2-min sliding window SE calculated based on NIRS (HbO, Hb, HbT, 9 DRE patients and 10 HC) and global MREG signal (10 DRE patients and 10 HC).

SENIRS		DRE patients	<	HC	p-value
		n = 9		n = 10	* < 0.05	** < 0.01
Time window [min]	0–2	1–3	2–4	3–5	4–6	5–7	6–8	7–9	8–10
Band
0–5 Hz						Hb*
0.009–0.08 Hz									Hb*
0.08–5 Hz	Hb*, HbT*	Hb*	Hb*	Hb*	Hb*				HbO*
0.12–0.4 Hz			HbT*	HbT*	HbO**, Hb*, HbT**		Hb*		HbT*
0.7–1.6 Hz
SEGlobal MREG		DRE patients n = 10		HC n = 10
0–5 Hz									DRE > HC*
0.009–0.08 Hz	DRE > HC*						DRE < HC*
0.08–5 Hz
0.12–0.4 Hz						DRE > HC*	DRE > HC*
0.7–1.6 Hz

HbO, oxygenated; Hb, deoxygenated; HbT, total hemoglobin. *: p < 0.05, **: p < 0.01.

Fig. 4

Group-differences between DRE patients and HC in different SE measures presented for 0.12–0.4 Hz, RFP. Orange and yellow indicate a significant group-difference at p < 0.05(*), and p < 0.01 (**), respectively. A) Static SEEEG for separate channels (mean ± STD) and spatial representation of EEG channels with indicators for increases in DRE patients. B) Dynamic sliding-window results in SEMREG. C) Dynamic sliding-window results in SENIRS HbT.

SE of the global MREG signal

DRE patients showed significantly elevated SEMREG in FB in minutes 8 to 10 (p = .038) in global MREG FIX, in RFP in minutes 5 to 7 –and 6 to 8 (p = 0.028 and p = 0.01) and in VLFP in minutes 0 to 2 (p = 0.014) and minutes 6 to 8 (p = 0.044) (Table 4). SE of the whole 10-min global MREG signal without FIX-correction showed no significant difference between study groups.

Discussion

We compared SE between DRE patients and HC in simultaneously measured synchronous MREG, EEG and NIRS signals and detected spatially overlapping differences between the groups in SE of all three modalities in 0–5 Hz frequency range. Interestingly, SEEEG in frontal midline and SEMREG in anterior cingulate gyrus were increased in DRE patients. We found lower SE in RFP in DRE patients compared to HC in both EEG and NIRS signals, which were in line with previous evidence from spatial SEMREG in cingulate gyrus. Taken together, using multimodal imaging techniques, we identified spatiotemporal alteration in brain areas at different frequency ranges in DRE patients.

Entropy and spectral entropy in epilepsy – DC-EEG

Epilepsy is a disease, which tends to cause rhythmic activity in certain frequencies but also irregular bursts (Fisher et al., 2014b), making a measure of signal complexity ideal to study epileptiform activity. It is known that synchronous neuronal activity and recruitment of surrounding neurons are required in epileptic activity (Schwartzkroin, 1997) resulting increase in regularity. This reduced complexity of EEG has been demonstrated during inter-ictal state of patients with epilepsy compared to HC (Kim et al., 2002; Mirzaei et al., 2010) suiting to our SEEEG results. SE has been shown to reflect epileptic pathology in intracranial (0.5–100 Hz) (Aarabi et al., 2009) and scalp EEG recordings (0–60 Hz) (Mirzaei et al., 2010) by showing higher SE during epileptic seizures compared to healthy and inter-ictal states. Using Shannon's SE it has been shown that SE is elevated in idiopathic epilepsy compared to controls in 6–12 Hz frequency range in EEG (Urigüen et al., 2017). Additionally, lower SE in ictal state has been demonstrated in baseline alpha sub-band (8–15 Hz) EEG compared to both HC and inter-ictal state as well as in whole EEG (0–60 Hz) in inter-ictal versus healthy states (Mirzaei et al., 2010).

In the current study, both higher SEEEG but also higher STD values were detected in DRE patients, which were located in midline frontal and parieto-occipital regions. Higher SE has also been found earlier in EEG in higher frequencies (6–12 Hz) (Urigüen et al., 2017), but wide 0–60 Hz EEG measures SE was previously lower in inter-ictal versus healthy states (Mirzaei et al., 2010). Urigüen et al. found that SE values decreased with time from the last seizure (Urigüen et al., 2017), which can have an effect for our results too. Interestingly, the SEEEG was decreased in DRE patients in RFP especially on the forehead and also almost in all channels, which supports the corresponding low SE in NIRS and global MREG in RFP.

fMRI signal - MREG BOLD

Thalamus has been linked to seizure initiation, propagation and spread in temporal lobe epilepsy (Bertram et al., 2008; Guye et al., 2006; Li et al., 2014b; Spencer, 2002), supporting our high SE results in right thalamus. Recently, a worldwide epilepsy study showed lower structural volume in right thalamus, raising the importance of a thalamo-cortical pathway in epilepsy (Whelan et al., 2018). Thalamo-cortical circuits are thought to drive the cortex with strong connections to limbic system (Çavdar et al., 2008) as well as to regulate consciousness and sleep (Hale et al., 2016; John, 2001; Llinás et al., 1998). Recently, thalamic functional connectivity was found to be altered during sleep in idiopathic generalized epilepsy (Bagshaw et al., 2017). Furthermore, thalamic modulation may have inhibitory effects on anterior cingulate cortex seizures (Chang and Shyu, 2014). It also has been shown that functional connectivity and BOLD low frequency power are reduced by mild anesthesia in cingulate cortex (Greicius et al., 2008). Both, thalamus and cingulate cortex are linked to altered functions in epilepsy.

Moreover, SE was higher in DRE patients in frontal pole and inferior frontal gyrus. Frontal pole is still poorly understood (Semendeferi et al., 2001), nevertheless, together with the anterior cingulate cortex and thalamus, the frontal pole is a part of a corticolimbic circuit (Navratilova and Porreca, 2014), which has been suggested to play a role in respiratory control (Evans, 2010) and epileptic seizures in addition to corticothalamic circuit (Holmes, 2008). We found higher SE in inferior frontal gyrus, which is responsible for speech and language function (Lazar et al., 2000; Paulesu et al., 1997). It is known that focal epileptic activity induced by inter-ictal changes or seizure spread can influence to speech and language organization (Janszky et al., 2003).

Hemodynamic signal – NIRS

Numerous studies show that epileptic seizures are associated with an increase in regional cerebral blood flow (CBF). A portable and noninvasive technique capable of measuring CBF related parameters thus could be potentially utilized in diagnosis of epilepsy. NIRS has already proved its potential in measuring hemodynamic changes in the brain during epileptic seizures (Jeppesen et al., 2015; Obrig, 2014; Peng et al., 2014; Watanabe et al., 2000; Watanabe et al., 2002). In general, most of these studies have reported changes in hemodynamic parameters (HbO, HbR, HbT), following epileptic activity. However, both animal and human studies have found also “pre-seizure” changes, where hemodynamic responses are disturbed several minutes or seconds before the electroencephalographic seizure onset (Mäkiranta et al., 2005; Osharina et al., 2017; Zhang et al., 2014). In other words, neurovascular coupling is disturbed in epileptic brains even in the absence of detectable epileptic seizure activity.

In this study, no epileptiform activity was detected with the EEG and thus we were able to evaluate SE as a commensurable marker for multimodal data during inter-ictal state. To the best of our knowledge this is the first study of SE NIRS signal in epilepsy. We found changed SE from all multimodal measures synchronously. In general there are only a few studies based on NIRS entropy in brain diseases. Decreased complexity of brain signals was linked to mortality and brain injury in the study of preterm infants (Sortica da Costa et al., 2017). If decreased entropy is assumed to be associated with a pathological state (Goldberger et al., 2002), this is also supported by our findings with DRE patients.

Just like Liu et al. described decreased heart rhythm complexity in DRE (Liu et al., 2017), our results suggest that respiratory and cardiac variation measured in frontal cortex are decreased in DRE. In 10-min measurements, patients had lower SE in HbO and Hb in RFP and in Hb in CFP, which was especially seen in first 6 min in dynamic measure. Interestingly, significantly lower SE values, in all other bands than VLFP were observed in DRE patients compared to HC. Decreased complexity in physiological frequency ranges can thus be seen in multimodal brain signals of DRE patients.

Physiological considerations

Autonomic functions such as respiration and HRV are known to be altered in epilepsy (Ansakorpi et al., 2011; Lotufo et al., 2012; O'Regan and Brown, 2005; Tomson et al., 1998). Our lower SE findings in NIRS and EEG measured from brain especially in parasympathetic frequency range (0.12–0.4 Hz, RFP) are supported by earlier studies (Jansen and Lagae, 2010; Mukherjee et al., 2009; Sevcencu and Struijk, 2010). Recent article by Liu et al. calculated multi-scale entropy from ECG in DRE patients and showed similar results than our finding of lower entropy in DRE patients in parasympathetic frequency range (Liu et al., 2017). Although EEG is usually though to have neuronal origin, respiratory-related DC shifts can be seen in EEG (Nita et al., 2004). Furthermore in our study, DRE patients had higher SE in anterior cingulate gyrus, which has been shown to be an area of autonomic functions such as respiratory control (Critchley et al., 2003; Woodward et al., 2008). Taken together the EEG and NIRS data indicate that altered autonomic functions may also be altered within the brain physiological signal in DRE. The autonomic mechanisms may be linked to the control of recently discovered physiological pulsation mechanisms of the brain (Kiviniemi et al., 2016). Our findings strongly suggest abnormal respiration pulsation mechanism in the brain in epilepsy, although we did not see difference in direct respiration measure. Findings are supported by previous literature that venous return from the brain is dictated by intra-thoracic respiratory pressure pulsations (Schaller, 2004) and respiration affects for movement of cerebrospinal fluid (Klose et al., 2000; Yamada et al., 2013).

Considering that SE is an established measure of consciousness (bispectral index, depth of anesthesia), findings between these connections may be linked to known altered vigilance in epilepsy (Danielson et al., 2011; Poirel, 1987). Vigilance changes in drug-resistant, focal epilepsy at pre-ictal and ictal state have been found by calculating permutation entropy scalp EEG (Bruzzo et al., 2008). Interestingly, global MREG showed the following results: in RFP patients tended to have lower entropy until 4–6-min, after which SE increased compared to controls. Vice versa to RFP, SE in VLFP decreased after 6-min in patients compared to controls. Global fMRI signal has been linked to vigilance and arousal by correlating with electrophysiology (Liu et al., 2018; Wong et al., 2013), suiting the altered vigilance in epilepsy (Danielson et al., 2011; Poirel, 1987), and this relationship should be studied further to understand the function between vigilance and disease itself. Furthermore, sleep disturbances together with respiratory changes have been linked to epilepsy (Manni and Tartara, 2000; Oliveira et al., 2000).

Preprocessing considerations

We used FIX to clean various types of noise from the data by auto-classifying bad and good components of the signal (Griffanti et al., 2014; Salimi-Khorshidi et al., 2014). We got highly localized results between DRE patients and HC from MREG signal after FIX, which decreases the possibility that difference between groups was caused by noise. It is already shown that FIX is an efficient method on noise removal in traditional fMRI BOLD (Boyacioglu et al., 2015; Feis et al., 2015; Pruim et al., 2015; Wang and Li, 2015), and our results support its efficiency also in MREG data. Like other fMRI researchers (Liu et al., 2018), also we used parallel analysis without and with FIX focusing on FIX results.

Limitations

Because of the limited number of good quality 10-min data in modalities studied (simultaneous MREG, NIRS and EEG) we were unable to get age- and gender-matched data, ending up with more females with older age in the DRE group. Complexity of heart rate dynamics associated with parasympathetic activity has been found to be higher in women and high frequency HF power tends to decrease with age (Ryan et al., 1994). In a recent publication on entropy, impairment of heart rate complexity seemed to be related to epileptic process itself, and not to age, AED regimen, severity nor duration of epilepsy (Liu et al., 2017).

As no inter-ictal or ictal epileptiform discharges were detected during the measurements, our results cannot be straight provided as results of simultaneous epileptiform activity. However, epilepsy is the whole brain disease covering many physiological aspects (e.g autonomic nervous system) and our findings may thus be related to other subclinical phenomena in epilepsy such as inhibitory mechanisms in parasympathetic and/or sympathetic nervous system (Ansakorpi et al., 2011; Tomson et al., 1998). Furthermore, as parasympathetic vagal nerve stimulation reduces epileptic seizure activity (Ben-Menachem, 2002; Wheless et al., 2018; Yang et al., 2007), our findings of elevated SE may reflect mechanisms inhibiting epileptic activity instead of interictal epileptic activity. The found abnormality in three different modality supports the epileptic origin of the findings in our study.

Heterogeneity of the patient sample was observable, still sharing the similarities as drug-resistance and focal type seizures. Overall, epilepsy is a complex disease with variance of etiology, AEDs used, seizure frequency etc. and in this drug-resistant epilepsy population epilepsy is a wide range of brain diseases which all manifest themselves as epileptic seizures. However, structural and functional MRI studies have shown that seizure propagation pathways are often similar in focal epilepsies (Centeno and Carmichael, 2014; Laufs et al., 2011; Löscher and Ebert, 1996; Steriade, 2005), which supports our choice to study this group.

Unfortunately, NIRS data was not always obtainable from DRE patients and therefore in five cases the NIRS data was taken only from one of the four scans available in a patient. However, excluding the non-matching NIRS measurements, the results still showed decrease in SE in concordance to frontal EEG data, especially significant in the respiratory range during first six minutes. We did not see difference in SE between groups in respiration belt data, but unfortunately, we had only eight data of DRE patients and six of HC available.

Conclusions

Altered SE in DRE patients in cingulate gyrus of MREG and in RFP of NIRS and EEG together support the prior observation that brain pulsation activity in parasympathetic respiratory frequency range is altered in epilepsy. Spatial SEMREG was significantly higher in DRE patients compared to HC in thalamus, which supports alterations in thalamic connections in epilepsy. SE findings could be linked to abnormal parasympathetic neurovascular control in epilepsy.

The following are the supplementary data related to this article.

Supplementary table 1

Scan number by modality used in SE analysis for DRE patients.

Supplementary material

Conflict of interests

The authors declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Funding

This work was supported by Health and Biosciences Doctoral programme -grant (HH), Jane and Aatos Erkko Foundation grant (VKi), Academy of Finland and Aivosäätiö TERVA grant 314497 (VKi,TM), Academy of Finland Grants 123772 & 275342 (VKi), The SalWe Research Program for Mind and Body (Tekes—the Finnish Funding Agency for Technology and Innovation, Grant No. 1104/10) (VKi), Novo Nordisk Foundation NNF17OC0030124 (TM), Finnish Medical Foundation (VKi, TT), Finnish Neurological Foundation, KEVO grants from Oulu University hospital (VKi), Epilepsy Research Foundation (JK), Finnish Cultural Foundation, North Ostrobothnia Regional Fund (JK), Orion Research Foundation sr (JK, TT), Tauno Tönning Foundation (JK).