Participants

Three patients were implanted with intracranial electrodes covering the lateral parietal surface on the right (participant 1, 64 electrodes and participant 3, 114 electrodes) or left (participant 2, 112 electrodes; Fig. 1). All three patients were implanted solely for clinical reasons and gave their informed consent to participate in the study, which was approved by the Stanford University Institutional Review Board. All patients were right handed, with normal IQ, with age 37, 43 and 45 years. In all participants, the lateral parietal electrodes of interest were void of pathological activity. Participant 1 was female and participants 2 and 3 were male. We emphasize that the aim of the current study was not to explore the lateralization of response or the effect of gender in numerical cognition, and thus our participants were not balanced on the basis of hemisphere coverage or gender.

Experimental task and stimuli

In the experimental task, participants judged with self-paced responses the accuracy of addition equations or memory statements. In the ‘arithmetic condition’, simple addition equations comprises adding one single digit with one double-digit number (n=48). These equations were randomly mixed with the ‘non-arithmetic condition’, which contained the memory statements without any numerical content (n=~258) and 5-s long fixation or rest trials (n=~32). The task was administered twice in each participant (either the same day (P1) or consecutive days (P2 and P3). Sentences had similar apparent visual length to the arithmetic equations and were presented in the same font and on the same screen background.

Electrode localization

We reconstructed individual three-dimensional brain images by aligning preoperative high-resolution T1-weighted magnetic resonance images (MRIs) with postoperative computed tomography (CT) scan images obtained after the electrodes implantation with an accuracy of ~±5 mm as described in our previous publications29,59,60,61. In brief, postoperative CT images were aligned to preoperative structural T1-weighted MRI whole-brain scans separately for each subject. Anatomical MRI data was reoriented to ‘AC–PC space’ by manually identifying the AC (anterior commissure), the PC (posterior commissure) and a third point in the midsagittal plane. MRI data were then resampled to 1 mm isotropic voxels using a b-spline image interpolation algorithm from SPM5 (http://www.fil.ion.ucl.ac.uk/spm). Postoperative CT images were then aligned to the T1 MRI anatomy scans using a mutual information algorithm from SPM5. After CT–MRI alignment, electrodes were identified in the coregistered CT image slices and their centroid coordinates recorded (subject T1 headspace).

Intracranial EEG recordings

Patients were implanted with flexible grid subdural platinum electrodes (AdTech Medical Instruments Corp.). The centre-to-centre interelectrode spacing was 10 mm and the diameter of each electrode’s exposed surface was 2.3 mm. Electrophysiological data were recorded simultaneously with a Tucker Davis Technologies (TDT) system and by Nihon Kohden Technology (NK) system. TDT data were synced with the display of a stimulus-presenting laptop by TTL pulses from a custom-built photodiode and NK system synced with a real-time online video camera. The TDT data were sampled at 3,052 Hz, whereas the NK data were sampled at 500 (P3) or 1,000 Hz (P1 and P2). One of the subdural electrodes served as a reference channel and an external ear lobule electrode as the ground contact for both TDT and NK systems.

Preprocessing and spectral decomposition

Neurologist specialized in epileptology (J.P.) examined the EEG data in all electrodes and excluded those with epileptic activity from further analysis. Notch filter from EEGLAB package (http://sccn.ucsd.edu/eeglab) was used to remove 60 Hz noise and its second and third harmonics. All EEG channels were re-referenced to a common average reference. EEG data were downsampled into 436 samples per second using MATLAB software. We aligned TDT and NK EEG data by calculating their maximum cross-correlation. Signals from each electrode contact were decomposed in two different ways. To make ERSP, we decomposed each channel into 42 bandwidths by filtering the signal by 42 custom-made band-pass filters. The central frequency of each band-pass filter was chosen from 1 to 232 Hz range with equal log-distance spacing. Second, each channel was also decomposed in δ (1–4 Hz), θ (4–8 Hz), α (8–12 Hz), β (15–25 Hz), low-γ (30–55 Hz) and high-γ (70–110 Hz) range by using band-pass filters from EEGLAB package. Finally, we applied Hilbert transform on decomposed signals to obtain the instantaneous amplitude of each frequency range. In this study, we call the instantaneous amplitudes of HFB range as the HFB trace.

Event-related spectral perturbation

The instantaneous power (square of amplitude) of each frequency was aligned and then averaged with respect to the stimulus presentation or to participant’s response time for each condition of the experiment. The trial length varied due to self-paced design of the experiment, and the number of trials in the arithmetic condition differed from the non-arithmetic condition. We matched the number of trials in the non-arithmetic condition to the arithmetic condition by choosing a subset of trials to equalize the median of trial length of both conditions. The statistical significance of power change was established by normalizing the averaged power with the mean and s.d. of the distribution of the surrogates of the entire experiment38. Each surrogate signal was produced by preserving the amplitude while randomizing the phase of the signal. Each time-frequency point in the surrogate ERSP matrix had a Gaussian distribution (central limit theorem) and, hence, Z-scores represented the significance of the change in power.

Significance test and d′-values

The trace of power change for each lateral parietal electrode was averaged over 1 s window (from 0.5 to 1.5 s after stimulus presentation) for every trial to calculate the mean-power in six bands of δ (1–4 Hz), θ (4–8 Hz), α (8–12 Hz), β (15–25 Hz), low γ (30–55 Hz) and high-γ (70–110 Hz). We calculated the d′- and P-values for all electrodes and for all frequency bands. Next, we normalized the mean power for each band by mean and s.d. of a generated surrogate mean-power distribution described above (the same as ERSP) by fitting the histogram of values to a Gaussian distribution with equal mean and s.d. (central limit theorem guarantees the goodness of fit). We used this mean power to test the significance of power change for each band, condition and electrode. We evaluated the null hypothesis that the mean of arithmetic trials was not different from that of non-arithmetic trials by applying non-paired Student’s t-test. We then quantified the extent to which the distribution of mean power in arithmetic trials was different from that of non-arithmetic trials by calculating the d′-values. The d′-value measures the distance between the means of these two distributions (centroid) in units of s.d. of non-arithmetic distribution. We sorted electrodes based on their high-γ d′-values (d′ of HFB) and visualized in a mean-maximum colour scale on each patient’s three-dimensional brain as depicted in Fig. 1f.

Characterization of HFB peaks

To quantify the occurrence of HFB activations, we defined the so-called HFB ‘peaks’ as discrete events in the HFB traces where the HFB instantaneous amplitude increased above a chosen threshold for amplitude and was sustained beyond a chosen duration threshold. Therefore, each HFB peak was determined by two threshold parameters: amplitude and duration (Fig. 2). Before searching for the optimal amplitude and duration thresholds, we smoothed the HFB traces with a standard Gaussian function (2 × σ=500 ms) to capture the local trend of HFB traces. As in the experimental task we knew which trials were arithmetic and which were non-arithmetic trials, we aimed to find the optimal values for duration and amplitude by which we could accurately identify the arithmetic trials from non-arithmetic trials. We used the ROC curve to identify the parameter that produced the highest sensitivity (hits/hits+misses) and specificity (correct rejection/correct rejection+false alarm). Hits were defined as the presence of HFB peak during arithmetic trials, misses as the absence of HFB peak during arithmetic trials, correct rejections as the absence of HFB peak in non-arithmetic trials and false alarms as the presence of HFB peak during non-arithmetic trials. The shape of the ROC curve of an electrode is determined by its d′-value; an electrode with d′- close to zero produces an ROC curve close to diagonal line (Fig. 2). On the ROC curve, the optimal amplitude threshold is the closest point to the upper left corner (a point with the highest sensitivity and specificity). For a given electrode, we tried a range of duration values from 500 ms (the width of smoothing Gaussian function) to 1 s for a given amplitude value in a step-wise fashion (Fig. 2a). The ROC curve with the largest area under the curve area determined the optimal duration value. We used one of the two sessions (training) of the first experiment to find the optimal amplitude and duration thresholds. Next, we used these values to separate the arithmetic trials from non-arithmetic trials in the other session (testing) of the first experiment and in the natural condition.

We use the term ‘natural condition’ even though being implanted with intracranial electrodes in a hospital may be different than the real ‘every-day’ life. Nevertheless, the daily life setting in the hospital can be considered ‘natural’, because patients were engaged in natural social interactions with other people as they would have done in any real-life setting, and our observations were made outside the artificial controlled experimental condition generated within a lab and without any obtrusive manipulations.

Selectivity of HFB response

We defined selectivity by the difference between HFB peak occurrence frequency during arithmetic and non-arithmetic conditions. We normalized this difference by summation of these two occurrence frequencies. We defined selectivity index as (FOA—FON)/(FOA+FON), where FOA is the frequency occurrence of HFB peaks during arithmetic condition and FON is the frequency occurrence of HFB peaks during non-arithmetic condition. The selectivity is a value between zero (no selectivity) and one (maximum selectivity). For example, the selectivity index value of 0.25 is equivalent to 66% increase in the number of HFB peak frequency during arithmetic condition.

Labelling natural events

We had access to 10 (P1)- and 6 (P2 and P3)-min-long simultaneous video and intracranial EEG data from the patients when they interacted with their environment. We chose a random video segment using the following specific inclusion criteria: (1) the video should be from the same day of recording as the experimental data to avoid the plausible problem of mismatch in the signal quality between the experimental and natural files; (2) patient must interact with others and sounds must be audible (if there was no interaction between the patient and others, we could not judge the behavioural correlates of the HFB responses); and (3) lack of seizures or significant epileptic activity in the rest of the brain within 2 h before or after video file to avoid pre-ictal or post-ictal changes in the signal quality.

Two independent reviewers (Stanford University undergraduate volunteers) evaluated the behavioural content of the words in conversations captured on the videos by time-stamping the start and the end of each word uttered during conversations (see Supplementary Table S1 and Supplementary Movie 1). Each reviewer then labelled the words that either had, or not, any numerical content.

HFB peaks during numerical events in natural condition

We selected videos in which the participants were awake and communicating with other people in their hospital environment. We defined the co-occurrence as the ratio of the number of HFB peaks with behavioural events in which words with numerical content were exchanged to the total number of HFB peaks during natural condition. To show that this co-occurrence ratio was not due to random alignment, we made a surrogate video by displacing the behavioural data from the end of the video to the beginning and then recalculated the co-occurrence ratio. It is noteworthy that the surrogate video data is only temporally misaligned but otherwise identical to the original video data. We made a Gaussian distribution of the co-occurrence ratio by making 500 surrogate videos. The P-value of the comparison between real and surrogate natural events is the probability of HFB peaks and numerical events being aligned by chance. We also estimated the chance level by averaging of the obtained co-occurrence values from the surrogate distribution.

The sparseness of HFB peaks during natural condition

The sparseness of the response during natural events was calculated by segmenting ECoG recording by non-overlapping 500 ms windows (as long as the minimum HFB peak duration). The sparseness62 of the activity of an electrode is calculated by Sparseness={1−(∑gi/n)^2/∑(gi^2/n)}/1–1/n}, where n is the total number of bins and gi is a bin that contains a HFB peak. The maximum sparseness is one when just one of the bins contains a HFB peak and the rest of bins are empty, and it is at minimum (zero) when all bins contain a HFB peak.