The Gene Expression Omnibus (GEO) database (18) is a sub-database of the National Center for Biotechnology Information (NCBI) and is used to store gene expression datasets. We downloaded six datasets from the GEO database, GSE107499 (Controls = 44, UC = 75), GSE11223 (Controls = 73, UC = 129), GSE38713 (Controls = 13, UC = 30), GSE53306 (Controls = 12, UC = 28), GSE94648 (Controls = 22, UC = 25), and GSE87466 (Controls = 21, UC = 87) (19–23). The GSE107499 dataset was used for the analysis of this study. Five data sets (GSE11223, GSE38713, GSE53306, GSE94648, and GSE87466) were used for independent validation. The raw gene expression data from the six GEO datasets were analyzed and standardized using the robust multiarray average method.

We used the GEO database’s GEO2R tool to screen and visualize differentially expressed genes using | log FC | ≥ 1.0 and adjusted p value <0.05 as screening criteria. The GEO2R function is implemented based on the “limma” R package. The “limma” R package is a generalized linear model-based differential expression screening method that can obtain DEGs between different comparison groups and controls (24). Specifically, we obtain the gene expression profile data set, remove the genes with an expression value greater than 50%, then use the “voom” function to transform the data, further using the “ImFit” function to perform multiple The data were then transformed using the “voom” function, and further multiple regression using the “ImFit” function was performed to further compute moderated t-statistics, moderated F-statistics, and log-odds of differential expression by empirical Bayes moderation of the standard errors toward a common value to finally obtain the significance of differences for each gene. The Metascape database (25) is a biological database that allows enrichment analysis online. The Metascape database was then used for Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analysis and visualization of the differentially expressed genes in UC. Then, the correlation analysis and identification of DEGs in the CRGs of the GSE107499 dataset were performed using the “corrplot” R package based on Spearman’s statistical method and the “limma” R package.

CIBERSORT is an analytical tool for estimating gene expression profiles and using gene expression data to make relative estimates of the abundance of cell types in mixed cell populations. CIBERSORT is based on linear support vector regression (26). We estimated the relative abundance of 22 immune cell types in each sample of GSE107499 gene expression data using the CIBERSORT algorithm and the LM 22 feature matrix. For each sample, CIBERSORT calculates an inverse fold product p value using Monte Carlo sampling. Each sample’s 22 immune cell proportions added up to 1 in total. Only samples with a p-value of 0.05 or lower were regarded as precise immune cell fractions. The relationship between CRGs and immune cells linked with UC was next examined. Initially, using spearman statistical methods, correlation coefficients between the expression of CRGs and the relative fraction of immune cells were looked at. A p value less than 0.05 was then identified as a significant association using the Spearman correlation coefficient. At last, the “corrplot” R tool (version 0.92) was used to show the findings.

Based on data from associated copper death gene expression profiles, unsupervised cluster analysis (27) was carried out. Using a k-means algorithm with 1,000 iterations, the 75 UC samples were divided into various clusters. We determined the optimal number of clusters based on a combination of cumulative distribution function (CDF) curves, consistency matrices, and consistency clustering scores >0.8 and selected the maximum number of subtypes (k = 9) for analysis. Principal component analysis (PCA) (28), one of the most widely used algorithms for dimensionality reduction of data, was subsequently performed on the two clusters after clustering. Finally, the clustered groups were analyzed for differences and correlations in the CRGs.

With the aid of the analytical technique WGCNA (29), it is possible to analyze the gene expression profiles of several samples, classify genes with similar expression patterns, and look into the associations between certain traits or phenotypes and modules. Using the optimal soft threshold, a weighted proximity matrix was created, following which a topological overlap matrix was created. A unique color is assigned to every module. The module signature genes are a representation of the overall gene expression profile in each module. The importance of modules demonstrates the association between modules and illness states. GSVA in order to clarify the variations in the collection of enriched genes between various CRGs clusters. By contrasting the GSVA scores between the various CRGs clusters, many expression pathways and biological processes were discovered. Subsequently, the investigation of immune cell infiltration was then repeated for several CRGs clusters.

In recent years, more and more machine learning and deep learning methods have been widely used in the medical field with outstanding results. XGB, SVM, RF, GLM, and adaptive boosting (AdaBoost) are some of the most commonly used methods for machine learning, while artificial neural networks (ANN), multilayer perceptron (MLP), and fully neural network (FNN), often referred to simply as neural networks, are the foundation of deep learning. Machine learning favors the interpretability of the model, while deep learning is more concerned with the accuracy of the model. Machine learning is more applicable to tabular data with a relatively small number of variables, while deep learning methods are specifically designed for large data and large feature sets and are more applicable to images or other data with a large number of variables (30, 31). Both AdaBoost and XGB are built based on the boosting algorithm. AdaBoost locates the deficiencies of the model by boosting the weights of the error points, while XGB locates the deficiencies of the model by counting the gradient. Therefore, compared with AdaBoost, the XGB model can use more kinds of objective functions (32). In addition, the dataset selected for this study is from the GEO database, which has a small sample size. Therefore, four machine learning methods—XGB, SVM, RF, and GLM—are used for this experimental study. The four machine learning models XGB, SVM, RF, and GLM were built using the “caret” R package (version 6.0.91) to find differentially expressed genes shared by the UC and turquoise modules. XGB is a supervised model that enables thorough comparisons of classification error and model complexity. It is built on a set of gradient-enhanced augmented trees (33). A binary classification model called SVM transfers the feature vector of an instance to a set of spatial points. The optimum separating hyperplane that maximizes the positive and negative sample intervals on the training set will be found via SVM, as will the best partitioning hyperplane that can distinguish between two classes and maximize the interval (34). RF is an integrated machine learning technique that predicts classification or regression utilizing different independent decision trees by randomly selecting samples and attributes during node splitting (35). GLM, a development of the linear model, establishes the mathematical expectation of the response variable by linking a linear combination of predictor variables (36). The four machine learning models discussed above were interpreted using the “DALEX” package (version 2.4.0), which was also used to show the residual distribution and feature significance. The “pROC” R package (version 1.18.0) was used to visualize the ROC curve. In conclusion, the primary predictive genes associated with UC were found to be the top five significant variables, and the best machine learning models were found.

To determine the prevalence of UC, a nomogram model was developed. Each predictor is given a score, and the “total score” of all predictors is calculated by adding the individual scores of the predictors mentioned above. Using calibration curves and DCA, the Nomogram model’s prediction accuracy was calculated. In addition, five datasets (GSE11223, GSE38713, GSE53306, GSE94648, and GSE87466) were used for independent validation analyses. The ROC curves for these datasets were constructed and visualized using SPSS 25.0 for the prediction model, which was used to validate the accuracy of the prediction model in differentiating between non-UC and UC patients.

Figure 2 displays the study’s flow chart. First, differential expression analysis on the GSE107499 dataset was conducted to identify 849 DEGs, comprising 309 up-regulated genes and 540 down-regulated genes (Figure 3A). The top 30 up-regulated and down-regulated genes were shown separately (Figure 3B). DEGs were also analyzed for GO and KEGG enrichment, and the findings demonstrated that the genes with differential expression were primarily engaged in immune-related pathways, such as leukocyte migration, leukocyte chemotaxis, cell chemotaxis in GO, and the cytokine-cytokine receptor interaction pathway in KEGG (Figure 3C). To elucidate the role of CRGs in the development of UC, the expression profiles of 19 CRGs were systematically evaluated using the GSE107499 dataset (Figure 4A). A total of 13 CRGs associated with UC were identified. In UC samples, the expression of NLRP3 and CDKN2A was significantly upregulated compared to normal samples. Conversely, the expression of NFE2L2, ATP7B, FDX1, LIAS, DLD, DLAT, PDHA1, PDHB, DBT, GCSH, and DLST was significantly decreased in UC samples. Differential analysis (Figure 4B) and correlation analysis (Figure 4C) were performed, followed by the visualization of chromosomal positions (Figure 4D).

The CIBERSORT algorithm was used to elucidate whether immune system differences existed between the UC and non-UC groups and to visualize differences in the proportions of the 22 infiltrating immune cell types, showing significant differences in B cells memory, T cells CD4 memory activated and Mast cells activated in UC patients (Figures 5A,B). The correlation analysis of the 13 CRGs with 22 immune cells showed that ATP7B showed a strong positive correlation with Eosinophils and Plasma cells, and a strong negative correlation with B cells memory and Macrophages M1. Eosinophils, T cells CD4 memory resting, and plasma cells all had a substantial positive connection with FDX1, as did T cells CD4 memory activated and T cells follicular helper (Figure 5C). These findings could imply that genes associated to cuproptosis control the infiltration of these immune cells, which in turn controls the development of UC.

Seventy five UC samples were clustered using a trustworthy clustering technique based on the expression profiles of 19 CRGs in order to determine the CRGs expression patterns in UC samples. As soon as the k value was set to 2, the number of clusters was at its most stable (Figure 6A), and within a minimum range of 0.2–0.6 for the consensus indices, the CDF curve changed (Figure 6B). The difference between the two CDF curves is shown by the area under the CDF curve when k = 2–9 (Figure 6C). Additionally, only when k = 2 was the concordance score for each subtype larger than 0.85 (Figure 6D). In summary, we finally divided the 75 UC samples into two clusters, including cluster 1 (n = 40) and cluster 2 (n = 35). The results of a subsequent PCA analysis showed significant differences between these two clusters (Figure 6E). In addition, analysis of differences between the CRGs and the two clusters after clustering showed that there remained 11 CRGs (PDHA1, NLRP3, LIAS, ATP7B, FDX1, DLD, DLAT, NFE2L2, CDKN2A, DLST, and DBT) that were significantly different between the two clusters (Figures 7A,B).

For cluster 1 and cluster 2, we kept track of immune cell infiltration using the CIBERSORT algorithm. The findings demonstrated that plasma cells, CD4 memory resting T cells, and activated NK cells remained substantially different between the two clusters (Figures 8A,B). A GSVA analysis was then performed on both clusters. Functional enrichment results showed enhanced protein lysine 6 oxidase activity and negative regulation of extracellular matrix disassembly in cluster 1, while positive regulation of microtubule nucleation and carboxylic ester hydrolase activity were enhanced in cluster 2 (Figure 8C). In addition, pathway enrichment results showed that Glycosphingolipid biosynthesis, Glycosaminoglycan biosynthesis – chondroitin sulfate / dermatan sulfate and complement and coagulation cascades were enhanced in cluster 1, while Biosynthesis of unsaturated fatty acids, Pentose phosphate pathway and Arginine and proline metabolism were enhanced in cluster 2 (Figure 8D).

To locate key UC-associated gene modules, we built co-expression networks and modules for clusters 1 and 2 using the WGCNA algorithm. Co-expressed gene modules were discovered when the scale-free R2 was equal to 0.9 and the soft threshold was set at 16 (Figure 9A). The dynamic cutting method produced eight different colored co-expression modules, and a heat map of the topological overlap matrix is also displayed (Figures 9B–D). These eight color modules’ genes were then sequentially used to examine the similarity and proximity of the co-expression of the module’s clinical characteristics (cluster 1 and cluster 2). The turquoise module, which contains 283 genes, exhibited the highest correlation with cluster 2 (Figure 9E). Furthermore, significant relationships between the turquoise module genes and the chosen module genes were shown by correlation analysis (Figure 9F).

We identified a total of 65 intersecting genes for DEGs of UC and Cluster 2 (Figure 10A). To further identify these intersecting genes as having high diagnostic value, based on the 65 cluster-specific DEGs’ expression characteristics in the UC training cohort, we developed four well-established machine learning models: XGB, SVM, RF, and GLM. The “DALEX” package (version 2.4.0) was used to examine the four models and to show the residual distributions for each model in the test set. The machine learning models for XGB and RF showed very little residuals (Figures 10B,C). Following that, based on root mean square error, the top 15 significant feature variables for each model were ranked (Figure 10D). In addition, on the basis of 5-fold cross-validation, by developing ROC curves, in the test set, we evaluated the four machine learning algorithms’ discriminative performance. The AUC of XGB was 0.981, the AUC of RF was 0.967, the AUC of SVM was 0.965, and the AUC of GLM was 0.921 (Figure 10E). In summary, the XGB model had the highest accuracy. The top 5 of these genes were PLXDC1, WAS, CTSK, PLCE1, and LIMD2.

We constructed a column line plot to evaluate the predictive accuracy of the XGB model (Figure 11A), and the scores for each gene were summed to obtain a total score, which corresponded to the disease risk value. To assess the precision of the column plots, calibration curves and DCA curves were created, which showed high accuracy (Figures 11B,C). The AUC of the five independent validation datasets were GSE11223 (AUC = 0.987), GSE38713 (AUC = 0.815), GSE53306 (AUC = 0.946), GSE94648 (AUC = 0.809), and GSE87466 (AUC = 0.981), suggesting high accuracy (Figures 11D–H).