The study protocol was approved by the Institutional Review Board of Imam Abdulrahman bin Faisal University. All participants received comprehensive information about the study and voluntarily signed a written consent form in accordance with the approved protocol. Written informed consent was obtained from the individuals for the publication of any potentially identifiable images or data included in this article. The study included 16 participants (10 males and 6 females) with an age of 30.5 ± 9.7. Of these, 5 (all male) had red-green color vision deficiency, while the remaining 11 had normal color vision. All subjects had no history of neurological or ophthalmological diseases. Table 1 summarizes the demographic details of the study participants, including age, gender, and health condition (CVD or healthy).

To record the EEG, an OpenBCI Cyton Biosensing Board (OpenBCI, Brooklyn, NY, United States), as shown in Figure 1A, was used. The board was connected to a PC via a USB dongle. The board has 8 channels and a 32-bit processor, and all the data were sampled at 250 Hz. The EEG was recorded using a cap with dry comb electrodes, shown in Figure 2B, which were placed according to the international 10–20 electrode system. These electrodes were chosen for their compatibility with the OpenBCI board, as well as their ease of use. Four channels were used: O1, O2, Cz, and Pz/, as shown in Figure 2A. Selecting these specific regions due to their significant association with cortical vision in the brain (Norcia et al., 2015; Srinivasan et al., 2006). The reference and ground electrodes were placed on the right and left earlobes, respectively, using conductive gel to enhance conductivity and reduce impedance. The OpenBCI GUI, shown in Figure 1C, was used for data acquisition. The experiment interface was designed in Python (version 3.12.0), and data analysis was conducted using both Python and MATLAB (version R2023b; The MathWorks, Inc., Natick, MA, United States). Figure 2B shows the experimental set-up, with the subject wearing the EEG cap equipped with dry comb electrodes. The subject is exposed to a screen with the stimulation function, while another screen displays the measured data during the experiment.

A specially designed stimulation setup, incorporating a software interface, was used to elicit SSVEPs. The interface, as shown in Figure 3 consists of two squares flickering at different frequencies (15 Hz and 18 Hz) with a series of Ishihara plates positioned above each square. The EEG was recorded for the subjects during the task. Each subject was asked to sit in front of a display in a laboratory with normal lighting. Each subject then focused on a flickering square (15 Hz or 18 Hz) positioned below a pre-identified number to assess whether they selected the correct plate, which indicated if the subject had CVD. Each subject underwent six sessions, with each session consisting of 10 trials. During each trial, subjects were asked to focus on the flickering square for 10 s, followed by a 10 s rest interval. During each trial, subjects were instructed to refrain from blinking, moving, or talking. However, during rest, they were allowed to blink, move, and talk to ensure comfort. By considering the electrode connections, ensuring signal clarity, and conducting the experiment, the total duration of the experiment was approximately 30 min per subject. Figure 4 illustrates the step-by-step procedure followed in this research to achieve its objectives.

The raw EEG data were collected using the OpenBCI GUI software (version 5.2.2), saved as. txt files, and then imported into Python and MATLAB for processing. First, the EEG data for all subjects were filtered using a 4th-order Butterworth bandpass filter (BPF) with a cutoff frequency range of 5–30 Hz, followed by a notch filter to eliminate 60 Hz powerline noise, both applied in the OpenBCI GUI software. Additionally, a second BPF with a range of 5–50 Hz was applied in Python to ensure better signal representation. Next, the EEG data were re-referenced to the average of all channels to improve the SNR. The data were then segmented into 10-s trials for each session. After segmentation, artifact rejection was performed using Independent Component Analysis (ICA) to remove artifacts-related components, such as eye-blinks, from the EEG signals. After the ICA-based removal process, visual inspection and amplitude-based rejection were conducted to discard trials in which any EEG channel exceeded ±100 μV. After that the 10 trials for each session were averaged to improve the SNR.

The SSVEP data segmentation was performed using Python. First, during data recording, a trigger was utilized to separate the data into two distinct conditions as previously mentioned. When the subject focused on the stimulation, the trigger was set to 0; conversely, when the subject was not stimulated, the trigger was set to 1 in each session. The data was segmented into 10-s intervals, alternating between “trial” (stimulation) and “rest” (no stimulation) segments. This structured segmentation approach facilitated data organization for subsequent analysis. With 10 trial and 10 rest periods in each session, the averages of these segments were calculated separately for each session and used in the classification process.

Feature extraction was conducted using traditional target identification algorithms, including fast Fourier transform (FFT), power spectral density (PSD), and canonical correlation analysis (CCA). PSD, a technique widely employed for detecting SSVEPs, allows for the extraction of frequency information from EEG signals. PSD values were computed using the Welch method with a Hann window to mitigate sharp transition effects on frequency content, and the results were subsequently converting to µV²/Hz. The 10-s EEG data windows (2,500 samples) were divided into 19 segments (250 samples) with a 50% overlap. On the other hand, CCA, known for its rapidity and straightforward integration, was utilized to identify target frequencies. CCA facilitates the identification of the target stimulation frequency through correlation coefficient analysis by constructing reference signals using sine and cosine templates. These methods were utilized for their efficacy in handling multichannel EEG signals and optimizing electrode and time parameters for SSVEP analysis (Hakvoort, 2010; Scarpino et al., 1990; Ma et al., 2022).

For PSD, the equation for the periodogram of each segment is given by Welch’s Method (2024), as shown in Equation 1. P k = X k 2   M (1) where P k represents the periodogram at frequency k , X k 2 denotes the magnitude squared of the Fourier transform of the signal at frequency k , and M represents the total number of samples in the signal.

The Welch estimate of the power spectral density is then given by Welch’s Method (2024), as shown in Equation 2. S k = 1 N × ∑ P n k (2) where S k represents the Welch estimate of the PSD at frequency k , N represents the number of frames used in the Welch method. P n k represents the periodogram of the n t h frame at frequency k . This equation calculates the average of all the periodograms at frequency k .

CCA involves examining two sets of variables X , Y . Unlike other correlation methods, CCA is capable of detecting robust linear relationships between multidimensional variables that might go unnoticed due to the coordinate system used. CCA’s effectiveness lies in identifying pairs of linear transformations ( W x , W y for these variable sets. These transformations maximize the correlation between the coordinate systems when applied (Hakvoort, 2010; Wang et al., 2014).

The resulting projections of these transformations are termed canonical variates and can be expressed as shown in Equation 3 (Hakvoort, 2010). W x , W y = argmax c o r r X W x , Y W y (3) where W x and W y represent the linear transformations, and X and Y represent the sets of variables.

CCA was performed for each stimulation frequency between a set of EEG signals in X and a set of Y f of SSVEP responses, as shown in Equation 4 (Hakvoort, 2010). Y f = sin 2 π f t cos 2 π f t . . . sin 2 π H f t cos 2 π H f t (4) Where f represents the stimulation frequency, H represents the number of harmonics, and T represents the sampling period, and f s represents the sampling rate. While CCA produces several correlation coefficients, the largest one was considered. In this study, the first 3 harmonics were included in the CCA analysis to capture the strongest components of the SSVEP response, as higher harmonics generally have lower SNRs and introduce more noise.

To analyze the neural responses, a two-way repeated measures ANOVA was performed using SPSS Statistics (IBM, Armonk, NY, United States), with frequency (15 Hz and 18 Hz) as a within-subjects factor and participant groups (normal vs CVD) as a between-subjects factor (Park et al., 2009). The statistical significance level (α) for all analyses was set at p ≤ 0.05. Mauchly’s test of sphericity was applied to verify the assumption of sphericity for the within-subjects comparisons, and Greenhouse-Geisser or Huynh-Feldt corrections were applied when necessary to adjust for any violations.

For each subject, a feature matrix was generated from the EEG data, comprising PSD and CCA features. From each selected EEG channel (O1, O2, Cz, and Pz) and target frequencies (15 Hz and 18 Hz), a total of 16 features per subject were extracted. Specifically, this feature set included 8 features from PSD (4 channel × 2 frequencies) and 8 features from CCA (4 channel × 2 frequencies) to efficiently capture relevant characteristics of EEG signals for distinguishing between healthy and CVD subjects. Given the class distribution, the feature matrix contained a total of 176 features for healthy subjects (16 features × 11 subjects) and 80 features for CVD subjects (16 features × 5 subjects).

To address this imbalance and improve classification performance, Synthetic Minority Over-sampling Technique (SMOTE) was implemented to balance the class distribution in our data (Chawla et al., 2002). SMOTE is a widely used oversampling algorithm that proportionally increases the representation of minority classes, thereby balancing the number of samples in each class. In this study, ensuring that the classifier model is trained on an equal number of samples from each class is essential to avoid model bias toward the majority class. To achieve this balance, SMOTE was applied to augment the minority class data (CVD) with synthetic samples. SMOTE was implemented on the combined feature set, generating synthetic samples for the minority class by interpolating between existing PSD and CCA values, resulting in an equal number of samples for each class in the training set. The SMOTE algorithm generates these synthetic samples using the following (Equation 5). x n e w = x i + λ x j − x i (5) where x n e w is the new synthetic sample generated, x i is an existing sample from the minority class, x j is a randomly chosen nearest neighbor of x i within the minority class, and λ is a randomly selected value between 0 and 1 for each synthetic sample, which determines the position of x n e w along the line segment connecting x i and x j .

After balancing the dataset using SMOTE, we applied three distinct classifiers (Decision Tree (DT), Support Vector Machine (SVM), and K-Nearest Neighbors (KNN)) to assess their performance in differentiating CVD from healthy cases. The DT classifier creates models to classify data by representing decision logics. This algorithm takes the form of a tree-like structure, comprising multiple levels, with the top node referred to as the root node. The internal nodes of the tree correspond to tests conducted on input variables, and based on the outcomes of these tests, the algorithm branches out to the appropriate nodes. The leaf nodes of the tree represent the decision outcomes. Building decision trees depends on calculating the entropy based on the data requirements (Uddin et al., 2019; Decision Tree, 2024). Figure 7A shows an illustration of the DT algorithm.

The equation that represents the entropy of a single attribute is shown in Equation 6 (Britannica, 2024). E S = ∑ i = 1 c p i ⁡ log 2 ⁡ p i (6) Where E S represents the entropy function, p i   represents the probability of an event S .

The SVM classifier has the ability to classify data depending on the Kernel method, which utilizes a linear classifier for non-linear problems, making it capable of classifying both linear and non-linear data. It transforms each data point into a different space and identifies a hyperplane that effectively separates the data into different classes. By representing each data point as a coordinate in space, the SVM algorithm determines the optimal hyperplane to distinguish the classes (Uddin et al., 2019; Boateng et al., 2020). Figure 7B shows an illustration of the SVM algorithm.

The equation that represents the SVM classifier is given by Equation 7. K x i , x i ′ = 1 + ∑ j = 1 p x i j x i j ′ ∧ d (7) Where K represents the Kernel function, x i , x i ′ represents the inner products of the training observations, and p represents the degree of the polynomial kernel.

The KNN classifier depends on the number of nearest neighbors, denoted as K . By choosing different values for K , the classification results for a given sample object will vary. The Euclidean distance function is used in the measurement process to assign the case to the nearest neighbor class (Uddin et al., 2019; Boateng et al., 2020). Figure 7C shows an illustration of the KNN algorithm, where the star represents the new object that is classified as black when K = 3 , while it is classified as red when K = 5 .

The equation that represents the Euclidean distance metrics used in the KNN classifier is given by Equation 8 (Boateng et al., 2020). D x , y = ∑ i = 1 N x i 2 _   y i 2 (8) where D x , y represents the distance between the points x and y , and N represents features numbers.

Moreover, the classification conducted in this study was evaluated using various methods. Initially, classification accuracy is a common method to evaluate the performance of classifying EEG signals in BCI applications. We used k-fold cross-validation to train and test extracted features for all classifiers, where k was set to 5. In this approach, the dataset was divided into 5 partitions, with 4 partitions iteratively used for training and the remaining partition used to test the model’s performance. The accuracy can be calculated using Equation 9. A c c u r a c y = N u m b e r   o f   c o r r e c t   S S V E P   r e s p o n s e s   a s s i g n e d   c l a s s e s   T o t a l   n u m b e r   o f   S S V E P   r e s p o n s e s × 100 (9)

To analyze the classification accuracies, - a two-way repeated measures ANOVA was performed using SPSS Statistics (IBM, Armonk, NY, United States), with classifier type (SVM, DT, KNN) and feature set (PSD, CCA, PSD + CCA) as within-subjects factors. The statistical significance level α for all analyses was set by p ≤ 0.05 .

The confusion matrix provides a comprehensive evaluation of the classification model’s performance by offering essential data to calculate precision, recall, and the F1-score, which collectively provide a complete picture of the model’s performance. The matrix consists of different situations (illustrated in Table 3), which are True Positive (TP), True Negative (TN), False Positive (FP), and False Negative (FN) for the three classes.

The F1 score can be calculated by first computing the precision and recall using Equations 10–12. P r e c e s i o n = T P   T P + F P (10) R e c a l l = T P   T P + F N (11) F 1 − s c o r e = 2 × P r e c e s i o n × R e c a l l   P r e c e s i o n + R e c a l l (12)

It is worth noting that precision is also referred to as the positive predictive value, recall is referred to as sensitivity or true positive rate, while specificity is referred to as the true negative rate. Also, sensitivity is associated with TP and FN, while specificity is associated with TN and FP.

SSVEP responses were observed throughout the experiment. Figure 8 illustrates the topographical maps of SSVEP amplitudes under three different conditions: rest, 15 Hz stimulation, and 18 Hz stimulation. These maps display the distribution of the brain activity across the scalp recorded from eight electrodes (O1, O2, Pz, C3, Cz, C4, Fp1, Fp2) placed according to 10–20 system. During the rest condition, the SSVEP activity is predominantly low across most of the scalp, particularly in parietal and occipital regions. However, the central region exhibits a higher amplitude, as indicated by the red coloration in the topographical map. This observation suggests that the brain maintains a certain level of baseline activity even without external visual stimuli. This persistent activity could be attributed to several factors, including alpha wave generation, spontaneous brain processes, and the influence of Default Mode Network (DMN) (Knyazev et al., 2011).

At 15 Hz stimulation, the topographical map reveals a notable increase in SSVEP amplitude, particularly in the occipital region (O1 and O2). This finding aligns with expectations, given that these areas of the brain are responsible for visual processing. The heightened amplitude in these areas indicates that the 15 Hz visual stimulus generates a strong neural response in the visual cortex. Furthermore, the map shows that while the activity is concentrated in the occipital region, it also extends into adjacent parietal areas, suggesting a broad activation of the brain’s visual processing networks when exposed to the 15 Hz frequency. As SSVEPs are primarily generated in the occipital lobe, a relative reduction in activity in the frontal and central regions is observed, as depicted in Figure 8. Similarly, during 18 Hz stimulation, a clear increase in SSVEP amplitude, particularly in the occipital region (O1 and O2), was observed. These observations suggest that the 15 Hz and 18 Hz frequencies effectively elicit strong neural responses and robustly stimulate the visual cortex, making them ideal choices for SSVEP-based CVD diagnostic.

The SSVEP responses were analyzed using four electrodes (O1, O2, Pz, Cz) which exhibited the highest activity during stimulation and are particularly significant for SSVEP analysis. Figure 9 shows the PSD analysis of SSVEP responses for two subjects: one with normal vision and another with CVD. The central portion of the figure presents the visual stimuli, with the target frequency set at 18 Hz. The graphs on either side illustrate the PSD results for the occipital electrodes (O1, O2) and central electrodes (Pz, Cz), providing a comprehensive view of the neural responses recorded from these specific electrode locations. This detailed analysis allows for a deeper understanding of how individuals with normal vision and CVD exhibit distinct brain activity patterns in response to visual stimuli, particularly in relation to the target frequencies presented. For the subject with normal vision, the PSD graph on the left reveals a distinct peak around 18 Hz, particularly evident in the occipital regions (O1 and O2). This robust neural response signifies that the subject’s brain is effectively synchronized with the 18 Hz stimulus, aligning precisely with the intended target frequency for focus and attention. In contrast, for the subject with CVD, the PSD graph on the right displays a peak response around 15 Hz instead of the target 18 Hz. This indicates that the subject focused on the 15 Hz, despite the 18 Hz stimulus being the target. This misalignment underscores the difficulties individuals with CVD encounter in accurately processing visual stimuli, shedding light on the distinct challenges they face in cognitive processing and perception. A distinct peak at 10 Hz (non-target) was observed for both stimulation frequencies but becomes more apparent when the subject focuses on 18 Hz. Moreover, upon comparing the average PSD ratios between normal vision and CVD subjects, the differences in their responses are evident in Figure 10. Ratios greater than 1 signify higher PSD in normal vision subjects, whereas ratios below 1 indicate higher PSD in CVD subjects. The PSD ratio is calculated as the average PSD in normal vision subjects divided by that in CVD subjects. In the 14–16 Hz range, normal vision subjects exhibit higher PSD for the 15 Hz target, whereas CVD subjects display higher PSD for the 18 Hz target. Conversely, in the 17–19 Hz range, normal vision subjects demonstrate higher PSD for the 18 Hz target, while CVD subjects showcase higher PSD for the 15 Hz target.

To examine the observed differences in neural responses, a two-way repeated-measures ANOVA was performed using color vision status (normal vs CVD) and frequency region (14–16 Hz vs 17–19 Hz) as factors. The results revealed a significant interaction effect between color vision status and frequency region on the average PSD ratio F 1.14 = 5.76 , p = 0.021 as shown in Figure 10. This finding indicates that the effect of color vision status on the PSD ratio differs depending on the frequency region being analyzed. Bonferroni-corrected post hoc comparisons further showed that normal vision subjects exhibited significantly higher PSD ratios at the 15 Hz target in the 14–16 Hz range and at the 18 Hz target in the 17–19 Hz range compared to CVD subjects, reflecting a potential difference in visual perception at these frequencies.

The classification accuracy is shown in Figure 11, offering a detailed insight into the performance of three distinct classifiers: DT, SVM, and KNN. These classifiers are evaluated across two pivotal target frequencies, which are 15 Hz and 18 Hz, providing a nuanced understanding of their efficacy at different frequency levels. Moreover, the evaluation is conducted using three feature sets: PSD, CCA, and a combined feature set comprising both PSD and CCA (PSD + CCA). This approach allows for a thorough examination of how these classifiers perform when leveraging distinct feature sets, shedding light on the interplay between feature selection and classification accuracy. This detailed analysis provides valuable insights into the optimal configurations for achieving high accuracy in this classification task. Each classifier is represented by three bars corresponding to the different feature sets, with error bars indicating the standard deviation.

In both conditions at 15 Hz and 18 Hz, the classifiers consistently achieve high accuracy levels, exceeding the 75% threshold. This indicates the effectiveness of the classifiers in accurately classifying data based on the chosen feature sets. At 15 Hz, depicted in Figure 11A, the SVM classifier demonstrates the highest accuracy, reaching 93.75% ± 1.5% with the combined feature set (PSD + CCA), notably outperforming the individual feature sets (PSD and CCA). A similar trend is observed with the DT classifier, which achieves an accuracy of 87.5% ± 1.3% with PSD + CCA, slightly below that of SVM. For the KNN classifier at 15 Hz, both the PSD and PSD + CCA feature sets yield the same accuracy of 81.25%, indicating limited performance for KNN with feature set combination.-When analyzing the results at 18 Hz (Figure 11B), it is evident that the SVM classifier outperforms both DT and KNN across most of feature sets, achieving an accuracy of 87.5% ± 1.3% with PSD + CCA feature set. This consistent superior performance of SVM at 18 Hz suggests its effectiveness in accurately classifying data at this specific frequency. Additionally, DT and KNN classifiers exhibit similar performance, at around 81.25% with the same feature set, though SVM maintains a clear accuracy advantage - when using PSD and PSD + CCA feature set. Overall, the highest accuracy is achieved by the SVM classifier using the combined feature set (PSD + CCA) for both frequencies.-The exceptional performance of SVM when combining multiple feature sets can be attributed to its proficiency in handling high-dimensional feature spaces effectively, a characteristic for which SVM is renowned (Cervantes et al., 2020). A two-way repeated-measures ANOVA was conducted with classifier type (SVM, DT, KNN) and feature set (PSD, CCA, PSD + CCA) as within-subjects factors. The analysis revealed no significant main effect of classifier type on classification accuracy F 2 , 28 = 1.34 , p = 0.28 , indicating that accuracy levels did not differ significantly among the classifiers. Similarly, no significant main effect of feature set was observed F 2 , 28 = 0.92 , p = 0.41 ), suggesting that the choice of feature set (PSD, CCA, or PSD + CCA) did not significantly impact classification accuracy. Moreover, based on the confusion matrices presented in Figure 12, the SVM classifier shows more accurate predictions and fewer misclassifications, particularly when using the combined feature set (PSD + CCA). This is shown by the high rates of correct classification, with SVM achieving 100% accuracy for CVD cases and 90.9% for normal cases using PSD + CCA, demonstrating its effectiveness in distinguishing between the two classes.