Recursive independent component analysis (ICA)-decomposition of ictal EEG to select the best ictal component for EEG source imaging

Journal Pre-proofs Recursive independent component analysis-decomposition of ictal EEG to select the best ictal component for EEG source imaging Mohammad Ashfak Habib, Fatimah Ibrahim, Mas S. Mohktar, Shahrul Bahyah Kamaruzzaman, Kheng Seang Lim PII: DOI: Reference:

S1388-2457(19)31362-8 https://doi.org/10.1016/j.clinph.2019.11.058 CLINPH 2009073

To appear in:

Clinical Neurophysiology

Accepted Date:

30 November 2019

Please cite this article as: Ashfak Habib, M., Ibrahim, F., Mohktar, M.S., Bahyah Kamaruzzaman, S., Seang Lim, K., Recursive independent component analysis-decomposition of ictal EEG to select the best ictal component for EEG source imaging, Clinical Neurophysiology (2019), doi: https://doi.org/10.1016/j.clinph.2019.11.058

This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

© 2019 Published by Elsevier B.V. on behalf of International Federation of Clinical Neurophysiology.

Recursive independent component analysis-decomposition of ictal EEG to select the best ictal component for EEG source imaging

Mohammad Ashfak Habib1,2,3, Fatimah Ibrahim1,2,*, Mas S. Mohktar1,2, Shahrul Bahyah Kamaruzzaman2,4 and Kheng Seang Lim2,4

1

Department of Biomedical Engineering, Faculty of Engineering, University of Malaya,

50603 Kuala Lumpur, Malaysia 2

Centre for Innovation in Medical Engineering (CIME), Faculty of Engineering, University

of Malaya, 50603 Kuala Lumpur, Malaysia 3

Department of Computer Science & Engineering, Chittagong University of Engineering &

Technology, Chittagong-4349, Bangladesh 4

Department of Medicine, Faculty of Medicine, University of Malaya, 50603 Kuala Lumpur,

Malaysia

*Corresponding

Author:

Professor Ir. Dr. Fatimah Ibrahim Head (Centre for Innovation in Medical Engineering) Department of Biomedical Engineering, Faculty of Engineering, University of Malaya 50603 Kuala Lumpur, Malaysia. Tel/Fax: +60379676878 Email: [email protected]

1

ABSTRACT Objective: This study aimed to present a new ictal component selection technique, named as recursive ICA-decomposition for ictal component selection (RIDICS), for potential application in epileptogenic zone localization. Methods: The proposed technique decomposes ictal EEG recursively, eliminates a few unwanted components in every recursive cycle, and finally selects the most significant ictal component. Back-projected EEG, regenerated from that component, was used for source estimation. Fifty sets of simulated EEGs and 24 seizures in 8 patients were analyzed. Dipole sources of simulated-EEGs were compared with a known dipole location whereas epileptogenic zones of the seizures were compared with their corresponding sites of successful surgery. The RIDICS technique was compared with a conventional technique. Results: The RIDICS technique estimated the dipole sources at an average distance of 12.86 mm from the original dipole location, shorter than the distances obtained using the conventional technique. Epileptogenic zones of the patients, determined by the RIDICS technique, were highly concordant with the sites of surgery with a concordance rate of 83.33%. Conclusions: Results show that the RIDICS technique can be a promising quantitative technique for ictal component selection. Significance: Properly selected ictal component gives good approximation of epileptogenic zone, which eventually leads to successful epilepsy surgery.

Keywords: Focal epilepsy, Epileptogenic zone, Ictal EEG, Source localization, Surgical, Independent component analysis.

2

HIGHLIGHTS: 

Recursive decomposition of ictal EEG and elimination of unwanted independent components (ICs) help to find the best ictal IC.



Quantitative feature is useful to reduce the dependency on visual inspection for ictal IC selection.



Cortical source of ictal rhythm can be estimated more concordantly from the best ictal component.

3

1. Introduction Epilepsy is one of the common chronic neurological disorders, which affect people of all ages. Around 65 million people worldwide have epilepsy (Ngugi et al. , 2010, Thurman et al. , 2011, Moshe et al. , 2015), and approximately 4.5% of all patients with epilepsy are potential candidates for surgical resection (Engel, 1993, Rosenow et al. , 2001). The success of focal epilepsy surgery depends largely on proper localization of epileptogenic zone. Commonly used diagnostic techniques for epileptogenic zone localization are: magnetic resonance imaging (MRI), functional MRI (fMRI), positron emission tomography (PET), single-photon emission computed tomography (SPECT), scalp electroencephalography (EEG), and invasive EEG (iEEG). Among all these techniques scalp EEG is most commonly used for epilepsy evaluation because of its high temporal resolution, suitability for long-term monitoring and low cost. Traditionally, EEG recordings are inspected visually. Simple visual inspection is not adequate for precise localization of epileptogenic zone because of the volume conduction effects of a series of barriers: scalp, skull, cerebrospinal fluid (CSF), etc. (Nunez et al. , 2006, He et al. , 2013). Electroencephalography source imaging (ESI) is a comparatively new model-based computational technique, which can consider volume conduction effects while localizing the possible cortical sources of EEG (Kaiboriboon et al. , 2012). Both ictal and interictal EEG can be used independently for ESI-based cortical source estimation. Ictal EEG is believed to be more reliable than interictal EEG for localizing the epileptogenic zone (Jayakar et al. , 1991), but ictal-EEG-based source imaging (i.e. ictal ESI) is comparatively difficult mainly because of the low Signal to Noise Ratio (SNR) of ictal EEG (Kaiboriboon et al. , 2012, Yang et al. , 2011). Ictal EEG measures cortical seizure discharges superposed with various artifacts, external noises and other background brain oscillations (Urrestarazu et al. , 2004).

4

The major obstacle in ictal ESI is to eliminate these unwanted portions of the ictal EEG before using it for source estimation. Up to now, only few studies (Assaf et al. , 1997, 1999, Lantz et al. , 1999, Merlet et al. , 2001, Boon et al. , 2002, Nam et al. , 2002, Beniczky et al. , 2006, Iriarte et al. , 2006, Leal et al. , 2006, Ding et al. , 2007, Jung et al. , 2009, Holmes et al. , 2010, Koessler et al. , 2010, Yang et al. , 2011, Lu et al. , 2012a, Lu et al. , 2012b, Kovac et al. , 2014, Habib et al. , 2016) have used ictal EEG for cortical source estimation through ESI. Almost all of those studies, except two (Nam et al. , 2002, Iriarte et al. , 2006), used digital bandpass filter in the EEG acquisition step and/or any other step(s) of their proposed technique to improve the SNR by minimizing high frequency noises, DC linear trends, movement and muscle artifacts. Even though the effective uses of filters are well-verified, they are inadequate for noise elimination because of their inability to discriminate between artifact and brain waves (Nam et al. , 2002, Iriarte et al. , 2003). Averaging of the ictal-events is another commonly used technique for the improvement of SNR (Assaf et al. , 1997, 1999, Merlet et al. , 2001, Beniczky et al. , 2006, Kovac et al. , 2014, Habib et al. , 2016), but such averaging can eliminate valuable cortical information and thus can lead to inaccuracies in ESI calculations (Chitoku et al. , 2003). A variety of procedures were proposed to correct ocular artifacts, but those were not very popular in ictal ESI studies due to their failure to deal with other noises. Rather than use those procedures, sometimes visible artifacts are screened through visual inspection (Holmes et al. , 2010) in addition to bandpass filters. Another well-known technique for signal and noise subspaces separation is the principal component analysis (PCA). None of the ictal ESI studies applied PCA because of its unrealistic assumption of orthogonality between neural activities and artifacts (Delorme et al. , 2004, Romero et al. , 2008).

5

Independent component analysis (ICA) is comparatively a more powerful statistical technique for decomposing a mixed dataset into statistically independent components. In the context of EEG signals, ICA can decompose each set of EEG into a series of spatially fixed and temporally independent components (Yang et al. , 2011). The ICA-based ictal ESI studies (Nam et al. , 2002, Iriarte et al. , 2006, Leal et al. , 2006, Jung et al. , 2009, Yang et al. , 2011, Lu et al. , 2012a) used ICA for EEG denoising and thus extracted the EEG signals of interest. Independent components (ICs) from extracerebral origin (e.g. muscle artifacts, eye movements, and power line noise, etc.) were removed for denoising (McMenamin et al. , 2010a), whereas the components that display ictal nature (i.e. ictal components) were extracted for epileptogenic zone localization. Proper identification of the ictal components is the major challenge for ICA-based ictal ESI techniques. Iriarte et al. (2006) depended only on visual inspection of the component time-courses and their corresponding reconstructed EEGs for selecting the ictal component(s). In addition to visual inspection, Nam et al. (2002) reviewed autocorrelogram for ictal component selection. They also conducted a quantitative analysis and used power spectral density (PSD) for ictal component selection. Components having higher proportion of power in the theta-band (exactly 2-10 Hz) were selected for source estimation. Jung et al. (2009) initially excluded the ICs of extracerebral origin by visually inspecting the scalp voltage topography and activation spectra (i.e. PSD). They also excluded the ICs with a residual variance of more than 20% before selecting the ictal components through visual inspection. Three other ICA-based ictal ESI studies (Leal et al. , 2006, Yang et al. , 2011, Lu et al. , 2012a) utilized time-frequency representation (TFR) in addition to the above mentioned tools and techniques for selecting the ictal components. They also used bootstrap statistical method and surrogate data, with further assistance of visual inspection, to evaluate the significance of the spectral changes in a TFR map or to assess the correlation between two different TFR maps.

6

The above review indicates that the existing ictal component selection techniques are highly dependent on visual inspection and are therefore susceptible to experimenter’s judgment or bias. In this context, the aim of this study was to develop a quantitative technique to select the best ictal component (BIC) for EEG source imaging with potential application in presurgical evaluation of refractory focal epilepsy. A recursive approach and a unique quantitative feature of ictal EEG decomposition were used for the BIC selection. The cortical sources of both simulated and real ictal EEG were estimated from their corresponding BICs. Dependency on manual visual inspection was evaluated against a conventional visual inspection-based (VI) technique. The estimated dipole sources of the simulated ictal EEGs were validated using their known dipole source, whereas the known sites of successful surgeries were used for the validation of the epileptogenic foci of the seizures. Since RIDICS technique uses a recursive approach, its practicality was evaluated by examining its execution time.

2. Materials and methods The RIDICS technique aims to find the BIC, which has the best correspondence with the rhythmic ictal discharges of the epileptogenic cortical zone. In order to do that, the RIDICS technique uses a unique quantitative feature of ictal EEG decomposition in a number of recursive steps. The background and the mathematical formulations of the quantitative feature are presented in the next subsection. 2.1 A quantitative feature of ictal EEG decomposition An ICA algorithm can decompose the recorded ictal EEG signals (X) into a series of IC time-courses (U) according to the following relation: 𝑿 = 𝑨𝑼

(1)

7

where A is the “mixing matrix”. Its columns contain the relative weights with which the corresponding ICs are projected to the EEG electrodes. If an n-channel EEG record (having T samples) is considered as the mixture of m independent components, then X, A and U represent n×T, n×m and m×T matrices respectively. Each IC (e.g. the pth IC) can be backprojected into an n-channel EEG (𝑿𝑝 ). Therefore, 𝑿𝑝 is a portion of 𝑿 and 𝑿 is the sum of all such back-projected EEGs:

𝑿=∑

𝑚

𝑿𝑝

(2)

𝑝=1

According to Eq. (1) the matrix 𝑿𝑝 is the product of 𝒂𝑝 (the pth column of A) and 𝒖𝑝 (the pth row of U) i.e. 𝑿𝑝 = 𝒂𝑝 𝒖𝑝

(3)

The Eq. (3) implies that the EEG signal of the qth electrode in 𝑿𝑝 can be expressed as the product of the qth element of ap and the IC time-course up: 𝒙𝑝𝑞 = 𝑎𝑞𝑝 𝒖𝑝

(4)

The strength of the ictal rhythm frequency f in 𝒙𝑝𝑞 can be estimated, using discrete Fourier transform, as follows: ℱ𝑘+1 (𝒙𝑝𝑞 ) = ∑

𝑇−1 𝑡=0

𝑝 𝑥𝑞(𝑡+1) 𝑒 −2𝜋𝑖𝑘𝑡/𝑇

(5)

where ℱ𝑘+1 (𝒙𝑝𝑞 ) represents the (k+1)th discrete value, which holds the strength of f in the signal 𝒙𝑝𝑞 . The symbols k and t are the indices of the discrete values in frequency domain and time domain respectively and both vary from 0 to T-1. If the sampling frequency of the signal be 𝑓𝑠 , then k can be obtained from the relation, 𝑘 = 𝑓𝑇/𝑓𝑠 . Moreover, i represents the imaginary unit (√−1). Substituting Eq. (4) into (5) yields the following expression: 8

𝑝

−2𝜋𝑖𝑘𝑡/𝑇 ℱ𝑘+1 (𝒙𝑞 ) = 𝑎𝑞𝑝 ∑𝑇−1 𝑡=0 𝑢𝑝(𝑡+1) 𝑒

∴ ℱ𝑘+1 (𝒙𝑝𝑞 ) = 𝑎𝑞𝑝 ℱ𝑘+1 (𝒖𝑝 )

(6)

If ℱ𝑘+1 (𝑿𝑝 ) be the column matrix that holds all the discrete values of ℱ𝑘+1 (𝒙𝑝𝑞 ) for all the electrodes (𝑞 = 1, 2, … , 𝑛) in 𝑿𝑝 , then according to Eq. (3) and (6): ℱ𝑘+1 (𝑿𝑝 ) = 𝒂𝑝 ℱ𝑘+1 (𝒖𝑝 )

(7)

Similarly, the strength of f in all the electrodes in 𝑿 can be represented by ℱ𝑘+1 (𝑿). Since the Fourier transform is a linear transform, ℱ𝑘+1 (𝑿) can be expressed, on the basis of Eq. (2), as follows:

ℱ𝑘+1 (𝑿) = ∑

𝑚

ℱ𝑘+1 (𝑿𝒑 )

𝑝=1

(8)

If the pth independent component be the BIC, then 𝑿𝒑 contain the highest amount of ffrequency activities with respect to the other portions of 𝑿. If the strength of f is negligible in all the portions of 𝑿, other than 𝑿𝒑 , then Eq. (8) can be rewritten as: ℱ𝑘+1 (𝑿) ≈ ℱ𝑘+1 (𝑿𝑝 )

(9)

Substituting Eq. (7) into (9) implies: ℱ𝑘+1 (𝑿) ≈ 𝒂𝑝 ℱ𝑘+1 (𝒖𝑝 )

(10)

Equation (9) and Eq. (10) show that the strength of ictal rhythm in the BIC-projected EEG has the strongest correlation with that of the recorded ictal EEG. This quantitative feature has been used in each recursive step of the RIDICS technique. The column matrices ℱ𝑘+1 (𝑿) and ℱ𝑘+1 (𝑿𝑝 ) have been termed as input rhythm magnitude (Y) and back-projected rhythm magnitude (Yp) respectively for the ease of explanation. Moreover, the term

9

component rhythm magnitude (yp) has been used for representing the discrete value ℱ𝑘+1 (𝒖𝑝 ). All the bold-face symbols including Y and Yp represent arrays while the normal font symbols represent single value. 2.2 Recursive ICA-decomposition for ictal component selection An EEG is a mixture of electrical signals derived from unlimited number of independent sources (both cerebral and extracerebral), but an ICA algorithm can decompose it into a limited number of ICs, which are mainly associated with the sources of stronger signals. Each weaker component of the mixture is either summed into a single IC or distributed into multiple ICs. Therefore, the ICA-decomposition of high-noise ictal EEG may distribute the ictal signal, generated from a single cortical region, into more than one IC (Onton et al. , 2006). The RIDICS technique intends to accumulate the majority (or all) of such distributed portions of an ictal signal into a single IC (i.e. the BIC) by using recursive ICA-decomposition and noise reduction in each recursive cycle. The major steps of the RIDICS technique are depicted in the simplified flowchart of Fig. 1. The recursion base case starts with the EEG epoch that contains mostly (or totally) the post-onset EEG. If pre-onset EEG is included, its length should not exceed one-third of the total epoch-length. Multiple segments of pre-onset and post-onset EEG can be concatenated for avoiding excessive artifacts in the input EEG epoch. The input EEG epoch passes through a wideband bandpass (1-45 Hz) filter, which removes the signals from the extracerebral origins. Before decomposing the filtered EEG epoch into ICs, the input rhythm magnitude Y is computed. The electrode with the highest corresponding value in Y is identified and defined as the most significant electrode (MSE). The ICA-decomposition step decomposes the filtered EEG epoch into n number of ICs, where n is also the number of electrodes in the input EEG. In the next step, the back-projected rhythm magnitude Yp is computed for each of all the ICs (i.e. for p = 1, 2, …, n) by using the expression of Eq. (6). The zero-lag cross-

10

correlation (Z) between Y and each of Yp matrices are estimated in the succeeding step. At the end of the base case the IC with the highest Z value (i.e. Zmax) is tentatively identified as the BIC. This base case is repeated for all other recursive cycles. If there is no improvement in the Zmax value with respect to that of the previous recursive cycle, the decomposition results are discarded. Otherwise, the results are stored and used for further processing. Both conditions lead to an IC elimination step that eliminates one unwanted IC from the stored decomposition results and regenerates back-projected EEG. The MSE of the regenerated EEG is passed through a narrowband bandpass filter. The resultant EEG is used as the input dataset for the next recursive cycle until any of the stopping criteria defined below is satisfied. The recursion stops when any of the following two conditions (stopping criteria) is met: (1) the Zmax value doesn’t increase for (n-1) consecutive recursive cycles, where n is the total number of components. (2) the value of Zmax exceeds the 75% of the autocorrelation value of Y. Whenever a stopping criterion is met, the regenerated dataset, for the next recursive cycle, is discarded and the last stored decomposition results are retrieved with their corresponding Z values. Finally, the IC with the highest Z value is considered as the BIC and an EEG is regenerated from that IC for the epileptogenic zone estimation. In the IC elimination step, the unwanted ICs are selected based on their Z values. The IC with the lowest Z value is discarded in every recursive cycle that has increment in the Zmax value. If the Zmax value decreases for r consecutive recursive cycles (where r = 1, 2, …, n-1), then in every rth cycle r components with comparatively lower Z values are discarded. In order to improve the degraded rank of the regenerated dataset matrix only the signal of MSE is filtered instead of conventional all channel filtration. The bandwidth of the narrowband bandpass filter is chosen to be ~5Hz (f-2 Hz to f+2 Hz) around the ictal rhythm frequency f.

11

The RIDICS technique was implemented in MATLAB R2013a (MathWorks, Inc., Natick, Massachusetts, USA) and confirmed with both simulated and real EEG data. 2.3 Simulated ictal-EEG generation Simulated ictal EEG data were generated by using a fixed dipole in a four-shell spherical head model consisting of scalp, skull, CSF and brain (Berg et al. , 1994). Outer radii of the compartments were chosen to be 85, 79, 72, and 71 mm respectively while their conductivities were chosen to be 0.33, 0.0042, 1, and 0.33 mho/m respectively (MacKinnon et al. , 2000, Tenke et al. , 2015). A dipole source (Fig. 2a) with 250 nAm moment and rhythmic sinusoidal 6 Hz ictal pattern (Fig. 2b), which is common in medial temporal lobe epilepsy (Nam et al. , 2002, De Vos et al. , 2007), was placed in the right temporal lobe. In Cartesian coordinates, the exact location of the dipole was: x = 58.65, y = 16.575 and z = – 3.91. Scalp potentials were computed at 33 electrode-locations, arranged according to 10-10 system, by multiplying the source wave of the dipole by the amplification factors obtained from the forward solution. The time-courses of the scalp potentials were stored with a sampling frequency of 500 Hz. In case of n-channel data, at least α×n2 data points are needed to obtain reliable ICAdecompositions, where α is the multiplication factor with a recommended value of 20 (Onton et al. , 2006, Groppe et al. , 2009). Therefore, 22,000 data points were generated for each set of simulated EEG. Similar to other ictal ESI studies (Jung et al. , 2009, Yang et al. , 2011, Lu et al. , 2012a), each 44 seconds EEG dataset was consisted of both pre-onset (12 seconds) and post-onset (32 seconds) EEG. The initial simulated dataset E was noise-free and therefore the pre-onset background activities were represented by zero potential. Fifty other sets of EEG were generated in a similar fashion with the same dipole source but adding different levels of coherent noises with E. The noise was coherent in the sense that there was quite a high correlation between signal amplitudes from electrodes that

12

were close to each other. The level of noise (root mean square value) was varied from 1 µV to 50 µV with a regular interval of 1 µV for generating the simulated datasets E1, E2 … and E50 respectively. The EEG responses around the ictal onset (1 second before and 2 seconds after the onset), recorded on the MSEs (i.e. the T8 channels), are shown in Fig. 2b. The MSE waveforms from three different datasets (E10, E25 and E50) with three different levels (low, medium, and high) of noise, along with their filtered (1 to 11 Hz bandpass) versions are included in this figure. It illustrates the changes of E with different levels of added noises. All the signals of Fig. 2b are in the same scale. Simulated EEG datasets were generated by using the BESA Simulator Version 1.1.0 (BESA GmbH, Germany). Since the simulator cannot produce more than 2000 samples of ictal EEG at a time without adding any interval, eight different epochs of simulated ictal EEG were concatenated (Yang et al. , 2011, Lu et al. , 2012a) for producing the 32 seconds of post-onset EEG. The SNR was calculated for the MSE of each generated ictal EEG dataset. The MSE of E was considered as the desired signal S for SNR calculation. The SNR of the jth dataset Ej was computed with the equation adopted from (Raz et al. , 1988). ∑𝑇𝑡=1(𝑆(𝑡))2 𝑆𝑁𝑅𝑗 = 𝑇 ∑𝑡=1(𝑆𝑗 (𝑡) − 𝑆(𝑡))2 where Sj is the MSE of Ej and j = 1, 2, …, 50. The data point 𝑡 = 1, 2, … , 𝑇 and the total data points T = 22,000 for all the simulated datasets. The SNR values that were computed for all the simulated datasets are illustrated in Fig. 2c. 2.4 Evaluation using simulated data An extension of the infomax algorithm of Bell et al. (1995), as implemented in the runica function of the EEGLAB toolbox (Delorme et al. , 2004), was used for all the ICAdecompositions of this study. Each n-channel EEG dataset was decomposed into n ICs. Therefore, 33 ICs were obtained from the ICA-decomposition of each simulated dataset. The 13

RIDICS technique identified one of those ICs as the BIC after one or more ICAdecomposition steps. The RIDICS technique is explained in the section 2.2. All the required Fourier transforms were applied on the 32 seconds post-onset EEG. The 12 seconds pre-onset EEG was included in the input datasets in order to fulfill the requirement of the VI technique that was used for the evaluation. The VI technique was adopted mainly from the study of Iriarte et al. (2006). Following that study, the IC time-courses and their back projected EEGs were inspected visually for identifying the components with likely ictal activities. Moreover, scalp voltage topography (Nam et al. , 2002, Jung et al. , 2009), PSD curve (Nam et al. , 2002, Jung et al. , 2009), and TFR map (Leal et al. , 2006, Yang et al. , 2011, Lu et al. , 2012a) were inspected visually for better identification of the ictal components. Instead of selecting a single ictal component, as done in the RIDICS technique, the existing studies (Nam et al. , 2002, Iriarte et al. , 2006, Jung et al. , 2009) selected one or multiple (usually 1 to 5) ictal components depending on the outcomes of the ICA-decompositions. Therefore, the VI technique identified one or three ictal components for each input EEG dataset. The selected ictal components as well as the BICs, obtained from the RIDICS technique, were later back-projected into EEG datasets. Thus, for each of the two techniques, 50 backprojected datasets (EB1, EB2 … EB50) were regenerated from their corresponding simulated datasets (E1, E2 … E50) and were used for evaluation. The performance of both the techniques were quantitatively evaluated against two references, namely the noiseless simulated ictal EEG dataset (E) and the known location of the current dipole (D) that caused the simulated ictal activities. Zero-lag cross-correlations between each of the back-projected datasets and E were computed as a measure of similarity. Higher zero-lag cross-correlation was considered as the better identification of the ictal component(s). Besides, an equivalent current dipole source and its Euclidean distance from D were estimated for each of the regenerated EEGs. Shorter distance indicates better estimation

14

of EEG source and therefore was considered as better identification of ictal component(s). The DIPFIT2 plug-in of EEGLAB toolbox was used for dipole source estimation from single-component-projected EEGs. An EEG that is back-projected from multiple ICs has different dipole sources at different time points (He et al. , 2018). Therefore, for each of such EEG datasets, the DIPFIT2 plug-in was applied on a set of event related potentials (ERPs) for estimating the equivalent current dipole source. The electrode potentials were extracted from an average map of 5 to 10 ictal events (each of 150 milliseconds) and from the time point where the peak global field power (GFP) occurred (Habib et al. , 2016). The 4 shell spherical head model that was used for simulated-EEG generation was also used for all the dipole estimations. 2.5 Patients and data acquisition This study was conducted with the use of anonymized retrospective data, and the study protocols were approved by the University of Malaya Research Ethics Committee, Malaysia. Thirty-four patients with drug-resistant focal epilepsy who underwent presurgical evaluation including MRI brain and long-term video-EEG monitoring were selected primarily from the database of the University of Malaya Medical Centre. Eight patients who underwent surgical resections with Engel’s class-I surgical outcome were included in this study. The demographics and the diagnosis of the selected patients are summarized in Table 1. Long-term video-EEGs were recorded with standard EEG setup of 21 electrodes according to the international 10-20 system with additional 2 electrodes at T1 and T2, using NicoletOne EEG/LTM system (Natus Medical Inc., Pleasanton, California, USA). Electrode impedances were kept below 10 kΩ, and sampling rate was set to 500 Hz (250 Hz for patient6 and patient-8). Seizures (range 2 - 4) were recorded in all patients and all the seizures were included in this study.

15

Anatomical details of every patient’s head were obtained from their corresponding MRI brain scans. Those details were used for generating patient-specific head models. Structural MRI scans were performed with a Signa HDxt 3.0T scanner (General Electric Healthcare, Wauwatosa, Wisconsin, USA). 2.6 Validation using real data An EEG epoch was extracted from each of the seizures for ICA-decomposition. Similar to the simulated datasets, each epoch consisted of both pre-onset and post-onset EEG and at least 20×n2 data points. Most of the EEG epochs were extracted from continuous time points. Other EEG epochs were obtained by concatenating multiple EEG segments of multiple time segments. It was to avoid excessive movement artifacts, abnormal signals from disconnected channels, etc. The ictal rhythm frequency f of an EEG epoch was identified from the vicinity of the ictal onset time. Since, both frequency and cortical source of an ictal rhythm may change with time, careful selection of f is very important for the success of the RIDICS technique. In this study, f was identified by analyzing the PSD curves of a 2 to 12 seconds segment of each EEG epoch. The EEG segments that were closer to the ictal onset and contained comparatively stable ictal rhythms with fewer artifacts were chosen. An example of such epoch selection is shown in Fig. 3. It is for the first seizure of patient 3. The RIDICS technique identified a BIC for each seizure by using the same procedure as described for the simulated datasets. An EEG dataset was regenerated for each of the seizures by backprojecting the corresponding BIC. No other ictal component selection technique was implemented on real datasets for comparing the results with those of the RIDICS technique. It is due to the fact that none of the existing approaches are suited for the strong assumption of selecting a single ictal component (i.e. the BIC). An epileptogenic zone was estimated from each of those regenerated EEG datasets and by following the procedure explained in (Habib et al. , 2016). Subject specific realistic

16

head models were generated from the corresponding MRI scans by using the BrainSuite software tool (Shattuck et al. , 2002). The source model, channel locations and the forward solution techniques were used in the same way as described in (Habib et al. , 2016). The linear distributed inverse model sLORETA was used for the inverse solutions. Epoching and averaging were not necessary for the single-component-projected EEGs, because the cortical source of such EEG is temporally independent and spatially fixed. Therefore, estimated locations of EEG sources were fixed with respect to time. All the estimated epileptogenic zones were compared with their corresponding sites of cortical surgery. Each ESI result was considered to be concordant with the surgery-site if both ESI focus and surgery site were located in the same lobe. Concordance rate was computed as the percentage of concordant results in the total number of considered results. Moreover an ESI result was considered as lateralized if both the estimated zone and the surgery-site were found in the same hemisphere. 2.7 Further validation and execution time analysis The RIDICS-outcomes with the simulated datasets were also visually inspected for further investigating their correctness. It was verified that the most appropriate IC with the most significant ictal rhythm were selected as the BIC. The time-courses of the BICs, their scalp voltage topographic maps, PSD curves and TFR maps were visually inspected for the purpose. The RIDICS technique recursively decomposed each EEG dataset into a number of ICs and discarded one or few of them in every recursive cycle. The objective was to eliminate the non-ictal components from the EEG datasets and thus facilitate the identification of the BICs. The gradual improvements of the EEG datasets and their decomposed ICs due to the recursive IC eliminations were investigated visually. The ICA-decomposition process itself is a resource and time consuming process and RIDICS technique uses ICA-decomposition recursively. Therefore, the practicality of the

17

RIDICS technique is highly dependent on its execution time. One may claim that the RIDICS technique consumes unacceptable amount of time. In order to encounter/oppose such claim, the required time for executing the RIDICS technique was measured. Since the ICAdecomposition step consumes the highest amount of time with respect to the other steps, the CPU execution time was recorded for every ICA-decomposition. In comparison with the real datasets, the simulated datasets contained higher number of data-points and thus consumed more time for ICA-decomposition. Therefore the practicality of the RIDICS technique was assessed on the basis of its execution times for simulated datasets.

3. Results Fig. 4 illustrates the Euclidean distances between D and each of the equivalent current dipoles that were estimated from the simulated EEG datasets. Shorter distance usually represents better estimation of dipole source and thus indicates better identification of ictal component(s). Therefore, Fig. 4 indicates that the RIDICS technique gives better ESI results in comparison with the VI technique. The dipoles obtained from the RIDICS-based ESI were closer to D for most of the (35 out of 50) datasets. Fig. 4 also indicates that the RIDICS technique works well for both low-noise and high-noise datasets. The outcomes of the RIDICS technique were as good as that of the VI technique for low-noise datasets. Slight differences in the dipole distance were due the nondeterministic nature of ICA-decomposition and for using different ICA-decomposition-outcomes for different techniques. Both the techniques identified a single ictal component for each of the low noise datasets. Moreover, the RIDICS technique did not require more than one recursive cycle to identify the BICs. These similarities with the VI technique indicate that the RIDICS technique can identify the BIC properly from low-noise dataset by using its unique quantitative feature (as described in section 2.1). In case of high-noise datasets, Fig. 4 illustrates that the outcomes of the RIDICS 18

technique were significantly better than those of the VI technique. Instead of single component, the VI technique identified multiple ictal components for most of the (E22, E26, E32, E33, E36, E37, E38, E40, E41, E42, E43, E45, E46, E47, E48, E49, and E50) high-noise datasets. Single ictal component was extracted from each of the remaining high-noise datasets and was used for dipole source estimation. The average of all the dipole distances, obtained for both low-noise and high-noise datasets, was calculated. Those average dipole distances were 12.86 mm (minimum 0.29 mm, maximum 45.08 mm), and 24.31 mm (minimum 0.39 mm, maximum 120.53 mm) for the RIDICS and the VI technique respectively. Shorter average dipole distance encourages the use of the RIDICS technique for the ESI solutions. Fig. 5 illustrates the similarities between the noiseless EEG dataset E and each of the regenerated datasets that were back-projected from the selected ictal component(s). The cross-correlation value is considered as the similarity measuring metric. Higher crosscorrelation between two signals was considered as the reflection of their higher similarity. Since E contained no noise, higher level of similarity to E indicates lower level of noise in the regenerated EEG. Source estimation from such noise-reduced EEG was expected to give better source localization. This expectation is supported by the dipole distances that are shown in Fig. 4. Moreover, higher cross correlation value also indicates better selection of ictal-component(s). According to Fig. 5, the RIDICS-technique gives highly correlated backprojected EEGs for both low-noise and high-noise datasets. On the contrary, in case of the VI technique, the average level of cross-correlations was comparatively low for the high-noise datasets. Out of 50 datasets, 31 produced higher cross-correlation values for using the RIDICS technique. Moreover the RIDICS technique maintained the level of cross-correlation more consistently than that of the VI technique. These findings are also in line with the findings of Fig. 4.

19

In order to validate with real ictal EEG, 24 seizures of 8 subjects were analyzed. The RIDICS technique identified a BIC for each of the seizures and a set of back-projected EEG was regenerated from the BIC. An epileptogenic zone was estimated for every set of regenerated EEGs. Table 2 presents those estimated cortical sources along with their corresponding sites of surgery. The concordance, as well as the discordance, of the estimated epileptogenic zones and their corresponding surgery-sites were evaluated and summarized in Table 3. It shows that the estimated epileptogenic zones were found concordant for 20 out of 24 seizures (concordance rate 83.33%). The cortical sources of three seizures (Seizure 1, Seizure 2 and Seizure 3) in patient 4 were not found concordant. Those cortical sources were found in the occipital, temporal or temporo-occipital region instead of the parieto-occipital region. Although those cortical sources were close to their corresponding site of surgery, those were not considered as concordant while calculating the concordance rate. Moreover, the estimated epileptogenic zone for the seizure 4 of patient 3 was also found discordant. Although some results were found discordant, the overall concordance rate of the RIDICS technique is admissible compared with that of the existing relevant techniques. Simulation-based evaluation results (Fig. 4 and Fig. 5) showed that the RIDICS technique performed well for both low-noise datasets (E1 to E25) and high-noise datasets (E26 to E50). Further investigations were carried out visually for verifying the correctness of those results. Fig. 6 and Fig. 7 illustrate the IC time-courses, topographic maps, and TFR maps for two different datasets, E23 and E41, respectively. These two datasets were selected purposefully. The first dataset (E23) corresponds to the worst performance of the RIDICS technique under low-noise condition. Amongst all the low-noise datasets, the highest dipole distance was obtained for E23. On the contrary, the VI technique provided better dipole distance, as well as better cross-correlation value, for this dataset. The other dataset (E41) was chosen because the RIDICS technique performed significantly better compared to the VI

20

technique for this dataset despite its high noise content. The time-courses with magenta color represent the selected ictal components including the BICs. The topographic maps of those selected ICs are presented on top of the corresponding IC time-courses. The TFR maps are included in Fig. 7 to justify the cause of selecting the ICs having no clearly-visible ictal rhythm. The IC time-courses that were used by the RIDICS and the VI technique are shown in Fig. 6a and Fig. 6b respectively. The third component of each set of these time-courses was considered as the ictal component. Visible ictal rhythm and definable topographic map support the correctness of these ictal components. Poor ictal rhythm of Fig. 6a, compared to that of Fig. 6b, indicates the reason behind the slightly inferior results of the RIDICS technique for E23. This dissimilarity of ictal rhythm occurred due to the nondeterministic nature of ICA decomposition. The RIDICS and the VI technique used two different decomposition-outcomes, which have different levels of accuracy in separating the ictal components. It is visible (in Fig.6) that the ictal component of Fig. 6a was separated less accurately compared to that of Fig. 6b. Longer dipole distance and lower cross-correlation values were obtained due to this poorly separated ictal component. Since the RIDICS technique identified the correct ictal component from poorly decomposed ICs (shown in Fig. 6a), it has the potential to work successfully with better decomposition-outcomes (shown in Fig. 6b) as well. Moreover, the success of any ICA-based ESI technique depends on proper ICA-decomposition. The IC time-courses of Fig. 7a illustrate that the RIDICS technique identified the BIC correctly for the high-noise dataset E41 and the topographic map also supports the correctness of the BIC selection. The ictal rhythm is clearly visible in the RIDICS-selected BIC, but the rhythm is not that clear in the ICs that were selected by the VI technique. Therefore it was difficult to become sure about the correctness of the VI-selected ICs by inspecting the IC

21

time-courses only. It was also difficult to make any positive decision about those ICs from the topographic maps. The selection was made based on the TFR maps. Strong existence of 6 Hz ictal rhythm frequency in the TFR map was considered for selecting the ICs. The TFR maps of three selected ICs (IC-5, IC-15 and IC-24) are illustrated in Fig. 7c. The time-courses of Fig. 7a were obtained after four recursive cycles of the RIDICS technique and in every cycle a small amount of noise was eliminated. Therefore after four recursive cycles the ictal rhythm became clearer in the RIDICS selected BIC. On the contrary, the VI technique selected the ICs from the initial decomposition outcomes. As a result, the ictal rhythm was contaminated by other signals and was hardly visible in the IC time-courses. These results provide a strong justification for the better performance of the RIDICS technique with the higher levels of noises in the datasets. The RIDICS technique eliminates noise in every recursive cycle. An example of such noise reduction is shown in Fig. 8. The epoch that was extracted from the first seizure of patient 3 (as shown in Fig. 3) is illustrated in Fig. 8a. The selected epoch of Fig. 3 is 100 times magnified in Fig. 8a. The RIDICS technique eliminated noise from this epoch in 21 recursive cycles and the outcome is shown in Fig. 8b. This EEG was regenerated BIC and was used for localizing the corresponding epileptogenic focus of Table 2. Both the epochs of Fig. 8 are in the same scale. The execution times that were recorded for the RIDICS technique are presented in Fig. 9. It discards the possible claim of excessive time requirement for the recursive operations of the RIDICS technique. The average time for analyzing each patient's dataset was 801.65 seconds (i.e. 13.36 minutes), which is quite acceptable with respect to other diagnostic modalities (e.g. SPECT). The average time required for each ICA decomposition step was 108.92 seconds while the average number of recursive-cycle required for analyzing each dataset was 7.36. The worst case (dataset) required 33 recursive cycles and took 4067.11

22

seconds (67.79 min) and that time is also within the acceptable limit. All the simulations and time computations were performed in a laptop computer with Intel Core i7-4700MQ processor, 32GB DDR3L RAM, and Windows 8.1 operating system.

4. Discussion Proper selection of IC provides better estimation of epileptogenic zone, but the IC selection step of the existing ICA-based ESI solutions are highly dependent on the manual visual inspections. In this study, a quantitative IC selection technique (RIDICS) has been proposed not only for reducing the high dependency on visual inspection but also for extracting less noisy ictal components for better estimation of epileptogenic zone. The RIDICS technique has been evaluated for both simulated and real ictal EEG datasets. The results demonstrate the practicality and usefulness of the RIDICS technique in the context of ictal ESI for presurgical epilepsy evaluation. The evaluation study with the simulated datasets showed that the RIDICS technique performed better for both the performance measuring metrics and for both (low and high) levels of noise. In comparison with the VI technique, RIDICS technique provides more accurate ictal component for better estimation of dipole source. The back-projected EEGs, obtained from the BICs, showed better similarities with the noise-free EEG. Actually, these two performance measuring metrics are correlated and they support each other. Improved results for both the measuring metrics support the correctness of the BICs, which were obtained by using the RIDICS technique. Although the VI technique performed satisfactorily with low-noise datasets, it could not identify the appropriate ictal components from highnoise datasets. It is due to the limitations of the manual visual inspections. On the other hand the main reason behind the better performance of the RIDICS technique is the recursive rejection of noise (or unwanted IC) in its every recursive cycle. Other studies (Vorobyov et 23

al. , 2002, Onton et al. , 2006, Delorme et al. , 2007, McMenamin et al. , 2010b) also support the usefulness of such IC rejection for EEG noise elimination. The RIDICS technique selects a single ictal component (i.e. the BIC) from each EEG dataset but the existing techniques, including the VI technique, do not always follow this strong assumption. Therefore, the comparison results between the RIDICS and the VI technique may not always indicate the superiority of one technique over the other. Those results actually support the usability of the RIDICS technique. One major motivation behind the comparative analysis was to evaluate the dependency of the RIDICS technique on manual visual inspection. Both the techniques are equally dependent on visual inspection for epoch selection from raw EEG, but for rest of the procedure RIDICS has excessively lower dependency on visual inspection. The VI technique needs rigorous visual inspection (on IC time courses, PSD curves, TFR maps, topographic maps, etc.) to select one or few ictal components. On the contrary, the RIDICS technique depends on manual inspection for identifying ictal rhythm and ictal rhythm frequency from EEG epoch and PSD curve respectively. No visual inspection is required to estimate the cortical source of a BIC, but multiple ictal events should be visually identified and then averaged for estimating a cortical source from the outcomes of the VI technique. Similarly, the RIDICS technique requires less visual inspection compared to all other existing techniques that select more than one ictal components for source estimation. The validation results for the patients' data analysis also support the strength of the RIDICS technique. Significantly high concordance rate was achieved due to the same reason (recursive noise elimination) as discussed for the simulated-dataset analysis. Moreover, almost all the estimated cortical sources (23 out of 24) were lateralized correctly in the same hemisphere as the surgery site. According to Table 3, four epileptogenic foci that were estimated from four different seizures (1 seizure of patient 3 and 3 seizures of patient 4) were

24

found discordant. The Fz channel contact became corrupted with high level of line noise immediately after the onset of the third seizure of patient 3. Therefore, the RIDICS technique could not perform properly for that seizure. Although the results for the first three seizures of patient-4 are discordant, those sources are lateralized in the same hemisphere as the surgery site. Besides concordance rate, lateralizing power was also considered as a performance measuring metric in various related studies (Nam et al. , 2002, Holmes et al. , 2010, Kovac et al. , 2014, Englot et al. , 2015, Habib et al. , 2016). According to that metric the performance of the RIDICS technique is highly appreciable. Straightforward inverse solutions (e.g. sLORETA) estimate the cortical sources for every discrete EEG time point (Ray et al. , 2007, Holmes et al. , 2010). The cortical source for an EEG time point is considered as a reliable epileptogenic zone if the source location remains unchanged (or stable) for at least few consecutive time points (Koessler et al. , 2010, Habib et al. , 2016). Indeed, various efforts were made to improve the stability of the estimated cortical sources, such as averaging of epileptic EEG events (Merlet et al. , 2001, Kovac et al. , 2014, Habib et al. , 2016), dynamic spatiotemporal imaging by EEG decomposition in sensor space and recombination in source space (Yang et al. , 2011, Lu et al. , 2012a), and so on. In this study, each epileptogenic zone was estimated from the EEG that was regenerated from a single selected ictal component (i.e. BIC). Since each BIC is temporally independent and spatially fixed in nature, their back projected EEGs also preserve those features. Consequently, each estimated cortical source from such EEG remained stable at the same site for the total duration of the EEG. This is an advantage of RIDICS-based ESI over other ESI techniques, which usually recombine multiple ictal components for EEG regeneration and source estimation. Although using one ictal component (BIC) for source estimation solves the instability problem, using multiple ictal components may provide a scope for analyzing ictal

25

propagation. This study suggests the potential use of one ictal component (BIC) for ESI and other ICA-based ESI studies (Nam et al. , 2002, Iriarte et al. , 2006) support such single component-based ESI. If analyzing the ictal propagation is important, then RIDICS technique can be easily modified for selecting multiple ICs based on their Z values. The time-courses of Fig. 6 highlighted an important problem, the stochastic behavior, of the ICA-decomposition algorithm. It means that multiple runs of the ICA-decomposition algorithm on the same dataset can produce slightly different results (Himberg et al. , 2004, Soldati et al. , 2013). The order of the components, their time-courses and their topographic maps may vary from run to run. Since ICA-decomposition results play a significant role in the accurate EEG source estimation, the VI technique was also applied on the ICs that were found from the first recursive cycle of RIDICS. The RIDICS technique outperformed the VI technique for that analysis as well. It proves that the RIDICS technique performed better not because of using any favorable ICA-decomposition outcome but due to its algorithmic strength. The RIDICS technique considered rhythmic ictal discharges and used the ictal rhythm frequency for selecting the BIC. Although the ictal EEG patterns are mostly rhythmic, nonstationary patterns or even no distinct pattern can also be found in the ictal EEG (Ebersole et al. , 1996). The RIDICS technique may not work with such ictal EEGs. Improvement of the RIDICS technique for source estimation from ictal EEG with non-periodic ictal patterns would be of importance for future studies.

5. Conclusions This study shows that the proposed RIDICS technique improved the accuracy in localizing the epileptogenic zones. This is likely due to its ability to correctly identify the best ictal component. This is supported by the lower average dipole distance for using the RIDICS 26

technique. High concordance rate with the resected region in RIDICS-based ictal ESI validates the clinical utility of the RIDICS technique. Implementation as part of conventional clinical setup and acceptable average time for recursive decomposition indicates the practicality of the RIDICS technique. Therefore, this study forms a basis for the future studies to utilize the RIDICS technique for presurgical evaluation in those with refractory focal epilepsy.

Acknowledgements This study is a part of the Malaysian Elders Longitudinal Research (MELoR) project which is fully funded by the High Impact Research (HIR) grant (UM.C/625/1/HIR/MOHE/ASH/02) form University of Malaya and Ministry of Higher Education (MOHE), Malaysia.

Conflict of Interest Statement None of the authors have potential conflicts of interest to be disclosed.

27

References Assaf BA, Ebersole JS. Continuous source imaging of scalp ictal rhythms in temporal lobe epilepsy. Epilepsia. 1997;38:1114-23. doi: 10.1111/j.1528-1157.1997.tb01201.x Assaf BA, Ebersole JS. Visual and quantitative ictal EEG predictors of outcome after temporal lobectomy. Epilepsia. 1999;40:52-61. doi: 10.1111/j.15281157.1999.tb01988.x Bell AJ, Sejnowski TJ. An information-maximization approach to blind separation and blind deconvolution. Neural Comput. 1995;7:1129-59. doi: 10.1162/neco.1995.7.6.1129 Beniczky S, Oturai PS, Alving J, Sabers A, Herning M, Fabricius M. Source analysis of epileptic discharges using multiple signal classification analysis. Neuroreport. 2006;17:1283-7. doi: 10.1097/01.wnr.0000230517.93714.f6 Berg P, Scherg M. A fast method for forward computation of multiple-shell spherical head models. Clin Neurophysiol. 1994;90:58-64. doi: 10.1016/0013-4694(94)90113-9 Boon P, D’havé M, Vanrumste B, Van Hoey G, Vonck K, Van Walleghem P, et al. Ictal source localization in presurgical patients with refractory epilepsy. J Clin Neurophysiol. 2002;19:461-8. doi: 10.1097/00004691-200210000-00009 Chitoku S, Otsubo H, Ichimura T, Saigusa T, Ochi A, Shirasawa A, et al. Characteristics of dipoles in clustered individual spikes and averaged spikes. Brain Dev-Jpn. 2003;25:1421. doi: 10.1016/s0387-7604(02)00104-3 De Vos M, Vergult A, De Lathauwer L, De Clercq W, Van Huffel S, Dupont P, et al. Canonical decomposition of ictal scalp EEG reliably detects the seizure onset zone. Neuroimage. 2007;37:844-54. doi: 10.1016/j.neuroimage.2007.04.041 Delorme A, Makeig S. EEGLAB: an open source toolbox for analysis of single-trial EEG dynamics including independent component analysis. J Neurosci Meth. 2004;134:9-21. doi: 10.1016/j.jneumeth.2003.10.009 Delorme A, Sejnowski T, Makeig S. Enhanced detection of artifacts in EEG data using higher-order statistics and independent component analysis. Neuroimage. 2007;34:1443-9. doi: 10.1016/j.neuroimage.2006.11.004 Ding L, Worrell GA, Lagerlund TD, He B. Ictal source analysis: localization and imaging of causal interactions in humans. Neuroimage. 2007;34:575-86. doi: 10.1016/j.neuroimage.2006.09.042 Ebersole JS, Pacia SV. Localization of Temporal Lobe Foci by Ictal EEG Patterns. Epilepsia. 1996;37:386-99. doi: 10.1111/j.1528-1157.1996.tb00577.x Engel J, Jr. Update on surgical treatment of the epilepsies. Summary of the Second International Palm Desert Conference on the Surgical Treatment of the Epilepsies (1992). Neurology. 1993;43:1612-7. doi: 10.1212/WNL.43.8.1612 Englot DJ, Nagarajan SS, Imber BS, Raygor KP, Honma SM, Mizuiri D, et al. Epileptogenic zone localization using magnetoencephalography predicts seizure freedom in epilepsy surgery. Epilepsia. 2015;56:949-58. doi: 10.1111/epi.13002 Groppe DM, Makeig S, Kutas M. Identifying reliable independent components via split-half comparisons. Neuroimage. 2009;45:1199-211. doi: 10.1016/j.neuroimage.2008.12.038 Habib MA, Ibrahim F, Mohktar MS, Kamaruzzaman SB, Rahmat K, Lim KS. Ictal EEG source imaging for presurgical evaluation of refractory focal epilepsy. World Neurosurg. 2016;88:576-85. doi: 10.1016/j.wneu.2015.10.096 He B, Ding L. Electrophysiological Mapping and Neuroimaging. In: He B, editor. Neural Engineering. Boston, MA: Springer US; 2013. p. 499-543. doi: 10.1007/978-1-46145227-0_12

28

He B, Sohrabpour A, Brown E, Liu Z. Electrophysiological source imaging: A noninvasive window to brain dynamics. Annu Rev Biomed Eng. 2018;20:171-196. doi: 10.1146/annurev-bioeng-062117-120853 Himberg J, Hyvärinen A, Esposito F. Validating the independent components of neuroimaging time series via clustering and visualization. Neuroimage. 2004;22:121422. doi: 10.1016/j.neuroimage.2004.03.027 Holmes MD, Tucker DM, Quiring JM, Hakimian S, Miller JW, Ojemann JG. Comparing noninvasive dense array and intracranial electroencephalography for localization of seizures. Neurosurgery. 2010;66:354-62. doi: 10.1227/01.NEU.0000363721.06177.07 Iriarte J, Urrestarazu E, Artieda J, Valencia M, LeVan P, Viteri C, et al. Independent component analysis in the study of focal seizures. J Clin Neurophysiol. 2006;23:551-8. doi: 10.1097/01.wnp.0000236579.08698.23 Iriarte J, Urrestarazu E, Valencia M, Alegre M, Malanda A, Viteri C, et al. Independent component analysis as a tool to eliminate artifacts in EEG: a quantitative study. J Clin Neurophysiol. 2003;20:249-57. doi: 10.1097/00004691-200307000-00004 Jayakar P, Duchowny M, Resnick TJ, Alvarez LA. Localization of seizure foci: pitfalls and caveats. J Clin Neurophysiol. 1991;8:414-31. Jung K-Y, Kang J-K, Kim JH, Im C-H, Kim KH, Jung H-K. Spatiotemporospectral characteristics of scalp ictal EEG in mesial temporal lobe epilepsy with hippocampal sclerosis. Brain Res. 2009;1287:206-19. doi: 10.1016/j.brainres.2009.06.071 Kaiboriboon K, Lüders HO, Hamaneh M, Turnbull J, Lhatoo SD. EEG source imaging in epilepsy—practicalities and pitfalls. Nat Rev Neurol. 2012;8:498-507. doi: 10.1038/nrneurol.2012.150 Koessler L, Benar C, Maillard L, Badier J-M, Vignal JP, Bartolomei F, et al. Source localization of ictal epileptic activity investigated by high resolution EEG and validated by SEEG. Neuroimage. 2010;51:642-53. doi: 10.1016/j.neuroimage.2010.02.067 Kovac S, Chaudhary UJ, Rodionov R, Mantoan L, Scott CA, Lemieux L, et al. Ictal EEG source imaging in frontal lobe epilepsy leads to improved lateralization compared with visual analysis. J Clin Neurophysiol. 2014;31:10-20. doi: 10.1097/WNP.0000000000000022 Lantz G, Michel C, Seeck M, Blanke O, Landis T, Rosén I. Frequency domain EEG source localization of ictal epileptiform activity in patients with partial complex epilepsy of temporal lobe origin. Clin Neurophysiol. 1999;110:176-84. doi: 10.1016/S00134694(98)00117-5 Leal AJ, Dias AI, Vieira JP. Analysis of the EEG dynamics of epileptic activity in gelastic seizures using decomposition in independent components. Clin Neurophysiol. 2006;117:1595-601. doi: 10.1016/j.clinph.2006.03.020 Lu Y, Yang L, Worrell GA, Brinkmann B, Nelson C, He B. Dynamic imaging of seizure activity in pediatric epilepsy patients. Clin Neurophysiol. 2012a;123:2122-9. doi: 10.1016/j.clinph.2012.04.021 Lu Y, Yang L, Worrell GA, He B. Seizure source imaging by means of FINE spatio-temporal dipole localization and directed transfer function in partial epilepsy patients. Clin Neurophysiol. 2012b;123:1275-83. doi: 10.1016/j.clinph.2011.11.007 MacKinnon CD, Verrier MC, Tatton WG. Motor cortical potentials precede long-latency EMG activity evoked by imposed displacements of the human wrist. Exp Brain Res. 2000;131:477-90. doi: 10.1007/s002219900317 McMenamin BW, Shackman AJ, Maxwell JS, Bachhuber DR, Koppenhaver AM, Greischar LL, et al. Validation of ICA-based myogenic artifact correction for scalp and sourcelocalized EEG. Neuroimage. 2010a;49:2416-32. doi: 10.1016/j.neuroimage.2009.10.010 29

McMenamin BW, Shackman AJ, Maxwell JS, Bachhuber DRW, Koppenhaver AM, Greischar LL, et al. Validation of ICA-based myogenic artifact correction for scalp and source-localized EEG. Neuroimage. 2010b;49:2416-32. doi: 10.1016/j.neuroimage.2009.10.010 Merlet I, Gotman J. Dipole modeling of scalp electroencephalogram epileptic discharges: correlation with intracerebral fields. Clin Neurophysiol. 2001;112:414-30. doi: 10.1016/S1388-2457(01)00458-8 Moshe SL, Perucca E, Ryvlin P, Tomson T. Epilepsy: new advances. Lancet. 2015;385:88498. doi: 10.1016/S0140-6736(14)60456-6 Nam H, Yim TG, Han SK, Oh JB, Lee SK. Independent component analysis of ictal EEG in medial temporal lobe epilepsy. Epilepsia. 2002;43:160-4. doi: 10.1046/j.15281157.2002.23501.x Ngugi AK, Bottomley C, Kleinschmidt I, Sander JW, Newton CR. Estimation of the burden of active and life-time epilepsy: a meta-analytic approach. Epilepsia. 2010;51:883-90. doi: 10.1111/j.1528-1167.2009.02481.x Nunez PL, Srinivasan R. Electric fields of the brain: the neurophysics of EEG: Oxford University Press, USA; 2006. Onton J, Westerfield M, Townsend J, Makeig S. Imaging human EEG dynamics using independent component analysis. Neurosci Biobehav R. 2006;30:808-22. doi: 10.1016/j.neubiorev.2006.06.007 Rankine L, Stevenson N, Mesbah M, Boashash B. A quantitative comparison of nonparametric time-frequency representations. Signal Processing Conference, 2005 13th European: IEEE; 2005. p. 1-4. Ray A, Tao JX, Hawes-Ebersole SM, Ebersole JS. Localizing value of scalp EEG spikes: a simultaneous scalp and intracranial study. Clin Neurophysiol. 2007;118:69-79. doi: 10.1016/j.clinph.2006.09.010 Raz J, Turetsky B, Fein G. Confidence intervals for the signal-to-noise ratio when a signal embedded in noise is observed over repeated trials. IEEE T Bio-Med Eng. 1988;35:646-9. doi: 10.1109/10.4598 Reyes B, Charleston-Villalobos S, Gonzalez-Camarena R, Aljama-Corrales T. Timefrequency representations for second heart sound analysis. Engineering in Medicine and Biology Society, 2008 EMBS 2008 30th Annual International Conference of the IEEE: IEEE; 2008. p. 3616-9. doi: 10.1109/IEMBS.2008.4649989 Romero S, Mañanas MA, Barbanoj MJ. A comparative study of automatic techniques for ocular artifact reduction in spontaneous EEG signals based on clinical target variables: a simulation case. Comput Biol Med. 2008;38:348-60. doi: 10.1016/j.compbiomed.2007.12.001 Rosenow F, Luders H. Presurgical evaluation of epilepsy. Brain. 2001;124:1683-700. doi: 10.1093/brain/124.9.1683 Shattuck DW, Leahy RM. BrainSuite: an automated cortical surface identification tool. Med Image Anal. 2002;6:129-42. doi: 10.1016/S1361-8415(02)00054-3 Soldati N, Calhoun VD, Bruzzone L, Jovicich J. The use of a priori information in ICA-based techniques for real-time fMRI: an evaluation of static/dynamic and spatial/temporal characteristics. Front Hum Neurosci. 2013;7:64. doi: 10.3389/fnhum.2013.00064 Tenke CE, Kayser J. Surface Laplacians (SL) and phase properties of EEG rhythms: Simulated generators in a volume-conduction model. Int J Psychophysiol. 2015;97:28598. doi: 10.1016/j.ijpsycho.2015.05.008 Thurman DJ, Beghi E, Begley CE, Berg AT, Buchhalter JR, Ding D, et al. Standards for epidemiologic studies and surveillance of epilepsy. Epilepsia. 2011;52 Suppl 7:2-26. doi: 10.1111/j.1528-1167.2011.03121.x 30

Urrestarazu E, Iriarte J, Alegre M, Valencia M, Viteri C, Artieda J. Independent component analysis removing artifacts in ictal recordings. Epilepsia. 2004;45:1071-8. doi: 10.1111/j.0013-9580.2004.12104.x Vorobyov S, Cichocki A. Blind noise reduction for multisensory signals using ICA and subspace filtering, with application to EEG analysis. Biol Cybern. 2002;86:293-303. doi: 10.1007/s00422-001-0298-6 Yang L, Wilke C, Brinkmann B, Worrell GA, He B. Dynamic imaging of ictal oscillations using non-invasive high-resolution EEG. Neuroimage. 2011;56:1908-17. doi: 10.1016/j.neuroimage.2011.03.043

31

Figure Legends Fig. 1: Simplified flow chart of RIDICS technique. ICA: Independent component analysis; RIDICS: recursive ICA-decomposition for ictal component selection.

Fig. 2: (a) Three different views of the location and orientation of the dipole. (b) Changes of MSE signals with the increase of added noises. (c) The SNRs (in decibel) of the MSE signals of the simulated EEG datasets. MSE: most significant electrode; SNR: signal to noise ratio.

Fig. 3: First seizure of patient 3. Dotted vertical line indicates the seizure onset time. Inner rectangle represents the epoch that was selected for the RIDICS technique. ICA: Independent component analysis; RIDICS: recursive ICA-decomposition for ictal component selection.

Fig. 4: Euclidean distance between the original current dipole location D and the estimated current dipole location DBIC, modeled from each BIC. BIC: best ictal component.

Fig. 5: Zero-lag cross-correlation values that were estimated between noiseless dataset E and each back-projected dataset EBIC regenerated from the BIC of each simulated ictal EEG dataset. BIC: best ictal component.

Fig. 6: Independent Component time-courses of E23 dataset and the topographic maps of the ictal components (highlighted with Magenta) obtained from (a) the RIDICS technique, and (b) the VI technique. ICA: Independent component analysis; RIDICS: recursive ICA-decomposition for ictal component selection; VI: visual inspection-based.

Fig. 7: Independent Component time-courses of E41 dataset and the topographic maps of the ictal components (highlighted with Magenta) obtained from (a) the RIDICS technique, and (b) the VI technique. (c) Shows the TFR maps for the selected ictal components in (b). 32

ICA: Independent component analysis; RIDICS: recursive ICA-decomposition for ictal component selection; VI: visual inspection-based; TFR: time-frequency representation.

Fig. 8: Noise elimination from ictal EEG by using the RIDICS technique. (a) the 100x view of the selected EEG epoch of Fig. 3, and (b) Ictal EEG, obtained from the EEG epoch of Fig. 8a. ICA: Independent component analysis; RIDICS: recursive ICA-decomposition for ictal component selection.

Fig. 9: Execution time required for every ICA-decomposition. ICA: Independent component analysis.

33

ABSTRACT Objective: This study aimed to present a new ictal component selection technique, named as recursive ICA-decomposition for ictal component selection (RIDICS), for potential application in epileptogenic zone localization. Methods: The proposed technique decomposes ictal EEG recursively, eliminates a few unwanted components in every recursive cycle, and finally selects the most significant ictal component. Back-projected EEG, regenerated from that component, was used for source estimation. Fifty sets of simulated EEGs and 24 seizures in 8 patients were analyzed. Dipole sources of simulated-EEGs were compared with a known dipole location whereas epileptogenic zones of the seizures were compared with their corresponding sites of successful surgery. The RIDICS technique was compared with a conventional technique. Results: The RIDICS technique estimated the dipole sources at an average distance of 12.86 mm from the original dipole location, shorter than the distances obtained using the conventional technique. Epileptogenic zones of the patients, determined by the RIDICS technique, were highly concordant with the sites of surgery with a concordance rate of 83.33%. Conclusions: Results show that the RIDICS technique can be a promising quantitative technique for ictal component selection. Significance: Properly selected ictal component gives good approximation of epileptogenic zone, which eventually leads to successful epilepsy surgery.

34

Start Noisy EEG Wideband bandpass filtering

ICA decomposition 𝒀𝑝 -matrix calculation for all ICs Estimation of cross-correlation (Z) between 𝒀 and each of all 𝒀𝑝

Recursion Base Case

𝒀-matrix calculation

Identification of the tentative BIC having the highest Z value (Zmax) Zmax Increased?

No

Recursionoutcomes rejection

Yes Recursionoutcomes preservatio

Unwanted IC elimination from stored outcomes Narrowband bandpass filtering of the MSE Stopping criteria Yesmet?

No

Discard obtained-EEG Retrieval of preserved recursion-outcomes Backprojected EEG End

35

36

37

Table 1: Demographics and diagnosis results of the patients. Patient No. Sex Age at Onset Age at Evaluation Diagnosis 1 M 13 28 Left temporal cortical dysplasia 2 M 12 28 Left MTS 3 M 14 20 Right OLE (lesion negative) 4 M 13 25 Left TLE (lesion negative) 5 M 5 38 Right MTS with bilateral IEDs 6 F 9 42 Left MTS 7 M 12 43 Right MTS 8 M 15 48 Left MTS F, female; IED, interictal epileptiform discharge; M, male; MTS, mesial temporal sclerosis, OLE, occipital lobe epilepsy; TLE, temporal lobe epilepsy

38

Table 2: Actual sites of surgery and ESI estimated epileptogenic foci for 24 seizures of 8 patients. Red cortical area denotes epileptic focus. Patient 1

2

3

Surgery

Seizure1

Seizure2

Seizure3

Seizure4

Left middle temporal cortectomy

-

Left antero-medial temporal lobectomy

-

Right occipital cortectomy

4

Left parietooccipital cortectomy

5

Right anteromedial temporal lobectomy

-

-

6

Left antero-medial temporal lobectomy

7

Right anterior temporal lobectomy

-

-

Left selective amygdalohippocampectomy

-

-

8

39

Table 3: Concordance or discordance of ESI results with surgery sites. Patient

Seizure1

Seizure2

Seizure3

Seizure4

1

√

√

√

-

2

√

√

√

-

3

√

√

X

-

4

X

X

X

√

5

√

√

√

-

6

√

√

√

√

7

√

√

-

-

8

√

√

-

-

√, Concordant; X, Discordant; -, Not Available

40