Introduction

SZ is a severe chronic mental disorder that occurs most occurs in adolescence or early adulthood1, and clinical symptoms are categorized into three main groups, including positive symptoms, negative symptoms, and general psychopathological symptoms. Negative symptoms earlier than other symptoms are considered core symptoms in SZ, and are associated with poor functional outcomes2,3. Therefore, it is crucial to effectively assess the effectiveness of early treatment of negative symptoms of SZ, to find individual prognostic biomarkers of SZ, and to guide the use of appropriate treatment in the clinic.

Numerous studies4,5,6,7 have attested to the disparities in brain imaging metrics between patients with SZ and healthy HC. These differentiating metrics are pivotal for uncovering individualized prognostic indicators for SZ. A multitude of research approaches zero in on the functional activities of specific brain regions. For instance, techniques such as regional homogeneity (ReHo), amplitude of low-frequency fluctuations (ALFF), and fractional amplitude of low-frequency fluctuations (fALFF)8,9,10 are dedicated to probing the functional activity patterns of each voxel or brain region. ReHo gauges the temporal congruence between a voxel and its neighboring voxels, thereby obliquely reflecting alterations in the synchronization of local neuronal activity. Multiple investigations11,12 have revealed that the ReHo values in the prefrontal cortex are correlated with clinical symptoms and can be harnessed to anticipate treatment responses, laying a theoretical groundwork for the application of machine learning in this ___domain. Guo et al.13 further corroborated this by demonstrating that, through the use of the support vector regression (SVR) method, they could accurately predict treatment responses in SZ patients based on abnormal ReHo values detected in the sensorimotor network and the right putamen. ALFF directly mirrors the intensity of spontaneous blood oxygen level-dependent (BOLD) signal activity within a region, offering valuable insights into local spontaneous neuronal activity. fALFF further mitigates the impact of physiological noise. It has been observed14,15 that regions like the postcentral gyrus and thalamus are linked to therapeutic efficacy in patients. Studies employing SVR analysis have indicated that a substantial decline in fALFF levels just one week after treatment in SZ patients is highly and positively associated with improvements in positive symptoms after eight weeks of treatment. Another analogous study discovered that fluctuations in fALFF within the right caudate nucleus region could serve as an effective biomarker for forecasting treatment responses in SZ patients16,17. These studies enable the identification of region-specific changes within the brain, facilitating a thorough exploration of the functional characteristics of localized brain regions. However, they fall short when it comes to providing a profound comprehension of the intricate interactions occurring between different brain regions.

The other mainstream approach lies in integrating the functional activities between voxels or brain regions. Functional connectivity (FC)18 gauges the co-activity of various brain regions by computing the signal correlations among them, which in turn reveals whether there are any irregularities in the synchronization between these regions. Currently, a significant portion of research on predicting patients’ responses to drugs related to abnormal FC focuses primarily on the striatum and hippocampus19,20. Similar investigations21,22 have also uncovered that enhanced connectivity between the two hemispheres of the brain is associated with more favorable treatment responses. Specifically, abnormalities in inter-hemispheric connectivity, along with certain abnormal values, can be used to forecast how patients with SZ will respond to olanzapine. These research methodologies can offer valuable insights into the functional connectivity across different brain regions and contribute to a better understanding of the structure and function of the entire brain network. However, they cannot highlight individual brain regions throughout the brain in a large amount of redundant information.

The methods mentioned above rely solely on a single feature extracted from resting-state functional magnetic resonance imaging (rs-fMRI). SZ, however, is an extremely complex disorder, and a solitary feature is far from adequate to comprehensively elucidate its multifaceted characteristics. Consequently, there are distinct limitations when it comes to classification accuracy and the prediction of treatment responses. A wealth of research has corroborated that if a model can integrate crucial information from multiple features simultaneously, it will contribute to enhancing prediction accuracy23. Nevertheless, feature fusion may usher in copious amounts of irrelevant information, which in turn escalates the complexity of the model. Hence, feature reduction techniques are essential for minimizing noise and errors, to optimize the model’s performance. Recursive Feature Elimination (RFE) is more result-driven. Once it takes in all the features, it gauges feature importance according to the training outcomes of the model. This enables it to efficiently select features that are highly relevant to the predicted target variable while weeding out those that have only a weak correlation with the target variable.

Current research has relatively neglected the use of multi-feature fusion, along with the integration of pre-treatment and post-treatment information, for treatment response prediction. In this study, we developed a multi-feature fusion recursive feature elimination random forest approach for SZ classification and treatment response prediction. Firstly, we computed three key features: regional homogeneity, fractional amplitude of low-frequency fluctuations, and functional connectivity. Subsequently, a two-sample t-test was implemented to filter out the interfering data elements. We then utilized the recursive feature elimination random forest (RFE-RF) method to conduct SZ classification. Moreover, the RR of the PANSS was adopted to forecast the treatment response of individual patients. Through this process, we were able to identify the neuroimaging biomarkers relevant to SZ classification and treatment response prediction.

The main contributions of this study can be summarized as follows:

  1. (1)

    A multi-feature fusion RFE-RF was proposed for schizophrenia classification and treatment response prediction.

  2. (2)

    The abnormalities in the visual and default mode networks are the potential neuroimaging biomarkers in distinguishing SZ from HC.

  3. (3)

    The abnormalities in the visual network, sensorimotor network, and right superior frontal gyrus are the important biomarkers of short-term treatment response of negative symptoms in SZ patients.

  4. (4)

    The abnormalities in the visual and default mode networks are important biomarkers of the short-term treatment response of positive symptoms.

Materials and methods

Figure 1 is the flow chart, and it includes data acquisition and preprocessing, feature extraction, statistical analysis, recursive feature elimination and random forest, and classification and prediction.

Fig. 1
figure 1

Flowchart of this study.

Participants

190 participants were included in this study. We collected 104 SZ patients from the Second Affiliated Hospital of Xinxiang Medical University in China and 86 HC from nearby communities. The inclusion criteria were set as follows: (1) 18–55 years-old, (2) Han nationality, (3) right-handedness, (4) normal IQ, (5) fulfillment of the clinical diagnostic criteria for schizophrenia in the DSM-IV, and (6) no history of psychiatric disorders in the HC and relatives in the third generation. Exclusion criteria were set as follows: (1) history of head trauma, (2) having neurological disorders, (3) alcohol or drug dependence, (4) suffering from unstable physical conditions, (5) breastfeeding or pregnancy status, and (6) contraindications to MRI. Symptoms and neurocognitive function of SZ were assessed by the PANSS and MCCB. The treatment response was evaluated by RR, and the RR can be calculated as follows

$$RR=\frac{{PANS{S_{Bassline}} - PANS{S_{8weeks}}}}{{PANS{S_{Bassline}}}}$$
(1)

Patients were treated with second-generation antipsychotics for eight weeks before MRI examination, and the selection of dosage and medication was performed by psychiatrists. This study was approved by the Ethics Committee of the Second Affiliated Hospital of Xinxiang Medical University and performed in accordance with the tenets of the Declaration of Helsinki, all patients and control subjects have provided written informed consent.

Data acquisition and preprocessing

We obtained the raw resting-state functional magnetic resonance imaging (rs-fMRI) data using a 3.0T MR scanner (Siemens, Verio) with an echo planar imaging sequence sensitive to blood oxygen level-dependent (BOLD) contrast. The imaging parameters were set as follows: repetition time of 2000 ms, echo time of 30 ms, flip angle of 90°, matrix size of 64 × 64, axial slice resolution of 3.4 × 3.4 mm2, slice thickness of 4 mm, and a 0.6-mm gap between slices. The data collection lasted for 8 min, yielding 240 time points.

The rs-fMRI raw data were preprocessed via the Data Processing Assistant for Resting-State fMRI (DPABI). This preprocessing involved multiple steps, including time point removal, slice timing correction, head motion correction, normalization, spatial smoothing, drift de-linearization, covariate regression, and filtering. At the start of MRI data collection, the scanner requires time to stabilize the magnetic field gradient, while the subjects need to adapt to the scanning environment and reach a resting state. It is commonly acknowledged that discarding the first 10 time points can enhance data quality and stability. Hence, we removed these initial time points. Participants with head movements exceeding 3 mm or rotations greater than 3° during scanning were excluded from the study. To register different subjects to the same standardized space, we utilized the echo planar imaging (EPI) template. Regarding spatial smoothing, an appropriately sized smoothing kernel can eliminate noise and boost the signal-to-noise ratio. However, if the kernel is too large, it may cause signal loss, and if it’s too small, residual noise will remain. A reasonable guideline for choosing the smoothing kernel size is twice the full-width half-maximum (FWHM) of the voxel size24. Given that the voxel size in this experiment was 3, we employed a Gaussian smoothing kernel with an FWHM of 6 mm for smoothing. Finally, we filtered the data of different subjects to mitigate both low-frequency drift and high-frequency noise.

Feature extraction

We divided the whole brain into 116 regions (90 cortical and subcortical non-cerebellar regions and 26 cerebellar regions) using the standardized Automatic Anatomical labeling (AAL) atlas. The fALFF, REHO, and FC were calculated, and Z-Score normalization was performed for these features.

fALFF has been demonstrated to stably correspond to local glucose metabolism in the brain, thus reflecting the brain’s energy expenditure25. fALFF can be calculated as follows. The bold signals of the whole brain voxel were converted into a frequency-___domain power spectrum by the Fourier transform, and the square root of each frequency-___domain power spectrum was calculated to obtain the average square root of the low-frequency (0.01–0.08 Hz) frequency-___domain power spectrum. fALFF is the low-frequency ALFF value divided by the full-band (0.01–0.25 Hz) power spectrum.

$${\text{fALF}}{{\text{F}}_{\text{i}}}=\frac{{\sum\nolimits_{{\text{k}}} {\sqrt {{X_{\text{i}}}{{\left[ {\text{k}} \right]}^2}} } }}{{\frac{1}{N}\sum\nolimits_{{{\text{m}}=0}}^{{N - 1}} {\sqrt {{X_{\text{i}}}{{\left[ {\text{m}} \right]}^{\text{2}}}} } }}{\text{,k}} \in \left[ {0.01\frac{{2\pi }}{{{\text{f}}_{S}^{2}}}{\text{,}}0.08\frac{{2\pi }}{{{\text{f}}_{S}^{2}}}} \right]$$
(2)

ReHo represents the local consistency of the time series of neighboring voxels, which is gauged by Kendall’s Coefficient of Concordance (KCC). Abnormalities in ReHo mirror disruptions in the local coordinated neuronal activity. Indeed, numerous studies have already detected ReHo abnormalities in SZ26,27. The KCC was calculated as follows.

$$W=\frac{{{{\sum {({R_i})} }^2} - n{{{\text{(}}\overline {R} )}^2}}}{{\frac{1}{{12}}{K^2}({n^3} - n)}}$$
(3)

Where n represents the number of time points, \({R_i}\) stands for means the total number of ranks of K voxel points at the i-th time point, \(\overline {R}\) is the mean value of R, and K is the total number of specific voxels and surrounding voxels.

FC serves as an indicator to reveal whether there are any anomalies in the synergistic activities or correlations among different regions in SZ. It can be calculated as follows. The time series of all voxels in each brain region were averaged. Pearson correlation coefficients were adopted to construct the functional connectivity matrix for individual subjects. Pearson correlation (PC) can be calculated by the linear correlation between ROIs as follows.

$$W_{{ij}} = \frac{{\sum\limits_{1}^{n} {(x_{i} - \bar{x}_{i} )^{T} (x_{j} - \bar{x}_{j} )} }}{{\sqrt {\sum\limits_{{i = 1}}^{n} {(x_{i} - \bar{x}_{i} )^{T} (x_{i} - \bar{x}_{i} )} } \sqrt {\sum\limits_{{j = 1}}^{n} {(x_{j} - \bar{x}_{j} )^{T} (x_{j} - \bar{x}_{j} )} } }}$$
(4)

Where \({x_i}\)and \({x_j}\)are the bold signal of any two ROIs, \({\bar {x}_i}\) and \({\bar {x}_j}\) are the average of \({x_i}\) and \({x_j}\).

Statistical analysis

Demographic and clinical data were analyzed using SPSS (Statistical Product and Service Solutions) 26.0 software using chi-square tests and two-sample t-tests (where applicable). Gender, age, and years of education were used as covariates at baseline, and two-sample t-tests were used for between-group comparisons of SZ and HC. Using Gaussian Random Field (GRF) theory based (voxel significance: p < 0.001, clustering significance: p < 0.05) for multiple comparisons to obtain the differential brain regions.

Recursive feature elimination and Random Forest

Recursive Feature Elimination (RFE) constructs a model that assigns importance ratings to each feature within the dataset and then ranks the features according to these ratings. During the RFE process, the model is initially trained on the complete set of original features. Subsequently, the least important feature, as determined by the model, is removed. With the remaining features, a new model is then trained. In each iteration of this process, RFE eliminates the least significant feature based on the feature weights. This iterative removal continues until only one feature remains. Eventually, we rank the features in terms of accuracy and select the subset that yields the highest accuracy, thereby identifying the optimal number of features. The importance ranking is computed using the Python feature importance function. In this function, every time a node splits during model construction, the contribution of each feature to impurity or the Mean Squared Error (MSE) is calculated. These individual contributions are then aggregated to obtain the final relative importance score for all features, which serves as the feature weights. This iterative procedure effectively filters out the most influential features from a large number of candidates, mitigating the impacts of redundancy, noise, and high data dimensionality, and enhancing the model’s generalization performance, accuracy, and predictive capabilities.

Random Forest (RF)28, an ensemble learning method, constructs and integrates multiple decision tree models. In the tree-building process, a bootstrap sampling strategy is adopted to randomly extract multiple sub-datasets from the original feature set. Each sub-dataset is then used to independently train a decision tree. This approach reduces the risk of over-fitting, thereby improving the model’s generalization ability. Additionally, it exhibits good robustness against noise. Random Forest offers high accuracy and strong robustness, conferring significant advantages in both classification and regression tasks. As a result, it has emerged as one of the most widely utilized machine learning algorithms in practical applications.

Classification and prediction

We used RFE-RF for SZ classification and treatment response prediction in this study. In SZ classification, we used the differential brain regions of SZ and HC as the original feature set, adopted a leave-one-out cross-validation strategy, and assessed the performance of the classifier models using accuracy, specificity, and sensitivity. These metrics can be calculated as follow.

$$ACC=\frac{{TP+TN}}{{TP+FP+TN+FN}}$$
(5)
$$SPE=\frac{{TN}}{{TN+FP}}$$
(6)
$$SEN=Recall=\frac{{TP}}{{TP+FN}}$$
(7)

Where TP is true positive, FN is false negative, TN is true negative, FP is and false positive. In this process, we identified the optimal number of features and reported the importance rankings of the top ten features for SZ classification.

In the treatment response prediction for SZ, we focused on abnormal changes in brain regions before and after patient treatment, using pre- and post-treatment information on abnormal brain regions. We employed the Recursive Feature Elimination with Random Forest (RFE-RF) method, combined with leave-one-out cross-validation, to comprehensively assess the performance of our model. The mean squared error (MSE) was selected as the key performance metric, while the RR of the PANSS total score, negative symptom score, and positive symptom score served as the prediction targets. During the feature selection process, we iteratively removed one feature at a time until only a single feature remained. This iterative approach enabled us to pinpoint the feature subset that corresponded to the lowest MSE value, thereby optimizing the model’s predictive ability. With the positive symptom score as one of the central prediction targets, we were able to determine the optimal number of features for the model. In this process, we identified the optimal number of features and again used the feature importance function to report the importance ranking of the top ten features for which the PANSS negative score RR and positive score RR were most relevant to predictive performance.

Result

Subject characteristics

This study consisted of 190 participants. Following strict quality control procedures and the application of exclusion criteria, 88 patients with SZ and 81 HC were ultimately incorporated into the research. The demographic and clinical profiles of the subjects are presented in Table 1. There was no significant difference in age and gender between SZ and HC, and there was a significant difference in years of education (p < 0.001). After eight weeks of treatment period, patients in SZ demonstrated notable improvements in the PANSS total score, negative symptom scores, and positive symptom scores when compared to their baseline scores. These improvements were all highly statistically significant, with P-values less than 0.001, as elaborated in Table 2.

Table 1 Demographic data and clinical characteristics of all subjects.
Table 2 Changes in PANSS scores before and after 8 weeks of treatment.

Multi-feature with statistical differences

For each participant, a multi-feature vector with a dimension of 6902 × 1 was constructed by concatenating several components. Specifically, it included 6670 FC features, 116 ReHo features, and 116 fALFF features. After conducting a two-sample t-test followed by false discovery rate (FDR) correction with a significance level of P < 0.001, 693 functional connection features, 4 ReHo features, and 3 fALFF features exhibited statistical differences, as shown in Fig. 2. Among the top ten brain regions showing statistically significant differences in FC, eight brain regions belonged to the visual network. Regarding the regions with significant differences in ReHo, they were the left inferior temporal gyrus, the precuneus, and the right inferior parietal marginal angular gyrus. As for the fALFF, the brain regions demonstrating statistical differences were the right superior frontal gyrus (orbital part), the right lenticular nucleus (putamen), and the left caudate nucleus.

Fig. 2
figure 2

(A) Corrected functional connectivity; numbers 1-110 are AAL template brain region numbers; (B) corrected ReHo; ITG.L is left inferior temporal gyrus; PCUN is Precuneus; IPL.R is right inferior parietal marginal angular gyrus; (C) corrected fALFF; ORBsup.R is Right Superior frontal gyrus(orbital part); PUT.R is Right Lenticular nucleus(putamen); CAU.L is Left Caudate nucleus.

Classification results

In the feature selection process, we used RFE-RF and pinpointed 181 optimal classification features. These features enabled us to obtain the best training outcomes, with an accuracy of 98.2%, a sensitivity of 97.7%, and a specificity of 98.8%. Moreover, in the SZ classification task, we achieved the best classification results, boasting an accuracy of 91.7%, a sensitivity of 90.9%, and a specificity of 92.6%. In the confusion matrix (Fig. 3), the top-left and bottom-right corners denote the quantity of correctly classified samples, while the top-right and bottom-left corners represent the number of misclassified samples. Currently, SZ lacks specific biomarkers for diagnosis, leading to a certain misdiagnosis rate in clinical practice. When compared with recent experiments29,30, both our false positive rate and false negative rate have declined. This reduction indicates that our RFE-RF approach holds practical value, pending further clinical validation.

Fig. 3
figure 3

perd: 0 indicates that the predicted label is health control; perd: 1 indicates that the predicted label is schizophrenia; True: 1 indicates that the predicted label is schizophrenia; True: 0 indicates that the predicted label is health control.

We further presented the importance rankings of the top ten most discriminative brain regions for SZ classification, Notably, over half of these regions exhibit intra-network connectivity abnormalities within the visual and default mode networks, as detailed in Table 3. Prior research has indicated that abnormalities in visual processing in SZ are caused by dysfunctions in multiple nodes of the visual neural network. These abnormal visual processing mechanisms have been demonstrated to contribute to various manifestations of cognitive dysfunction, as documented in reference31. Abnormalities in the default mode network have also been recurrently reported. Such findings solidify the correlation between default mode network dysfunction and both the positive and negative symptoms of SZ, as supported by reference32.

Table 3 Top 10 most discriminating categorized features.

Treatment response prediction results

We have also clarified the optimal features number for the RR prediction of the PANSS total score, negative symptom score, and positive symptom score. The optimal feature combination for the positive symptom’s RR included only one feature for the ReHo of the left Precuneus, and all the other features were functional connections. It demonstrated that the ReHo of the left Precuneus has a unique role in predicting the RR for the PANSS positive symptoms. The optimal number of the predicted features and MSE are shown in Table 4. Using the best features to predict the treatment response, the predicted RR values of the PANSS total score, negative symptoms, and positive symptoms were all positively correlated with the actual RR values (Fig. 4).

Fig. 4
figure 4

(A) Predicted and actual RR of individual-level PANSS total scale scores were positively correlated (r = 0.982, p < 0.001) (B) Predicted and actual RR of individual-level PANSS negative symptom subscale scores were positively correlated (r = 0.985, p < 0.001) (C) Predicted and actual RR of individual-level PANSS positive symptom subscale scores were positively correlated (r = 0.984, p < 0.001). RR were positively correlated (r = 0.984, p < 0.001).

Table 4 Characterization categories screened and predicted results.

Biomarkers for treatment response prediction

We have also identified the optimal features for predicting the RR of the PANSS negative and positive symptoms after eight weeks of treatment (Fig. 5). As shown in Table 5, the optimal features for predicting the RR of the PANSS negative symptom scores predominantly involved the visual network, the sensorimotor network, and the superior frontal gyrus in default mode network. Regarding the prediction of the RR of the PANSS positive symptom scores, as presented in Table 6, the principal optimal features mainly involved the default mode network and the visual network.

Table 5 The ten most important features of RR for negative symptom scoring.
Table 6 The ten most important features of RR for positive symptom scoring.
Fig. 5
figure 5

(A) FC selected for predicting the RR of the PANSS total scores. (B) FC selected for predicting the RR of the PANSS Negative scores. (C) FC selected for predicting the RR of the PANSS Positive scores.

Discussion

In this study, we discovered that, in contrast to HC, SZ patients predominantly exhibited abnormal intra-network connectivity within the default mode network and the visual network. We utilized these abnormal brain features to distinguish SZ from HC. We successfully predicted the short-term treatment response at the individual level using the information on abnormal brain regions both before and after treatment. This research also pinpointed the neuroimaging biomarkers correlated with positive and negative symptoms during the treatment of SZ patients. Once clinically validated, these biomarkers had significant potential for application in the precise diagnosis and targeted treatment of SZ.

SZ is an intricate psychiatric disorder, entailing a wide array of cognitive, emotional, behavioral, and perceptual abnormalities. The visual network plays a pivotal role in the manifestation of these SZ-related abnormalities33,34,35. The visual network is tasked with processing visual information, facilitating our understanding of the ever-changing external world. It further enables us to concentrate visual attention, filter out irrelevant visual cues, prioritize and process more critical visual data, and store visual memories. The visual network is interconnected with motor control regions, emotional processing centers, and social cognitive areas, among others. Any impairment to this network can disrupt the transmission and integration of brain information, rendering it arduous to interpret information from multiple dimensions, including external perception, executive function, and emotion36. This study revealed that eight of the top ten brain regions, which showed statistically significant differences in functional connectivity between SZ and HC, were part of the visual network. In the classification of SZ and HC, a substantial portion of the top ten brain region features were associated with the visual network and the default mode network. This further underscores the prominent contribution of visual network abnormalities in SZ.

The abnormalities in the visual networks have also been widely observed in previous studies, and there has been substantial evidence of structural and functional abnormalities in the visual network. These abnormalities include structural alterations in the visual network disturbances in functional connectivity37,38,39, and changes in vision information processing40,41. The previous studies42,43 have also identified abnormalities within the default mode networks in SZ patients. Collectively, these findings corroborate our observations regarding the abnormalities in the visual network and the default mode network among SZ patients.

Negative symptoms constitute the core symptoms of SZ, and the main manifestations encompass loss of motivation and pleasure, impoverished thinking, and emotional apathy. These symptoms are a crucial determinant of the clinical and functional prognosis for SZ patients. Hence, treating negative symptoms has an important application value44. In this study, we discovered that the abnormalities in the visual network, sensorimotor network, and the right superior frontal gyrus played a substantial role in predicting the treatment response to negative symptoms in SZ patients after short-term drug-treatment. The sensorimotor network plays a central role in the detection and processing of sensory inputs as well as in the preparation and execution of motor functions. It was the first resting-state brain network to be identified and proposed in fMRI studies45. Multiple investigations have revealed that abnormal indicators of the sensorimotor network and visual network are correlated with negative symptoms in patients46,47,48. Similarly, the right superior frontal gyrus has been repeatedly associated with negative symptoms in SZ49,50. This corroborates our finding that the abnormalities in the visual network, the sensorimotor network, and the right superior frontal gyrus are the neural mechanisms for negative symptoms.

Positive symptoms also form an important part of the symptoms of SZ. Hallucinations, delusions, and abnormal behavior are the most prominent features. The default mode network is known to be activated during the rest state and is associated with individual mental processes51. Many studies have confirmed that the default mode network and visual network are associated with positive symptoms52,53,54. In this study, we specifically selected the ReHo values of the precuneus, which serves as a fundamental node hub within the default mode network. The abnormalities of the precuneus are also widely recognized as playing a role in SZ55,56. This further supports the concept that SZ results from the interaction of two or more brain regions rather than a lesion confined to a single brain region57.

K-fold cross-validation is a prevalent performance assessment approach within the realm of machine learning. This method entails partitioning the dataset into k subsets. Subsequently, k–1 of these subsets are utilized for training the model, while the remaining single subset serves as the test set. This entire procedure is iterated k times, ensuring that each subset gets the opportunity to act as the test set precisely once. By doing so, this technique yields more consistent and trustworthy results for evaluating the performance of a model. Leave-one-out cross-validation represents a special instance of k-fold cross-validation. In this case, every individual data point is successively designated as the test set, with all the other data points being used for training. Given the relatively limited scale of our dataset, we decided to employ leave-one-out cross-validation. This choice was made to optimize the utilization of the available data. However, this method requires higher computational costs. Hence, in the future, it will be necessary to explore more efficient validation methods.

There are some limitations in our research. Firstly, future studies should include more samples to enhance the generalizability of the model. Secondly, our study solely relied on rs-fMRI data, so we cannot determine whether the functional changes are due to disease-related structural alterations in the brain, and we have missed potential information provided by sMRI changes. We plan to incorporate structural brain data in future research. Finally, many studies58 have suggested that demographic factors, such as gender, age, and years of education, are associated with SZ. In our study, the variable of years of education demonstrated significant differences. Although it has been incorporated as a covariate in the analysis, these demographic factors should be independently incorporated into the model in future research initiatives.

Conclusion

In this study, we developed a multi-feature fusion recursive feature elimination random forest for SZ classification and treatment response prediction. We also identified the associated neuroimaging biomarkers. The abnormalities within the visual network and the default mode network are the potential neuroimaging biomarkers in distinguishing SZ from HC. To predict the individual treatment response of patients, we utilized the RR score derived from the PANSS. Through this approach, we were able to pinpoint the biomarkers conducive to treatment response prediction. Specifically, the visual network, the sensorimotor network, and the right superior frontal gyrus were found to be predictive of the short-term drug-treatment response regarding negative symptoms in SZ patients. Meanwhile, the abnormalities within the visual network and the default mode network could forecast the short-term drug-treatment response related to positive symptoms in SZ. These biomarkers have application value in the accurate diagnosis and treatment of SZ after clinical validation.