RNA-seq expression data from a total of 548 patients with UCEC and corresponding clinical information were obtained from the TCGA database (https://portal.gdc.cancer.gov/projects/TCGA-UCEC). LncRNAs were annotated by human gene annotation files (GRCh38.p12), which were downloaded from the Ensembl database (https://asia.ensembl.org/index.html). We excluded UCEC patients from TCGA dataset according to the following criteria: 1) patients with incomplete clinical information; 2) patients with overall survival (OS) time or follow-up time less than 30 days. Finally, 524 UCEC patients were selected in our study.

Among the tissues of these patients, 23 paired UCEC tissues and adjacent non-cancerous tissues were selected to screen out DEGs by using the edgeR package (version 3.32.1) (11). 1) Filter the data to remove genes with low counts. Genes that have count-per-million (CPM) values above one in at least two samples were kept; 2) Normalize the data. Trimmed Mean of M-values (TMM) normalization method was used for the data analysis; 3) Explore the data and estimate the dispersion; 4) Investigate DEGs. The false discovery rate (FDR) was applied for multiple testing correction of raw P values through the Benjamini-Hochberg method. A |log2 fold change (FC)| > 2 and an FDR < 0.05 were set as the threshold for identifying DEGs. Besides, the samples were clustered based on gene expression data by performing principal component analysis (PCA). Differentially expressed mRNAs (DEmRNAs) and lncRNAs (DElncRNAs) were used for further analysis.

WGCNA was used to construct the co-expression network of DElncRNAs and DEmRNAs in UCEC according to the following main steps: 1) Select the weighting coefficient, β; 2) Transform the gene expression profiles into an adjacency matrix; 3) Use the adjacency matrix to define a separate measure of similarity, the Topological Overlap Matrix(TOM); 4) Perform hierarchical clustering for TOM-based dissimilarity (dissTOM) to obtain the hierarchical clustering tree; 5) Identify the modules via dynamic branch cutting methods; 6) Calculate the module eigengene (ME) of each module; 7) Identify the modules that were most strongly related to the clinical traits (12, 13). Next, gene ontology (GO) annotation and Kyoto encyclopedia of genes and genomes (KEGG) pathway enrichment were performed using the clusterProfiler R package (version 3.18.1) (14). Significant enriched functions and pathways were visualized by ggplot2 R package (version 3.3.5).

Two methods were applied to screen the hub genes of the network. In the first method, the network screening function based on GS (representing the correlation between the gene and a given clinical trait) and MM (representing the correlation between the gene and a given module) was used to directly identify hub genes. A cut-off criteria (|MM| > 0.8 and |GS| > 0.4) was employed to obtain key genes with high connectivity in the clinically significant modules. In the second method, a PPI network was built based on the online STRING database (version 11.5) and analyzed using Cytoscape (version 3.8.0) (15). Common genes that scored in the top 30 by all five methods in CytoHubba were selected as key nodes of UCEC, and MCODE in Cytoscape was used to perform module analysis (16, 17). Overlapping genes identified by both WGCNA and PPI analyses were designated as potential hub genes for further validation and analysis (18, 19).

Based on the median expression value of the hub genes, samples of TCGA UCEC were divided into low- and high-expression phenotypes. Then, GSEA (version 1.52.1) was performed to detect which KEGG pathways were enriched (20). Terms with FDR < 0.05 were visualized by ggplot2 R package to investigate potential functions of the hub genes.

The lncRNA-miRNA-mRNA ceRNA and lncRNA-TF-mRNA networks in the co-expression modules were constructed based on several online databases (Starbase, miRWalk, LncTarD, RNA Interactome and LncMAP) and visualized in Cytoscape (version 3.8.0) (21–25).

To identify prognostic hub genes by integrating the clinical information of patients suffering from UCEC in TCGA, Kaplan-Meier curve analysis of the samples with hub genes was conducted with the survival R package (version 3.2-10), and significant differences in survival were determined with the log-rank test. P values <0.05 were regarded as significant.

The human UCEC cell lines Ishikawa and KLE were purchased from Procell Life Science &Technology Co., Ltd. and cultured in Dulbecco’s modified Eagle’s medium (DMEM) and DMEM/F12 medium respectively, supplemented with 10% fetal bovine serum and 1% penicillin/streptomycin and maintained with 5% CO₂ at 37°C in a humidified incubator.

Total cellular RNAs were extracted from Ishikawa and KLE by using TRIzol™ Reagent (Thermo Fisher Scientific Inc., United States) according to the manufacturer’s protocol. cDNA was synthesized using MMLV Reverse Transcriptase (Anhui Toneker Biotechnology Co., Ltd., China). The qRT-PCR analysis on mRNA and lncRNA was performed using 2 × PCR SYBR Green Mix buffer (Anhui Toneker Biotechnology Co., Ltd., China). The PCR process ran 40 cycles for 15s at 95°C and for 1 min at 60°C in Applied Biosystems 7500 Real-Time PCR System (Thermo Fisher Scientific Inc., United States). The 2[-Delta Delta C(T)] was adopted to calculate relative expression with GAPDH as an internal control (26). The primers used were listed in Table 1 .

Three oligonucleotides specific for short hairpin RNAs (shRNAs) against DUXAP8 were inserted into the lentiviral expression vector (pLVX-shRNA1) respectively. Lentiviral vectors were then cotransfected with packaging plasmids psPAX2 and pMD2G into HEK-293T cells. The viral supernatant was harvested 48 and 72 h after transfection, filtered through a 0.45 μm millipore filter, aliquoted and stored at -80°C for subsequent use. Stable cell lines were obtained by selection with puromycin. The interference sequences targeting DUXAP8 were listed in Table 2 .

Cells were seeded in 96-well plates at a density of 4×10³ cells/well in triplicate and incubated at 37°C. Cell viability was determined by a Cell Counting Kit-8 (CCK-8) kit (Dojindo Laboratories, Japan) according to the manufacturer’s instructions. The absorbance at 450 nm was measured after incubation at 37°C for 2 h with microtiter plate reader (BioTek, United States), and cell viability was calculated.

Alisertib (MLN8237) was purchased from Selleck, dissolved in dimethyl sulfoxide (DMSO) and kept frozen at -20°C. Cells were dispensed into 96-well plates at a density of 4×10³ cells in 100 μL of complete medium with different concentrations of Alisertib (0.01, 0.1, 1, 10, 50, 100 μM) for 72 h. The final concentration of DMSO in medium was kept at 0.2% for all groups including control samples (no Alisertib, DMSO only). CCK-8 assays were then performed to detect cell viability.

Cells were seeded into 6-well plates at a density of 5×10² cells, followed by incubation at 37 °C for 1-3 weeks. Then, colonies were fixed with glutaraldehyde (6.0% v/v) and stained with crystal violet (0.5% w/v). The colony was defined to consist of at least 50 cells and colonies were quantified by using ImageJ image analysis software (27).

Cells were seeded at a density of 3×10⁵ cells/mL and treated with the Alisertib (1μM) for 48 h. For cell cycle distribution analysis, cells were harvested and fixed with 70% cold ethanol at -20°C overnight. The fixed cells were then resuspended and stained with propidium iodide (PI)/RNase staining buffer (BD Biosciences, United States) according to the manufacturer’s instructions. The cells were analyzed by flow cytometry (BD Biosciences, United States) and the proportions of cells at each phase of the cell cycle were determined and compared using the ModFit LT software. For cell apoptosis analysis, cells were harvested and stained with Annexin V-FITC/PI apoptosis kit (MULTISCIENCES, China) following the manufacturer’s instructions and the data were analyzed with Flowjo software. Annexin V⁺/PI^– cells were considered in early apoptosis, Annexin V⁺/PI⁺ populations denoted late apoptotic/necrotic cells (28).

Statistical analyses are based on RStudio software (version 3.6.3) and GraphPad Prism (version 8.0). The experimental data are presented as the means ± standard deviation (SD). The independent Student’s t-test was used to compare equivalent variables with normal distribution between two groups. Statistical differences between groups with unequal sample sizes and different variances were assessed using unpaired t-tests with Welch correction (Welch t-tests) (29). Multiple group comparisons were performed with one-way analysis of variance (ANOVA) followed by Tukey’s multiple comparison test. P < 0.05 was considered statistically significant. All data generated or analyzed during the current study are included in the published article and Supplementary Materials .

The schematic overview of this study was shown in Figure 1 . The RNA-seq expression dataset and corresponding clinical information for patients with UCEC were obtained from the TCGA UCEC cohort. Once the raw data was normalized, DEG analysis was performed between 23 paired UCEC tissues and adjacent non-cancerous tissues. In total, 2569 DEmRNA strands (1295 upregulated and 1274 downregulated) and 1457 DElncRNA strands (733 upregulated and 724 downregulated) were identified at the threshold of FDR < 0.05 and |log2FC| > 2 and selected for further analysis. Volcano plots and heatmaps were plotted to show the distribution of DEmRNAs ( Figure 2A ) and DElncRNAs ( Figure 2B ) that were differentially expressed in UCEC tissues in comparison to nontumorous tissues. It was also demonstrated by PCA that the 23 tumors were clustered separately from their paired adjacent non-cancerous tissues ( Figures 2C, D ).

In the next stage, the co-expression network of all DEGs in UCEC, including 2569 DEmRNA strands and 1647 DElncRNA strands, was constructed by using the WGCNA R package. When the weighting coefficient β was 3, the independence degree was > 0.8, indicating that the scale free topology of the network could be appropriately assessed by the power value according to the WGCNA algorithm ( Figures 3A, B ). Therefore, β = 3 was selected as a soft threshold to construct the co-expression modules. The gene expression profiles were transformed into the adjacency matrix, and dissTOM for DEGs was then obtained ( Figure 3C ). 8 modules were detected by the method of dynamic tree cutting and merging similar modules ( Figure 3D ).

The Pearson correlation coefficients for the MEs of all modules and the clinical information were calculated to identify which modules were related to the clinical traits. As shown in Figure 4A , the Blue module (R = 0.50, P = 1e-34) and Brown module (R = 0.49, P = 2e-33) were significantly positively associated with neoplasm histologic grade. Moreover, groups of correlated eigengenes were identified through hierarchical clustering dendrogram and heatmap construction, which indicated that the Blue and Brown modules were most significantly associated with neoplasm histologic grade ( Figure 4B ), suggesting that genes in the two modules might have potential roles in tumor differentiation and progression. Therefore, the two modules mentioned above were selected as the clinically significant module for further analysis. GO term enrichment analysis demonstrated that genes in the two modules were significantly enriched in the following biological processes which are ‘Mitotic nuclear division’, ‘Nuclear division’, ‘Chromosome segregation’, ‘Microtubule cytoskeleton organization involved in mitosis’ and ‘Nuclear chromosome segregation’ ( Figure 4C and Table S1 ). In the KEGG pathway enrichment analysis, the top 10 enriched pathways were ‘Cell cycle’, ‘Oocyte meiosis’, ‘cellular senescence’, ‘Cytokine-cytokine receptor interaction’, ‘Viral protein interaction with cytokine and cytokine receptor’, ‘Chemokine signaling pathway’, ‘Human T-cell leukemia virus 1 infection’, ‘PPAR signaling pathway’, ‘Progesterone-mediated oocyte maturation’ and ‘P53 signaling pathway’ ( Figure 4D and Table S2 ).

The highly connected genes of the Blue and Brown modules were investigated as potential key factors related to tumor differentiation. Based on the cut-off criteria (|MM| > 0.8 and |GS| > 0.4), 58 genes with high connectivity in the clinically significant modules were identified (30, 31) ( Figure 5A ). Additionally, we also constructed a PPI network for all genes in the two modules based on STRING database and Cytoscape software, which was composed of 532 nodes and 9911 edges. 20 common genes that scored in the top 30 by all five methods in CytoHubba were selected as key nodes of UCEC in PPI analysis ( Figure 5B ). MCODE in Cytoscape was used to perform module analysis. We found that all of the 20 common genes were in module 1, which was the fairly significant module (MCODE score = 105.433) in all modules and potentially played an important role in UCEC progression ( Figure 5C and Table S3 ).

As shown in Figure 5D , 10 genes identified by both WGCNA and PPI analyses were designated as potential candidates for further validation and analysis (AURKA, BIRC5, BUB1, CCNA2, CCNB1, CDCA8, DLGAP5, KIF2C, NCAPG and TPX2). Based on the TCGA data, the expression levels of the 10 genes could be used to significantly distinguish well-differentiated UCEC (neoplasm histologic grade G1-G2) from poorly differentiated UCEC (neoplasm histologic grade G3) ( Figure 6A ) (P < 0.0001). The protein expressions of most of these genes were higher in UCEC tissues compared to non-tumor tissues, according to the Human Protein Atlas database (32) ( Figure 6B ), except that no protein expression was detected for KIF2C gene and there was no expression data for BUB1 gene in this database. Moreover, Kaplan-Meier curve analysis showed that 6 of the 10 genes were significantly negatively related to the 5-year overall survival (OS) of patients with UCEC (AURKA, BUB1, CDCA8, DLGAP5, KIF2C and TPX2) ( Figure 7 ). Therefore, these 6 genes were selected as representative genes for further investigation. Then, GSEA was performed to explore the possible pathogenesis of the 6 genes in UCEC. Five common gene sets, ‘Cell cycle’, ‘DNA replication’, ‘Mismatch repair’, ‘Homologous recombination’ and ‘Oocyte meiosis’, were enriched in the sample group with high gene expression (FDR < 0.05; Figure 8 ). Overall, these gene sets were tightly associated with cell proliferation, indicating a vital role of the 6 genes in regulating the cell cycle and proliferation of UCEC.

Besides, to explore the molecular mechanism of UCEC-related lncRNA, lncRNA-miRNA-mRNA and lncRNA-TF-mRNA networks were constructed from several online databases (Starbase, miRWalk, LncTarD, RNA Interactome and LncMAP), according to the DElncRNAs and DEmRNAs from the co-expression Blue and Brown modules. The lncRNA-miRNA-mRNA ceRNA networks consisted of 4 DElncRNA strands, 4 DEmiRNA strands, and 8 DEmRNA strands ( Table S4 ). The lncRNA-TF-mRNA networks consisted of 16 DElncRNA strands, 37 DEmRNA strands, and 20 transcription factors ( Figure 9A and Table S5 ). 3 lncRNA strands in the constructed networks were identified as being significantly associated with decreased 5-year OS (AC015849.16, DUXAP8 and DGCR5), with higher expression levels in UCEC tissues compared to non-tumorous tissues ( Figure 9B ). Moreover, the expressions of AC015849.16, DUXAP8 and DGCR5 were significantly upregulated in poorly differentiated UCEC tissues compared with well-differentiated UCEC tissues ( Figure 9C ) (P < 0.0001).

To confirm the reliability of the hub genes from WGCNA and PPI, lncRNA-miRNA-mRNA ceRNA, and lncRNA-TF-mRNA network, we verified the expression patterns of these genes in two endometrial cancer cell lines, including Ishikawa (histological grade 1; G1) and KLE (histological grade 3; G3) (33). All of these hub genes presented higher expression in KLE than in Ishikawa cells ( Figure 10A ), which were consistent with the results obtained for the TCGA UCEC cohort.

Among these hub genes, AURKA was the most significant gene associated with the survival rate in UCEC patients. Therefore, AURKA was selected as a representative gene for further verification. We evaluated the biological impact of inhibiting AURKA with Alisertib (MLN8237), an investigational, oral, and selective inhibitor of AURKA in UCEC cells in vitro. The cell viability (3 days), in response to Alisertib treatment, was determined by CCK-8 assay. It was discovered that Alisertib decreased cell viability in a dose-dependent manner ( Figure 10B ). To further assess the effects of AURKA on cell cycle distribution, Ishikawa cells were treated with Alisertib (1 μM) for 24h and the DNA content was then measured by flow cytometry. As shown in Figure 10C , incubation of cells with Alisertib resulted in a marked G2/M phase arrest from 14.0% at basal level to 84.6%. Moreover, treatment with Alisertib also increased the total proportion of apoptotic cells (early + late apoptosis). Compared to the control group (0.01% DMSO), there was 2.3-fold increase in total apoptotic cells when Ishikawa cells were incubated with Alisertib (1 μM) ( Figure 10D ). In the long-term clonogenic survival assay, after a single overnight treatment with Alisertib (1 μM), a significant reduction in the number of viable cell colonies in UCEC cells was demonstrated ( Figure 10E ). Collectively, these results indicate that inhibition of AURKA using Alisertib can effectively suppress cancer cell viability, modulate cell cycle distribution and induce apoptosis in UCEC.

We also selected DUXAP8, a strand of pseudogene-derived lncRNA that has been shown to be associated with the occurrence and development of several kinds of tumors, as a representative lncRNA for further verification. Little knowledge is currently available on the role of DUXAP8 in endometrial cancer. Therefore, we assessed the effect of DUXAP8 on UCEC cells proliferation and colony formation in vitro. DUXAP8-shRNA/GFP lenti-virus (versus control virus) was used to knock down DUXAP8 in human UCEC cells ( Figure 10F ). Colony formation assay results showed that clonogenic survival was inhibited following the downregulation of DUXAP8 ( Figure 10G ). In addition, the growth curves detected by CCK-8 showed that DUXAP8 knockdown significantly impaired UCEC cells growth ( Figure 10H ).