Human subjects

This study was performed in epilepsy patients undergoing surgical evaluation with iEEG monitoring at New York University Medical Center (NYUMC). The protocol was approved by the NYUMC Institutional Review Board and the Clinical Trials Registration number was NCT02263274 (www.clinicaltrials.gov). Subjects were eligible according to pre-established criteria, including: (1) age over 18 years old; (2) undergoing invasive monitoring for seizure localization for epilepsy surgery; and (3) ability to provide informed consent or have a legal guardian who could consent. Exclusion criteria included (1) significant cognitive impairment (IQ < 70), (2) facial or forehead skin breakdown that would interfere with surface electrode placement, (3) contraindication to MRI, (4) known adhesive allergy, (5) space occupying lesion, and (6) subjects who had an electrographic seizure for 1 h prior to the stimulation procedure. All patients (n = 17) or their caregivers provided informed consent. Subjects were enrolled between December 2013 and May 2017. A table listing subject characteristics is included in Supplementary Table 1.

Sleep staging

Stimulation was performed after patients had entered at least 5 min of continuous NREM sleep, during a daytime nap or nocturnal sleep. Initial sleep staging was performed by visual online analysis, for the presence of a slow-wave and spindle activity as detected in the real-time iEEG seen at the bedside by a physician board-certified in clinical neurophysiology, as well as by direct clinical observation. As part of standard practice at NYU, an extradural lead is customarily placed near the vertex of the craniotomy to aid with spindle detection. When the patient aroused or drifted into a lighter stage of sleep, stimulation was stopped. The iEEG segments were later confirmed by a second board-certified physician to be consistent with N2 and N3 sleep. During offline processing, we selected segments of NREM sleep for analysis by comparing raw spectrograms of sleep (after artifact subtraction), to demonstrate that the depth of NREM sleep is similar across testing conditions (Supplementary Fig. 6).

iEEG recordings

iEEG was recorded from implanted subdural platinum-iridium electrodes embedded in silastic sheets (2.3 mm diameter contacts, 10 mm center–center spacing, Ad-Tech Medical Instrument, Racine, WI) or depth electrodes (1.1 mm diameter, 5–10 mm center–center spacing). The decision to implant, placement or recording electrodes, and the duration of invasive monitoring were determined solely on clinical grounds and without reference to this study. Electrodes were arranged as grid arrays (8 × 8 contacts, 10 or 5 mm center-to-center spacing), linear strips (1 × 4 to 12 contacts), or depth electrodes (1 × 8 or 12 contacts), or some combination thereof. Subdural electrodes covered extensive portions of lateral and medial frontal, parietal, occipital, and temporal cortex of the left and/or right hemisphere.

Within 24 h after surgical implantation of electrodes, patients underwent a post-operative brain MRI to confirm subdural electrode placement. Electrode were localized and mapped onto the pre-implant and post-implant MRI using geometric models of the electrode strips/girds and the cortical surface72.

Here, we present an efficient method to accurately localize intracranial electrode arrays based on pre-implantation and post-implantation MR images that incorporates array geometry and the individual’s cortical surface.

Clinical (macroelectrode) recording equipment

Recordings from grid, strip, and depth electrode arrays were made using a NicoletOne C64 clinical amplifier (Natus Neurologics, Middleton, WI), bandpass filtered from 0.16–250 Hz and digitized at 512 Hz. ECoG signals were referenced to a two-contact subdural strip facing toward the skull near the craniotomy site. A similar two-contact strip screwed to the skull was used for the instrument ground.

NeuroConn DC stimulator

The DC-STIMULATOR PLUS (NeuroConn, Germany) is a CE-certified medical device for conducting noninvasive TES in humans. The stimulator is a micro-processor-controlled constant current source, which continuously monitors electrode impedance, and detect insufficient contact with the skin. The device is battery powered, and therefore electrically isolated from the clinical recording electrodes and equipment.

Low-frequency tACS

We performed 0.75 and 1 Hz sinusoidal tACS on 13 epileptic patients with implanted subdural and depth electrodes. Seven subjects were stimulated during NREM sleep (four daytime nap; three nocturnal sleep) and six subjects were stimulated during waking rest, eyes closed. Patients were over 18 years old and fluent in English. Subjects were excluded if they had frequent (>2) electroclinical seizures in the 24 h preceding stimulation. Patient characteristics and electrode coverage are summarized in Supplementary Table 1 and Fig. 1.

All subjects tolerated scalp stimulation. All subjects who were stimulated during night-time sleep (N = 3, subjects 11–13) and most subjects during an afternoon nap (N = 4, subjects 7–10) were able to sleep through trials at stimulation intensities between 0.5 and 2 mA. One subject (subject 13) woke from sleep and reported an itching sensation during one stimulation block with 2.5 mA current intensity. For the tACS experiments, we recorded and analyzed from 1700 electrodes without artifacts (mean of 131 electrodes per subject, example electrode in Fig. 3b). There were no complications from stimulation, and no induced electrographic seizures. One patient (subject 7) had a typical electroclinical seizure during stimulation. Because this patient had frequent spontaneous seizures, it was determined by the patient’s epileptologist that stimulation was unlikely to have caused the seizure.

Furthermore, we enrolled one subject (S14) who had a bilateral subdural strip and depth survey, to perform a precise replication experiment of prior protocols. This patient had multiple target clinical events captured, which were non-epileptic in nature. He did not have any interictal or ictal activity captured during 1 week of monitoring, even while medications were being withdrawn. In other words, the patient did not demonstrate any epilepsy-related pathophysiology.

We reviewed the hour of iEEG recording prior to stimulation to exclude recent seizures. We performed a pre-stimulation clinical assessment (including assessment of the stimulation skin site and neurologic examination). A physician (AL) was present at the bedside during the entire procedure to monitor for safety. The patient’s iEEG recording was monitored in real time at the bedside during stimulation for seizures.

For patients S1–S13, two stimulating electrodes were placed medially over the frontal and occipital poles (2 cm × 2 cm rubber electrodes) for patients S1–S10. In patients S10–S13 the frontal electrode was offset from midline by 3 cm (S10 and S13 left frontal; S11 and S12 right frontal) to minimize distance from stimulating electrodes to recording electrodes. In one patient S10, the stimulation electrodes were placed at three additional locations to examine the effects of electrode placement on extracerebral and intracerebral current spread. Subjects were covered with a nickel-cadmium shroud to reduce environmental artifact during recording, and other sources of environmental noise (60 Hz) were minimized in the patient area. The stimulation protocol used the NeuroConn DC Stimulator Plus (NeuroConn, Germany), with a biphasic sinusoidal waveform at 0.75 and 1 Hz, at variable intensities between 0.3 and 2 mA, for 10 s (cycles) to determine the peak intensity at which amplifier saturation occurred. Thereafter, subjects were stimulated with TES at 0.75 and 1 Hz, at variable intensities up to the peak intensity, for a duration between 5 and 10 min (10 min for A1–A6 and 5 min for S1–S4). Up to four blocks of stimulation were applied, until subjects woke up. The more than ten-fold increase of the subdurally recorded iEEG amplitude, compared to the EEG signal73, allowed for simultaneous recording and stimulation (up to saturation levels of the amplifiers). Stimulation was immediately stopped in the event of an electrographic seizure (S7). A repeat clinical assessment (including assessment of stimulation skin site and neurologic examination) was performed after stimulation.

For S14, who was enrolled to perform a replication experiment of prior protocols, we selected a surgical patient who had a bilateral strip and depth survey. There were two windows that were cut into the patient’s surgical bandage to allow electrode placement at the F3/F4 positions. Stimulation electrodes (8 mm Ag/Cl ring type) were applied bilaterally, with anodes at F3/F4 and cathodes on each mastoid. To test for acute effects, we utilized a stimulation protocol using trapezoidal waveform (0.33 s ramp up/0.33 s steady state/0.33 s ramp down/0.33 s zero current), 0.75 Hz, 0–0.26 mA, for 8 cycles on/8 cycles off (10.66 s ON/10.66 s OFF), for 80 cycles. Stimulation was started after the first 5 min of NREM sleep.

Low-frequency noise-burst auditory stimulation

Auditory stimulation consisted of 50 ms pink-noise bursts (1/f spectrum) repeated at a rate of 1 of 0.75 Hz to four patients (subjects A1–A4). Sound was delivered via flat-profile headphones (Bedphones, Millwood, NY), which were placed on the patient’s ears. Placement of the headphones required access under the surgical bandage and was administered by an epilepsy physician (AL), with patient verification of correct positioning. Acoustic pulses (50 ms, pink noise, 5 ms on/off ramps) were digitally generated and delivered via a laptop placed at the bedside. The sound level of the stimulation was manually adjusted for each patient to maximize comfort (e.g., ability to sleep with sound playing in the background). The volume was recorded and estimates of the sound level presented to each subject were assessed via an ear simulator (KEMAR Head and Torso simulator, Knowles Research, coupled to a B&K type 3134 Pressure microphone and B&K type 2230 Sound Level Meter, Bruel & Kjaer, Denmark). Resultant peak sound pressure level estimates for each subject are 72 dB SPL (subject A1), 68 dB SPL (subject A2), 46 dB SPL (subject A3), 70 dB SPL (subject A4).

Similar to the procedure used during tACS stimulation, acoustic stimulation was presented at repetition rates of 0.75 and 1 Hz, in blocks of 5 min during NREM sleep. For each subject, we collected a block at each stimulation rate during NREM sleep, which was visually confirmed offline. The order of the presentation blocks was randomized across subjects. An awake control condition was performed for each patient, to verify the presence of acoustic ERPs. Additional time-locked TTL triggers were generated for each stimulus presentation and recorded by the EEG amplifier’s DC input to aid in offline analysis.

iEEG data preprocessing

All electrodes were inspected for signal quality by plotting spectrogram, raw voltage, and the power spectrum. We recorded from a total of 1700 electrodes for this analysis. Electrodes were discarded based on high 60 Hz noise (likely due to poor contact impedance), amplifier saturation (clipping), or poor removal of the tACS artifact (due to non-stationarity of stimulation artefact, typically resulting from patient movement). Example of artefact-free recording electrodes are shown in Supplementary Figs. 2 and 3. The electrodes remaining for each subject during tACS were: 122/126 (S1), 112/126 (S2), 117/126 (S3), 78/84 (S4), 111/124 (S5), 155/254 (S6), 80/128 (S7), 116/122 (S8), 89/100 (S9, 1 Hz tACS), 93/100 (S9, 0.75 Hz), 83/124 (S10), 49/98 (S11), 79/102 (S12), 88/188 (S13, 1 Hz tACS), 93/188 (S13, 0.75 Hz tACS). During endogenous sleep the electrodes remaining for each subject during night 1 were: 123/126 (S7), 103/122 (S8), 96/100 (S9), 91/124 (S10), 70/98 (S11), 94/102 (S12), 170/188 (S13). During night 2: 121/126 (S7), 113/122 (S8), 89/100 (S9), 96/124 (S10), 82/98 (S11), 93/102 (S12), 164/188 (S13).

Measurement and modeling of electric fields

During tACS, the current alternates in directionality between two stimulating electrodes. This alternation results in a sinusoidal signal that can be used to determine the magnitude of the stimulation voltages. Magnitude was estimated by fitting a sinusoid to the signal at each electrode location and estimating the magnitude of the fitted signal. The output of this processing was plotted and manually inspected electrode by electrode.

The measured voltage in each location is then used to derive the projected electric field, by subtracting potential values between adjacent electrode pairs and dividing by their distance, resulting in V/m. The adjacent electrode was defined as the closest electrode within a 10 mm vicinity for cortical electrodes on the same grid array and linear strip, and 5 mm vicinity for depth electrodes on the same strip, to reflect the different inter-electrode distance. It is important to realize that this only captures a fraction of the field magnitude at any given location as the field orientation may not be parallel to the direction of two neighboring electrodes. The distant stimulating electrodes were expected to generate the strongest field intensities on the cortical surface directly under the scalp electrodes40. However, recording electrodes lay predominantly orthogonal to field direction (parallel to the cortical surface). Thus, the measured field projections will not capture maximal intensities, except in the rare circumstance that a depth array is precisely underneath one of the stimulating electrode and oriented toward a second, distant stimulating electrode.

The computational models were built following our previous work74. Briefly, the MRI for each patient was automatically segmented by the New Segment toolbox75 in Statistical Parametric Mapping 8 (Wellcome Trust Centre for Neuroimaging, London, UK) in Matlab (R2013a, MathWorks, Natick, MA). Segmentation errors were corrected first by a customized Matlab script74 and then by hand in ScanIP software (v4.2, Simpleware Ltd., Exeter, UK). The field of view of the clinical MRI scans was extended down to the neck by co-registering a standard head74, and pasting the lower portion of the standard head to the model. The 2 × 2 cm stimulation electrodes were positioned on the model using CAD software. For each patient, a finite element model was generated from the segmentation data and then the electric potential distribution was computed assuming 1 mA current through the stimulation electrodes. Tissue conductivities were adjusted to minimize the mean-square difference between predicted and measured potentials. With these calibrated models, we then computed electric fields throughout the brain. Electric potentials of model and measurements corresponded closely, with correlation values of r = 0.95 ± 0.04 (mean ± Std across patients, N = 1545 electrodes across 10 subjects). Electric field is the spatial derivative of these potentials. They are estimated as the difference in electric potential between neighboring electrodes, divided by the distance. This is the electric field projected on the orientation of the electrode pair40.

Phase-amplitude coupling

PAC measures the degree to which the amplitude of a high-frequency oscillations, A HF (t), is aligned with the phase of a lower frequency, ϕ HF (t). We were interested in the interaction between the amplitude of spindle activity band at 14 Hz with the phase of endogenous slow oscillations at 1 Hz (or 0.75 Hz) during sleep as well as entrainment to the applied stimulation (tACS and acoustic). This section refers to 14 Hz activity, but the identical analysis was done for power amplitudes in the theta, alpha, and gamma bands (see next paragraph). To measure entrainment of spindle oscillations during endogenous sleep, we compared spindle power against the phase of the endogenous slow oscillation activity (ϕ LFE (t), low-frequency endogenous). To measure entrainment to tACS, we used the phase of the electrical stimulation artifact (ϕ LFS (t), low-frequency stimulation; Supplementary Fig. 1). To obtain the phase during acoustic stimulation we used the delay from the onset of each noise burst. To obtain the phase during tACS and remove the stimulation artifact, we first modeled the 1 Hz artifact as a linear superposition of sines and cosines at multiples of a base frequency (harmonics up to 40 Hz) by fitting the amplitude of each sine/cosine and the base frequency. An example of this fitting procedure is shown in Supplementary Fig. 4 with the top row indicating the raw signals and the bottom row showing the signals after the fitted harmonic artifact has been subtracted. The resulting harmonic fit captures the 1 Hz stimulation artifact including any harmonic distortion that may have resulted from amplifier nonlinearities. We calculated the stimulation phase from the harmonic fit (Supplementary Fig. 1D). To obtain the phase during endogenous sleep, first we applied a complex-valued Morlet wavelet filter centered at 1 Hz (or 0.75 Hz) with a bandwidth of 1 Hz (or 0.75 Hz; in humans, the center frequency of slow-wave activity is often assumed to be <1 Hz)48. The instantaneous phase can be directly extracted from the complex-valued filtered signal. Peak and trough of the slow oscillation are indicated by 0° and 180°, respectively, which represent the cortical upstate and downstate as discussed in the main text.

To obtain the instantaneous amplitude of the high-frequency rhythm during tACS we subtract the harmonic fit (Supplementary Fig. 4) and filtered the residual iEEG signal with a complex-valued Morlet bandpass for spindle activity (fc = 14 Hz, bandwidth = 7 Hz), for alpha activity (fc = 10 Hz, 5 Hz bandwidth) and for gamma activity (fc = 90 Hz, bandwidth = 40 Hz, Supplementary Fig. 1C). The instantaneous HF power was obtained by taking the absolute value of the filtered signal (Supplementary Fig. 1E). Patient movement and other artifacts resulted in outliers in the HF amplitude estimates. We removed these by removing cycles with excessive power as follows. For each electrode and each LF cycle we compute the mean power of the HF band (square amplitude averaged over one LF cycle). Cycles are removed as outliers if their mean power exceeds two times the interquartile range across all cycles in an electrode. To obtain the amplitude of the high-frequency rhythm during endogenous sleep and auditory stimulation we used the same procedures starting with the raw iEEG signal (no harmonic fit is needed). Outlier rejection was done as before based on the mean power in the high-frequency band. In addition, for endogenous slow wave we excluded cycles for which the slow oscillation amplitude was below 50 µV. For the analysis of the acoustic stimulation we used the same procedures as with the endogenous sleep except that the phase was defined based on the time since the onset of the noise burst, and cycle duration determined from the average inter-stimulus interval (TTL pulse) for each subject.

PAC is measured here using the modulation index r, which is defined as the absolute value of the time average, r = | < z(t)> |, of the complex-valued quantity, \(z(t) = {A_{HF}}(t){e^{i{\phi _{LF}}(t)}}\) (Supplementary Fig. 1F). Time average < z(t) > is computed over cycles and time within a cycle. If z(t) does not have a radially symmetric distribution this will cause the time average r (modulation index) to be different from zero. This can be a result of (1) the amplitude A HF (t) is consistently higher at a certain phase, or (2) ϕ LF (t) is not uniformly distributed in time. The phase is uniformly distributed in the case of tACS and acoustic stimulation. However, during endogenous sleep, phase is extracted from the slow oscillations, which is not sinusoidal resulting in non-uniform phase distributions. Consequently, we applied a histogram equalization to the phase distribution, ensuring that non-zero modulation index is only a result of modulated amplitude coupling, which we confirmed by testing for significance using surrogate data with constant HF amplitude.

Randomized surrogate data to estimate statistical significance of PAC

Significance was determined by randomizing the phase of each LF cycle and thus creating the distribution of the modulation index r, under the null hypothesis of no PAC. Phase randomization makes no assumptions on the distribution of the data, except that the phase is uniform, which has been addressed in the paragraph above. For each slow oscillation cycle the low-frequency phase, ϕ LF (t), was incremented by a random value uniformly distributed between 0° and 180°. The phase was shifted by the same random phase for all electrodes but independently for different cycles. The randomizing procedure was repeated 10,000 times and the p-value was measured as the rate of a random phase having a modulation index higher than the original data. The minimum numerical p-value possible, given the number of randomizations, was 10−4 (1/number of shuffles). All the p-values were corrected for multiple comparisons across electrodes using FDR correction76 with q < 0.05. No correction was performed across segments (as we use segments to determine how reliable a potential entrainment is over longer periods of time) or frequency bands (as these were planned comparison). Uncorrected p-values for electrodes that had p < 0.05 in at least one stimulation block are shown in Supplementary Figs. 7–9 for tACS, endogenous slow-wave (no stimulation), and acoustic stimulation, respectively.

Note that we examine the acute effects of modulation within a single cycle, thus permitting many single “trials” (N = 300 cycles for 1 Hz and N = 225 for 0.75 Hz within each 5-min block, repeated 2–6 times for each stimulation intensity; see number of blocks tested for each subject in Supplementary Table 1). This large number of trials permits, in principle, detection of small changes in power within a cycle, on the order of a few percent (e.g., assuming independent noise: sqrt(1/300) = 5.6%). See Supplementary Note 3 for an extended discussion on statistical power.

To aid comparison with endogenous sleep and tACS stimulation condition, we analyzed entrainment for different sleep duration segments. The longer the recordings (approximately 20 min), the more accurate the mean vector strength technique is at measuring entrainment (approximately 50% of electrodes entrained). For shorter durations of endogenous sleep (5 min of iEEG) comparable to the duration of tACS blocks, there is still modulation of spindle activity but in only 10% of the electrodes.

Spindle event detection

Spindle detection follows existing methods50. Briefly, the signal is bandpass filtered in the spindle band as above, and instantaneous power is obtained as the square amplitude of the filtered signal. For each channel a detection threshold is defined as six times the median of the instantaneous power. For segments that cross this threshold a lower threshold (1 standard deviation of the power) is applied to detect onset and offset of the spindle event. Only spindles whose duration was between 0.2 and 2 s were considered for further analysis50, 77. To avoid false positive detections due to patient movements or high-frequency artifacts spindle events that coincide with increases in broadband power were discarded. Events exhibiting broadband power increases (p > 0.10, comparing the maximum broadband power in that event vs. the distribution of broadband power during the whole recording) were excluded.

Slow-wave detection during auditory stimulation

Additional analyses were performed on the acoustic stimulation data to assess the physiological origin of the PAC entrainment. First, a slow-wave detection algorithm was applied during stimulation segments to assess whether underlying slow-wave activity was altered during stimulation. In contrast to tACS, this is possible for auditory stimulation as there are no electric stimulation artefacts in the iEEG signal. Following previous literature61, 78, slow-wave detection consisted of bandpass filtering the waveform (fc = 0.75 Hz, bandwidth = 0.75 Hz) and using a zero-crossing algorithm to identify events where two subsequent negative zero-crossing (e.g., from positive to negative) were within the range of 0.5–2 s (2–0.5 Hz, bandwidth = 1.5 Hz). Subsequently, once a slow-wave event was identified, the down-state was identified as the minimum voltage within this event. Second, time-frequency response (TFR) functions were calculated to assess whether stimulation trials consisted of physiological sleep spindles. TFRs were computed using a 6-cycle Morlet wavelets between 5 and 25 Hz and averaging the resulting spectrograms (locked to stimulation onset) across all noise bursts (Fig. 6). Significance between conditions (sleep vs. baseline) was assessed via paired sample t-test (one-tailed, p < 0.05) compared across all stimulation trials (n = 300).

Analysis of changes in power before and after stimulation

Qualitative comparison between the NREM during the two nights without stimulation and the stimulation period show no evident changes in the power spectrum (Supplementary Fig. 6). To determine if there were significant changes in power at the different frequency bands before vs. after stimulation we used the Chronux toolbox79 (http://chronux.org/; version 2.12). Briefly, for each electrode, 30 s immediately after each 5-min stimulation block were compared against 30 s preceding the first block. Differences in power in the slow oscillation (0.5–1 Hz) and spindle bands (10 and 14 Hz) are shown in Supplementary Fig. 12. Statistical significance of power changes is computed with the Chronux toolbox for the fast spindle band and slow oscillation band and FDR corrected (q < 0.05).

Entrainment with electric stimulation assessed after stimulation block

In subject S14 we collected enough stimulation blocks to evaluate slow-wave entrainment in the intervals immediately following tACS. For this, slow-wave entrainment was tested following previous reports37. Briefly, the 10-s stimulation-free interval immediately after stimulation was fit with a sinusoid (0.5–1.5 Hz) and phase coherence was calculated across all trials. Statistical significance was calculated using a circular test statistic (Rayleigh test for non-uniformity). In addition, to test for after-effects, we utilized the same trapezoidal waveform, frequency, and intensity in five separate 5-min blocks during NREM sleep (5 min ON/5 min OFF). We compared these sessions against two control periods of endogenous NREM sleep for spindle (10, 14 Hz) entrainment. The power in the stimulation-free intervals was calculated for slow oscillations (0.5–1 Hz), slow spindle activity (8–12 Hz), and fast spindle activity (12–15 Hz). This quantity was then compared to the power during a different night, when the subject was in an equal sleep state. Slow oscillations and slow spindle activity were averaged across locations close to Fz. And fast spindle activity was determined across parietal locations.

Code availability

The code used to generate the main findings of the current study are available from the corresponding author upon reasonable request.

Data availability

The data sets generated during and/or analyzed during the current study are available from the corresponding author upon reasonable request.