Bayesian Optimization of Machine Learning Classification of Resting-State EEG Microstates in Schizophrenia: A Proof-of-Concept Preliminary Study Based on Secondary Analysis

Resting-state electroencephalography (EEG) microstates reflect sub-second, quasi-stable states of brain activity. Several studies have reported alterations of microstate features in patients with schizophrenia (SZ). Based on these findings, it has been suggested that microstates may represent neurophysiological biomarkers for the classification of SZ. To explore this possibility, machine learning approaches can be employed. Bayesian optimization is a machine learning approach that selects the best-fitted machine learning model with tuned hyperparameters from existing models to improve the classification. In this proof-of-concept preliminary study based on secondary analysis, 20 microstate features were extracted from 14 SZ patients and 14 healthy controls’ EEG signals. These parameters were then ranked as predictors based on their importance, and an optimized machine learning approach was applied to evaluate the performance of the classification. SZ patients had altered microstate features compared to healthy controls. Furthermore, Bayesian optimization outperformed conventional multivariate analyses and showed the highest accuracy (90.93%), AUC (0.90), sensitivity (91.37%), and specificity (90.48%), with reliable results using just six microstate predictors. Altogether, in this proof-of-concept study, we showed that machine learning with Bayesian optimization can be utilized to characterize EEG microstate alterations and contribute to the classification of SZ patients.


Introduction
Schizophrenia (SZ) is a psychiatric disorder whose neurobiological underpinnings are still largely unknown. One of the most widely used techniques in SZ research is electroencephalography (EEG), which measures electrical activity in the brain. EEG signals reflect the oscillation and large-scale synchronization of underlying neural populations [1] and therefore can be used to investigate the synchrony and dynamics of neural circuits in SZ patients [2,3]. In patients with SZ, EEG abnormalities have been reported from studies measuring event-related potential paradigms, including face processing and mismatch negativity [4,5], sleep EEG characteristics [6,7], graph analysis measures of resting state EEG [8], and EEG power spectra, which are the most commonly computed features [9][10][11]. Statistical features, such as mean, skewness, and kurtosis were also investigated in SZ [12,13]. However, most of these EEG-based features have failed to take advantage of the EEG millisecond resolution [14][15][16][17][18].
Microstates are topographic maps of resting state brain activity that exploit the EEG temporal resolution. Specifically, Lehmann et   We identified four microstates (i.e., A, B, C, and D, Figure 2), which reflected the activity of the EEG channels with almost 70% of the total topographic variance, and labeled the topography at each GFP peak as one of these microstates. Five categories of features were computed for these microstates (Figures 3 and S1): (1) the average number of times per second that each microstate occurred during the EEG recording (occurrence); (2) the average amount of time each microstate lasted after it occurred (duration); (3) the percentage of total recording time in which each microstate was dominant (coverage); [19] (4) the average GFP during microstate dominance (mean GFP) [28]; and (5) the average

Microstate Analysis
The entire EEG signal can be represented by a small set of topographic maps that alternate at discrete intervals when EEG is viewed as a topography of electric poten-tials evolving over time [46]. These topographic maps are named microstates. These images/topographies do not gradually merge or overlap in time; instead, a single map dominates for roughly 80-120 milliseconds before abruptly switching to another map (i.e., microstates are quasi-stable periods of single maps) [19]. To extract the relevant maps, the root-mean-squared potential differences at all N electrodes (i.e., V i (t)) from the mean of instantaneous potentials across electrodes (i.e., V mean (t)) were computed to characterize the EEG global field power (GFP) [47].
Topographies of the local maxima of the GFP curve were identified based on the above equation. We then used a modified K-means clustering approach and global explained variance (GEV) criteria [27,48,49] by setting the re-run parameter to 20 times, the convergence threshold at 10 −6 , and the maximum number of iterations to 1000 [33,44].
We identified four microstates (i.e., A, B, C, and D, Figure 2), which reflected the activity of the EEG channels with almost 70% of the total topographic variance, and labeled the topography at each GFP peak as one of these microstates. Five categories of features were computed for these microstates (Figures 3 and S1): (1) the average number of times per second that each microstate occurred during the EEG recording (occurrence); (2) the average amount of time each microstate lasted after it occurred (duration); (3) the percentage of total recording time in which each microstate was dominant (coverage); [19] (4) the average GFP during microstate dominance (mean GFP) [28]; and (5) the average correlation between each labeled GFP peak map (i.e., A) and the corresponding microstate template (microstate map correlations) [34]. All microstate features are summarized in Table 1.   Four normalized microstates (i.e., A, B, C, and D) of resting-state EEG recordings were obtained for patients diagnosed with schizophrenia and healthy control subjects; SZ: patients diagnosed with schizophrenia; HC: healthy control subjects.

Figure 3.
Results of the predictor importance scores for twenty microstate features extracted from resting state EEG recordings of SZ patients and HC subjects.

Predictors' Ranking and Classification
We ranked the 20 features/predictors described above using chi-square tests ("fsc-chi2") embedded in the Statistics and Machine Learning Toolbox of MATLAB software (MathWorks, Inc., Natick, MA, USA, 2020) ( Figure 3). Then, an optimized machine-learning approach was applied to classify the combination of ranked microstate features for SZ classification by adding each feature incrementally. Of note, machine learning classification can be subject dependent or subject independent. In the subject-independent classification, the model is tested using unseen subject's data, (i.e., the data that are used to test the model are not included in the training phase), while here in the subject-dependent classification, training and testing sets are randomly split, and EEG-segment data from the same subject are included in both sets [31]. Because our cohort of SCZ and HC subjects was rather small, subject-dependent classification was employed. In each run, 80% of the data was utilized in the training process with a 5-fold cross-validation method to avoid overfitting, and the remaining data (20%) were used to test for the best-fitted trained model. Histograms of subject-based features used for this study to show the within-group subject variabilities of microstate features were also computed (see Figures S5 and S6). Accuracy, area under the curve (AUC), sensitivity, and specificity were used to assess classifier performance. Furthermore, we used the leave-one-out validation method by excluding the data of one individual per group each time to examine the subject-independent design in addition to the current approach. The trained models on the remaining dataset were tested on these two individuals, who each time were left out. Our results are comparable in both designs (output measures in the range of about 90% for subject independent vs. subject dependent) in Table S3. The classifier algorithms were selected from a pool of existing classifiers based on prior EEG studies that highlighted the important application role of ML as a real-time health monitoring system for stroke prognostics [50], access ischemic stroke-derived cortical impairment [51] and other biomedical engineering works [50][51][52][53][54][55]. The Shapley additive explanations (SHAP) or Shapley values of features were computed using the "shapley" function embedded in the Statistics and Machine Learning Toolbox of MATLAB2022b, which explains the deviation of the prediction for the query points from the average prediction, due to the feature ( Figure 4). For each query point, the sum of the Shapley values for all features corresponds to the total deviation of the prediction from the average [56,57]. The details of the machine learning optimization analysis are described below.
subject variabilities of microstate features were also computed (see Figures S5 and S6). Accuracy, area under the curve (AUC), sensitivity, and specificity were used to assess classifier performance. Furthermore, we used the leave-one-out validation method by excluding the data of one individual per group each time to examine the subject-independent design in addition to the current approach. The trained models on the remaining dataset were tested on these two individuals, who each time were left out. Our results are comparable in both designs (output measures in the range of about 90% for subject independent vs. subject dependent) in Table S3. The classifier algorithms were selected from a pool of existing classifiers based on prior EEG studies that highlighted the important application role of ML as a real-time health monitoring system for stroke prognostics [50], access ischemic stroke-derived cortical impairment [51] and other biomedical engineering works [50][51][52][53][54][55]. The Shapley additive explanations (SHAP) or Shapley values of features were computed using the "shapley" function embedded in the Statistics and Machine Learning Toolbox of MATLAB2022b, which explains the deviation of the prediction for the query points from the average prediction, due to the feature ( Figure 4). For each query point, the sum of the Shapley values for all features corresponds to the total deviation of the prediction from the average [56,57]. The details of the machine learning optimization analysis are described below.

Bayesian Optimization of Classification
To perform Bayesian optimization of machine learning classification, we employed an algorithm comprising 2 principal steps (Equations (S2) and (S3) in Table S1), where DATA 1:t−1 = {x n , y n } t − 1 n = 1 defined the training dataset with the t − 1 observation of an unknown function (Table S1 and Figure S2). To automatically choose the machine learning algorithm with tailored hyperparameters, the Statistics and Machine Learning Toolbox TM (MATLAB and Release 2020b, The MathWorks, Inc., Natick, MA, USA) was utilized with the "fitcauto" function [37,38,56]. A multi-TreeBagger model of the objective function was included in the Bayesian optimization approach of "fitcauto". The objective function of this model differed from the Gaussian process model implemented by other machine learning toolbox functions using Bayesian optimization, and the next point to be examined was determined by an acquisition function (i.e., expected improvement). The output of the "fitcauto" algorithm was the point with the lowest objective function value among the points assessed during the optimization. This method automatically chose the best machine learning method for training data from among the most applicable machine learning methods (e.g., discriminant analysis . When the optimization process was completed, "fitcauto" returned the trained model for the entire train dataset to perform classification [57].

Results
Four microstates were identified: A, B, C, and D that exhibited right-frontal leftposterior, left-frontal right-posterior, midline frontal-occipital, and midline frontal topographies respectively. Microstate scalp topographies of patients diagnosed with SZ and healthy control subjects are shown in Figure 2. We also computed and compared twenty microstate features, including occurrence, duration, coverage, and microstate correlation maps between SZ patients and HC subjects. We found that most of the predictors were significantly different between groups after Bonferroni correction for multiple comparisons (α < 0.008, Table 2, Figure S1). We then ranked microstate parameters, which indicated that the features from microstates C and D had the highest predictor importance scores (Figures 3 and 4). Specifically, Ocurrence_C, Coverage_C, MsMC_C, Duration_C, MsMC_B, and Coverage_D were the highest ranked parameters. Furthermore, to compare the performance obtained using ranked microstate features in classifying patients diagnosed with SZ and HC subjects, we implemented models that incrementally considered ranked features (e.g., model 1 included only the Occurrence_C feature, and model 2 considered Occurance_C and Coverage_C as input) The number of features that were fed to the machine learning approach based on the ranked order is presented in Figure 3. Furthermore, the contributions of the value of the feature to the difference between the actual prediction and the mean prediction is estimated as Shapley values in Figure 4. The x-axis indicates the variable name, and the value next to them is the mean SHAP value. On the y-axis is the SHAP value that indicates how much the change in features can positively or negatively affect the probability of prediction.
Classification performance using the optimized machine learning approach relative to the quadratic SVM [27,58] showed that the highest output measure results were obtained when using 19 ranked features as the input of the optimized ML algorithm [ACC = 90.93%, AUC = 0.90%, sensitivity = 91.37%, specificity = 90.48%] (Tables 3 and S2, Figures 5, S3 and S4). Furthermore, comparable results in terms of accuracy, sensitivity, and ACU were obtained with the optimized ML algorithm by using the first six ranked features. The best-fitted model selected by the optimized algorithm using 19 features was SVM with Gaussian kernel and 'Ensemble' when using the 6 most important features, while quadratic SVM was used in the recent study on 19 microstate features.   (1) shows the best-fitted model results acquired for using six ranked features. Number (2) shows the highest output measures when using 19 input features.

Discussion
We employed an optimized machine learning approach to microstate measures and examined their potential for SZ classification. By applying the Bayesian optimized machine learning approach to ranked microstate measures, we were able to discern resting EEG recordings of SZ from HC subjects with high sensitivity, specificity, and accuracy ( Table 3, Figures 5 and S2). We also established that with only six features, we could efficiently classify SZ using microstate analyses (Figures 3-5 and S3). Overall, findings from this proof-of-concept study show that optimized ML applied to microstate features could contribute to the identification of patients with SZ relative to HC subjects.
In line with previous studies [18,27,28,30,59], we found four microstates (i.e., A, B, C, and D) that had a similar topography in SZ and HC groups. These four microstates explained the global topographic variance and have been suggested to represent distinct functions of the brain [25,27,48]. Specifically, microstates A and B have been associated with the processing of different sensory modalities and with the mental visualization of the situation [60], whereas types D and C have been implicated in attention regulation and default mode functionality, respectively [27]. Although the types and topographies of these microstates were similar between SZ and HC subjects, individual microstate features differed across groups. For instance, the occurrence of microstate B was significantly increased in SZ vs. HC groups, in line with another study that also reported an association of this altered microstate parameter with the positive symptoms of patients with SZ [30]. We also found the mean GFP was significant in SZ vs. HC individuals across all four types  (1) shows the best-fitted model results acquired for using six ranked features. Number (2) shows the highest output measures when using 19 input features.

Discussion
We employed an optimized machine learning approach to microstate measures and examined their potential for SZ classification. By applying the Bayesian optimized machine learning approach to ranked microstate measures, we were able to discern resting EEG recordings of SZ from HC subjects with high sensitivity, specificity, and accuracy (Table 3, Figures 5 and S2). We also established that with only six features, we could efficiently classify SZ using microstate analyses (Figures 3-5 and S3). Overall, findings from this proofof-concept study show that optimized ML applied to microstate features could contribute to the identification of patients with SZ relative to HC subjects.
In line with previous studies [18,27,28,30,59], we found four microstates (i.e., A, B, C, and D) that had a similar topography in SZ and HC groups. These four microstates explained the global topographic variance and have been suggested to represent distinct functions of the brain [25,27,48]. Specifically, microstates A and B have been associated with the processing of different sensory modalities and with the mental visualization of the situation [60], whereas types D and C have been implicated in attention regulation and default mode functionality, respectively [27]. Although the types and topographies of these microstates were similar between SZ and HC subjects, individual microstate features differed across groups. For instance, the occurrence of microstate B was significantly increased in SZ vs. HC groups, in line with another study that also reported an association of this altered microstate parameter with the positive symptoms of patients with SZ [30]. We also found the mean GFP was significant in SZ vs. HC individuals across all four types of microstates, likely indicating that SZ patients have a higher level of synchronization (i.e., increased cortical power) during resting EEG recordings [61]. The MsCM, a new microstate feature that was calculated for this study showed higher values in HC vs. SZ across types A, B, and C, thus suggesting that GFP topographies are more consistently repeated (i.e., less variability in topographic patterns) in HC relative to SZ patients. Besides MsCM, Type C had reduced occurrence and coverage in SZ vs. HC individuals. These findings are in line with results from a recent meta-analyses of microstate research [18,30]. Of note, microstate C has been linked to the functionality of the saliency network, including the anterior cingulate, inferior frontal gyrus, and insula [20], and aberrant activity in the salience network has been consistently reported in SZ [62][63][64]. Thus, alterations in microstate C parameters further point to dysfunction in this network in SZ. Here, microstate Type D showed an increase in occurrence, duration, and coverage, while previous studies reported a decrease in these parameters [18,27]. Although these discrepancies could be related to methodological differences as well as medication status [15,65], other microstate studies reported an increase in these Type D characteristics, consistent with our findings [30,66]. Microstate D features have been linked to flexible aspects of attention because of their association with the frontoparietal attention network [20]. Studies indicated that impairments of microstates of class D in SZ are associated with deficits in context update, attentional processes, and executive control, which are often observed in these patients [18,67].
Since we observed alterations in most of the microstate parameters computed, we wanted to assess whether some alterations were more relevant than others in differentiating individuals with SZ from HC. We therefore computed the predictor importance score for each of these 20 parameters and found that features from microstate Type C and D were ranked the highest. A reduction in the occurrence and the coverage of microstate Type C were the two top-ranked features, while coverage and duration of microstate D, both of which were increased in SZ vs. HC, were also highly ranked. Given the implication of types C and D in default mode and dorsal attention respectively, these findings suggest that alterations in these domains may be more relevant for SZ classification [20,22]. At the same time, some microstate B features were ranked high as well, including the MsMC that was computed in this study for the first time, and therefore these parameters should be considered and assessed in future studies of SZ. In contrast, Type A microstate features were ranked lower in our study, a finding in agreement with previous work, showing that these features were relatively intact in patients with SZ vs. HC [27].
Multivariate analysis based on machine learning algorithms provides an opportunity to understand the SZ classification by analyzing many features simultaneously. A handful of studies have utilized microstate features and have tested their accuracy and precision to classify SZ with multivariate patterns and have suggested the efficacy of this approach [14,27]. Building on this body of evidence, in this study, we used this machine learning approach to evaluate multiple, ranked microstate features at the same time. As such, we created a more generalized cross-validated model with increased AUC, sensitivity, and specificity. In particular, the optimized machine learning approach employed here was able to achieve greater than 90% efficacy based on the Gaussian SVM using 19 features and greater than 88% accuracy with the ensemble as the best fitted using just 6 features. Importantly, a recent study that used the same dataset for SZ classification reported highest output measures as [Acc = 75.64%" Sensitivity = 71.93%, and Specificity = 75.50%] with quadratic SVM [27]. Several factors, including the number of EEG trials, the type and number of microstate features (i.e., MsMC), and the feature selection method (i.e., feature importance score calculation) may have contributed to this difference. The accuracy, specificity, and sensitivity scores reported here are comparable with some deep learning approaches that use EEG microstates for SZ classification [68][69][70][71][72][73]. The fact that findings obtained from our optimized ML approach were comparable to deep learning methods in terms of performance but outperformed in computation and time of processing (i.e., our method [order of minutes] vs. [order of hours]) potentially provide a more rapid, efficient way to achieve an optimized SZ classification.
This study has several limitations that should be addressed in future studies. For example, although the number of epochs was large enough for applying the optimized classification approach, the sample size of the SZ and HC groups was rather small; thus, we decided to apply subject-dependent machine learning approach classification on 5-second segments of data. We also chose this approach to make our results more comparable to a previous study using the same dataset [27]. Compared to the subject-dependent method, the subject-independent method may offer greater generalizability regarding learning the HC vs. SZ labels rather than from the individual's signature. Therefore, future work on larger groups of SZ patients is needed to confirm the findings from this proof-ofconcept study on a larger dataset, and also by applying a subject-independent classification. Nonetheless, in this study, we run a leave-one analysis and found that the main findings did not change (Table S3). Additionally, even though age and gender were matched, obtaining enough data to reflect the broader range of ages in both genders is necessary for the generalization of the trained models. Thus, to increase classification accuracy and develop an accurate model reflecting the general population, the classification performance of microstate characteristics should be examined using data from larger cohorts encompassing the lifespan. This will contribute to establish how specific features of individuals with SZ vs. HC are captured by the machine learning classification method presented here, in line with a personalized medicine approach. Furthermore, in the present study, four microstates were identified which encompassed at least 70% of the global explained variance (GEV) ( Figure S7) [14,15,45,49]. While we found that GEV did not significantly change when the number of microstates increased, this could still affect classification performance. Future work should, therefore, also assess whether more than four microstates are identified and whether a different number of archetypes may affect the classification of SZ. Relatedly, in the present study, the characteristics of Type C and D microstates were among the highest ranked features, thus indicating that these parameters may be more reliable diagnostic features than Type A and B features in schizophrenia classification, which eventually could have relevant implication in the day-to-day clinical psychiatry practice [20,22]. Of note, each patient enrolled in this study underwent a medication washout period of at least seven days before the EEG recordings were performed. Nonetheless, future work should confirm these findings in medication-naïve patients and/or more thoroughly assess the possible impact of antipsychotic medications on EEG microstate parameters.

Conclusions
By employing for the first time an optimized machine learning approach on microstate measures of resting EEG recordings we achieved higher accuracy, sensitivity, and specificity of SZ patients compared to conventional classification methods, even with just six microstate predictors. Furthermore, our results showed that ranking microstate features was critical to optimize this process. Further studies should confirm and extend these findings on datasets involving larger cohorts of SZ patients. Eventually, the novel machine learning approach employed here may help establish EEG microstates as neurophysiological biomarkers that contribute to the classification of SZ.
Supplementary Materials: The following supporting information can be downloaded at: https: //www.mdpi.com/article/10.3390/brainsci12111497/s1, Table S1: Bayesian optimization workflow for the training dataset; Table S2: Hyperparameters of the machine learning models; Table S3: Hyperparameters of the machine learning models for subject independent design; Figure S1: The distribution of twenty microstate features that were extracted in this study; Figure S2: Flowchart of input data and Bayesian optimization of ML approach; Figure S3: Results of AUC computed for each model; Figure S4: ROC curves for the best-fitted models; Figure S5: Subject-based histogram of HC microstate features; Figure S6: Subject-based histogram of SZ microstate features; Figure S7: Global explained variance (GEV) rate for choosing different numbers of archetypes map.