In this study, experimental data were derived from the public emotional EEG dataset SEED-IV [32]. SEED-IV is an advanced version of the original SEED, and contains 62-channel EEG signals for four emotions recorded in15 subjects. Scalp voltages were recorded using the ESI NeuroScan System (SynAmps RT, Compumedics Neuroscan Inc., Charlotte, NC, USA) at a 200 Hz sampling rate. The SEED-IV dataset used audio-visual stimuli to elicit four target emotions: neutral, sad, fear, and happy. Subjects were instructed to watch the video clips and experience the corresponding emotions during each trial. Each subject completed three separate sessions, and the stimuli for the three sessions were completely different. A session comprised a total of 24 trials, encompassing six trials per emotion. The duration of each video clip was approximately two minutes, and each trial started with a 5 s hint and ended with a 45 s self-assessment and relaxation. To facilitate data analysis, we organized the four trials related to each of the four different emotions into a block. Hence each session consisted of six blocks.

Generally, unprocessed EEG signals contain electrical noises induced by facial muscle activity and eye movements. EEG signals were preprocessed with the open-source EEGLAB toolbox [33]. First, unnecessary channels were eliminated, and channel locations were edited. The CB1 and CB2 channels were removed due to their nonconformity with the international 10–20 system. Subsequently, the reference electrode standardization technique (REST) [34], a widely recommended reference method, was applied to re-reference EEG signals to the point at infinity. An finite impulse response (FIR) filter was then utilized for bandpass filtering in the 1–70 Hz range, and a 50-Hz notch filter was employed to eliminate power-line interference. Next, independent component analysis (ICA) was executed. Finally, ICLabel [35] and ADJUST (Automatic EEG artifact Detection based on the Joint Use of Spatial and Temporal features) [36] were utilized to remove independent components associated with blinks, eye movements, and electromyography (EMG). Wavelet packet transform was used to decompose the EEG signal of each channel into six sub-bands: Delta (1–4 Hz), Theta (4–8 Hz), Alpha (8–13 Hz), Beta (13–30 Hz), Low Gamma (30–45 Hz), and High Gamma (45–70 Hz). Following this, the 60-channel EEG signals of each frequency band were segmented into 5-second time windows to calculate the phase-amplitude coupling values. The window length of 5 s was determined by considering the maximum inclusion of training samples while capturing genuine phase-amplitude coupling values effectively. A study on the segmentation of EEG signal for classification was also referenced [37]. The Supplementary Material provided further evidence that a 5-second time window was reasonable.

The mean vector length (MVL) [38] is a regular method to quantify the coupling strength between low-frequency signal phase and high-frequency signal amplitude, as it is more sensitive to the coupling strength than other methods [39]. For the original EEG signal filtered into low- and high- frequency signals denoted as $\mathrm{Low}(t)$ and $\mathrm{High}(t)$, Hilbert transform was used to estimate both the instantaneous phase of the low-frequency signal ${\varphi }_{\text{Low}}(t)$ and the instantaneous amplitude of the high-frequency signal $A_{\text{High}}(t)$. The mean vector length is calculated using Eqn. 1:

(1)$$\mathrm{MVL}(f_{\text{Low}},f_{\text{High}})=\left|\frac{1}{T}\sum _{t=1}^T{A}_{\text{High}}(t){\mathrm{e}}^{\mathrm{j}{\varphi }_{\text{Low}}(t)}\right|$$

Where $t$ represents time, and $T$ is the total number of time points. Higher MVL values generally indicate stronger phase-amplitude coupling.

PCB is calculated by averaging the complex vectors of phase angles as shown in Eqn. 2:

(2)$$\mathrm{PCB}=\frac{1}{T}\sum _{t=1}^T{e}^{i\varphi f_{\text{Low}}(t)}$$

Driel et al. [31] introduced a modified version of the regular MVL method, which can uniform the distribution of phase angles so as to debias MVL method. The debiased mean vector length (DMVL) is calculated using Eqn. 3.

(3)$$\mathrm{DMVL}(f_{\text{Low}},f_{\text{High}})=\left|\frac{1}{T}\sum _{t=1}^T{A}_{\text{High}}(t)\left(e^{i\varphi f_{\text{Low}}(t)}-\mathrm{PCB}\right)\right|$$

PCB needs to be linearly subtracted from each phase angle before evaluating the phase-amplitude coupling.

In this study, the phase-amplitude coupling values were calculated using the regular MVL method (hereafter referred to as PAC) and modified version of the regular MVL method (hereafter referred to as DPAC).

For each channel, the PAC and DPAC values were calculated using formulas (1) and (3), respectively, for each 5-s time window. These values were then averaged across the time windows within each emotion for 15 subjects. The Wilcoxon signed rank test was employed to compare the differences among the four emotions for PAC and DPAC, as mentioned in the main text or figure legends. The 95% confidence level threshold was set, and false discovery rate (FDR) was used to adjust the p-value for multiple comparisons. “*” represents p 0.05 and “**” represents p 0.01. The phases of low-frequency signals were randomly shuffled for 1000 iterations and a permutation test was conducted. Statistical Program for Social Sciences (IBM SPSS Statistics 27, IBM Corp., Armonk, NY, USA) was used for statistical analysis.

To further assess the effect of PCB on emotional EEG phase amplitude coupling, PAC and DPAC were used as features for emotion recognition, and their performance was compared. Support Vector Machine (SVM) is a binary linear classifier that distinguishes classes by finding hyperplanes that maximize the margins between two classes [40]. The LIBSVM toolbox (https://www.csie.ntu.edu.tw/~cjlin/libsvm/) is an integrated software that implements support vector classification, and supports multi-class classification [41]. SVM with a linear kernel function and a penalty factor C = 10 was applied to emotion recognition.

When evaluating performance, it is important to note that if the EEG signal of the entire session is randomly split into K folds, both the training and test sets may contain segments from the same trial. Information could be leaked from the training set to the test set due to the correlated data distribution of the segments from the same trial. This may result in an inflated recognition performance on the test set. To evaluate the classification results more reasonably, leave-one-block-out cross-validation method was adopted in this study, where one block refers to four trials of four different emotions.

The EEG signals possess a high dimensional feature space owing to the number of channels. Some features might have a small contribution to the classification task. To mitigate the risk of overfitting and to develop a simpler classification model, the Minimal-Redundancy-Maximal-Relevance (mRMR) [42] algorithm was utilized to rank features according to their importance to the four-class classification tasks. The most important M features were sequentially selected for emotion recognition, recording the classification accuracy for each frequency band pair.

The comodulogram, which depicts the strength of phase-amplitude coupling, is typically used to visualize the interaction between different frequency bands. The horizontal axis represents the phase frequency, while the vertical axis represents the amplitude frequency, and the color corresponds to the coupling strength.

To examine the coupling strength between low-frequency phase and high-frequency amplitude across all frequency bands during different emotional states, phases were extracted from 14 sub-bands within the frequency range of 1–29 Hz (with a bandwidth of 2 Hz and a step size of 2 Hz), and amplitudes from 17 sub-bands within the frequency range of 1–69 Hz (with a bandwidth of 4 Hz and step size of 4 Hz), as determined by the related work of comodulogram [21]. PAC and DPAC were calculated separately for all frequency pairs in each channel. The average comodulograms of four emotions were obtained by averaging the PAC and DPAC across 15 subjects within each emotion. From 60 channels, 60 groups of comodulograms were obtained.

As the frontal lobe is associated with emotions [43], comodulograms of channels in the frontal lobe were averaged. This process was then repeated for DPAC. Fig. 1 depicts the coupling patterns of neutral, sad, fear, and happy emotions. Specifically, the coupling patterns of neutral and happy emotions exhibited similarities, mainly manifesting in the phase frequency range of 1–5 Hz and the amplitude frequency range of 10–30 Hz. On the other hand, the two negative emotions, sad and fear emotions, displayed similar patterns, particularly in the phase frequency range of 1–5 Hz and the amplitude frequency range of 30–69 Hz. As seen from the bottom panel of Fig. 1, the regions of colored in red, which representing the stronger coupling in corresponding frequency ranges, exhibited smaller size when compared to those in the top panel of Fig. 1. This discrepancy can be attributed to the weaker coupling strength after removing PCB.

Fig. 1.

Average emotion comodulograms for channels in the frontal lobe (N = 15 subjects). The top panel shows PAC and the bottom panel shows DPAC. The title of each column gives the emotion. Color bar corresponds to the strength of coupling. Note the similar patterns between PAC and DPAC, and with removal of the phase clustering bias leading to weaker coupling strength. PAC, phase-amplitude coupling; DPAC, debiased phase-amplitude coupling.

As per the emotion comodulograms mentioned earlier, emotional EEG phase-amplitude coupling mainly occurs between six pairs of frequency band: Delta and Beta (D-B), Delta and Low Gamma (D-LG), Delta and High Gamma (D-HG), Theta and Low Gamma (T-LG), Theta and High Gamma (T-HG), and Alpha and High Gamma (A-HG). To further compare the phase-amplitude coupling of these six frequency band pairs, the PAC and DPAC were calculated for all EEG signal channels for the four emotions, and then averaged across channels. As shown in Fig. 2, DPAC and PAC had similar patterns for the four emotions in the six frequency band pairs. However, the Wilcoxon sign ranks test demonstrated that the DPAC values were significantly smaller than those of PAC (Z = 3.408, p = 0.001).

Fig. 2.

Statistical differences among emotions (neutral, blue; sad, orange; fear, yellow; happy, purple) in the average of all channels (N = 15 subjects, **p 0.01, *p 0.05, Wilcoxon signed rank tests). (A) PAC. (B) DPAC. Error bars represent standard error of mean (S.E.M.). A-HG, Alpha and High Gamma; D-B, Delta and Beta; D-LG, Delta and Low Gamma; D-HG, Delta and High Gamma; T-LG, Theta and Low Gamma; T-HG, Theta and High Gamma.

Furthermore, Wilcoxon signed rank tests showed that DPAC and PAC differed significantly between some emotions, as indicated by the significance markers shown in Fig. 2. In particular, after removing PCB, DPAC between fear and happy emotions (Z = 2.360, p = 0.018), as well as neutral and happy emotions (Z = 2.521, p = 0.012), in the D-B frequency band pair both showed significant differences.

As EEG signals have numerous channels the spatial distribution of phase-amplitude coupling strength on the scalp surface can be observed through brain topographic maps. With 60 EEG channels, average brain topographic maps were obtained across the 15 subjects. Figs. 3,4 present the brain topographic maps of PAC and DPAC for six frequency band pairs and four emotions, respectively. Each row denotes an emotion and each column represents a frequency band pair. PAC and DPAC are linearly mapped to the normal scale [0,1].

Fig. 3.

Topographic maps of PAC for 60 channels (N = 15 subjects). The color bars represent the phase-amplitude coupling value, which is linearly mapped to the normal scale [0,1].

Fig. 4.

Topographic maps of DPAC for 60 channels (N = 15 subjects). The color bars represent the phase-amplitude coupling value, which is linearly mapped to the normal scale [0,1].

Figs. 3,4 reveal that PAC and DPAC are strong in both the prefrontal and temporal lobes, and also exist in the occipital lobe. The brain topographic maps of PAC and DPAC also appear to exhibit notable symmetry between the left and right hemispheres.

A comparison of Figs. 3,4 showed that the topographic maps of phase-amplitude coupling exhibit greater differences among the four emotions and the six frequency band pairs both prior to and post the removal of PCB. A statistical analysis of the observed phenomena was also conducted. The difference values of 60 channels between any two frequency band pairs were calculated under each emotion. Similarly, the difference values of 60 channels between any two emotions were calculated under each frequency band pair. For all these difference values, the difference between PAC and DPAC were significant by the Wilcoxon signed rank test (Z = 6.7359, p 0.001). The spatial distributions of phase-amplitude coupling values in D-B and D-LG frequency band pairs were most affected by PCB. Fig. 4 also illustrates that the asymmetry patterns of the left and right hemispheres become more obvious in the D-B frequency band pair.

Since the frontal lobe is frequently reported associated with emotional EEG in previous literature [44], PAC and DPAC of each subject’s frontal lobe channels (F8, F6, F4, F2, FZ, F1, F3, F5, F7) were averaged. Wilcoxon signed rank tests revealed the differences among the four emotions for six frequency band pairs. The significant test results are depicted in Fig. 5.

Fig. 5.

Statistical differences among emotions in the frontal lobe (N = 15 subjects, **p 0.01, *p 0.05, Wilcoxon signed rank tests). (A) PAC. (B) DPAC. The horizontal axis denotes six pairs of the four emotions, namely, neutral-sad (N-S), neutral-fear (N-F), neutral-happy (N-H), sad-fear (S-F), sad-happy (S-H), and fear-happy (F-H), while the vertical axis represents the six frequency band pairs. The color bars represent the p-value.

Results showed that removing PCB enlarged the difference in coupling strength between fear and happy emotions significantly, but seemingly reduced differences in coupling strength between neutral and fear emotions. Notably, a statistically significant difference existed in EEG phase-amplitude coupling in the frontal lobe between neutral and happy emotions, while the p-values are almost equal whether or not PCB was removed. The six p-values in the third column of PAC and DPAC were compared by the Wilcoxon signed rank test and the difference was not significant. This implies that the difference in EEG phase-amplitude coupling in the frontal lobe between neutral and happy emotions was not influenced by PCB (Z = 0.105, p = 0.916). After removing PCB, DPAC between fear and happy emotions (Z = 2.439, p = 0.015), as well as neutral and happy emotions, in the D-B frequency band pair showed significant differences (Z = 2.296, p = 0.022). Similar results were also obtained in Section 3.2.

The emotion recognition outcomes depicted in Fig. 6 present the average classification accuracy of 15 subjects as a function of feature numbers on the horizontal axis. The results indicated that with or without removal of the phase clustering bias, each frequency band pair exhibited a pattern of increasing initially and then decreasing classification accuracy as the feature numbers increased, which could be attributed to classifier overfitting due to the high number of feature dimensions.

Fig. 6.

Emotion recognition performance of PAC and DPAC with different feature numbers. (A) PAC. (B) DPAC. The horizontal axis denotes the number of features, the vertical axis represents the average classification accuracy of 15 subjects.

The optimal number of features with the best classification accuracy of PAC vs. DPAC were 26 vs. 24 (A-HG), 18 vs. 16 (D-B), 14 vs. 12 (D-LG), 14 vs. 14 (D-HG), 22 vs. 18 (T-LG), and 26 vs. 24 (T-HG). Their corresponding accuracy is shown in Table 1. The performance in the four-class emotion recognition task is above the chance level (25%) for each subject. The Wilcoxon signed rank test revealed that the number of features with the best classification accuracy of PAC and DPAC were significantly different (Z = 2.121, p = 0.034), while the difference between their best classification accuracy was not statistically significant (Z = 0.210, p = 0.833). Furthermore, the top-ranked features for both PAC and DPAC were mainly from the prefrontal and lateral temporal lobes, especially FP2 and T7 channels, a phenomenon irrelevant to the subject. Table 1 also showed that D-LG and T-HG achieved higher classification accuracy among the six band pairs of PAC. A-HG and T-LG achieved higher classification accuracy among the six frequency band pairs of DPAC.