TY - JOUR AU - Classen, Joseph AB - Abstract Abnormal phase-amplitude coupling between β and broadband-γ activities has been identified in recordings from the cortex or scalp of patients with Parkinson’s disease. While enhanced phase-amplitude coupling has been proposed as a biomarker of Parkinson’s disease, the neuronal mechanisms underlying the abnormal coupling and its relationship to motor impairments in Parkinson’s disease remain unclear. To address these issues, we performed an in-depth analysis of high-density EEG recordings at rest in 19 patients with Parkinson’s disease and 20 age- and sex-matched healthy control subjects. EEG signals were projected onto the individual cortical surfaces using source reconstruction techniques and separated into spatiotemporal components using independent component analysis. Compared to healthy controls, phase-amplitude coupling of Parkinson’s disease patients was enhanced in dorsolateral prefrontal cortex, premotor cortex, primary motor cortex and somatosensory cortex, the difference being statistically significant in the hemisphere contralateral to the clinically more affected side. β and γ signals involved in generating abnormal phase-amplitude coupling were not strictly phase-phase coupled, ruling out that phase-amplitude coupling merely reflects the abnormal activity of a single oscillator in a recurrent network. We found important differences for couplings between the β and γ signals from identical components as opposed to those from different components (originating from distinct spatial locations). While both couplings were abnormally enhanced in patients, only the latter were correlated with clinical motor severity as indexed by part III of the Movement Disorder Society Unified Parkinson’s Disease Rating Scale. Correlations with parkinsonian motor symptoms of such inter-component couplings were found in premotor, primary motor and somatosensory cortex, but not in dorsolateral prefrontal cortex, suggesting motor domain specificity. The topography of phase-amplitude coupling demonstrated profound differences in patients compared to controls. These findings suggest, first, that enhanced phase-amplitude coupling in Parkinson’s disease patients originates from the coupling between distinct neural networks in several brain regions involved in motor control. Because these regions included the somatosensory cortex, abnormal phase-amplitude coupling is not exclusively tied to the hyperdirect tract connecting cortical regions monosynaptically with the subthalamic nucleus. Second, only the coupling between β and γ signals from different components appears to have pathophysiological significance, suggesting that therapeutic approaches breaking the abnormal lateral coupling between neuronal circuits may be more promising than targeting phase-amplitude coupling per se. Parkinson’s disease, phase-amplitude coupling, source analysis, spatiotemporal characteristics Introduction Over the past decade, exaggerated phase-amplitude coupling (PAC) between β and broadband-γ activities has been detected in subthalamic nucleus (STN), motor cortex, and interactions between the two regions of patients with Parkinson’s disease via subdural electrocorticography (ECoG) and local field potential (LFP) recordings (de Hemptinne et al., 2013; Kondylis et al., 2016). It has also been shown that enhanced PAC, as well as parkinsonian bradykinesia, can be reduced by dopaminergic therapy (Miller et al., 2019) or through STN deep brain stimulation (DBS) (De Hemptinne et al., 2015; van Wijk et al., 2016; Malekmohammadi et al., 2018), suggesting that abnormally enhanced PAC might represent a neural circuit abnormality that is close to the causative pathophysiological mechanism underlying parkinsonian bradykinesia. Recently, two EEG studies have shown that even non-invasive techniques can be used to differentiate between patients and healthy individuals, and that the state of drug-induced improvement of mobility in patients with Parkinson’s disease can be differentiated from the akinetic state by means of PAC (Swann et al., 2015; Miller et al., 2019). However, two issues related to abnormal PAC cannot yet be satisfactorily addressed using data recorded at EEG sensors: (i) the spatial localization of the brain regions from where the pathological coupling originates is incompletely explored. EEG sensor signals are poor in localizing the brain areas generating the underlying neuronal activity, because each signal represents a mixture of differently weighted contributions from numerous brain areas. This issue is important as the spatial distribution of abnormal PAC informs about its possible origin within anatomically segregated loops of the basal-ganglia-thalamocortical circuit (BGTC), which may relate to distinct motor and non-motor functions (DeLong and Wichmann, 2010); and (ii) the origin of PAC-involved β and γ activity at a local scale and its relationship with the pathophysiology of parkinsonian motor impairment remain unclear. Both slow and fast activities could originate from a single oscillator, in which slow and fast rhythms are nested and tightly phase-coupled (Gast et al., 2020). Alternatively, slow and fast activities may originate from relatively independent neural circuits co-located in the same brain region, in which β and γ networks are only weakly coupled. In this scenario, the β-signal would only modulate the gain of the γ network, resulting in PAC without phase-phase coupling (Schroeder and Lakatos, 2009). However, although the spatial distributions of these two oscillators would be different, it is not a priori clear that they can be distinguished within the methodological limits of scalp EEG. Each of these different scenarios may have important implications for the mechanisms of PAC and its pathophysiological significance. EEG source localization techniques, based on detailed and individual biophysical head models, offer the possibility to reconstruct the original source signals and thereby reveal the spatiotemporal architecture of the underlying brain activity. Here we used advanced source analysis to more precisely locate PAC in distinct regions of individual brains. To differentiate the fine spatial structure of the involved oscillators within each region, we used independent component analysis (ICA). ICA separates spatiotemporal independent subnetworks linearly mixed in the original signals, and thus provides a unique spatial signature for every separated subnetwork. These techniques helped us to investigate which spatiotemporal properties of regional PAC may be most relevant for the parkinsonian motor dysfunctions. Our findings provide new evidence on the pathophysiological mechanisms of PAC in Parkinson’s disease. Obtaining such information by non-invasive EEG may have implications for the development of neurostimulation techniques. Materials and methods The study protocol was approved by the local Ethics Committee (reference number: 147/18-ek). Participants Twenty-one patients with Parkinson’s disease, according to the current diagnostic criteria (Postuma et al., 2015), were recruited from the outpatient clinic of the Department of Neurology, University of Leipzig Medical Center. Motor function of all patients was assessed before the experiments, using the International Parkinson and Movement Disorders Society Unified Parkinson’s Disease rating scale part III (MDS-UPDRS III) (Goetz et al., 2008). Twenty-three age- and sex-matched healthy control subjects were recruited through various media. All participants were right-handed as confirmed by the Edinburgh Handedness Inventory (Oldfield, 1971). Written informed consent was obtained from all study participants. Exclusion criteria and information on individual structural MRI for all participants are described in the Supplementary material. Two patients were excluded from the cohort because of serious tremor activity interfering with the EEG recording in one patient and failure to obtain an MRI in the other. One healthy subject was excluded from the analysis because of an abnormal finding in MRI. EEG signal recording and preprocessing Signals were recorded using a 64-channel EEG system (eegoTMmylab, ANT Neuro) with 24-bits resolution and sampled at 2000 Hz. Vertical electrooculography and bipolar EMG of the first dorsal interosseous (FDI) muscle were also recorded. In patients, EMG was recorded from the hand side that was prominently affected by the disease, as indicated by the bradykinesia MDS-UPDRS III hemibody scores. In healthy controls, the side of FDI-EMG recordings was pseudorandomly chosen to eventually match the respective subsample sizes of patients. During the experiment, patients were in a practically defined ‘OFF medication’ state, after overnight withdrawal (at least 12 h) of parkinsonian medication. Prior to the EEG recording, the positions of the EEG electrodes and fiducial markers were acquired by a 3D optical digitization system (EEG Pinpoint, Localite). Recordings comprised 5 min of rest followed by several movement tasks. Only the data from the resting period were analysed for this report. For resting state recordings, subjects were asked to relax and fixate on a white cross displayed at the centre of a black computer screen 80 cm in front of them. EEG signals were processed using the EEGLAB toolbox (Delorme and Makeig, 2004) and custom MATLAB scripts. Signals were first detrended and then high-pass filtered at 0.5 Hz (eegfilt.m in EEGLAB with FIR default settings). Channels with poor data quality were excluded either based on visually detected artefacts (long-term large spiking-like activities) or if their power spectra failed to follow the canonical 1/f pattern (Miller et al., 2009). Components containing eye movement artefacts, channel noise, 50 Hz line noise, and heart beating were detected using ICA and subsequently removed. Artefacts from transitory muscle activity in EMG or EEG signals, or other environmental noise, were visually detected and marked in the continuous raw data. If more than 50% of data were contaminated by artefacts, the whole dataset was excluded from later analysis. Through this procedure, two healthy control subjects were excluded from the analysis. There were no significant differences between patients and controls regarding the data volume after removing all artefacts (P = 0.474). Afterwards, all data were re-referenced to the average of all electrodes. Region of interest-based source analysis The procedure for source analysis is displayed in Fig. 1A. It mainly relies on functions from the Fieldtrip toolbox (Oostenveld et al., 2011), EEGLAB toolbox, and custom MATLAB scripts. First, we constructed a lead-field matrix for each subject, which maps potential source activities on the cortex to electrical potentials at the EEG sensors (Supplementary material). A spatial filter was then constructed by using a least-square minimum variance beamformer (Van Veen and Buckley, 1988) based on the covariance matrix of artefacts-removed EEG signals and on the lead-field matrix. Figure 1 Open in new tabDownload slide Workflow for analysis on the source level. (A) Workflow of region of interest (ROI)-based EEG source analysis. We first used individual MRI images and electrode positions to create forward models, which describe the geometry and electrical properties of the head tissues and map the electrical brain activity to the electrical potentials at the EEG electrodes. The model is displayed as a brain outline (head model in light pink and source model in dark pink) with sensors (green dots) (top left). Using the HCP atlas information, we created the region of interest lead-field matrix. After preprocessing of EEG sensor signals artefacts-removed and artefacts-marked EEG datasets were prepared for further analysis. A data covariance matrix was computed using artefacts-removed EEG signals. We applied beamformer techniques to project artefacts-removed EEG signals from sensor to source space. PCA-ICA was applied to artefacts-removed source signals of each region of interest to obtain the ICA weighting matrices. We then implemented PCA-ICA weighting matrices on artefacts-marked EEG signals to obtain artefacts-marked time series data of ICA source components. (B) For each subject, we computed single KL-MI values (average KL-MI values in the phase frequency range from 13 to 30 Hz and amplitude frequency range from 50 to 150 Hz) to create a pairwise PAC matrix among ICA components. Each column represents the KL-MI values when each component provided the β-phases. Each row represents the KL-MI values when each component provided the γ-amplitudes. Figure 1 Open in new tabDownload slide Workflow for analysis on the source level. (A) Workflow of region of interest (ROI)-based EEG source analysis. We first used individual MRI images and electrode positions to create forward models, which describe the geometry and electrical properties of the head tissues and map the electrical brain activity to the electrical potentials at the EEG electrodes. The model is displayed as a brain outline (head model in light pink and source model in dark pink) with sensors (green dots) (top left). Using the HCP atlas information, we created the region of interest lead-field matrix. After preprocessing of EEG sensor signals artefacts-removed and artefacts-marked EEG datasets were prepared for further analysis. A data covariance matrix was computed using artefacts-removed EEG signals. We applied beamformer techniques to project artefacts-removed EEG signals from sensor to source space. PCA-ICA was applied to artefacts-removed source signals of each region of interest to obtain the ICA weighting matrices. We then implemented PCA-ICA weighting matrices on artefacts-marked EEG signals to obtain artefacts-marked time series data of ICA source components. (B) For each subject, we computed single KL-MI values (average KL-MI values in the phase frequency range from 13 to 30 Hz and amplitude frequency range from 50 to 150 Hz) to create a pairwise PAC matrix among ICA components. Each column represents the KL-MI values when each component provided the β-phases. Each row represents the KL-MI values when each component provided the γ-amplitudes. Source signals were then selected corresponding to particular cortical regions of interest defined by the Human Connectome Project atlas, which is based on multi-modal cortex parcellation (Glasser et al., 2016). In the atlas, 180 cortical areas are defined. Because abnormal cortical PAC of Parkinson’s disease has been previously reported in recordings from C3 and C4 EEG electrodes (Swann et al., 2015; Miller et al., 2019), which only broadly reflect neuronal activity originating from sensorimotor brain regions, we aimed at further spatial specification of PAC particularly focusing on the regions involved in motor control. Thus, we kept a fine resolution in the regions ‘somatosensory and motor cortex’ and ‘paracentral lobular and mid cingulate cortex’. Specifically, as described by Glasser et al. (2016), the region ‘somatosensory and motor cortex’ contained three (sub)regions, namely primary motor cortex (M1), primary somatosensory cortex and primary somatosensory complex, and the region ‘paracentral lobular and mid cingulate cortex’ contained three (sub)regions, namely cingulate motor area, supplementary motor area, and area 5. For all other parts of the brain we retained the coarser resolution of 20 regions as defined by Glasser et al. (2016) to avoid statistical power issues. As a result, 26 regions of interest were defined covering the whole brain. For the signals in each region of interest, we derived spatiotemporal patterns by using a combination of principal component analysis (PCA) and ICA (Jonmohamadi et al., 2014; Jonmohamadi and Jones, 2016). Specifically, PCA was applied to reduce the noise by keeping the main components, which can explain 95% of the data variance. To detect statistically independent stationary sources, we applied ICA to the PCA-truncated source signals for each region of interest separately (PCA and ICA were computed by runica.m in EEGLAB with the rank of 95% data variance). Then, we applied the beamformer weights and PCA-ICA weights calculated from the artefact-removed EEG signals to the artefact-marked signals, to generate continuous, breakpoint-free data, required for the PAC analysis. The sizes of the 26 regions of interest, as indexed by the average numbers of dipoles across subjects, and the average numbers of ICA components for patients and controls in all 26 regions of interest, are presented in Supplementary Table 1. None of the 26 regions of interest showed significant differences in the number of ICA components between patients and controls (two-tailed Wilcoxon rank-sum test, P > 0.05). Phase-amplitude coupling calculation We calculated PAC both on sensor and source signals by means of the Kullback-Leibler-based modulation index (KL-MI) (Tort et al., 2008). We computed KL-MI values in phase-amplitude signal pairs in which the frequencies of phase signals range from 4 to 40 Hz in 2-Hz steps, and that of amplitude signals from 32 to 200 Hz in 4-Hz steps. The details for PAC calculation are described in the Supplementary material. In sensor space, we calculated KL-MI values for each sensor signal separately. To enhance the spatial specificity of the sensor signals, we applied the current source density procedure (Supplementary material) before computing PAC. At the source level, we calculated KL-MI values in ICA component pairs within each region of interest, that is, we not only calculated the KL-MI values for each component separately, but also investigated the coupling between different components within the same region of interest. Next, a single KL-MI value of each sensor or component pair was computed by averaging the KL-MI values over the a priori defined phase frequency range of 13–30 Hz, and the a priori defined amplitude frequency range of 50–150 Hz (Swann et al., 2015). This yielded, for a region of interest with n components, n × n KL-MI values, as a pairwise PAC matrix shown in Fig. 1B. Finally, for comparison between groups, we computed a single PAC value for each region of interest of each subject. This was obtained by averaging across the n × n KL-MI values. Considering the different contributions of the ICA components to the original source signals, they were weighted by the percentage of variance (pvar) (Supplementary material) for each component pair before averaging. Power spectral density We used the MATLAB function ‘pwelch’ to calculate the power spectral density (PSD). The window was set to Hamming window in 2000 ms with 50% overlap. We then computed the normalized log PSD values for each ICA component. Waveform analysis for phase-locking estimation When estimating the strength of phase-phase coupling between β and γ activities, the common method is computing the n:m phase-locking value in which the phases of β and γ activities and the ratio between phases are directly defined (Ermentrout, 1981; Langdon et al., 2011). However, the direct phase definition of the broadband-γ activity is not available as it is composed of many different frequencies with different phases. Thus, to test whether some of the γ-phases are indeed locked to the β-phase, we defined the phase-phase coupling between β and broadband-γ activities as correlation of the γ signal with β-phase across β-cycles. Therefore, we applied the method of waveform analysis as described below. To examine whether the PAC-involved γ-activities were also phase-phase coupled to β-activities, we generated averaged waveforms for γ-activities (50–150 Hz, z-score normalized) time-locked to the β-phase (13–30 Hz, z-score normalized, then Hilbert transformed) zero crossings. We also averaged the envelopes of γ-activities (from Hilbert transformation), time-locked to the β-phase. For calculating the averaged waveforms, we extracted epochs of 160 ms (a period of around two cycles at 13 Hz) from both signals centred at the β-peak. To avoid overlapping, we used only every second epoch. If all phases of the broadband-γ activity were identical in each epoch (which would require phase reset), the amplitude of averaged γ oscillations would be identical to the average (over β-cycles) of the envelopes. Otherwise, it would be reduced. The maximum amplitudes of the averaged waveforms in a β-cycle were identified as the phase-locked amplitude (PLA) to β-phase. The mean PLA of component pairs represented the PLA for each subject. The PLAs of γ oscillations were defined as PLAosci, while those of γ envelopes were defined as PLAamp. By comparing PLAosci and PLAamp within subjects, and between groups, we examined to what degree abnormal PAC in Parkinson’s disease is accompanied by phase-phase coupling. The contributions of ICA components to phase-amplitude coupling For each PAC value of a component pair, one ICA component provides the β-phase and another (or the same) component provides the γ-amplitude. The overall contribution of a particular component to the β-phase information or γ-amplitude information in PAC of a region of interest can then be quantified by summing all PAC values that involve that component as a β-contributor or as a γ-contributor, respectively. For each region of interest, the values can be arranged in a matrix, as illustrated in Fig. 1B. In that matrix, the sum of each column represents the cumulative contribution of β-phases by the respective component, while the sum of each row represents the cumulative contribution of γ-amplitudes by the respective component. To obtain the relative contributions of a component to the β-phases and γ-amplitudes, respectively, these values were normalized to the sum of the entire PAC matrix of the region of interest. Finally, to estimate the inequality of the relative contributions of ICA components to β-phases and γ-amplitudes of the observed PAC in a region of interest, a Gini coefficient (Farris, 2010) was calculated in each case. A high Gini coefficient indicates that most contribution to the PAC comes from a few components, while a low Gini coefficient means that there are no such prominent contributors. By comparing the Gini coefficients between two groups, we circumvented the problem that ICA components differ between subjects and this approach enabled us to investigate whether in patients with Parkinson’s disease some subnetworks in the same region of interest stand out among others with regard to their contribution to the β-phases and γ-amplitudes involved in the observed PAC. Phase-amplitude coupling weighted ICA topography Each ICA source component defines an independent network with a unique spatiotemporal pattern (Chen et al., 2013) on the individual cortical surface. For each ICA component, we computed the topography as the spatial distribution of the absolute values of the ICA weights. To characterize the topographic distributions of the sources contributing to β-phases and coupled γ-amplitudes for each subject, we computed the average ICA topographies (as ‘β-topography’ and ‘γ-topography’, respectively) by averaging the topographies of all ICA components after weighting them according to their relative contributions to β-phases and γ-amplitudes, respectively. To compare the ICA topographies of the left or right hemisphere across subjects, they were mapped and morphed to the left hemisphere of a template source model from ‘fsaverage_sym’. Absolute Spearman correlation coefficients were calculated for each subject to compare the spatial similarity between the β- and γ-topography, and also between subjects to estimate the spatial similarity of the topographies across subjects. Statistical analysis We defined ‘contralateral’ hemisphere for patients as the hemisphere contralateral to the clinically more affected hemibody, and for controls as the hemisphere contralateral to the selected hand side for EMG recording (see above). For statistical comparison, we applied non-parametric statistical tests. Specifically, the two-tailed Wilcoxon ‘rank-sum’ test was used for comparisons between groups and the two-tailed Wilcoxon ‘sign-rank’ test for comparisons within subjects. For correlation analysis, we applied Spearman correlation analysis to investigate the relationship between PAC and other features including PSD, PLAs, and clinical scores. To deal with multiple comparisons, we applied false discovery rate (FDR) correction to the P-values (FDR P-values < 0.05). Data availability Personal data are protected by data privacy statements signed by all subjects. The data can be made available upon specific request taking into account the opinion of the local data privacy board. Results Characteristics of participants After exclusion of five participants (see ‘Materials and methods’ section), data from 19 patients (six females, age: 60.9 ± 10.8 years) and 20 healthy controls (eight females, age: 62.6 ± 7.9 years) were included in the analysis. Clinical characteristics of the patients are detailed in Table 1. The ‘contralateral’ hemisphere was the right hemisphere in 10 patients and the left hemisphere in nine patients corresponding to their more affected body side. Accordingly, in controls, the ‘contralateral’ hemisphere was the right hemisphere in 10 subjects and the left hemisphere in 10 subjects. There was no significant age difference between patients and controls (P = 0.527). On average, patients had been diagnosed with Parkinson’s disease 5.3 ± 4.5 years before the current investigation. Table 1 Characteristics of Parkinson’s disease patients ID . Sex . Age, years . Disease duration, years . Clinically more affected body side/ MDS-UPDRS III hemibody scores . Total MDS-UPDRS III (medication OFF) . l-DOPA equivalent dose, mg/day . 01 Male 57 4 Right/7 15 310 02 Female 46 1 Left/5 6 210 03 Male 54 1 Left/10 22 400 04 Female 65 3 Left/6 10 152 05 Male 55 1 Left/8 12 500 06 Male 49 4 Right/10 30 735 07 Male 67 12 Right/12 20 683 08 Female 74 11 Left/11 24 400 09 Male 75 3 Left/10 30 525 10 Male 76 2 Right/5 12 100 11 Female 78 4 Left/13 32 300 12 Female 58 6 Left/7 19 955 13 Male 71 12 Left/11 35 1395 14 Male 56 2 Right/6 11 355 15 Female 57 3 Left/15 26 930 16 Male 63 6 Right/10 20 400 17 Male 61 17 Right/12 29 1185 18 Male 38 2 Right/13 21 520 19 Male 57 6 Right/8 14 930 ID . Sex . Age, years . Disease duration, years . Clinically more affected body side/ MDS-UPDRS III hemibody scores . Total MDS-UPDRS III (medication OFF) . l-DOPA equivalent dose, mg/day . 01 Male 57 4 Right/7 15 310 02 Female 46 1 Left/5 6 210 03 Male 54 1 Left/10 22 400 04 Female 65 3 Left/6 10 152 05 Male 55 1 Left/8 12 500 06 Male 49 4 Right/10 30 735 07 Male 67 12 Right/12 20 683 08 Female 74 11 Left/11 24 400 09 Male 75 3 Left/10 30 525 10 Male 76 2 Right/5 12 100 11 Female 78 4 Left/13 32 300 12 Female 58 6 Left/7 19 955 13 Male 71 12 Left/11 35 1395 14 Male 56 2 Right/6 11 355 15 Female 57 3 Left/15 26 930 16 Male 63 6 Right/10 20 400 17 Male 61 17 Right/12 29 1185 18 Male 38 2 Right/13 21 520 19 Male 57 6 Right/8 14 930 Open in new tab Table 1 Characteristics of Parkinson’s disease patients ID . Sex . Age, years . Disease duration, years . Clinically more affected body side/ MDS-UPDRS III hemibody scores . Total MDS-UPDRS III (medication OFF) . l-DOPA equivalent dose, mg/day . 01 Male 57 4 Right/7 15 310 02 Female 46 1 Left/5 6 210 03 Male 54 1 Left/10 22 400 04 Female 65 3 Left/6 10 152 05 Male 55 1 Left/8 12 500 06 Male 49 4 Right/10 30 735 07 Male 67 12 Right/12 20 683 08 Female 74 11 Left/11 24 400 09 Male 75 3 Left/10 30 525 10 Male 76 2 Right/5 12 100 11 Female 78 4 Left/13 32 300 12 Female 58 6 Left/7 19 955 13 Male 71 12 Left/11 35 1395 14 Male 56 2 Right/6 11 355 15 Female 57 3 Left/15 26 930 16 Male 63 6 Right/10 20 400 17 Male 61 17 Right/12 29 1185 18 Male 38 2 Right/13 21 520 19 Male 57 6 Right/8 14 930 ID . Sex . Age, years . Disease duration, years . Clinically more affected body side/ MDS-UPDRS III hemibody scores . Total MDS-UPDRS III (medication OFF) . l-DOPA equivalent dose, mg/day . 01 Male 57 4 Right/7 15 310 02 Female 46 1 Left/5 6 210 03 Male 54 1 Left/10 22 400 04 Female 65 3 Left/6 10 152 05 Male 55 1 Left/8 12 500 06 Male 49 4 Right/10 30 735 07 Male 67 12 Right/12 20 683 08 Female 74 11 Left/11 24 400 09 Male 75 3 Left/10 30 525 10 Male 76 2 Right/5 12 100 11 Female 78 4 Left/13 32 300 12 Female 58 6 Left/7 19 955 13 Male 71 12 Left/11 35 1395 14 Male 56 2 Right/6 11 355 15 Female 57 3 Left/15 26 930 16 Male 63 6 Right/10 20 400 17 Male 61 17 Right/12 29 1185 18 Male 38 2 Right/13 21 520 19 Male 57 6 Right/8 14 930 Open in new tab Phase-amplitude coupling is enhanced in multiple cortical sources β-γ PAC as computed by the KL-MI method in 26 brain regions (see ‘Materials and methods’ section) showed a markedly inhomogeneous distribution across the brain, which was similar between patients and controls. However, comparison between patients and controls revealed specific regional differences as illustrated in Fig. 2. First, we analysed the differences in PAC by pooling the homologous regions of both hemispheres (Fig. 2A). In patients, β-γ PAC was enhanced, relative to controls, in six regions of interest comprising inferior frontal cortex (IFC, FDR P = 0.034, Z = 2.80), dorsolateral prefrontal cortex (DLPFC, FDR P = 0.034, Z = 2.88), premotor cortex (PMC, FDR P = 0.023, Z = 3.33), M1 (FDR P = 0.037, Z = 2.63), primary somatosensory cortex [Brodmann area (BA)3, FDR P = 0.037, Z = 2.68], and primary somatosensory complex (BA1,2, FDR P = 0.034, Z = 2.80). Z-statistic was highest for PMC. Because of the asymmetry in clinical parkinsonian motor symptoms, we also separately considered the hemispheres located contralaterally and ipsilaterally to the clinically more affected hemibody. In the contralateral hemisphere of patients, five regions of interest including DLPFC (FDR P = 0.035, Z = 2.99), PMC (FDR P = 0.035, Z = 3.25), M1 (FDR P = 0.035, Z = 2.80), BA3 (FDR P = 0.035, Z = 2.77), BA1,2 (FDR P = 0.035, Z = 2.74) showed higher PAC than in control subjects (Fig. 2B). The differences did not reach statistical significance in the ipsilateral hemisphere although we did not explicitly find significant differences of PAC between contralateral and ipsilateral hemispheres in either patients or controls. The comodulograms of the median KL-MI values in all five regions shows enhanced PAC in patients (Fig. 2B). PAC from temporal regions of interest (i.e. lateral temporal cortex, FDR P = 0.172, Z = 1.62; medial temporal cortex, FDR P = 0.412, Z = 0.96) did not differ between patients and controls, suggesting that the enhanced PAC in patients was unlikely to be driven by muscle artefacts from pericranial muscles. As we did not find any differences in the number of ICA components between patients and controls in any of the 26 regions of interest, the PAC differences between patients and controls are unlikely to be caused by differences in the number of components. Moreover, we applied a Pearson correlation test to estimate the correlation in each region of interest between the number of components and PAC values for 39 subjects. It showed no significant correlation in any of the 26 regions of interest (FDR P > 0.2). Figure 2 Open in new tabDownload slide PAC observed at source and sensor levels. (A) Analysis of PAC distribution across 26 regions of interest of the whole brain anatomy when averaging PAC from two hemispheres. Left: Regions marked as red showed significant differences between patients and controls after FDR correction [inferior frontal cortex (IFC), DLPFC, PMC, primary motor cortex (M1), primary somatosensory cortex (BA3), primary somatosensory complex (BA1,2)]. Right: Box plot showing the PAC of the six regions of interest that presented significant differences between patients and controls (marked by asterisks). (B) Analysis of PAC distribution across 26 regions of interest of the whole brain anatomy in the contralateral hemisphere. Top left: Five regions marked as red showed a significant difference between patients and controls after FDR correction. Top right: Box plot showing the PAC of the five regions of interest that presented significant difference (marked by asterisks). Bottom: Group comodulograms showing the median of single KL-MI across subjects in each group in the five regions of interest. (C) Analysis of PAC on sensor level. Left: Electrodes F3, F4; FC3, FC4; C3, C4; CP3, CP4, which are related to frontal and motor-sensory areas, were selected for comparison to the results at source level. Right: Box plot showing the PAC extracted from the four electrodes on the contralateral hemisphere of both patients and control subjects. None of the differences reached statistical significance. To refer collectively to the electrodes relating to the patients’ most affected hemispheres we used a notation with a forward slash (e.g. F3/4), because the most affected hemispheres varied between patients. For instance, F3/4 refers either to F3 or F4, depending on the patient’s dominantly affected side. Figure 2 Open in new tabDownload slide PAC observed at source and sensor levels. (A) Analysis of PAC distribution across 26 regions of interest of the whole brain anatomy when averaging PAC from two hemispheres. Left: Regions marked as red showed significant differences between patients and controls after FDR correction [inferior frontal cortex (IFC), DLPFC, PMC, primary motor cortex (M1), primary somatosensory cortex (BA3), primary somatosensory complex (BA1,2)]. Right: Box plot showing the PAC of the six regions of interest that presented significant differences between patients and controls (marked by asterisks). (B) Analysis of PAC distribution across 26 regions of interest of the whole brain anatomy in the contralateral hemisphere. Top left: Five regions marked as red showed a significant difference between patients and controls after FDR correction. Top right: Box plot showing the PAC of the five regions of interest that presented significant difference (marked by asterisks). Bottom: Group comodulograms showing the median of single KL-MI across subjects in each group in the five regions of interest. (C) Analysis of PAC on sensor level. Left: Electrodes F3, F4; FC3, FC4; C3, C4; CP3, CP4, which are related to frontal and motor-sensory areas, were selected for comparison to the results at source level. Right: Box plot showing the PAC extracted from the four electrodes on the contralateral hemisphere of both patients and control subjects. None of the differences reached statistical significance. To refer collectively to the electrodes relating to the patients’ most affected hemispheres we used a notation with a forward slash (e.g. F3/4), because the most affected hemispheres varied between patients. For instance, F3/4 refers either to F3 or F4, depending on the patient’s dominantly affected side. Previous studies reported enhanced PAC in Parkinson’s disease patients when PAC was directly computed from signals of C3 and C4 electrodes (Swann et al., 2015) located over motor regions. However, in our study, patients and controls did not differ when analysis was performed on signals recorded by C3/C4 electrodes (contralateral: P = 0.126; ipsilateral: P = 0.267), or on signals from other electrodes, including F3/F4 (FDR, contralateral: P = 0.438; ipsilateral: P = 0.693), FC3/FC4 (FDR, contralateral: P = 0.178; ipsilateral: P = 0.747), CP3/CP4 (FDR, contralateral: P = 0.178; ipsilateral: P = 0.693), which are related to frontal and sensorimotor areas (Fig. 2C). Even when sensors from both hemispheres were collectively considered (FDR, F3,4 P = 0.553, FC3,4 P = 0.177, C3,4 P = 0.177, CP3,4 P = 0.177), no difference was present between patients and controls. We examined the possibility that enhanced PAC may be driven by altered spectral power in either the β or γ band in patients. In M1 (Fig. 3A), β power differed significantly between patients and controls in contralateral (FDR, P = 0.040) and ipsilateral hemispheres (FDR, P = 0.019). Furthermore, γ power differed significantly in the contralateral hemisphere (FDR, P = 0.011), but not in the ipsilateral one (FDR, P = 0.093). The other four regions of interest showed similar results (Supplementary Fig. 1). Figure 3 Open in new tabDownload slide Relationship between power spectral density (PSD) and PAC in M1. (A) Left: Box plots of spectral power in the β frequency range (13–30 Hz) in patients and control subjects. Right: Box plots of spectral power in the γ frequency range (50–150 Hz). Significant differences between patients and controls are marked by asterisks. (B) Left: Non-significant correlation between power peak frequency at β range and peak PAC frequency for phase. Right: Non-significant correlation between power peak frequency at γ range and peak PAC frequency for amplitude. (C) Three examples of comodulograms of PAC and the related β and γ power from patients. Left: Phase and amplitude frequency of the largest PAC are related to the β and γ power peak. Middle: Lack of PAC even when there is a clearly discernible β power peak. Right: The phase frequency and the amplitude frequency of the largest PAC are not related to the β and γ power peaks. Figure 3 Open in new tabDownload slide Relationship between power spectral density (PSD) and PAC in M1. (A) Left: Box plots of spectral power in the β frequency range (13–30 Hz) in patients and control subjects. Right: Box plots of spectral power in the γ frequency range (50–150 Hz). Significant differences between patients and controls are marked by asterisks. (B) Left: Non-significant correlation between power peak frequency at β range and peak PAC frequency for phase. Right: Non-significant correlation between power peak frequency at γ range and peak PAC frequency for amplitude. (C) Three examples of comodulograms of PAC and the related β and γ power from patients. Left: Phase and amplitude frequency of the largest PAC are related to the β and γ power peak. Middle: Lack of PAC even when there is a clearly discernible β power peak. Right: The phase frequency and the amplitude frequency of the largest PAC are not related to the β and γ power peaks. We calculated the Spearman correlation between PAC and β and γ power for each region of interest of each hemisphere (Supplementary Table 2). We found no relationship between β power and PAC in any of the five regions of interest. While γ power was not significantly correlated with PAC in either DLPFC, PMC or M1, both metrics were correlated in BA3 and BA1,2. As power was not correlated with PAC in any of the three regions exhibiting the largest PAC difference between patients and controls (DLPFC, PMC or M1), the enhancement of power is unlikely to underlie the enhancement of PAC in patients. Moreover, we compared the frequencies at which the largest PAC value was noted with the frequencies where the largest β or γ power was detected. A correlation test applied on each signal pair in M1 showed no correlation between the peak PAC frequency and the power peak frequency in either the β or γ band (Fig. 3B). Similar results were obtained in the other four regions of interest (data not shown). Three examples from the dataset of patients are displayed in Fig. 3C. Together, the above findings indicate that the detected PAC enhancement in patients was unlikely to be driven by differences in spectral power. The preceding analysis has dealt with the macroscopic distribution of abnormal PAC. In the following section, we investigate the mechanisms of enhanced PAC in Parkinson’s disease at a mesoscopic scale. We perform detailed analyses regarding the composition and origin of the signals abnormally coupled in patients. To enhance the clarity of presentation and to enable comparison with previous reports on PAC recorded from ECoG of M1, we first present results of M1. Findings from the other regions (DLPFC, PMC, BA3 and BA1,2) harbouring abnormal PAC are also reported. Neural circuits contributing to abnormal phase-amplitude coupling While PAC indicates coupling between γ-amplitude and β-phase, its presence does not imply nor exclude the possibility of additional phase-phase coupling. Therefore, we addressed this issue separately. The phase-phase coupling between β and γ activity can be distinguished when γ activity is averaged across epochs centred at a specific β-phase. This idea is illustrated in Fig. 4A, which shows two cases. In both cases the amplitude of high-frequency oscillation (HFO) is coupled to the phase of the low-frequency oscillation (LFO), but only in one case was there additional phase-phase coupling between HFO and LFO. The average waveform of HFO in the non-phase locked case results in a relatively flat waveform compared with that of the envelope of HFO, while in another case the HFO was retained. Therefore, the largest amplitude of the averaged waveform, termed PLA, can be used for estimating the phase relationship between PAC-involved activities. The results in M1 revealed that the PLAamp was larger than PLAosci both in the contralateral (Fig. 4B, left) and ipsilateral hemisphere for each subject, which rules out that the γ activities were completely phase-phase coupled to β-phase in both patients (FDR P < 0.001) and controls (FDR P < 0.001). In addition, in the contralateral hemisphere (Fig. 4B, left), we found no significant difference of PLAosci between patients and controls (FDR P = 0.399), but a significant difference in PLAamp (FDR P = 0.005). In the ipsilateral hemisphere, we found no significant difference between the two groups in either PLAosci (FDR, P = 0.189) or PLAamp (FDR, P = 0.479). Correlation analysis revealed no significant correlation between PAC and the PLAosci (Fig. 4B, middle) in patients (FDR P = 0.278) and controls (FDR P = 0.456) in the contralateral hemisphere, but a significant positive correlation between PAC and the PLAamp (Fig. 4B, right), both in patients (FDR P < 0.001, R = 0.80) and controls (FDR P = 0.016, R = 0.58). The results for DLPFC and PMC were similar to M1, but slightly different for somatosensory areas (Fig. 4C). However, even in somatosensory areas, the PLAosci values were still lower than PLAamp values (FDR P < 0.001). These findings demonstrate that the PAC-involved β and γ activities are not completely phase-phase coupled, and the PAC between β and γ activities that are not strictly phase-phase coupled becomes stronger in the pathological state. Figure 4 Open in new tabDownload slide Waveform estimation of PLA. (A) Simulation of phase-locked and non-phase locked conditions on PAC-involved oscillations. Left: Two oscillatory signals in which the amplitude of the high-frequency oscillation (HFO, red) is coupled to the peak of the low-frequency oscillation (LFO, black). In the top case, the phase relation between both signals varies over time, whereas in the bottom case, the HFO is phase-phase coupled to the LFO. Middle: The varying phase relation (top) results in a relatively flat average waveform (red) of the LFO, when compared to the average of the instantaneous envelope of HFO (blue). In contrast, in the bottom part, the HFO is retained when it is phase-phase coupled to the LFO. Right: In both cases the KL-MI values are identical. (B) Experimental results on PLA in M1 on the contralateral hemisphere. Left: Paired plot showing that PLAosci (PLA of γ oscillations) were always lower than PLAamp (PLA of γ envelopes) within subjects. Furthermore, PLAamp were significantly larger in patients than in control subjects. Significant differences are marked by asterisks. Middle: No significant correlation of PLAosci with mean PAC in patients and controls. Right: Significant correlation of PLAamp with mean PAC in patients (red line) and controls (black line). (C) Experimental results on PLA in BA3 in the contralateral hemisphere. Left: Paired plot showing that PLAosci was always lower than PLAamp in patients and controls. Furthermore, PLAamp was significantly larger in patients than in controls. Significant differences between patients and controls are marked by asterisks. Middle: Significant correlation of PLAosci with mean PAC in patients (red line) but not in control subjects. Right: Significant correlation of PLAamp with mean PAC in patients (red line) and controls (black line). Lines are only displayed if the correlation was statistically significant. Figure 4 Open in new tabDownload slide Waveform estimation of PLA. (A) Simulation of phase-locked and non-phase locked conditions on PAC-involved oscillations. Left: Two oscillatory signals in which the amplitude of the high-frequency oscillation (HFO, red) is coupled to the peak of the low-frequency oscillation (LFO, black). In the top case, the phase relation between both signals varies over time, whereas in the bottom case, the HFO is phase-phase coupled to the LFO. Middle: The varying phase relation (top) results in a relatively flat average waveform (red) of the LFO, when compared to the average of the instantaneous envelope of HFO (blue). In contrast, in the bottom part, the HFO is retained when it is phase-phase coupled to the LFO. Right: In both cases the KL-MI values are identical. (B) Experimental results on PLA in M1 on the contralateral hemisphere. Left: Paired plot showing that PLAosci (PLA of γ oscillations) were always lower than PLAamp (PLA of γ envelopes) within subjects. Furthermore, PLAamp were significantly larger in patients than in control subjects. Significant differences are marked by asterisks. Middle: No significant correlation of PLAosci with mean PAC in patients and controls. Right: Significant correlation of PLAamp with mean PAC in patients (red line) and controls (black line). (C) Experimental results on PLA in BA3 in the contralateral hemisphere. Left: Paired plot showing that PLAosci was always lower than PLAamp in patients and controls. Furthermore, PLAamp was significantly larger in patients than in controls. Significant differences between patients and controls are marked by asterisks. Middle: Significant correlation of PLAosci with mean PAC in patients (red line) but not in control subjects. Right: Significant correlation of PLAamp with mean PAC in patients (red line) and controls (black line). Lines are only displayed if the correlation was statistically significant. We also characterized abnormal PAC further by investigating whether the enhanced PAC of patients could be caused by β and γ activities from different ICA components. We segregated the KL-MI values of phase-amplitude pairs originating from two different components (PACinter) and those from the same components (PACiden). Compared with controls, in the contralateral hemisphere of patients (Fig. 5A), significantly enhanced PACiden was found in DLPFC (FDR P = 0.018), PMC (FDR P = 0.006), M1 (FDR P = 0.016), BA1,2 (FDR P = 0.026), and BA3 (FDR P = 0.074), where the difference was marginally significant. In addition, we found a significant difference of PACinter between patients and controls in all five regions of interest (FDR, DLPFC P = 0.001, PMC P = 0.001, M1 P = 0.001, BA3 P = 0.001, BA1,2 P = 0.001). Although the KL-MI values were relatively higher when phase-amplitude pairs were from identical components, we also found the coupling between different components to be enhanced in patients compared with controls. This result indicates that the PAC-involved β and γ activities could be generated from different subnetworks, which provides support for our conclusion that abnormally enhanced coupling comprises distinct subnetworks in at least five brain regions of patients. Figure 5 Open in new tabDownload slide PAC computed fromidentical or different ICA components. (A) Box plots show PAC differences for all regions of interest between patients and controls from interactions between different (left) or identical ICA components (right), in the contralateral hemisphere. Significant differences between patients and controls are marked by asterisks. (B) Correlation between PAC and clinical scores in the five regions of interest (black lines designate statistically significant correlations). Left: Correlation between mean PAC on the overall pairwise matrix with MDS-UPDRS III hemibody scores. Middle: Correlation between mean PAC from the interaction between different components (PACinter) with MDS-UPDRS III hemibody scores. Right: Correlation between mean PAC from identical components (PACiden) with MDS-UPDRS III hemibody scores. Figure 5 Open in new tabDownload slide PAC computed fromidentical or different ICA components. (A) Box plots show PAC differences for all regions of interest between patients and controls from interactions between different (left) or identical ICA components (right), in the contralateral hemisphere. Significant differences between patients and controls are marked by asterisks. (B) Correlation between PAC and clinical scores in the five regions of interest (black lines designate statistically significant correlations). Left: Correlation between mean PAC on the overall pairwise matrix with MDS-UPDRS III hemibody scores. Middle: Correlation between mean PAC from the interaction between different components (PACinter) with MDS-UPDRS III hemibody scores. Right: Correlation between mean PAC from identical components (PACiden) with MDS-UPDRS III hemibody scores. To examine the relationship between PAC and the severity of motor impairment of patients with Parkinson’s disease, we computed the correlation between the PAC values in contralateral and ipsilateral hemispheres with total MDS-UPDRS III hemibody scores. The correlations between PAC and clinical scores in all five regions of interest are shown in Fig. 5B, and the correlation coefficients and P-values are presented in Table 2. We found the overall PAC (no segregation between component pairs within a given region) to be significantly correlated with MDS-UPDRS III hemibody scores in M1, but not in the other four regions of interest. We also computed the correlation between clinical scores and PACinter as well as PACiden. The PACinter showed a significant correlation with UPDRS hemibody scores in PMC, M1, BA3, and BA1,2. However, we did not find a significant correlation between clinical scores and PACiden in any of these regions of interest. Moreover, we did not find any significant correlation between PAC values from DLPFC and clinical scores. Table 2 Correlation between PAC and MDS-UPDRS III hemibody scores ROI . Overall PAC . PACinter . PACiden . DLPFC P = 0.491, R = 0.12 P = 0.407, R = 0.21 P = 0.940, R = 0.01 PMC P = 0.088, R = 0.28 P = 0.004*, R = 0.49 P = 0.707, R = 0.15 M1 P = 0.027*, R = 0.36 P = 0.007*, R = 0.46 P = 0.135, R = 0.30 BA3 P = 0.169, R = 0.22 P = 0.033*, R = 0.39 P = 0.499, R = 0.19 BA1,2 P = 0.065, R = 0.30 P = 0.003*, R = 0.49 P = 0.252, R = 0.25 ROI . Overall PAC . PACinter . PACiden . DLPFC P = 0.491, R = 0.12 P = 0.407, R = 0.21 P = 0.940, R = 0.01 PMC P = 0.088, R = 0.28 P = 0.004*, R = 0.49 P = 0.707, R = 0.15 M1 P = 0.027*, R = 0.36 P = 0.007*, R = 0.46 P = 0.135, R = 0.30 BA3 P = 0.169, R = 0.22 P = 0.033*, R = 0.39 P = 0.499, R = 0.19 BA1,2 P = 0.065, R = 0.30 P = 0.003*, R = 0.49 P = 0.252, R = 0.25 P-values were corrected (for each region of interest, ROI). *P-values (FDR P < 0.05). PACiden = phase-amplitude coupling from the same components; PACinter = phase-amplitude coupling between different components. Open in new tab Table 2 Correlation between PAC and MDS-UPDRS III hemibody scores ROI . Overall PAC . PACinter . PACiden . DLPFC P = 0.491, R = 0.12 P = 0.407, R = 0.21 P = 0.940, R = 0.01 PMC P = 0.088, R = 0.28 P = 0.004*, R = 0.49 P = 0.707, R = 0.15 M1 P = 0.027*, R = 0.36 P = 0.007*, R = 0.46 P = 0.135, R = 0.30 BA3 P = 0.169, R = 0.22 P = 0.033*, R = 0.39 P = 0.499, R = 0.19 BA1,2 P = 0.065, R = 0.30 P = 0.003*, R = 0.49 P = 0.252, R = 0.25 ROI . Overall PAC . PACinter . PACiden . DLPFC P = 0.491, R = 0.12 P = 0.407, R = 0.21 P = 0.940, R = 0.01 PMC P = 0.088, R = 0.28 P = 0.004*, R = 0.49 P = 0.707, R = 0.15 M1 P = 0.027*, R = 0.36 P = 0.007*, R = 0.46 P = 0.135, R = 0.30 BA3 P = 0.169, R = 0.22 P = 0.033*, R = 0.39 P = 0.499, R = 0.19 BA1,2 P = 0.065, R = 0.30 P = 0.003*, R = 0.49 P = 0.252, R = 0.25 P-values were corrected (for each region of interest, ROI). *P-values (FDR P < 0.05). PACiden = phase-amplitude coupling from the same components; PACinter = phase-amplitude coupling between different components. Open in new tab Multiple abnormal couplings between spatially distributed different sources Since PACinter was correlated with the clinical severity, we asked whether the enhanced PACinter in patients resulted from the abnormal enhancement of all or just a few dominant components. As the ICA components are generally different between subjects, this question cannot be examined directly. However, the pathological dominance of a subset of components in patients can be inferred indirectly by comparing the distribution of relative contributions of components to PAC between groups. An abnormal contribution of only a subset of components in patients would be likely, if the relative contributions of components were more unevenly distributed in patients than in controls. By contrast, if there were no outstanding contributors to the PAC enhancement in patients, the contributions of ICA components to the β-phases or γ-amplitudes of PAC would likely be similarly distributed in patients and control subjects. The evenness/unevenness of the distribution can be estimated by the Gini coefficient, a metric widely used in social statistics (see ‘Materials and methods’ section). Since the enhanced PAC in patients was stronger in the contralateral hemisphere, this analysis was done on this hemisphere. Gini coefficients in M1 (Fig. 6A, left) were significantly different between patients and controls suggesting that the contributions of different ICA components to the β-phases (FDR, P = 0.005) and γ-amplitudes (FDR, P = 0.044) depended on the disease state. Results were similar in the other four regions of interest (Supplementary Fig. 2). This finding indicates that, compared with relatively equivalent contributions of subnetworks to PAC in control subjects, a few subnetworks in prefrontal and sensorimotor areas of patients contribute β-phases or γ-amplitudes to exaggerated PAC. Besides, both in patients and controls, we found significant correlations between PAC strengths and Gini coefficients for β-phases and γ-amplitudes of PAC in M1 (Fig. 6A, middle and right) and in the other four regions of interest (Supplementary Table 3). Those findings indicate that the enlarged PAC was related to enhanced contributions from a few subnetworks. Figure 6 Open in new tabDownload slide Correlation analysis for PAC weighted topography estimation in the contralateral hemisphere. (A) Left: Box plot of Gini coefficients in M1. Comparison between patients and controls. Middle and right: Gini coefficients of contributions to β-phases and γ-amplitudes were correlated with PAC values. (B) Examples of weighted average β-topographies and the corresponding γ-topographies in each of the five regions of interest. ICA weights were z-score normalized. The spatial patterns are similar but not identical between β- and γ-topographies in each region of interest. (C) Top: Histogram of spatial similarity coefficients as indicated by absolute Spearman correlation coefficients of M1 showing raw correlation values between β- and γ-topographies for patients and controls (the dashed lines represent the median coefficients of each group). Bottom: Significant correlation between β-γ topographic similarity and the Gini coefficient of component contributions to β-phases and γ-amplitudes (red and black lines indicate statistically significant correlations). (D) Histogram of Spearman correlation coefficient of weighted average topographies in M1 across subjects. Note the large variability of similarities of topographies across subjects. Black dashed lines represent the median coefficient of each distribution, and show low similarity among subjects. Topo = topography. Figure 6 Open in new tabDownload slide Correlation analysis for PAC weighted topography estimation in the contralateral hemisphere. (A) Left: Box plot of Gini coefficients in M1. Comparison between patients and controls. Middle and right: Gini coefficients of contributions to β-phases and γ-amplitudes were correlated with PAC values. (B) Examples of weighted average β-topographies and the corresponding γ-topographies in each of the five regions of interest. ICA weights were z-score normalized. The spatial patterns are similar but not identical between β- and γ-topographies in each region of interest. (C) Top: Histogram of spatial similarity coefficients as indicated by absolute Spearman correlation coefficients of M1 showing raw correlation values between β- and γ-topographies for patients and controls (the dashed lines represent the median coefficients of each group). Bottom: Significant correlation between β-γ topographic similarity and the Gini coefficient of component contributions to β-phases and γ-amplitudes (red and black lines indicate statistically significant correlations). (D) Histogram of Spearman correlation coefficient of weighted average topographies in M1 across subjects. Note the large variability of similarities of topographies across subjects. Black dashed lines represent the median coefficient of each distribution, and show low similarity among subjects. Topo = topography. Finally, we examined the question of how the networks contributing to β-phases and γ-amplitudes involved in PACinter were spatially organized on a local scale. We studied the ‘β-topography’ and ‘γ-topography’ for each subject in the contralateral hemisphere as introduced in the ‘Materials and methods’ section. The β-topography and γ-topography in a local region represent the averages of the spatial distributions of the ICA components, after weighting them by their respective contributions to the PAC in terms of β-phases and γ-amplitudes. Figure 6B illustrates examples of β-topographies and related γ-topographies of the five regions of interest. The example shows a similar but not identical spatial pattern between both topographies within each region of interest. For statistical estimation of the spatial similarity between β- and γ-topographies within each subject in M1, we computed the absolute Spearman coefficient for each subject, which is displayed in the histogram (Fig. 6C). This showed that the β- and γ-topographies for each individual were highly spatially correlated, both in patients (R = 0.96 ± 0.03) and in control subjects (R = 0.98 ± 0.02). However, there was a strong tendency for the β-γ topographic similarities of patients to be reduced compared to controls (Wilcoxon ranksum test, P = 0.051). This observation indicates that abnormally enhanced interactions became more prevalent between a few spatially distinct subnetworks in patients than in control subjects. The lower spatial similarity within patients was also observed in BA1,2 (P = 0.026), but not in the other three regions of interest (Supplementary Fig. 3). Moreover, the β-γ topographic similarities in M1 were negatively correlated with the Gini coefficients for component contributions to γ-amplitudes (FDR P < 0.001, R = −0.64), as well as to β-phases (FDR P < 0.001, R = −0.58) (Fig. 6C). The other four regions of interest showed similar results (Supplementary Fig. 4). Furthermore, when estimating the topographic similarity across subjects, the median values of absolute Spearman coefficient histogram (Fig. 6D) were low in the similarities of β-topographies among patients (median R = 0.30) and controls (median R = 0.25), and γ-topographies among patients (median R = 0.30) and controls (median R = 0.24). The other four regions of interest showed similar results. These findings indicate that there was no consistent spatial pattern of β-γ PAC across subjects. Discussion Using source localization techniques, the present study has comprehensively characterized enhanced PAC in scalp EEG recorded from patients with Parkinson’s disease, providing novel insights into abnormal brain synchronization on different anatomical scales. The present study has localized abnormal PAC in DLPFC, PMC, M1 and somatosensory areas (including BA3, BA1,2), and thus to more regions of the brain than previously known. Direct ECoG recordings have revealed abnormal PAC in patients with Parkinson’s disease in signals from subdural electrodes overlying the precentral gyrus, thus probably reflecting activity in PMC/M1 (de Hemptinne et al., 2013). Although ECoG has superior spatial resolution compared to EEG (Buzsáki et al., 2012), its spatial coverage is limited by the small size of the electrode strip that can be placed through a burr hole. Enhanced PAC in patients with Parkinson’s disease has also been observed in signals from EEG electrodes overlying the sensorimotor region (Swann et al., 2015; Miller et al., 2019). Using advanced source localization algorithms, the present findings appear to be the first to provide non-invasive evidence of abnormal PAC located in the DLPFC and somatosensory areas. The contribution of somatosensory area is noteworthy as recordings were done at rest, which implies minimal activity of reafferent signals from contracting muscles. Furthermore, it has been hypothesized that the ‘hyperdirect’ cortico-subthalamic pathway is essential for generating PAC (de Hemptinne et al., 2013). Somatosensory areas do not monosynaptically project to STN, which receives its monosynaptic excitatory input mainly from neurons located in M1 as well as from the supplementary motor area (Nambu et al., 1996) and prefrontal regions (Haynes and Haber, 2013; Bruni et al., 2018). Therefore, the presence of PAC in the somatosensory area indicates that abnormal PAC at rest is unlikely to exclusively involve brain regions connected by the hyperdirect tract, but may also comprise other BGTC circuits or cortico-cortical connections. A similar conclusion has been reached based on recordings from globus pallidus internus (GPi) (Connolly et al., 2015; Malekmohammadi et al., 2018). Although the spatial resolution of scalp EEG is limited even with advanced source localization, our findings strongly suggest that abnormal PAC in Parkinson’s disease, apart from M1, is present in a variety of cortical regions involved in motor control, even at rest. PAC was clearly abnormally enhanced in DLPFC of patients with Parkinson’s disease. However, its magnitude was not correlated to MDS-UPDRS III. As MDS-UPDRS III exclusively captures motor symptoms of Parkinson’s disease, the lack of correlation in DLPFC, on the one hand, is readily explained by the fact that DLPFC is predominantly involved in associative BGTC circuits, which are implicated in executive and other cognitive dysfunctions rather than motor functions (DeLong and Wichmann, 2015; Magrinelli et al., 2016). On the other hand, as executive dysfunction is part of the clinical characteristics in Parkinson’s disease, enhanced PAC in DLPFC suggests that abnormal brain synchronization may also have a role in behavioural domains outside of motor control (Oswal et al., 2013). As the DLPFC has direct connections to the basal ganglia network and the functional activity in DLPFC can be altered through STN DBS (Boertien et al., 2011; DeLong and Wichmann, 2015), it will be interesting to see whether PAC may be correlated with cognitive symptoms in Parkinson’s disease and how responsive it is to different clinical interventions. Previous reports have examined the relationship between the magnitude of PAC and clinical severity of parkinsonian bradykinesia. Correlation between PAC and clinical scores has been noted in LFP recordings from human STN (López-Azcárate et al., 2010; van Wijk et al., 2016; Ozturk et al., 2020), GPi (Tsiokos et al., 2017), or in ECoG from PMC/M1 (Malekmohammadi et al., 2018). In our study, correlations between PAC and clinical motor scores were noted in all regions with enhanced PAC, except DLPFC. While previous work has demonstrated a correlation of therapy-induced changes in PAC derived from EEG with therapy-induced changes in clinical severity (Miller et al., 2019), our study demonstrates a direct relationship between the strength of native PAC computed from scalp EEG and clinical severity scores in the OFF medication state. This observation not only enhances the validity of our findings, but may additionally underline the pathophysiological significance of enhanced PAC in cortical sources. PAC was statistically significantly enhanced only in the clinically more affected hemisphere, which agrees with previous reports (Shreve et al., 2017). In animal models of Parkinson’s disease, PAC did not appear before advanced parkinsonian stages (Connolly et al., 2015; Devergnas et al., 2019). As previous ECoG have been performed in patients undergoing DBS (de Hemptinne et al., 2013), advanced stages of Parkinson’s disease are implicated. In contrast, the majority of patients in the present study were at relatively early stages of the disease. Therefore, PAC may have been expressed more weakly in the present patients and hence may have been too low to be detectable reliably in the less severely affected hemisphere. This could also explain why we failed to find the PAC abnormality in EEG sensor signals. As we applied commonly used procedures for recording, preprocessing, and PAC calculation, which were similar to those used in previous studies (Swann et al., 2015; Miller et al., 2019), it seems unlikely that the differences between the results obtained in previous studies and in the present study could result from methodological details. Rather, we suggest that for PAC to be detectable in sensor signals (Miller et al., 2019), patients may need to be in more advanced stages of the disease. The fact that application of realistic head modelling and inverse solutions were required for PAC detection from EEG signals suggests that advanced signal processing technology may be required if PAC should be used more widely in future non-invasive closed-loop paradigms. Our results emphasize that the transformation of EEG sensor signals to source signals may be important to enhance the sensitivity of the biomarker. Moreover, a more precise location of abnormal PAC could help to find a more effective neuromodulation target for specific symptoms. Our analyses also provide insight into the potential mechanism of enhanced PAC in patients with Parkinson’s disease on a mesoscopic scale. Using waveform analysis, we ruled out that the PAC-involved β and γ activities of patients are completely phase-phase coupled. This finding specifically renders the possibility unlikely that PAC arises because a recurrent (cortical) network, which physiologically oscillates in the γ range (Ray et al., 2008), is strongly driven by the output of an enhanced (subcortical) β-oscillator (Spiegler et al., 2011). Our findings suggest, instead, that PAC arises to a significant extent through the activity of two physiologically distinct oscillators which likely correspond to two spatially distinct cortical substrates. This conclusion also has implications for another hypothesis about PAC, which was developed through time-domain analysis of β oscillations in ECoG (Cole et al., 2017) or scalp EEG recordings (Jackson et al., 2019). According to this hypothesis, the non-sinusoidal wave shape of the β oscillations is a consequence of over-synchronization of local firing of action potentials in the cortex (De Hemptinne et al., 2015; Voytek and Knight, 2015). In principle, PAC could result from excessively synchronized neural activity in the BGTC loop (Cole et al., 2017) and may still be correlated with clinical motor impairment. Nonetheless, our results suggest that attenuation of PAC by decorrelating excessive synchrony through DBS (Voytek and Knight, 2015) might fail to break up PAC across spatially separated circuits, which might represent a pathophysiologically more relevant target. Our conclusion, assuming the existence of spatially and physiologically distinct oscillators contributing to PAC, was substantiated by further investigation of ICA components. We found enhanced PAC when phase-amplitude pairs were composed of β and γ signals from identical, but also from distinct components. Importantly, PAC was only correlated with the clinical severity scores when computed between different components. Because ICA is not suited to distinguish between temporally independent components if they originate from exactly the same cortical patch, pathological coupling across distinct components implies coupling between different cortical columns. However, despite the fact that the electric activity was generated in the cortex, our findings do not contradict a model in which one or both of the cortical oscillators are driven by activity from subcortical nuclei. The nature of the coupling could even be routed entirely in distinct subcortical projections driving independent cortical oscillations (Pasquereau and Turner, 2011; Shimamoto et al., 2013; Belluscio et al., 2014). We found that Gini coefficients were higher in patients than in controls. This finding indicates that some components contributed much more to PAC in Parkinson’s disease than others. This result provides additional support for the assertion that synchronization between the activities of spatially distinct neuronal circuits was of greater pathophysiological relevance than that between the activities of the same topographical origin. Our conclusion was further corroborated when the spatial organization of ICA components was explicitly tested. In M1 and BA1,2, where β- and γ-topographies were generally similar for each individual, this similarity was reduced in patients. Moreover, we found that the higher the Gini coefficients, the greater the spatial discrepancy between the source distributions of β and coupled-γ activities. The implication is that the abnormal enhancement of PAC in Parkinson’s disease reflects the abnormal dominance of a relatively small number of coupled subnetworks. Animal experiments have shown that the loss of dopamine impairs the directionality and hierarchical organization of normal β and γ propagation through different BGTC pathways (West et al., 2018). Against this background, we suggest that in the parkinsonian state the dopamine depletion may disrupt the regular operation of the BGTC loops by reinforcing and altering the connectivity between subnetworks so that the physiological segregation of certain feedback loops in basal ganglia circuits is lost. This may then lead to more widespread coupling between diverse circuits. Coupling of circuits with similar intrinsic frequencies will be evident as abnormal synchronization and reduction of dimensionality in the EEG signals, whereas coupling of circuits with different frequencies will be evident as enhanced PAC. PAC-weighted topographies in patients did not display a uniform parkinsonian pattern, but differed considerably interindividually. This finding indicates a pathophysiological heterogeneity which may reflect different clinical motor phenotypes as well as the variability of the disease stage. This further implies that pathophysiologically based future treatment of Parkinson’s disease by non-invasive brain stimulation may need to be highly individualized and dynamically adjusted. Because all analyses were done on recordings while the patients were at rest, our findings do not allow us to establish a direct link between physiological abnormalities and the clinical phenotype. As certain properties of PAC derived from EEG were correlated with the parkinsonian motor symptoms, it appears likely that abnormal PAC, even if derived from the cortex, is more than a marker for the parkinsonian state. However, future studies will have to explore abnormal cross-frequency coupling during motor behaviour and to investigate the effect of treatment interventions such as dopamine replacement therapy and DBS on cross-frequency coupling. In the present study, patients with marked rest tremor were excluded and all recordings were done in patients in early to moderate disease stages. Therefore, the generalizability of our findings to tremulous or more severe parkinsonian phenotypes remains unclear. Larger variability of clinical severity and a larger number of patients may allow more robust assessments of the nature of the correlation between pathological synchronization and individual items of the clinical phenotype. Acknowledgements We thank all patients and healthy volunteers who participated in our study. We also thank Alhuda Dabbagh for assisting in data collection, and Juanli Zhang for the discussion of analysis methods. Funding This work was supported by the CortExplorer program (P1140048) of the Hertie Foundation (Gemeinnützige HERTIE-Stiftung) awarded to J.C. R.G. is recipient of a scholarship of the International Max Planck Research School on Neuroscience of Communication (IMPRS Neurocom).. Competing interests The authors report no competing interests. Supplementary material Supplementary material is available at Brain online. References Belluscio MA , Escande MV , Keifman E , Riquelme LA , Murer MG , Zold CL. Oscillations in the basal ganglia in Parkinson's disease: role of the striatum . Basal Ganglia 2014 ; 3 : 203 – 12 . Google Scholar Crossref Search ADS WorldCat Boertien T , Zrinzo L , Kahan J , Jahanshahi M , Hariz M , Mancini L , et al. Functional imaging of subthalamic nucleus deep brain stimulation in Parkinson's disease . Mov Disord 2011 ; 26 : 1835 – 43 . Google Scholar Crossref Search ADS PubMed WorldCat Bruni S , Gerbella M , Bonini L , Borra E , Coudé G , Ferrari PF , et al. Cortical and subcortical connections of parietal and premotor nodes of the monkey hand mirror neuron network . Brain Struct Funct 2018 ; 223 : 1713 – 29 . Google Scholar PubMed OpenURL Placeholder Text WorldCat Buzsáki G , Anastassiou CA , Koch C. The origin of extracellular fields and currents—EEG, ECoG, LFP and spikes . Nat Rev Neurosci 2012 ; 13 : 407 – 20 . Google Scholar Crossref Search ADS PubMed WorldCat Chen JL , Ros T , Gruzelier JH. Dynamic changes of ICA‐derived EEG functional connectivity in the resting state . Hum Brain Mapp 2013 ; 34 : 852 – 68 . Google Scholar Crossref Search ADS PubMed WorldCat Cole SR , van der Meij R , Peterson EJ , de Hemptinne C , Starr PA , Voytek B. Nonsinusoidal beta oscillations reflect cortical pathophysiology in Parkinson's disease . J Neurosci 2017 ; 37 : 4830 – 40 . Google Scholar Crossref Search ADS PubMed WorldCat Connolly AT , Jensen AL , Bello EM , Netoff TI , Baker KB , Johnson MD , et al. Modulations in oscillatory frequency and coupling in globus pallidus with increasing parkinsonian severity . J Neurosci 2015 ; 35 : 6231 – 40 . Google Scholar Crossref Search ADS PubMed WorldCat de Hemptinne C , Ryapolova-Webb ES , Air EL , Garcia PA , Miller KJ , Ojemann JG , et al. Exaggerated phase–amplitude coupling in the primary motor cortex in Parkinson disease . Proc Natl Acad Sci USA 2013 ; 110 : 4780 – 5 . Google Scholar Crossref Search ADS PubMed WorldCat De Hemptinne C , Swann NC , Ostrem JL , Ryapolova-Webb ES , San Luciano M , Galifianakis NB , et al. Therapeutic deep brain stimulation reduces cortical phase-amplitude coupling in Parkinson's disease . Nat Neurosci 2015 ; 18 : 779 – 86 . Google Scholar Crossref Search ADS PubMed WorldCat DeLong M , Wichmann T. Changing views of basal ganglia circuits and circuit disorders . Clin Eeg Neurosci 2010 ; 41 : 61 – 7 . Google Scholar Crossref Search ADS PubMed WorldCat DeLong MR , Wichmann T. Basal ganglia circuits as targets for neuromodulation in Parkinson disease . JAMA Neurol 2015 ; 72 : 1354 – 60 . Google Scholar Crossref Search ADS PubMed WorldCat Delorme A , Makeig S. EEGLAB: an open source toolbox for analysis of single-trial EEG dynamics including independent component analysis . J Neurosci Methods 2004 ; 134 : 9 – 21 . Google Scholar Crossref Search ADS PubMed WorldCat Devergnas A , Caiola M , Pittard D , Wichmann T. Cortical phase–amplitude coupling in a progressive model of parkinsonism in nonhuman primates . Cereb Cortex 2019 ; 29 : 167 – 77 . Google Scholar Crossref Search ADS PubMed WorldCat Ermentrout GB. n: m Phase-locking of weakly coupled oscillators . J Math Biol 1981 ; 12 : 327 – 42 . Google Scholar Crossref Search ADS WorldCat Farris FA. The Gini index and measures of inequality . Am Math Mon 2010 ; 117 : 851 – 64 . Google Scholar Crossref Search ADS WorldCat Gast R , Schmidt H , Knosche TR. A mean-field description of bursting dynamics in spiking neural networks with short-term adaptation . Neural Comput 2020 ; 32 : 1615 – 34 . Google Scholar Crossref Search ADS PubMed WorldCat Glasser MF , Coalson TS , Robinson EC , Hacker CD , Harwell J , Yacoub E , et al. A multi-modal parcellation of human cerebral cortex . Nature 2016 ; 536 : 171 – 97 . Google Scholar Crossref Search ADS PubMed WorldCat Goetz CG , Tilley BC , Shaftman SR , Stebbins GT , Fahn S , Martinez‐Martin P , et al. Movement Disorder Society‐sponsored revision of the Unified Parkinson's Disease Rating Scale (MDS‐UPDRS): scale presentation and clinimetric testing results . Mov Disord 2008 ; 23 : 2129 – 70 . Google Scholar Crossref Search ADS PubMed WorldCat Haynes WI , Haber SN. The organization of prefrontal-subthalamic inputs in primates provides an anatomical substrate for both functional specificity and integration: implications for Basal Ganglia models and deep brain stimulation . J Neurosci 2013 ; 33 : 4804 – 14 . Google Scholar Crossref Search ADS PubMed WorldCat Jackson N , Cole SR , Voytek B , Swann NC. Characteristics of waveform shape in Parkinson’s disease detected with scalp electroencephalography . eNeuro 2019 ; 6 : ENEURO.0151-19.2019 . Google Scholar Crossref Search ADS PubMed WorldCat Jonmohamadi Y , Jones RD. Source-space ICA for MEG source imaging . J Neural Eng 2016 ; 13 : 016005 . Google Scholar Crossref Search ADS PubMed WorldCat Jonmohamadi Y , Poudel G , Innes C , Jones R. Source-space ICA for EEG source separation, localization, and time-course reconstruction . NeuroImage 2014 ; 101 : 720 – 37 . Google Scholar Crossref Search ADS PubMed WorldCat Kondylis ED , Randazzo MJ , Alhourani A , Lipski WJ , Wozny TA , Pandya Y , et al. Movement-related dynamics of cortical oscillations in Parkinson’s disease and essential tremor . Brain 2016 ; 139 : 2211 – 23 . Google Scholar Crossref Search ADS PubMed WorldCat Langdon AJ , Boonstra TW , Breakspear M. Multi-frequency phase locking in human somatosensory cortex . Prog Biophys Mol Biol 2011 ; 105 : 58 – 66 . Google Scholar Crossref Search ADS PubMed WorldCat López-Azcárate J , Tainta M , Rodríguez-Oroz MC , Valencia M , González R , Guridi J , et al. Coupling between beta and high-frequency activity in the human subthalamic nucleus may be a pathophysiological mechanism in Parkinson's disease . J Neurosci 2010 ; 30 : 6667 – 77 . Google Scholar Crossref Search ADS PubMed WorldCat Magrinelli F , Picelli A , Tocco P , Federico A , Roncari L , Smania N , et al. Pathophysiology of motor dysfunction in Parkinson’s disease as the rationale for drug treatment and rehabilitation . Parkinsons Dis 2016 ; 2016 : 1 – 18 . Google Scholar Crossref Search ADS WorldCat Malekmohammadi M , AuYong N , Ricks-Oddie J , Bordelon Y , Pouratian N. Pallidal deep brain stimulation modulates excessive cortical high β phase amplitude coupling in Parkinson disease . Brain Stimul 2018 ; 11 : 607 – 17 . Google Scholar Crossref Search ADS PubMed WorldCat Miller AM , Miocinovic S , Swann NC , Rajagopalan SS , Darevsky DM , Gilron R , et al. Effect of levodopa on electroencephalographic biomarkers of the parkinsonian state . J Neurophysiol 2019 ; 122 : 290 – 9 . Google Scholar Crossref Search ADS PubMed WorldCat Miller KJ , Sorensen LB , Ojemann JG , den Nijs M. Den Nijs M. Power-law scaling in the brain surface electric potential . PLoS Comput Biol 2009 ; 5 : e1000609 . Google Scholar Crossref Search ADS PubMed WorldCat Nambu A , Takada M , Inase M , Tokuno H. Dual somatotopical representations in the primate subthalamic nucleus: evidence for ordered but reversed body-map transformations from the primary motor cortex and the supplementary motor area . J Neurosci 1996 ; 16 : 2671 – 83 . Google Scholar Crossref Search ADS PubMed WorldCat Oldfield RC. The assessment and analysis of handedness: the Edinburgh inventory . Neuropsychologia 1971 ; 9 : 97 – 113 . Google Scholar Crossref Search ADS PubMed WorldCat Oostenveld R , Fries P , Maris E , Schoffelen J-M. FieldTrip: open source software for advanced analysis of MEG, EEG, and invasive electrophysiological data . Comput Intell Neurosci 2011 ; 2011 : 1 – 9 . Google Scholar Crossref Search ADS PubMed WorldCat Oswal A , Brown P , Litvak V. Synchronized neural oscillations and the pathophysiology of Parkinson's disease . Curr Opin Neurol 2013 ; 26 : 662 – 70 . Google Scholar Crossref Search ADS PubMed WorldCat Ozturk M , Abosch A , Francis D , Wu J , Jimenez‐Shahed J , Ince NF. Distinct subthalamic coupling in the ON state describes motor performance in Parkinson's disease . Mov Disord 2020 ; 35 : 91 – 100 . Google Scholar Crossref Search ADS PubMed WorldCat Pasquereau B , Turner RS. Primary motor cortex of the parkinsonian monkey: differential effects on the spontaneous activity of pyramidal tract-type neurons . Cerebral Cortex 2011 ; 21 : 1362 – 78 . Google Scholar Crossref Search ADS PubMed WorldCat Postuma RB , Berg D , Stern M , Poewe W , Olanow CW , Oertel W , et al. MDS clinical diagnostic criteria for Parkinson's disease . Mov Disord 2015 ; 30 : 1591 – 601 . Google Scholar Crossref Search ADS PubMed WorldCat Ray S , Crone NE , Niebur E , Franaszczuk PJ , Hsiao SS. Neural correlates of high-gamma oscillations (60–200 Hz) in macaque local field potentials and their potential implications in electrocorticography . J Neurosci 2008 ; 28 : 11526 – 36 . Google Scholar Crossref Search ADS PubMed WorldCat Schroeder CE , Lakatos P. Low-frequency neuronal oscillations as instruments of sensory selection . Trends Neurosci 2009 ; 32 : 9 – 18 . Google Scholar Crossref Search ADS PubMed WorldCat Shimamoto SA , Ryapolova-Webb ES , Ostrem JL , Galifianakis NB , Miller KJ , Starr PA. Subthalamic nucleus neurons are synchronized to primary motor cortex local field potentials in Parkinson's disease . J Neurosci 2013 ; 33 : 7220 – 33 . Google Scholar Crossref Search ADS PubMed WorldCat Shreve LA , Velisar A , Malekmohammadi M , Koop MM , Trager M , Quinn EJ , et al. Subthalamic oscillations and phase amplitude coupling are greater in the more affected hemisphere in Parkinson’s disease . Clin Neurophysiol 2017 ; 128 : 128 – 37 . Google Scholar Crossref Search ADS PubMed WorldCat Spiegler A , Knösche TR , Schwab K , Haueisen J , Atay FM. Modeling brain resonance phenomena using a neural mass model . PLoS Comput Biol 2011 ; 7 : e1002298 . Google Scholar Crossref Search ADS PubMed WorldCat Swann NC , de Hemptinne C , Aron AR , Ostrem JL , Knight RT , Starr PA. Elevated synchrony in Parkinson disease detected with electroencephalography . Ann Neurol 2015 ; 78 : 742 – 50 . Google Scholar Crossref Search ADS PubMed WorldCat Tort AB , Kramer MA , Thorn C , Gibson DJ , Kubota Y , Graybiel AM , et al. Dynamic cross-frequency couplings of local field potential oscillations in rat striatum and hippocampus during performance of a T-maze task . Proc Natl Acad Sci USA 2008 ; 105 : 20517 – 22 . Google Scholar Crossref Search ADS PubMed WorldCat Tsiokos C , Malekmohammadi M , AuYong N , Pouratian N. Pallidal low β-low γ phase-amplitude coupling inversely correlates with Parkinson disease symptoms . Clin Neurophysiol 2017 ; 128 : 2165 – 78 . Google Scholar Crossref Search ADS PubMed WorldCat Van Veen BD , Buckley KM. Beamforming: a versatile approach to spatial filtering . IEEE Assp Mag 1988 ; 5 : 4 – 24 . Google Scholar Crossref Search ADS WorldCat van Wijk BC , Beudel M , Jha A , Oswal A , Foltynie T , Hariz MI , et al. Subthalamic nucleus phase–amplitude coupling correlates with motor impairment in Parkinson’s disease . Clin Neurophysiol 2016 ; 127 : 2010 – 9 . Google Scholar Crossref Search ADS PubMed WorldCat Voytek B , Knight RT. Dynamic network communication as a unifying neural basis for cognition, development, aging, and disease . Biol Psychiatry 2015 ; 77 : 1089 – 97 . Google Scholar Crossref Search ADS PubMed WorldCat West TO , Berthouze L , Halliday DM , Litvak V , Sharott A , Magill PJ , et al. Propagation of beta/gamma rhythms in the cortico-basal ganglia circuits of the parkinsonian rat . J Neurophysiol 2018 ; 119 : 1608 – 28 . Google Scholar Crossref Search ADS PubMed WorldCat BGTC basal-ganglia-thalamocortical circuit DLPFC dorsolateral prefrontal cortex ECoG electrocorticography ICA independent component analysis KL-MI Kullback-Leibler-based modulation index M1 primary motor cortex; MDS-UPDRS III = Movement Disorders Society Unified Parkinson’s Disease rating scale part III PAC phase-amplitude coupling PLA phase-locked amplitude PMC premotor cortex © The Author(s) (2020). Published by Oxford University Press on behalf of the Guarantors of Brain. All rights reserved. For permissions, please email: journals.permissions@oup.com This article is published and distributed under the terms of the Oxford University Press, Standard Journals Publication Model (https://academic.oup.com/journals/pages/open_access/funder_policies/chorus/standard_publication_model) TI - Spatiotemporal features of β-γ phase-amplitude coupling in Parkinson’s disease derived from scalp EEG JF - Brain DO - 10.1093/brain/awaa400 DA - 2021-03-03 UR - https://www.deepdyve.com/lp/oxford-university-press/spatiotemporal-features-of-phase-amplitude-coupling-in-parkinson-s-LP76CD2yoW SP - 487 EP - 503 VL - 144 IS - 2 DP - DeepDyve ER -