Participants
Participants were recruited through a prospective, observational cohort study at four institutions across the United States as a part of the Developmental Synaptopathies Consortium (Clinical Trial NCT02461420): Icahn School of Medicine at Mount Sinai, University of Texas Southwestern, Rush University Medical Center, and Boston Children’s Hospital. Stanford University and the National Institute of Mental Health also participated in the overarching study, but because they only collected phenotyping data and did not collect EEG, participants recruited at those institutions are not included here. In total, 31 individuals with PMS and 17 TD individuals had EEG completed. Participants with PMS were included if they had pathogenic deletions or mutations of the SHANK3 gene; clinical reports were reviewed to confirm this information. Typically-developing individuals were matched at the group level with PMS participants on chronological age and sex. TD individuals were excluded if they had a diagnosis of any intellectual disability, ASD, or other learning, developmental, psychiatric, or neurological disorders as determined by parent report. All participants were 4 to 19 years of age (inclusive). Informed written consent was obtained from legal guardians and assent was obtained from participants when appropriate. Table 1 shows demographics for participants with adequate EEG data for inclusion (see below).
Phenotypic Data
To examine how our EEG measures related to developmental abilities and ASD phenotypes among individuals with PMS, the following assessments were conducted: the Vineland Adaptive Behavior Scales (Vineland II): Survey Interview Form (37), the Autism Diagnostic Observation Schedule, 2nd edition (ADOS-2) (38), the Autism Diagnostic Interview-Revised (39), the Autism Diagnostic Criteria Checklist from the Diagnostic and Statistical Manual of Mental Disorders, 5th edition (40), the Short Sensory Processing Profile (SSP) (41), and the Repetitive Behavior Scale-Revised (RBS-R) (42). A psychologist determined ASD diagnosis either on the basis of the study’s assessments or clinical experience when the participant was seen clinically on a regular basis. Additionally, to assess non-verbal cognitive ability, participants were either given the Mullen Scales of Early Learning (MSEL) (43) the Stanford Binet-5 (SB-5) (44), or the Differential Ability Scales, 2nd edition (DAS-II) (45). Data for a non-verbal intelligence quotient (NVIQ) was compiled depending on the test given: for participants given the SB-5, the sum of the non-verbal subtests were taken; for participants given the MSEL, the mean of the visual reception developmental quotient (visual reception age equivalent score / age in months) and fine motor developmental quotient (fine motor age equivalent score / age in months) was taken; for participants given the DAS-II, the non-verbal reasoning standard score was taken. Finally, for participants who had experienced seizures, a seizure history was collected.
EEG Acquisition / Processing
Continuous EEG was collected for up to 10 minutes. Participants viewed a silent movie of their choice during EEG recording as is common practice in individuals with neurodevelopmental disorders (47). EEG was recorded using either a 128-channel Hydrocel Geodesic Sensor Net or a 32 channel ActiveTwo Biosemi net. Data were sampled at either 512Hz or 1000Hz (all files were later resampled to 250Hz). Impedances were kept below the recommendations for the specific EEG system being used prior to recording. For a subset of individuals with PMS, continuous EEG was again collected approximately 12 weeks after the initial recording. These subsequent recordings were used in place of initial recordings if the initial recording did not meet data quality thresholds (n = 1); all other analyses were performed using the initial recording.
Files were processed using the Batch EEG Automated Processing Platform (BEAPP) (48). Within BEAPP, the Harvard Automated Preprocessing Pipeline for EEG (HAPPE), which was developed specifically to optimize preprocessing of developmental EEG data with potentially high levels of artifact and short recordings, was used to automate preprocessing and artifact minimization (49). Data were first filtered using a 1 Hz high-pass filter and a 100 Hz low-pass filter. Data were then downsampled to 250 Hz for optimal performance of the HAPPE pipeline. With the exception of Cz, which was used as a reference electrode in some sites’ systems, only electrodes in the international 10-20 system were included in this analysis (18 total) to allow standardization of analyses across net types. Epochs of signal with any channel’s amplitude >40 μV (the HAPPE default threshold, reflecting the reduced signal amplitude that results from wavelet-thresholding and independent components analysis in HAPPE) were removed. EEG recordings were removed from further analysis if they exceeded thresholds for HAPPE data quality as per (50) in one or more of the following output parameters: percent good channels, mean retained artifact probability, median retained artifact probability, percent of independent components rejected, and percent variance retained after artifact removal. Data were subsequently re-referenced using an average reference, and then segmented into 2 second windows for power and PAC analysis. For each participant, 150 segments (300 seconds of data) were randomly selected; files with fewer than 150 segments of data at this stage were not analyzed. Primary power and PAC metrics were then obtained using code added to the BEAPP software.
Power Analyses
Power was computed across frequencies using a three taper multitaper window (51). Power was then computed for a number of frequency bands: Delta [1 – 4Hz), Theta [4 – 8Hz), Alpha [8 – 12Hz), Beta [12 – 30Hz), and Gamma [30 – 55Hz). Total power was computed as all frequencies between [1 – 55Hz]. The power at each frequency band and the overall 1-55Hz range was computed by adding the power spectral density over the frequency range of interest.
To capture each frequency band’s relative contribution to total power, the relative power at each frequency band was computed as the power at each frequency band divided by the total power. Power values were then averaged across electrodes. Visual inspection of the power spectra, averaged across the occipital channels analyzed in this study (O1 and O2), was used to identify the peak alpha frequency of each participant.
PAC analysis
Modulation Index: To capture the presence of coupling, PAC was first quantified using the Modulation Index (MI) (52). Because the data are not time locked to any specific task, we focus on PAC in the alpha-gamma range, where prior studies have shown abnormalities in other neurodevelopmental disorders using resting or non-time-locked data (23,24). For each frequency pair, the raw signal in each segment was exported from MATLAB into Python and filtered into a range of alpha (8-12 Hz in 2 Hz steps) and gamma (here, 28-56 Hz to allow for division into 4 Hz steps) frequencies using code adapted from Dupré la Tour et al., 2017. Alpha frequencies were filtered using a constant bandwidth of 2 Hz, while gamma frequencies were filtered using an upper sideband variable bandwidth, so as to avoid including phase frequencies in the amplitude frequencies. In detail, for each gamma frequency, the lower passband cutoff was 2 Hz below the gamma frequency, and the upper passband cutoff was set as the alpha phase frequency plus the gamma amplitude frequency (24). For example, for the combination of 40 Hz gamma and 8 Hz alpha, the lower limit of the gamma amplitude filter was 38 Hz (40 – 2), while the upper limit of the filter was 48 Hz (40 + 8). Filtering at this step consisted of a zero-phase cosine-based filter to extract the real component, and then a sine-based filter to extract the imaginary component, resulting in a complex-valued output signal (Dupré la Tour et al., 2017). The alpha phase time series, or gamma amplitude time series, were obtained from this complex signal. The phases of the alpha signal were then binned into 18 20o intervals (-180o to 180o), and the mean of the amplitude of the gamma signal occurring within each phase bin was calculated. Mean gamma amplitude values in each phase bin were then normalized by dividing each bin value by the sum of all bin values. Data were then imported into MATLAB, where the amplitude of the gamma signal at each phase bin of the alpha signal was then averaged together across segments. The MIraw was then computed as the Kullback-Leibler divergence of the gamma amplitude distribution from a uniform distribution (52). We then employed a time-shift procedure to control for factors that may generate spurious phase-amplitude coupling. In detail, for each participant, 200 surrogate MI values (MIsurr) were generated by repeating the procedure after offsetting gamma amplitude from the alpha phase distribution by a randomized time shift between 0.1 to 1.9 seconds. A normalized MI (z-MI) was then computed as the z-score of the MIraw compared to the distribution of MIsurr values (54). The z-MI at each alpha and gamma frequency combination was then averaged to obtain a single overall alpha-gamma PAC value for each participant, at each electrode.
Phase Bias: Modulation Index captures the extent of coupling. We additionally set out to quantify whether gamma amplitude increased closer to the rising or falling phase of the alpha waveform, and to what degree. To do so, we employed a metric termed phase bias, drawing on measures of phase preference (32) and prior findings that this tends to show a bimodal distribution (i.e., with fast activity occurring maximally at either the negative or positive phases (corresponding to rising and falling phases, respectively) of the slower waveform (55). Specifically, we quantified the phase bias of the gamma amplitude to the positive phases of the alpha waveform; i.e, the relative change of gamma amplitude (gammaamp) during the positive phases (0 – 180o) of the alpha waveform. Thus, phase bias is calculated as (Σgammaamp in positive phases of the alpha waveform) / Σ(gammaamp in all phases of the alpha waveform) - 0.5. Importantly, here, a cosine-based filter was used to extract alpha phase; as a result, 0o corresponds to the peak of alpha, +90o corresponds to the falling zero crossing of alpha, 180o / -180o corresponds to the trough of alpha, and -90o corresponds to the rising zero crossing of alpha. Therefore, a phase bias >0 indicates gamma amplitude increases at the falling phase of the alpha waveform, and a phase bias <0 indicates gamma amplitude increases at the rising phase of the alpha waveform. Additionally, a larger distance from 0 (where gamma amplitude does not increase preferentially at either positive or negative phases of alpha) indicates stronger phase bias. The phase bias at each alpha frequency and gamma high frequency combination was then averaged to obtain a single overall alpha-gamma phase bias value.
Statistical Analysis
Group Comparisons: We first set out to test whether power or phase-amplitude coupling metrics differed between groups. Because most metrics were not normally distributed, all group comparisons were performed using a non-parametric test (independent samples Mann-Whitney U) unless otherwise specified. Relative power in each frequency band was compared between groups, and an independent samples t-test was used to compare peak alpha frequency between groups. To test whether overall PAC metrics differed in individuals with PMS as compared to typically developing individuals, group comparisons were first performed on z-MI and phase bias data averaged across all 10-20 electrodes. Subsequently, because PAC has been shown to differ between anterior and posterior scalp areas (33), these group comparisons were repeated after averaging PAC metrics across all anterior 10-20 electrodes (Fp1, Fp2, F3, F4, F7, F8, Fz) and then posterior 10-20 electrodes (P3, P4, P7, P8, Pz, O1, O2). Finally, these comparisons of overall, anterior, and posterior z-MI and phase bias were repeated between individuals with PMS diagnosed with ASD (N=11), and individuals with PMS diagnosed without ASD (N=14). Data were analyzed in SPSS (IBM Corp, 2016).
Clinical Associations: All associations were performed using linear regression analysis. Because PAC has been shown to change with age (33), we tested whether age was associated with PAC metrics among all participants. Additionally, to test whether the relationship between ln(z-MI) and age was different in individuals with PMS as compared to TD individuals, a regression was performed, with ln(z-MI) as the dependent variable, and age, group, and age by group included as independent variables. The association between alpha power and PAC metrics (averaged across all electrodes) was additionally examined. To test how PAC associated with behavioral phenotype in individuals with PMS, linear regression analysis was performed between PAC metrics (z-MI and phase bias) and the following measures: Vineland Adaptive Behavior Composite, Vineland Socialization Composite, ADOS comparison score, SSP, RBS-R, and NVIQ. Additional linear regressions were performed between z-MI and the 6 behavior sub-scales of the RBS-R (Restricted Interest, Sameness, Ritualistic, Compulsive, Self-Injurious, and Stereotypic). Because z-MI did not demonstrate a normal distribution, linear regressions were performed on the natural log transformation of z-MI; one negative z-MI value was not included in this analysis. Age was included as a control variable in all regressions.
Clinical Comparisons: In individuals with PMS, we tested whether PAC measures differed by a number of categorical clinical variables, including: sex, presence of an ASD diagnosis, presence of a seizure history (at least one seizure event experienced), and whether the participant has a SHANK3 mutation or deletion. All comparisons were computed using a Mann-Whitney U test. For all associations between EEG and clinical measures, a Benjamini-Hochberg correction was applied to power, ln(z-MI) and phase bias clinical correlations separately (FDR = .1).
Participants with Insufficient EEG data
In total, there were 33 individuals with PMS and 17 TD individuals in the study. After removal of participants who did not complete EEG (n = 2 with PMS), had insufficient data quality (n = 4 with PMS, n = 1 TD), or had insufficient data length (n = 1 with PMS and n = 1 with TD), 26 individuals with PMS and 15 TD individuals remained for further analysis. Compared to PMS participants included in this dataset, the 7 PMS participants excluded for unusable EEG data were more likely to be male (6/7). Otherwise, they were not significantly different in age (mean = 9.94 , SD = 4.58, p = 0.7587), they demonstrated similar prevalence of ASD diagnosis (4/7) and seizure history (2/7) and there were no significant differences in included vs. excluded participants with PMS on the Vineland Adaptive Behavior Composite, the Vineland Communication Composite, ADOS severity score, NVIQ, SSP, or RBS-R (p > .05). All tests were performed using an independent samples Mann-Whitney U test.