Integrated analysis of anatomical and electrophysiological human intracranial data.

Human intracranial electroencephalography (iEEG) recordings provide data with much greater spatiotemporal precision than is possible from data obtained using scalp EEG, magnetoencephalography (MEG), or functional MRI. Until recently, the fusion of anatomical data (MRI and computed tomography (CT) images) with electrophysiological data and their subsequent analysis have required the use of technologically and conceptually challenging combinations of software. Here, we describe a comprehensive protocol that enables complex raw human iEEG data to be converted into more readily comprehensible illustrative representations. The protocol uses an open-source toolbox for electrophysiological data analysis (FieldTrip). This allows iEEG researchers to build on a continuously growing body of scriptable and reproducible analysis methods that, over the past decade, have been developed and used by a large research community. In this protocol, we describe how to analyze complex iEEG datasets by providing an intuitive and rapid approach that can handle both neuroanatomical information and large electrophysiological datasets. We provide a worked example using an example dataset. We also explain how to automate the protocol and adjust the settings to enable analysis of iEEG datasets with other characteristics. The protocol can be implemented by a graduate student or postdoctoral fellow with minimal MATLAB experience and takes approximately an hour to execute, excluding the automated cortical surface extraction.

INTRODUCTION

Intracranial EEG (iEEG) allows simultaneous recordings from tens to hundreds of electrodes placed directly on the neocortex (electrocorticography, ECoG), or intracortically (stereoelectroencephalography, SEEG). In humans, the most common implementation of iEEG is when non-invasive techniques such as scalp-EEG and MRI do not provide sufficient information to guide surgery in medication refractory epilepsy patients. Each electrode reflects the activity of tens of thousands of neurons[1,2], and the recording and stimulation of these neural populations allow for identification of epileptogenic zones, as well as for mapping of functionally eloquent areas of human cortex to guide neurosurgery. The outcome of these procedures can be directly observed when the neural or behavioral response is straightforward such as speech arrest or muscle movement with direct stimulation[3]. Any more complex empirical study requires accurate knowledge of an electrode’s location in relation to the brain’s anatomy that is linked to the local electrophysiological signal. This integrated information is key to basic and clinical research work aimed at understanding human neural and cognitive processing[4,5].

Human iEEG analysis has traditionally relied on stand-alone and ad hoc workflows for the separate analysis of anatomical and functional aspects of the iEEG data, presenting researchers with a series of challenges to realize the full potential of this exceptionally powerful tool. To process the neuroanatomical data research labs are tasked with assembling software combinations for the conversion of file formats (e.g., DICOM to NIfTI using MRIConvert), coregistration of anatomical scans (e.g., CT to MRI using SPM[6], FSL[7], or AFNI[8]), localization of electrodes (e.g., BioImage Suite[9]), and the sorting and labeling of electrodes to match the format of the functional recording file (manually, or using custom software). This technological obstacle is receiving increasing attention in the form of more efficient workflows for localizing and visualizing electrodes[10-16], but no protocol exists that allows researchers to efficiently process the anatomical data within a single work environment, and seamlessly fuse with the electrophysiological data and its subsequent analysis. Ideally, in the light of scientific reproducibility[17], such a protocol should also provide complete transparency in the evolution of raw data into results and illustrative representations, allowing for a convenient and efficient exchange of data and workflows between researchers. These two components are particularly valuable in a growing field where the analysis of data is uniquely complex, but where the gold standard for that analysis is yet to be defined.

Here, we describe – at the implementation level – a comprehensive protocol to address the series of challenges associated with both the anatomical and functional aspects of human iEEG analysis. The protocol is directly integrated with the MATLAB-based open source FieldTrip toolbox (Box 1), offering the opportunity to readily and flexibly build on a continuously growing set of analysis techniques that have already been developed and employed by a large research community. The FieldTrip toolbox supports the data formats of most popular electrophysiological data acquisition systems and shares analysis code with other software packages such as SPM and EEGLAB[18]. In contrast to the host of proprietary programs currently available for the analysis of electrophysiological data, the central tenet of FieldTrip is to provide complete transparency in order to promote a deeper understanding of the analysis techniques and enhance the quality of the scientific work that depends on these techniques. Accordingly, all computer code is fully accessible and the well-defined data structures contain full provenance to facilitate sharing between researchers. Our aim is to utilize these open source features to advance the field of human iEEG by promoting interaction within and across methodologically contiguous research areas (e.g., non-invasive electrophysiology such as EEG or MEG).

All code of the protocol is directly integrated with, and freely available through FieldTrip[53]. This MATLAB-based open source toolbox offers advanced analysis methods for electrophysiological data, such as event-related averaging, frequency and time-frequency analysis, source modeling (for EEG and MEG), connectivity analysis, classification, real-time data processing, and (non)parametric statistical inference. The implementation as a toolbox allows users to perform elaborate and structured analyses of large data sets using the MATLAB command line and batch scripting. Tutorial documentation, answers to frequently asked questions, and example code are available online as a wiki: http://www.fieldtriptoolbox.org. The toolbox’s infrastructure allows users and developers to relatively easily extend the functionality and implement new algorithms. Over the past decade, the FieldTrip toolbox has grown to an estimated 5000 users.

To get started with FieldTrip, download the most recent version from its homepage or GitHub, and set up your MATLAB path.

FieldTrip functionalities, recognizable by an prefix, typically have a single output argument and one or two input arguments, the first input argument being configuration structure .

Here, input is processed by according to parameters specified in the fields, in this case applying a 1 Hz high-pass filter. Each function’s optional parameters are available in the respective function’s header (type ) and examples are given on the wiki.

The structure that holds the parameters to the algorithm at the present level is automatically appended to the output data structure, i.e. . Configuration structures used at previous levels are kept in , , and so on. This nesting of previous configurations allows for access to the full data provenance at any level of the analysis pipeline (see for visualizing the pipeline as a flowchart).

Application of the protocol

Our protocol is especially useful for studying human neural and cognitive processes with intracranial EEG. Human iEEG analysis is uniquely complex because it requires dealing with both neuroanatomical and large electrophysiological data sets. The scope of iEEG encompasses a wide range of basic and clinical research, varying from studies of higher-order cognition[19,20] to the localization and understanding of the sources and features of epileptogenic activity[21,22]. The methodological challenges that iEEG researchers face can be grouped into obstacles that are common to most empirical work and obstacles that are study-specific. This protocol aims to resolve the former, while providing adequate support and flexibility for the latter.

Advantages and limitations of the protocol

The main advantages of our protocol are that it (i) guides the researcher from the multitude of raw intracranial data files to integrated observations, in a fast and efficient way, (ii) is directly integrated with a comprehensive and open source hub for electrophysiological data analysis, (iii) can be readily adapted and automated, (iv) is completely transparent and (v) produces reproducible workflows and data sets that can be easily shared and generalized to other research modalities. The main limitation is that the MATLAB command line interface requires some basic programming knowledge, which may require more initial learning as compared to the execution of computer commands through a (black box) graphical user interface. However, the use of computer commands can be relatively easily mastered by virtue of using this protocol, paving the way for batch scripting in order to efficiently deal with repeated analyses within and across subjects and, ultimately, for a deeper understanding of the underlying algorithms.

Human intracranial data sets are approached from various angles and come in different shapes and sizes, so it is critical for a protocol to strike the right balance between efficiency and flexibility. This need is further amplified by the relatively unique nature of intracranial data, typically imposing greater demands on alternative options and strategies in the analysis than non-invasive data recorded with more standardized hard- and software in dedicated laboratory settings. Besides providing a quick guide to interpretable results, our protocol allows for easy switching between methods to accommodate different cases and situations. By changing a single parameter at execution, one can for instance readily apply a different fusion cost function or filter setting. Utilizing this versatility should not negatively impact continuation with the protocol. In fact, the full and automatic provenance, in combination with the systematic file naming, encourages adapting to the circumstances by alleviating potential concerns regarding oversight and reproducibility.

The spatiotemporal precision of intracranial EEG provides a unique window on neural processing. The size and dimensions of this window, however, may grow disproportionally large with certain types of analyses, complicating the overall interpretability of the data. Starting from the two dimensions of the raw neural signal (channels and time), a time-frequency analysis, for instance, results in 3 dimensions in the output (power as function of channel, time and frequency), whereas between-channel connectivity analysis expands the combinatorial space to 4 dimensions. Our protocol addresses this dimensionality issue and illustrates how the interactive manipulation of anatomically informed graphical representations of the neural data facilitates the inspection of the multi-dimensional outcome of an iEEG analysis, taking maximum advantage of the groundwork laid by the integrated processing of the anatomical and functional data.

Integration with FieldTrip

In addition to the complete transparency that comes with an open source toolbox, the integration with FieldTrip provides unique benefits to iEEG researchers by allowing them to build on algorithms for reading in raw data of various formats, data preprocessing, event-related potential analysis, spectral analysis, source modeling, connectivity analysis, classification, real-time data processing, and statistical inference. Applied to human iEEG data, these methods permit characterization of neural information flow with a level of detail inaccessible to noninvasive techniques. Additionally, invasive and non-invasive human electrophysiology can be directly overlaid using very similar analysis pipelines for an integrative perspective of neural processing, or a comparison of MEG/EEG source reconstruction methods with iEEG.

The open source development model allows for a relatively easy extension of the protocol. For instance, several techniques exist to compensate for electrode displacement due to the “brain shift” phenomenon explained below[11,12,23-30]. Given different strengths and weaknesses, these techniques may need to be evaluated on a case-by-case basis. FieldTrip’s modular architecture facilitates developers to incorporate new techniques and users to subsequently employ those techniques by virtue of changing a single parameter at function call. In a similar vein, the protocol can be extended to a number of exciting new research areas. These include single- and multiunit recordings, ‘NeuroGrid’ recordings[31], wireless ‘Neural Dust’ recordings[32], (deep) brain stimulation[33,34], and multimodal imaging[12]. Supported by a growing community of developers committed to the ongoing push to improve data analysis methods, we will coordinate with these new electrophysiological endeavors and continue sharing analysis code with other software packages.

Compatibility with FreeSurfer

The protocol is compatible with the freely available software package FreeSurfer[35]. Although optional, processing of the anatomical MRI with FreeSurfer (Step 6) offers several advantages for subsequent analysis and data interpretation. Processing the MRI with FreeSurfer results in the creation of a cortical mesh, consisting of approximately equally sized triangles that form a topological sphere for each of the cerebral hemispheres. This cortical mesh is particularly convenient for an anatomically realistic representation of the electrophysiological data on the neocortex (e.g., bottom center in Fig. 1). A smoothed version of the extracted cortical surface can be used in the compensation for electrode displacement due to brain shift (Step 22). Moreover, FreeSurfer automatically registers the subject’s brain to a template brain on the basis of its cortical gyrification pattern, an aspect of brain structure that remains difficult to accurately normalize using volume-based registration techniques due to its complexity and variability across subjects[36,37]. Our protocol uses the resulting surface registration maps to link electrode positions to their template homologues (Step 29). Finally, FreeSurfer-generated atlases are convenient for representations of neural and anatomical data for a single subject (Step 52), since they are defined in native subject space. Other supported atlases are defined in standardized (e.g., MNI) space and require the added step of transforming electrode positions to that space.

Figure 1

Overview of the procedure

The protocol is grounded in two parallel but interrelated workflows. The anatomical workflow minimally consists of the preprocessing and fusion of the anatomical images and electrode placement. The functional workflow encompasses the preprocessing of the neural recordings, but may also include follow-up activities such as event-related averaging, time-frequency and statistical analysis. The electrode placement activity offers the opportunity to directly link anatomical locations to electrode labels corresponding to the neural recordings, allowing for an early seamless integration of the two workflows to facilitate anatomically informed data exploration and visualization.

Human intracranial data

Anatomical images, typically MRI and CT scans, are used as part of the epilepsy diagnostic and surgical procedures. A pre-implant MRI shows the anatomy of the head including the brain and is used to identify structural abnormalities. An MRI is also instrumental in guiding SEEG electrode implantation subsequent to the clinical decision to record intracranially. A post-implant CT shows high-intensity objects such as the electrodes and skull but lacks details of brain anatomy. To obtain knowledge of an electrode’s location in relation to the brain’s anatomy, the two scans have to be fused.

Following fusion of the pre- and postoperative anatomical images, electrodes that have been surgically placed on the cortical surface occasionally appear “buried” within the cortical tissue, sometimes more than a centimeter deep[38-43]. This electrode displacement is typically due to “brain shift”, the inward sinking of the brain post-implant most commonly observed with electrocorticographic surface grid electrodes. The brain shift reflects tissue displacement, caused by the electrodes themselves, and by subdural fluid loss or accumulation. As noted, the displacement is most pronounced directly below a craniotomy and is usually minimal for implants solely involving burr holes[43]. It is important to account for this brain shift in order to accurately align electrode specific signals with the local cortical anatomy. Several labs have developed realignment techniques to compensate for electrode displacement due to brain shift, reducing localization error to under 3 mm when compared to intraoperative photographs[11,23-30]. Our protocol currently supports two of these techniques to project electrode grids back to the cortical surface while accounting for a grid’s shape and orientation[23,30].

Electrode localization can also be done using post-implant MRIs, although these are not commonly acquired in a clinical setting. These scans show the brain anatomy after electrode implantation, so brain shift is not an issue. In a T1-weighted MRI, electrodes appear dark, due to the magnetic susceptibility artifact. This is generally not an issue for recordings with depth electrodes (SEEG), where the electrodes are visible as dark voids in the higher intensity brain tissue. Electrode grids and strips (ECoG), on the other hand, are placed directly on the cortical surface. This complicates their identification, as the electrodes are surrounded by cerebral spinal fluid, which also appears dark on a T1 scan (but see [25,44-46] for workarounds). Depending on the availability of a post-implant MRI of sufficient quality that clearly shows the electrodes, the CT preprocessing and fusion Steps 9–15 may be left out, and electrode localization may be done on the post-implant MRI. However, if the post-implant MRI is of unsatisfactory quality regarding brain anatomy, for instance due to electrode induced MR signal distortion, we recommend fusing the post-implant MRI with the pre-implant MRI, as if it were a post-implant CT.

Neural recordings are typically part of the ongoing clinical monitoring and come in various file formats. Each data channel represents, as a function of time, the electric potential difference, obtained with either a bipolar or referential electrode scheme. That is, the electrodes are pairwise linked or referenced to a single, common electrode during acquisition. The latter montage has the benefit that the recordings can be easily re-montaged to a more preferred scheme in the offline analysis[47]. The markers or triggers for stimulus onset times and responses are typically recorded simultaneously in a dedicated channel, allowing for precise synchronization of experimental scenarios with the neural recording.

Overview of the procedure

The protocol is grounded in two parallel but interrelated workflows, as shown in Figure 1. The first workflow entails the processing of anatomical data. Its main activities constitute the preprocessing and fusion of the anatomical images, and electrode placement (Steps 1–19). Secondary activities that are also discussed include cortical surface extraction with FreeSurfer, brain shift compensation, spatial normalization, and anatomical labeling (Steps 6 and 20–33). Generally, the anatomical workflow aims to obtain estimates of the electrode locations in relation to the individual and atlas-based brain anatomy, which is a one-time procedure for each subject. The second workflow focuses on improving the signal-to-noise ratio and extracting the relevant features from the electrophysiological data, while preparing for subsequent analyses. It minimally encompasses the preprocessing of the neural recordings, but may also include follow-up activities such as time-frequency and single-subject or group-level statistical analysis (Steps 34–45). Generally, the specifics of the functional workflow depend ultimately on the clinical or research question at hand and contingencies in the experimental paradigm.

The two workflows become intrinsically connected for the first time during the electrode placement activity (Step 17), which offers the opportunity to directly link anatomical locations to electrode labels corresponding to the neural recordings. This activity involves an interactive electrode placement tool designed for efficient yet precise identification of electrodes in even the most challenging cases. The integration of the two workflows culminates in an interactive and anatomically informed data exploration tool and the ability to represent functional and epileptiform neural activity overlaid on cortical and subcortical surface models, in figure or video format (Steps 46–56).

Implementation and adaptation of the procedure

All implementations run on a single universal platform (MATLAB, except for FreeSurfer) to support relatively easily automated procedures for dealing with repeated analyses within and across subjects. We recommend that the user construct a single script for a single subject by incrementally copy-pasting code from this protocol into the MATLAB editor (Supplementary File 1), and evaluating segments of that script in the MATLAB command window. Once the script produces satisfactory results, it can be converted into a batch analysis by breaking it into separate components. By looping around the separate components for all subjects, the entire analysis pipeline for all subjects in a study can easily be executed and intermediate results can be saved and evaluated.

The whole batch can be documented and shared, or re-evaluated with different parameter settings as appropriate. By virtue of changing single parameters at a function call, one can for instance readily alternate between various fusion, localization, projection, normalization, filtering, re-montaging, and spectral estimation algorithms to accommodate different cases and situations. The output data structures are kept consistent across the different algorithms, and the parameters to the used algorithm are appended to allow for access to the full data provenance at any level of the analysis pipeline (Box 1).

Experimental design

The example iEEG data set was acquired at the Medical Center of the University of California, Irvine. The Office for the Protection of Human Subjects of the University of California, Berkeley, approved the study and the subject gave informed consent. The data set includes a pre-implant MRI, a post-implant CT, a post-implant MRI, and neural recordings from 96 ECoG and 56 SEEG electrodes that were implanted as part of the preparation for the epilepsy surgery (see Materials). The neural data were recorded in the context of an experiment that required the patient to press a button with the right hand when hearing a target tone. The original data set (after defacing the imaging data with ) and the processed results are available for download from ftp://ftp.fieldtriptoolbox.org/pub/fieldtrip/tutorial/SubjectUCI29.zip. Raw DICOM images and recording files are not shared to protect the subject’s identity.

We choose this iEEG data set for three reasons. First, it contains neural recordings from both cortical grid (ECoG) and stereotactically inserted depth electrodes (SEEG), requiring strategies for dealing with each type as well as their combination in the analysis. Second, the pre-implant MRI is not of the best quality (a contrast agent was used), electrodes of adjacent cortical grids have seemingly merged with one another in the post-implant CT, and there is significant electrode displacement due to a subdural hygroma contributing to brain shift. These issues reflect real world challenges in intracranial data analysis, allowing us to demonstrate the application of our protocol to non-ideal data. Finally, the experimental paradigm is simple enough to need no further explanation, yet requires performing all the fundamental steps underlying the analysis of intracranial data recorded using a more complex experimental paradigm[19,48]. We demonstrate the analysis of task-related high-frequency-band activity (~70 to 150 Hz), a prominent neural signature in intracranial data that has been associated with neuron population level firing rate[5,49-52]. Many other supported analyses such as event-related potential analysis, connectivity analysis, and statistical analysis have been described in detail elsewhere[53-55].

Pre-implant T1-weighted MRI (Magnetic Resonance Image, Siemens 3T TrioTim).

Post-implant CT (Computerized Tomography, Philips iCT 256).

Post-implant T1-weighted MRI (Magnetic Resonance Image, Siemens 1.5T Avanto). This scan is not used in the procedure but nevertheless included for completeness.

64-contact cortical grid with left parietal coverage (Integra, 8 × 8 layout, 10 mm inter-electrode spacing, labels have a LPG prefix)

32-contact cortical grid with left temporal coverage (Integra, 4 × 8 layout, 10 mm inter-electrode spacing, labels have a LTG prefix)

8-contact linear depth electrode targeting left amygdala (Ad-Tech, 5 mm inter-electrode spacing, labels have a LAM prefix)

8-contact linear depth electrode targeting left hippocampus head (Ad-Tech, 5 mm inter-electrode spacing, labels have a LHH prefix)

8-contact linear depth electrode targeting left hippocampus tail (Ad-Tech, 5 mm inter-electrode spacing, labels have a LTH prefix)

8-contact linear depth electrode targeting right amygdala (Ad-Tech, 5 mm inter-electrode spacing, labels have a RAM prefix)

8-contact linear depth electrode targeting right hippocampus head (Ad-Tech, 5 mm inter-electrode spacing, labels have a RHH prefix)

8-contact linear depth electrode targeting right hippocampus tail (Ad-Tech, 5 mm inter-electrode spacing, labels have a RTH prefix)

8-contact linear depth electrode targeting right occipital cortex (Ad-Tech, 5 mm inter-electrode spacing, labels have a ROC prefix)

All neural recordings were acquired using a Nihon Kohden recording system with a JE-120A amplifier (Nihon Kohden Corporation, Tokyo, Japan), analogfiltered above 0.01 Hz, and digitally sampled at 5 KHz

MATLAB environment (MathWorks, Natick, MA; installation and licensing through http://www.mathworks.com)

FieldTrip toolbox (Box 1, freely available at http://www.fieldtriptoolbox.org)

FreeSurfer software suite for cortical surface extraction (optional; freely available at http://www.freesurfer.net)

AFNI (*.head, *.brik)

Analyze (*.img, *.hdr)

ANT (*.mri)

DICOM (*.dcm, *.ima)

FreeSurfer (*.mgz, *.mgh)

MINC (*.mnc)

NIfTI (*.nii, *.nii.gz)

Anywave (*.ah5)

BCI2000 (*.dat)

BESA (*.besa)

Blackrock (*.nev, *.ns#)

Cambridge Electronic Design (*.smr)

European Data Format (*.edf)

GTec (*.mat, *.hdf5)

Micromed (*.trc)

Neuralynx (*.ncs, *.nse, *.nts, *.nst, *.ntt, *.nev)

Neuromag (*.fif)

Neuroscope (*.eeg, *.dat, *.xml)

Nihon Kohden (*.m00)

Plexon (*.ddt, *.nex, *.plx)

and various EEG, MEG, NIRS, and eye-tracker data formats

Specify the subject ID. This ID will be used in the file naming, in addition to 988information about the type of data (e.g., MRI, CT), the coordinate system it is in (e.g., ACPC, MNI), and the process(es) that were applied to it (e.g., f for fusion). For example, a CT scan that is aligned to the ACPC coordinate system and that has just been fused with the anatomical MRI is written out to file as subjID_CT_acpc_f.nii.

Import the anatomical MRI into the MATLAB workspace using . The MRI comes in the format of a single file with an .img or .nii extension, or a folder containing a series of files with a .dcm or .ima extension (DICOM; Supplementary File 2 may aid in the search and visualization of a DICOM series).

Determine the native orientation of the anatomical MRI’s left-right axis using (Box 2 and Supplementary Video 1).

CRITICAL STEP To correctly fuse the MRI and CT scans at a later step, accuracy in demarcating the right hemisphere landmark in the following step is important for avoiding an otherwise hard to detect flip of the scan’s left and right orientation.

Align the anatomical MRI to the ACPC coordinate system[56], a preferred convention for the FreeSurfer operation optionally used in a later step. In this coordinate system, the origin (coordinate [0,0,0]) is at the anterior commissure (AC), the Y-axis runs along the line between the anterior and posterior commissure (PC), and the Z-axis lies in the midline dividing the two cerebral hemispheres. Specify the anterior and posterior commissure, an interhemispheric location along the midline at the top of the brain, and a location in the brain’s right hemisphere. If the scan was found to have a left-to-right orientation in the previous step, the right hemisphere is identified as the hemisphere having larger values along the left-right axis. Vice versa, in a right-to-left system, the right hemisphere has smaller values along that axis than its left counterpart (Supplementary Video 2).

Write the preprocessed anatomical MRI out to file.

Coordinate systems define the orientation and units of the X-, Y-, and Z-axes of an anatomical volume in addition to an origin point along the brain’s midline (e.g., anterior commissure). Here we provide a guideline for determining the native coordinate system of the MRI and CT scans and, in particular, whether they have a left-to-right or a right-to-left orientation. Knowledge of the orientation of the left-right axis of the scan's native coordinate system provides the necessary context for demarcating the right hemisphere landmark in the succeeding alignment step. Although the interpretation of posterior-anterior and inferior-superior axes is straightforward from axial, coronal, or sagittal slices of the brain, differentiating left and right requires a three-dimensional context. To accomplish this, we recommend using , which depicts an anatomical volume as three intersecting, orthogonal slices and labels the X-, Y-, and Z-axes. This allows determining which of these three axes represents the left-right axis and, importantly, whether that axis has a left-to-right or a right-to-left orientation (Supplementary Video 1).

Visualize the coordinate system of the MRI or CT:

Determine which of the three axes, X, Y, or Z, runs through or along the left-right axis of the subject’s head. This axis is the left-right axis for this anatomical volume.

Determine the orientation of the left-right axis. If the values on the left-right axis increase to the right (indicated by a + sign), then the scan has a left-to-Right orientation. If the values on the left-right axis increase to the left, then the scan has a right-to-Left orientation.

Write down the orientation of the scan’s left-right axis.

Execute FreeSurfer’s functionality from the Linux or MacOS terminal (Windows via VirtualBox), or from the MATLAB command window as below. This set of commands will create a folder named ‘freesurfer’ in the subject directory, with subdirectories containing a multitude of FreeSurfer-generated files.

PAUSE POINT FreeSurfer’s fully automated segmentation and cortical extraction of the anatomical MRI currently may take up 10 hours or more. For tutorial purposes, the example data set contains the output from FreeSurfer, a folder named ‘freesurfer’, for continuation with the protocol.

Import the extracted cortical surfaces into the MATLAB workspace and examine their quality. Repeat the following code using to visualize the pial surface of the right hemisphere.

?TROUBLESHOOTING

Import the FreeSurfer-processed MRI into the MATLAB workspace for the purpose of fusing with the CT scan at a later step, and specify the coordinate system to which it was aligned in Step 4.

Import the anatomical CT into the MATLAB workspace. Similar to the MRI, the CT scan comes in the format of a single file with an .img or .nii extension, or a folder containing a series of files with a .dcm or .ima extension (Supplementary File 2 may aid in the search and visualization of a DICOM series).

In case this cannot be done on the basis of knowledge of the laterality of electrode implantation, determine the native orientation of the anatomical CT’s left-right axis using , similarly to how it was done with the anatomical MRI in Step 3 (Box 2 and Supplementary Video 1).

CRITICAL STEP To correctly fuse the MRI and CT scans at a later step, accuracy in demarcating the right and left preauricular landmark in the following step is important for avoiding an otherwise hard to detect flip of the scan’s left and right orientation.

Align the anatomical CT to the head surface coordinate system by specifying the nasion (at the root of the nose), left and right preauricular points (just in front of the ear canals), and an interhemispheric location along the midline at the top of the brain (Supplementary Video 3). The CT scan is initially aligned to the head surface coordinate system, given that the ACPC coordinate system used for the MRI relies on neuroanatomical landmarks that are not visible in the CT.

Automatically convert the CT’s coordinate system into an approximation of the ACPC coordinate system, the same system the anatomical MRI was aligned to.

Fuse the CT with the MRI, a necessary step to link the electrode locations in the anatomical CT to their corresponding locations in the anatomical MRI[57,58]. Given that both scans are from the same subject and their common denominator is the skull, a rigid body transformation suffices for their alignment under normal circumstances (the default technique when using the SPM-method in FieldTrip).

Carefully examine the interactive figure that is produced after the coregistration is completed, showing the MRI and fused CT superimposed. A successful fusion will show tight interlocking of CT-positive skull (in blue) and MRI-positive brain and skin tissue (in red).

CRITICAL STEP Accuracy of the fusion operation is important for correctly placing the electrodes in anatomical context in a following step.

?TROUBLESHOOTING

Write the MRI-fused anatomical CT out to file.

Import the header information from the recording file, if possible. By giving the electrode labels originating from the header as input to in the next step, the labels will appear as a to-do list during the interactive electrode placement activity. A second benefit is that the electrode locations can be directly assigned to labels collected from the recording file, obviating the need to sort and rename electrodes to match the electrophysiological data.

Localize the electrodes in the post-implant CT with , shown in Figure 2. Clicking an electrode label in the list will directly assign that label to the current crosshair location (Supplementary Video 4). Several in-app features facilitate efficient yet precise navigation of the anatomical image, such as a zoom mode, a magnet option that transports the crosshair to the nearest weighted maximum with subvoxel accuracy (or minimum in case of a post-implant MRI), and an interactive three-dimensional scatter figure that is linked to the two-dimensional volume representations. Furthermore, passing on the pre-implant MRI, , to allows toggling between CT and MRI views for the identification of specific electrodes based on their anatomical location. Generally, electrode #1 is the electrode farthest away from the craniotomy or burr hole in case of depths and single-row strips. Careful notes taken during surgery and recording are critical for determining the numbering of grid and multi-row strip electrodes.

Figure 2

Interactive electrode placement

Clicking an electrode label in the main panel on the left will directly assign that label to the current crosshair position in the CT scan. Several features facilitate precise navigation of the anatomical CT, such as a zoom mode, a magnet option that transports the crosshair to the nearest weighted maximum (or minimum in case of a post-implant MRI), and the interactive three-dimensional scatter figure shown on the right.

Examine whether the variables in resulting electrode structure match the recording parameters, e.g., the number of channels stored in the field. The electrode and channel positions are stored in the and fields, respectively. The field contains the original electrode positions. With exception of possible brain shift compensation, this field is not adjusted. The channel positions in the field are initially identical to the electrode positions but may be updated to accommodate offline adjustments in channel combinations, i.e. during re-montaging. For bipolar iEEG data, the best considered channel position is in between the two corresponding electrode positions. The chanpos field is used for overlaying the neural data on (sub-)cortical models during data visualization. The field is a matrix with the weight of each electrode into each channel, which at this stage merely is an identity matrix reflecting one-to-one mappings between electrodes and channels.

Save the resulting electrode information to file.

In case of “brain shift”, a displacement of brain tissue and electrodes post-implant, realignment of electrode grids to the preoperative cortical surface may be necessary. To prevent electrodes from being incorrectly placed in the nearby cortical sulci during back-projection, create a smooth hull around the cortical surface generated by FreeSurfer[59].

Save the hull to file.

Project the electrode grids to the surface hull of the implanted hemisphere. Given that different grids can move independently from one another and that the projection algorithm specified in considers the global electrode configuration of a grid[30], it is recommended to realign electrode grids individually by running separate realignment procedures for each grid. Here, we realign the electrodes of the left parietal grid followed by the electrodes of the left temporal grid (LPG and LTG respectively) and store the updated grid electrode information in a new variable together with the unaltered coordinates of the depth electrodes.

Visualize the cortex and electrodes together and examine whether they show expected behavior (Fig. 3).

Figure 3

Brain shift compensation

In some patients, compensation for electrode displacement due to brain shift after implantation may be necessary. In this particular case, a subdural hygroma at the top of the brain caused severe electrode displacement in a direction opposite to the more commonly observed inward shift (left). Realigning electrode grids to the cortical surface can compensate for electrode displacement due to brain shift (right). The thin black lines indicate each electrode’s path from its localized origin on the left to its projected location on the right.

CRITICAL STEP Accuracy of the realignment operation is important for correctly placing the electrodes in anatomical context in a following step.

?TROUBLESHOOTING

Save the updated electrode information to file.

To generalize the electrode coordinates to other brains or MNI-based neuroanatomical atlases in a later step, register the subject’s brain to the standard MNI brain. The volume-based registration technique considers the overall geometry of the brain[60] and can be used for the spatial normalization of all types of electrodes, whether depth or on the surface.

Use the resulting deformation parameters to obtain the electrode positions in standard MNI space.

Visualize the cortical mesh extracted from the standard MNI brain along with the spatially normalized electrodes and examine whether they show expected behavior (top right in Fig. 4).

Figure 4

Spatial normalization

On the left are the electrodes on the individual cortical sheet. The top right shows the electrodes on the standard MNI brain after volume-based registration. The bottom right shows the electrodes on FreeSurfer’s fsaverage brain after surface-based registration. Compared to volume-based registration, with surface-based registration the original grid geometry is no longer preserved as electrodes are moved from one brain to another according to the curvature pattern of the cortex.

CRITICAL STEP Accuracy of the spatial normalization step is important for correctly overlaying the electrode positions with a brain atlas in a following step.

?TROUBLESHOOTING

Save the normalized electrode information to file.

To generalize the electrode coordinates to other brains in a later step, map the electrodes onto FreeSurfer’s fsaverage brain. The surface-based registration technique solely considers the curvature patterns of the cortex[35] and thus can be used for the spatial normalization of electrodes located on or near the cortical surface. In the example case, this pertains to all electrodes of the left parietal and temporal grids.

Visualize FreeSurfer’s fsaverage brain along with the spatially normalized electrodes and examine whether they show expected behavior (bottom right in Fig. 4).

CRITICAL STEP Accuracy of the spatial normalization step is important for correctly overlaying the electrode positions with a brain atlas in a following step.

Save the normalized electrode information to file.

FieldTrip supports looking up the anatomical or functional labels corresponding to the electrodes in a number of atlases, including the AFNI Talairach Tournoux atlas[61], the AAL atlas[62], the BrainWeb data set[63], the JuBrain cytoarchitectonic atlas[64], the VTPM atlas[65], and the Brainnetome atlas[66], in addition to the subject-tailored Desikan-Killiany and Destrieux atlases produced by FreeSurfer[67,68]. With exception of the above FreeSurfer-based atlases, these atlases are in MNI coordinate space and require the electrodes to be spatially normalized (Step 25). First, import an atlas of interest, e.g., the AAL atlas, into the MATLAB workspace.

Look up the corresponding anatomical label of an electrode of interest, e.g., electrode LHH2, targeting the left hemisphere’s hippocampus. Supplementary File 3 represents a tool that automatically overlays all channels in an electrode structure with all of the above atlases and stores the resulting anatomical labels in an excel table (e.g., SubjectUCI29_electable.xlsx in the zip file).

?TROUBLESHOOTING

Define the trials, that is, the segments of data that will be used for further processing and analysis. This step produces a matrix containing for each segment the begin and end sample in the recording file. In the case of the example provided in the shared data, the segments of interest begin 400 ms before tone onset, are marked with a ‘4’ in the trigger channel, and end 900 ms thereafter.

Import the data segments of interest into the MATLAB workspace and filter the data for high-frequency and power line noise (see the documentation of for filtering options).

Examine whether the variables in the output data structure match the recording and preprocessing parameters, i.e. the sampling rate ( ), number of recording channels ( ), and segmentation into the experiment’s twenty-six trials ( , and their respective time axes in ).

Add the structure originating from the anatomical workflow and save the preprocessed electrophysiological data to file. The advantage of adding the electrode information at this stage is that it will be kept consistent with the neural data going forward, as when applying the same montage used for the neural recordings to the channel positions.

Inspect the neural recordings using and identify channels or segments of non-interest, for instance segments containing signal artifacts or (in this case) epileptiform activity. Mark the bad segments by drawing a box around the corrupted signal. Write down the labels of bad channels.

CRITICAL STEP Identifying bad channels is important for avoiding the contamination of other channels during re-montaging in Step 40.

Remove any bad segments marked in the above step.

Re-montage the cortical grids to a common average reference in order to remove noise that is shared across all channels. Box 3 provides a background on re-montaging. Bad channels noted in Step 38 can be excluded from this step by adding those channels to with a minus prefix. That is, if one were to exclude the LPG1 channel from the list of LPG and LTG channels.

Apply a bipolar montage to the depth electrodes. This can be done in a similar manner as in Step 40, but by selecting single channel labels for and . Alternatively, create a more elaborate scheme with (see the documentation of ). Here, we combine for each depth electrode shaft the 8 unipolar channels into 7 bipolar channels, using the weights defined in the 7×8 field. We also create new labels indicating the bipolar origin of the data, e.g., “RAM1-RAM2”, “RAM2-RAM3”, and so on. Note that because we added the structure to the data in Step 37, the same montage is automatically applied to the channel positions as well, with the resulting field containing the mean locations of all electrode pairs that comprise a bipolar channel.

Combine the data from both electrode types into one data structure for the ease of further processing.

Save the re-referenced data to file.

The recorded electrophysiological signals are a mixture of signal-of-interest and noise, both neural and non-neural. The main objective of the preprocessing of the neural recordings is to improve the signal-to-noise ratio of the data while optimally preparing it for follow-up analysis. Re-montaging to a different referencing scheme, also known as a montage, may aid in the removal of noise that is shared across multiple channels. The common average re-referencing technique, for instance, involves taking the average potential from all channels and subtracting this global noise estimate from the potential in each channel[47,69-72]. We demonstrated how to apply this technique to the cortical grid electrodes in our example case.

Depth electrodes, located inside the brain and using differently sized and shaped contact points, have a different sensitivity distribution and capture different types of activity and levels of noise[2]. There is currently no consensus on the preferred montage for depth-electrode recordings and, thus, what electrodes to use as references[73-76]. White matter signals may not be as silent as one would intuitively expect, and bipolar signals, despite being relatively clean, miss out on activity that had the same amplitude on the two consecutive electrodes prior to their combination[5,77]. Different options may need to be tested and evaluated per case, taking into account the purpose of any follow-up analysis[71]. For instance, see [55] for a discussion of connectivity analysis in relation to the referencing scheme.

Decompose the signal in time and frequency bins. The configuration options and determine the frequencies and time-points of interest, in this case from 5 to 200 Hz in steps of 5 Hz, and 300 ms prior to tone onset until 800 ms thereafter in steps of 10 ms.

Save the time-frequency data to file.

For an anatomically informed exploration of the multidimensional outcome of an analysis, create a layout based on the three-dimensional electrode locations. This layout is a symbolic representation in which the channels are projected on the two-dimensional medium offered by paper or a computer screen. The layout is complemented by an automatic outline of the cortical sheet that is specified in . The option allows selecting channels whose two-dimensional distances are used to determine the plotting box sizes in the following step.

Express the time-frequency representation of neural activity at each channel in terms of the relative change in activity from a baseline interval.

Visualize the time-frequency representations overlaid on the two-dimensional layout. The generated figure is interactive, so that selecting a group of channels will launch another figure representing the average time-frequency representation over those channels (Fig. 5). Selecting a certain frequency and time range in that time-frequency representation will launch yet another figure showing the topographical distribution of activity in the selected interval, and so on (Supplementary Video 5).

Figure 5

Interactive plotting

Fast browsing through various anatomically informed representations of the neural data can help address the multidimensionality of intracranial EEG data.

For an anatomically realistic representation of cortical activity, overlay a surface model of the neocortex with the spatial distribution of the high frequency-band activity. First, extract high-frequency-band activity during a time interval of interest.

Visualize the spatial distribution of high-frequency-band activity on a cortical mesh of the subject’s brain.

Add the electrodes to the figure (Fig. 6). By looping around Steps 49 to 51 while breaking down the time in terval of interest specified with in consecutive steps, it becomes feasible to observe the spatiotemporal dynamics of neural activity occurring in relation to known experimental structure and behavior (Supplementary Video 6). See for capturing and assembling time-lapse movies.

Figure 6

ECoG data representation

Task-induced high-frequency-band activity relative to a baseline interval, plotted on a cortical surface mesh of the subject’s brain.

For depth recordings, create an integrated representation of neural activity and anatomy by interpolating neural data from each bipolar channel in a spherical cloud, which can then be overlaid on a surface mesh of any deep brain structure. First, create a volumetric mask of the regions of interest (ROI). Here, we generate a mask for the right hippocampus and amygdala from the cortical parcellation and subcortical segmentation produced by FreeSurfer.

Create a triangulated and smoothed surface mesh on the basis of the volumetric masks.

Identify the subcortical electrodes of interest.

Interpolate the high-frequency-band activity in the bipolar channels on a spherical cloud around the channel positions, while overlaying the neural activity with the above mesh. By repeating the current step for neural data corresponding to consecutive time intervals, similarly to the process outlined in Step 51, it becomes feasible to create time-lapse movies of the spatiotemporal dynamics of deep-brain activity (Supplementary Video 7 shows the spatiotemporal evolution of epileptiform activity in a separate subject).

To create a more definitive image of the neural activity at particular positions, generate two-dimensional slices through the three-dimensional representations. This combination provides the most complete and integrated representation of neural and anatomical data (Fig. 7).

Figure 7

SEEG data representation

Task-induced high-frequency-band activity relative to a baseline interval, plotted as point clouds around a triangulated mesh of the subject’s amygdala and hippocampus in the right hemisphere. The two-dimensional planes on the right correspond to the slices in the image on the left.

Steps 2–5, Preprocessing of the anatomical MRI: ~2 min

Steps 6–8, Cortical surface extraction with FreeSurfer (optional): ~10 hrs

Steps 9–12, Preprocessing of the anatomical CT: ~2 min

Steps 13–15, Fusion of the CT with the MRI: ~3 min

Steps 16–19, Electrode placement: ~15 min

Steps 21–24, Brain shift compensation (optional): ~5 min

Steps 25–28, Volume-based registration (optional): ~2 min

Steps 29–31, Surface-based registration (optional): ~2 min

Steps 32–33, Anatomical labeling (optional): ~1 min

Steps 34–43, Preprocessing of the neural recordings: ~10 min

Steps 44–45, Time-frequency analysis (optional): ~2 min

Steps 46–48, Interactive plotting: ~3 min

Steps 49–51, ECoG data representation: ~1 min

Steps 52–56, SEEG data representation: ~2 min

Box 2, Coordinate system determination: ~1 min

ANTICIPATED RESULTS

Upon completion of the protocol, one should obtain an integrated representation of neural and anatomical data. The exact results depend ultimately on the clinical or research question at hand, contingencies in the experimental paradigm, and decisions made during the execution of the protocol. We demonstrated the analysis of spatiotemporal neural dynamics occurring in relation to known experimental structure and relatively simple behavior, namely the pressing of a button with the right hand when hearing a target tone (Fig. 5–7, Supplementary Video 6). However, with small adaptations of the protocol it is feasible to track the spatiotemporal evolution of epileptiform activity with high precision (Supplementary Video 7), or to perform group-level investigations of fine-grained decision-related neural dynamics in human orbitofrontal cortex[48]. A precise fusion of the anatomical images with the electrophysiological data is key to reproducible analyses and findings. Hence, it is important to examine the outcome of any critical step, as we have done in this protocol (e.g., Fig. 3 and 4).

TABLE 1

Troubleshooting table

Step	Problem	Possible reason	Solution
7	Unsatisfactory quality of cortical surfaces	Insufficient quality of the MRI	Repeat step 6 on another MRI or manually correct the topological defects (see FreeSurfer website)
14	Severe misalignment of CT and MRI	Failure of the automatic CT conversion in step 12	Directly align to the ACPC system in step 11 by virtue of educated guesses of the commissure locations
	Imperfect alignment of CT and MRI	A left-right flip of either MRI or CT	Re-examine the native orientations of the MRI and CT in steps 3 and 10, and redo the preprocessing of the affected scan
	Imperfect alignment of CT and MRI	MRI and CT contain different head anatomies	Repeat step 13 with a different cost function (type help ft_volumerealign)
17	Electrodes hard to identify in the 2D ortho plot	Cortical grid orientation not aligned with any of the 2D planes	Identify electrodes in the 3D scatter figure (tick the scatter checkbox)
23	Severe deformation of the electrode grid	Incorrect pairing of neighboring electrodes in space	Repeat step 22 with alternate settings (type help ft_electroderealign)
27	Unsatisfactory quality of the volume-based registration	Insufficient quality of the MRI	Repeat step 25 with an alternate cost function or template version (type help ft_volumenormalise)
33	No anatomical label found	No overlap of electrode position with any anatomical mask	Increase the search radius around the electrode by increasing cfg.maxqueryrange (type help ft_volumelookup)