The intersection of differential genes between UC and normal samples analyzed from the GSE38713, GSE87473, and GSE92415 datasets yielded 1,678 genes that were downregulated in UC and 1,605 genes upregulated in UC. Metascape is an analytical website integrated with functional analyses of GO and KEGG (12). We set the screening conditions as Min overlay = 3 and Min enrichment = 1.5 and used Metascape to predict the enrichment of differential gene intersection. Statistical significance was set at p-value <0.01 (13). We used the R package clusterProfiler to predict the different degrees of enrichment at various levels of expressed genes in the data.

The GEO database can be used as a common functional genome warehouse for high-throughput gene expression data, chips, and microarrays to qualitatively study the expression of disease-related genes (https://www.ncbi.nlm.nih.gov/geo/). We obtained the GSE38713, GSE16879, GSE87473, and GSE92415 datasets from the GEO database (14–16). The GSE38713 dataset, containing colon biopsy samples from 13 non-inflammatory controls and 30 UC patients, was conducted using Affymetrix Human Genome U133 Plus 2.0 Array. The basic analysis is contained in the manuscript published in the Gut journal with PMID 23135761. The GSE87473 (PMID: 29401083), containing colon mucosal biopsies from 21 healthy participants and 106 UC patients, was conducted using Affymetrix HT HG-U133+ PM Array Plate. The GSE92415 dataset, containing colonic mucosal biopsy samples from 21 healthy participants and 53 UC patients, was conducted using Affymetrix HT HG-U133+ PM Array Plate. The basic analysis is contained in the manuscript published in Gastroenterology with PMID 23735746. The GSE16879 dataset contains 48 UC samples and 12 normal samples. Total RNA was isolated from intestinal mucosal biopsies, labeled, and hybridized to Affymetrix Human Genome U133 Plus 2.0 Arrays. The basic analysis is contained in the manuscript published in PLoS One with PMID 19956723. The GSE16879 is a mucosal expression profile microarray in patients with inflammatory bowel disease before and after the first infliximab treatment (Table S3) (17). Our annotation information for the reference platform data converted the probe IDs into the corresponding gene symbols. We obtained the raw “CEL” files for the two datasets, called the AffyPLM package in R, and normalized them using the RMA algorithm (18). On differential expression analysis, p-values <0.05 were considered statistically significant.

STRING (version 11.5) was used to build a protein–protein interaction (PPI) network for the differential genes (19). Cytoscape (Version 3.9.1) was used for visual analysis of the interaction network between proteins (20). Molecular Complex Detection (MCODE) (Version 2.0) uses the inherent relationship between proteins in the network to find gene clusters and build a foundation for subsequent analysis (21). We used Cytoscape to read the PPI network files and then used the MCODE plugin to identify the key modules in the PPI network.

The ConsensusClusterPlus package in R was used to classify patients with quality differences in cuproptosis-related genes. The k value was selected and determined by defining the inflection point of the sum of errors.

Differential genes of cluster patterns were screened using the limma package. The expression profiles of the cluster patterns were obtained. The parameters were set at p <0.05.

The degree of each sample was obtained using the single-sample gene set enrichment analysis (ssGSEA) algorithm, and after downloading the gmt format gene-set data required for the analysis, the samples in GSE87473 were immunoscored using the R package GSVA, containing a total of 23 immune-related scores. The differential expression of immune checkpoint genes (PD-L1, IDO1, HAVCR2, PDCD1LG2, CD86, ICOS, TNFRSF9, and CTLA4) and their response to anti-TNF therapy were analyzed (22–24).

While the feeding and drinking habits of each mouse were documented, a daily and routine inspection of the appearance of each mouse (vitality and hair) was also performed. Body weight, stool consistency, and degree of intestinal bleeding were recorded, and the disease activity index (25) was calculated (Table S14). The colon was promptly removed when mice were sacrificed using the cervical dislocation technique on the 7^th day, and colon lengths were obtained. A portion of the colon samples was preserved at −80°C for further studies, while the remainder of the colon tissue was sliced (1 cm) and fixed in paraformaldehyde (4%) for histological testing (26).

The tissue sample was embedded in paraffin after fixing with paraformaldehyde (4%). The tissue was cut to a thickness of 5 mm and dyed with fuel. The slices were then examined under a light microscope for blind analysis. Table S15 displays the histological scores (26).

Freshly frozen tissues were maintained at −80°C. First, the RNA associated with the colon tissue proteins was isolated using TRIzol reagent and a Thermo reverse transcription kit. Second, using an ABI7500 fluorescence quantitative PCR instrument and PerfectStart^® Green qPCR SuperMix amplification reagent, PCR amplification and detection were performed with actin as the internal control. The expression levels were examined using the 2^−ΔΔCt approach. The primer sequences were as follows: DLAT-F: GAGGTGCTGTTGGTACGGAA, DLAT-R: ACGAGTTTGCTTCGGGAACT; DLD-F: AGTCGTGTGTACCGCTCCTT, DLD-R: CACTGTCACGTCAGCCTCAA; and PDHA1-F: GGATGGAGCTAAAGGCGGAT, PDHA1-R: TCCGTAGGGTTTAT GCCAGC.

Cells were lysed using RIPA buffer (Tris 20 mM, NaCl 150 mM, KCl 20 mM, MgCl₂ 1.5 mM, glycerol 10%, Triton X-100 1%, pH 7.5), and the protein concentration was measured by the BCA Protein Assay Kit (C503021; Sangon Biotech, Shanghai, China). Equal amounts of proteins were separated by 10%–12% sodium dodecyl sulfate-polyacrylamide gel electrophoresis and were transferred onto polyvinylidene fluoride membranes, and the target proteins were finally detected using standard Western blotting protocols and visualized using the Super Signal West Pico Plus Luminol/Enhancer Solution (UC279012; Thermo, Waltham, USA). β-Actin was used as the loading control. The primary antibodies used were listed as follows unless otherwise specified: anti-β-actin antibody (AC038; ABclonal, Wuhan, China), anti-DLD antibody (A5220; ABclonal), anti-PDHA1 antibody (A13687; ABclonal), anti-DLAT antibody (A8814; ABclonal). The secondary antibodies used were as follows unless otherwise specified: horseradish peroxidase (HRP)-conjugated goat anti-rabbit IgG (H + L) (AS014; ABclonal).

Most analyses used the R software, partly using the GraphPad Prism 8 software and the Sangerbox (27) website. In all statistical analyses, p <0.05 was significantly different. The p-value of the differential gene selection has been corrected (https://www.biocon.5octor.org/, https://www.graphpad.com/, http://sangerbox\2/login.html).

We performed a differential analysis of the expression between the normal and UC tissue samples in three datasets (GSE38713, GSE87473, and GSE92415) using the limma package (Figures 1A–C and Table S1). The upregulated and downregulated differentially expressed genes intersected in GSE38713, GSE87473, and GSE92415 (Table S2). We finally obtained 1,605 upregulated genes and 1,678 downregulated genes (Figures 1D, F). We used Metascape to predict the enrichment of the intersection of the differential genes (Figures 1E, G). Most of the upregulated genes were enriched in cytokine signaling in the immune system and positive regulation of cell motility, and the downregulated genes were enriched in the metabolism of lipids, generation of precursor metabolites, and energy and lipid biosynthetic processes.

The most critical modules in the PPI network were recognized using the MCODE plugin in Cytoscape. This article indicates the top six subnetworks in the networks constructed by the downregulated genes. The weight of these networks from left to right is 33.314, 14.067, 10.471, 8.227, 6.909, 6.727, 6.250, 6.000, 5.000, and 4.889 (Figure 1H and Figure S1). We then identified and marked three cuproptosis-associated genes in the top two ranked networks and analyzed them.

We analyzed the expression levels of the three cuproptosis-associated genes mentioned above in GSE38713, GSE92415, and GSE87473 (Figures 2A-C) based on the cuproptosis-associated genes studied by Tsvetkov et al. (26). To determine whether these cuproptosis-related genes affected disease progression in UC patients, we chose GSE87473 with a large sample size to assess the expression levels of cuproptosis-associated genes in UC and normal samples. The results indicate that seven genes (DLAT, CDKN2A, PDHA1, PDHB, DLD, GLS, and FDX1) were differentially expressed. We found high expression levels of GLS and CDKN2A in UC tissues (p < 0.01) and high expression levels of the other five genes in the normal tissues (p < 0.001) (Figures 2D, E). The results revealed that cuproptosis-associated genes indicated highly heterogeneous expression in normal and UC patients, suggesting that the expression of different cuproptosis-associated genes can affect disease progression in UC patients. Among them, 10 cuproptosis-associated gene expression profiles are presented in the Supplementary Material (Table S3). By analyzing the expression data from multiple ulcerative colitis, we found three cuproptosis-associated genes in the core subnetwork of the PPI network (Figures S2A–F). We finally selected the three hub genes DLD, PDHA1, and DLAT for the follow-up studies.

Using the ConsensusClusterPlus package in R, the cuproptosis pattern was classified according to the expression levels of three cuproptosis-associated genes and identified as k = 3 (Figures 2F–H). Finally, we identified three different modification patterns through unsupervised clustering: 21 cases of CuClusterA, 34 cases of CuClusterB, and 51 cases of CuClusterC (Table S4). Furthermore, PCA indicated that the UC samples could be completely distinguished (Figure 2I). After unsupervised clustering, we studied cuproptosis-related genes in different CuClusters (Figures 2J, K). The strongest expression of cuproptosis-associated genes was found in CuClusterB (p < 0.001).

Next, we compared the differences in immune cell infiltration among the three clusters using ssGSEA analysis. The results indicated that multiple immune cells, except dendritic cells, were significantly enriched in the three CuClusters (Figure 3A and Table S5). Correlations between individual CuClusters and immune infiltration were analyzed using the online website Sangerbox (Figure S3). We then explored the association between these three cuproptosis-associated genes and immune infiltration (Figure 3B). DLAT expression positively correlated with CD56 dim natural killer cells (Figure 3C; p < 0.001). DLD expression positively correlated with activated CD8T cells, activated CD4T cells, eosinophils, γ-δ T cells, type 2 T helper cells, CD56 bright natural killer cells, and type 17 T helper cells. However, the correlation was only significant for γ-δT cells (p < 0.01) and type 2 T helper cells (Figure 3D; p < 0.05). PDHA1 expression positively correlated with CD56 dim natural killer cells, CD56 bright natural killer cells, and type 17 T helper cells and significantly correlated with type 17 T helper cells (Figure 3E; p < 0.05). Moreover, we found that the three CuClusters, except IL15, were significantly differentially expressed in inflammatory factors (Figures 3F, G). Overall, cuproptosis-associated genes can affect the mucosal immune microenvironment of UC patients.

The number of differentially expressed genes in the three CuClusters was analyzed using the limma package, with a threshold of |logFC| >1 and p <0.05. CuClusterA–CuClusterB screened 2,616 differentially expressed genes (Table S6). CuClusterA–CuClusterC screened 2,658 differentially expressed genes (Table S7). CuClusterB–CuClusterC screened 122 differentially expressed genes (Table S8). Subsequently, the differential genes between the three CuClusters were intersected, and 10 differential genes were identified (Figure 4A). The pathway and enrichment information of the 10 genes were analyzed using the clusterProfiler R package, and two of the genes were not enriched for annotation (Figure 4B and Table S9). These genes were significantly enriched in extracellular matrix organization, extracellular structure organization, and external encapsulating structure organization (p < 0.05). In addition, these genes were also enriched in the collagen-containing extracellular matrix, platelet alpha granules, and basement membrane (p < 0.001). Among the cellular components, these genes were significantly enriched in heparin binding, glycosaminoglycan binding, and sulfur compound binding (p < 0.001). Differential genes were uploaded to the Sangerbox website for enrichment analysis (Figure 4C). Furthermore, enrichment in the KEGG pathway was observed for malaria, ECM receptor interaction, Staphylococcus aureus infection, and protein digestion and absorption (Figure 4D and Table S10; p < 0.05). After analyzing the expression levels of these 10 genes in GSE87473, we found that the expression of COLA5A3, CHN1, H19, THBS2, SELP, KIAA0125, and EGFL6 was significantly upregulated in UC patients (Figure 4E; p < 0.05). The expression levels of SLC9A2, FAM55A, and ANXA13 were significantly downregulated (Figure 4F; p < 0.05).

We classified UC patients into different genomic subtypes according to the differential expression levels of 10 genes associated with the cuproptosis phenotype. GeneClusters were well-distinguished (Figures 4G–J and Table S11). We found that these 10 genes were significantly different in GeneClusters (p < 0.05), with high expression of SLC9A2 and FAM55A in GeneCluster I and H19, EGFL6, CHN1, THBS2, KIAA0125, COL5A3, and SELP in GeneCluster II. The expression of ANXA13 was high in GeneCluster III (Figures 4K, L). Analysis of the results indicated that these 10 genes had the highest expression levels in GeneCluster II. On analyzing the relationship between CuClusters and GeneClusters using Sangerbox, we found that nearly all CuClusterA patients were from GeneCluster I and nearly all GeneCluster II patients were from CuClusterC (Figure 4M). These 10 genes have a strong association with cuproptosis-associated genes, suggesting that cuproptosis has a potential role in UC treatment.

To explore the impact of cuproptosis phenotype-associated hub genes in UC immune infiltration, we evaluated the immune infiltration relationship between GeneCluster I, GeneCluster II, and GeneCluster III. The results indicated that CD56 dim natural killer cells were more abundant in GeneCluster I. Type 17 helper T cells were more abundant in GeneCluster III. Other immune cells were more abundant in GeneCluster II (Figure 5A, p < 0.05). Correlations between individual GeneClusters and immune infiltration were analyzed using Sangerbox (Figure S4). Notably, the expression of eight important immune checkpoints significantly differed between the CuClusters and GeneClusters (Figures 5B, C). We found that CD86, ICOS, TNFRSF9, CTLA4, and LAG3 were expressed at their highest levels in CuClusterC (p < 0.05), and PDCD1LG2 was the most highly expressed gene in CuClusterA (p < 0.05). We also found that, except for PDCD1LG2, the other seven immune monitoring sites exhibited the highest expression in GeneCluster II (p < 0.05), suggesting that the hub genes associated with the cuproptosis phenotype can affect immunotherapy.

The current treatment for patients with refractory and severe UC is infliximab, usually in combination with steroids, which was approved after clinical trials on these patients. Studies have indicated that infliximab treatment improves clinical presentation and induces mucosal healing in UC patients (29–32). Therefore, we further validated the expression levels of the 10 cuproptosis phenotype-associated hub genes in GSE16879. However, the expression results for only eight hub genes were available in GSE16879. The expression levels of these eight hub genes were analyzed in responders before the first infliximab treatment and after the first infliximab treatment (Figures 5D, E and Table S12). Furthermore, we analyzed the expression of a single gene in responders before and after the first infliximab treatment (Figure 5F and Figure S5). Correlations between eight hub genes were analyzed in responders before and after the first infliximab treatment (Figure 5G). The results indicated that the correlation between KIAA0125 and SLC9A2 and the correlation between ANXA13 and SELP were both significantly enhanced, while the correlation between ANXA13 and KIAA0125 and the correlation between THBS2 and COLA53 were both significantly reduced after the first infliximab treatment compared with that before the first infliximab treatment. This suggests that the altered correlation between the genes is due to drug sensitivity. Based on the different expression levels of the eight hub genes, we performed a cluster analysis using the ConsensusClusterPlus R package, k = 2 (Figures 5H–J and Table S13). Notably, PCA indicated that the UC samples could be completely distinguished (Figure 5K). Furthermore, infliximab treatment was significantly different between the two gene clusters (Figures 5L, M). The hub genes (COLA53, CHN1, THBS2, SELP, KIAA0125, EGFL6) were significantly highly expressed in gene cluster α (p < 0.001), whereas SLC9A2 was more expressed in cluster β than in gene cluster α (p < 0.05). ANXA13 was not significantly differentially expressed in clusters α and β. Cluster α responded better to infliximab compared with cluster β (Figure 5N), suggesting that cuproptosis-associated hub gene expression can have implications in infliximab treatment.

In summary, the results suggested that high expression levels of some cuproptosis phenotype-associated hub genes in UC patients are sensitive to infliximab treatment. Analysis of the expression levels of cuproptosis phenotype-associated hub genes can facilitate the selection of clinical dosing.

The body weight of the mice in each group indicated different changes. We observed that the dextran sulfate sodium (DSS) group indicated a significant and continuous decrease on the 5th day of modeling (vs. control, p < 0.001) (Figure 6A). Notably, colon length is an essential index of the degree of colonic inflammation. Our investigations revealed that colon length in the DSS group was substantially shorter than that in the control group (Figure 6B). The mice were sacrificed on the 7th day. The mice in the treatment group weighed less than those in the control group (Figure 6C). We found that the colon/body weight was higher in DSS-treated mice than in the control group mice (Figure 6D). The DSS-treated mice indicated a reduced colon length (Figure 6E). The colon weight/length ratio increased after DSS administration (Figure 6F). The weight of the spleen in the DSS group was significantly higher than that in the control group (Figure 6G). An increase in the spleen/body weight ratio was also observed after the acute administration of DSS (Figure 6H). Hematoxylin and eosin-stained slides indicated that the epithelial cells of the mouse colon mucosa were structured and intact, and goblet cells were abundant, without atrophy, deformation, necrosis, or infiltration of inflammatory cells in the control group. In contrast, typical inflammatory symptoms in the acute phase, disappearance of normal pathological structures, destruction of crypt structures, obvious reduction or even disappearance of goblet cells, and large infiltration of inflammatory cells in the mucosa and submucosa were observed in the DSS group (Figure 7A). Furthermore, DAI and histological inflammation scores increased after acute administration of DSS (Figures 7B, C). We detected DLAT, DLD, and PDHA1 gene expression levels by using real-time quantitative polymerase chain reaction (qRT-PCR) and Western blot analysis (Figures 7D-H). Consistent with the results of bioinformatic analysis, DLAT, DLD, and PDHA1 were lower than those in the control group. In summary, these three genes play important roles in the development of UC.