Exploring potential therapeutic strategy for hepatocellular carcinoma and COVID-19 using bioinformatics analysis

Background Hepatocellular carcinoma (HCC) constitutes an important contributor to fatalities. Coronavirus disease 2019 (COVID-19) frequently presents with complications such as respiratory distress, systemic inflammatory responses, and damage to various organs. Several studies have investigated the relationship between COVID-19 and mortality in patients with liver cancer, but there are few research on the relationship between them. This study is to explore the correlation between the two diseases and drugs treating them. Methods The Gene Expression Omnibus (GEO) database provides gene datasets of COVID-19 patients and HCC patients. Through differential gene analysis and weighted gene co-expression network analysis, we determined 223 genes represented in HCC and COVID-19. We then used functional annotation, protein–protein interaction network construction, predictive model development and verification, prognostic value analysis, and miRNA–gene network construction. Besides, we created a drug–hub–gene network by predicting possible medications that interact with hub genes using the Drug–Gene Interaction Database (DGIdb). Ultimately, we applied immunohistochemistry to ascertain the hub genes expression. Results This study revealed that eight core genes (RRM2, TPX2, DTL, CDT1, TYMS, CDCA5, CDC25C, and HJURP) co-existed in both HCC and COVID-19 and were differentially expressed in both HCC and normal tissues. CDC25C, RRM2, CDCA5 , and HJURP had diagnostic value (AUC > 0.8) and prognostic value (adjusted P -value < 0.05). Genome enrichment analysis indicated that eight genes may function in liver cancer through engagement in the cell cycle, DNA replication, etc. In liver cancer samples, these genes were significantly and adversely associated with plasma cells while RRM2 was positively associated with neutrophil and NK cell activation and with dendritic cell resting. Using the miRNAnet database and DGIdb, 9 transcription factors, 7 miRNAs, and 51 drugs or molecular compounds were predicted to interact with the hub genes. Finally, RRM2 expression showed significant variation in clinical specimens, and analysis of the association of RRM2 with immunomodulators indicated that RRM2 was closely connected to MICB and CD276. Conclusions Our study revealed several metabolic genes related to HCC and COVID-19. Moreover, potential drugs related to central genes were predicted. These findings may provide new ideas for treating COVID-19 and HCC.


Background
The incidence of the COVID-19 triggered by the SARS-CoV-2 pathogen has escalated and presented a major threat to global public health, with a continually growing number of affected individuals [1].Hepatocellular carcinoma is a frequently diagnosed form of primary cancer.HCC is a fatal consequence of cirrhosis and can be caused by congenital hereditary alterations or acquired hepatopathy such as HBV/HCV infection, alcoholic liver disease, and aflatoxin contamination.The high incidence of HCC attributed to these risk factors.The lack of research on hepatocarcinogenesis and the failure to block or reverse malignant transformation have resulted in unsatisfactory outcomes for HCC patients [2].Several studies have found a correlation between SARS-CoV-2 and tumors.Additionally, the risk of SARS-CoV-2 infection and serious complications is about 2.31 times in pulmonary malignancy patients than in the general population [3].An Italian research represented that pulmonary malignancy patients with COVID-19 infection had a more fatality rate than those uninfected patients [4].An increasing number of researchers have indicated that COVID-19 shares similarities with potential risks of prostate cancer.Chakravarty et al. reported that COVID-19-associated proinflammatory disease in nearby tissues could trigger prostate cancer.According to their research, inflammation may be a driver of prostatic carcinoma malignant transformation [5].The association between COVID-19 infection and prostatic carcinoma is linked to a high level of TMPRSS2 expression [6].Stipp MC et al. reported a common pathogenic profile between breast cancer and COVID-19, represented by the potential for replication of inflammatory mediators and SARS-CoV-2 in metastatic cancer cells [7].
Due to the recent emergence of COVID-19, there are limited studies on the association between liver cancer and COVID-19, and even fewer studies have explored their common molecular mechanisms through advanced bioinformatics approaches.Numerous studies have suggested that systemic immune suppression brought on by chemotherapy or surgery may put cancer patients at higher risk of severe COVID-19 and worsen their prognosis [8,9].Following COVID-19 infection, Mallet et al. observed a noteworthy increase in 30-day mortality in patients with alcohol use disorders, alcohol-induced liver damage, cirrhosis, and HCC [10].Specifically, patients with HCC who are infected with SARS-CoV-2 are more likely than those without cancer to experience complications, ICU admission, and death [11].According to Martínez et al., SARS-CoV-2 infection was present in 48% of patients who succumbed to hepatic malignancies and bile duct malignancies.Additionally, 18.4% of HCC patients die shortly after diagnosis [12].Common pathogenic factors were also identified.VWF, a biomarker used to predict the progression of hepatocellular carcinoma, has been mentioned repeatedly [13,14].It is involved in the formation of microvascular thrombosis during the development of COVID-19 [15].Del Valle et al. found strong expression of IL-1β in COVID-19 patients [16].IL-1β is an inflammatory cytokine responsible for inducing PD-L1 expression, which further contributes to immune escape from HCC [17].TNF has been shown by Nakagawa H et al. to play a role in hepatocellular carcinoma pathogenesis [18].The levels of the proinflammatory cytokine IFN-γ are increased in hepatocellular carcinoma (HCC) patients.IFN-γ further mediates immune escape from HCC by increasing PD-L1 expression [19].Karki et al. found that the synergism of TNF-α and IFN-γ triggers SARS-CoV-2 infection [20].
In order to identify copathogenic genes and medications and investigate their pathogenesis, our goal was to perform molecular analysis on COVID-19 and HCC.Two gene datasets (GSE54236 and GSE177477) were obtained from the GEO database [21].The limma package and weighted gene co-expression network analysis were used to identify prominent genes associated with COVID-19 and hepatocellular carcinoma.In addition, we used the STRING website [22] and Cytoscape software [23] to create a PPI network.We examined the enabling modular genes and identified the hub genes.We also investigated the transcription factors associated with the chosen genes.Finally, we used DGIdb [24] to predict the drug compounds that interact with the pivotal genes.Finally, we conducted an immunohistochemical analysis of the hub genes.Based on these findings, we concluded that RRM2 may act as an immune checkpoint inhibitor in hepatocellular carcinoma.Our study provides new insights into the molecular mechanisms of hepatocellular carcinoma, which can be further utilized in the development of novel therapeutic strategies for diagnosis and treatment.

GEO dataset downloading and data preprocessing
We entered the terms "liver cancer" or "COVID-19" into the GEO database to select related gene expression datasets and then gained the liver cancer dataset GSE54236 and COVID-19 dataset GSE177477 from the GEO database and conducted differential gene analysis to produce volcano plots.Then, we screened the DEGs in the GSE54236 and GSE177477 datasets.(|log2FC|> 1, p < 0.05).

Weighted gene co-expression networks and module analysis
For DEGs in the HCC and COVID-19 datasets, coexpression networks with relative clinical distinctiveness were established using the WGCNA R package.The samples of the two datasets exhibited good clustering.The pickSoftThreshold function was employed to ascertain the soft threshold power and used for automated network construction.The results were obtained through the analysis of the topological overlap matrix, which produced module assignments labeled by color and module eigengenes (MEs).Ongoing efforts are made to evaluate the associations between MEs and clinical characteristics using Pearson correlation tests.The modules with a |ME|> 0.3 and an adjusted P-value < 0.05 were considered momentous in the interactions with clinical features [25].The four modular genes with the strongest positive or negative correlation with disease were screened using the Venn diagram overlap to obtain common DEGs between the four modules.These genes are potentially associated with the development of hepatocellular carcinoma and COVID-19.

Gene ontology and KEGG analysis of DEGs
To explain biological functions and further interactions, two DEGs related to COVID-19 and HCC were annotated with Gene Ontology (GO) analysis and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis.Gene Ontology, a bioinformatics tool, provides brief annotations of gene products based on their functions and biological pathways [26].The Kyoto Encyclopedia of Genes and Genomes pathway documents host genetic pathway data from different species [27].Adjusted P values < 0.05 were regarded as statistically significant.

Protein-protein interaction network construction
The Search Tool for the Retrieval of Interacting Genes serves to explore the relation between some proteins and establish a complex regulatory protein-protein interaction (PPI) network [28].Cytoscape serves to validate this network.The core functional modules were scrutinized utilizing Cytoscape's plug-in molecular compound detection technology, and the module genes' involvement for KEGG and GO analyses was then determined through a bioinformatics website.An interaction with a collective score above 0.4 was deemed statistically significant.

Selection and analysis of hub genes
We applied the Cytoscape cytoHubba plug-in for verifying the hub genes.We ascertained seven commonly used algorithms (MCC, MNC, Degree, Closeness, Radiality, EPC, DMNC) by cytoHubba to assess and select hub genes.After that, we gained the hub genes.Afterward, we created a co-expression network of these hub genes using GeneMAGIA [29], a reliable tool for identifying internal associations within gene sets.We examined the correlation between shared genes and pathways.

Differential expression profiling of hepatic cancer genes with the Human Protein Atlas
The expression levels of common hub genes in people who suffer from hepatic cancer are displayed using a boxplot.Differential expression analysis of common pivotal genes was carried out on normal and tumor tissues.Besides, we used the Human Protein Atlas (HPA) database [30] to gain immunochemical data for the relevant genes to ascertain their expression levels.This database facilitates the assessment of protein levels in both cancerous and non-carcinomatous tissues, alongside the assessment of OS in individuals with hepatocellular carcinoma.

Kaplan-Meier plotter
Kaplan-Meier plotter [31] was able to estimate the correspondence between the expression of all genes (mRNAs, miRNAs, proteins) and survival in 21 tumor kinds that include breast cancer, ovarian neoplasm, pulmonary malignancy, and stomach carcinoma.Resources of the database included GEO, EGA, and TCGA.We used the Kplot website to analyze the survival of common pivotal genes.This analysis evaluated the prophylactic value of several common pivotal genes (adjusted p-value < 0.05), all of which were statistically significant.With respect to the common situation, the change in the predictive proportion (endpoint incidence) between the two cohorts is proportional to the degree of division (divergence) between the two lines.

Prognostic gene identification and functional enrichment analysis
After identifying the central genes, we used the R package to produce nomogram tables to forecast disease progression and receiver operating characteristic curve (ROC) plots to determine whether the genes are valuable for disease diagnosis.Gene set enrichment analysis [32] calculates whether a priori defined genome indicates a statistically important and consistent distinction between two biological states.We performed GSEA on each gene using the clusterProfilter package and considered the first five pathways and the last five pathways of significant enrichment for characterization.In this study, gene enrichment analysis was used to elucidate significant functional and pathway differences.The first five and last five enrichment pathways were obtained as predefined gene sets for enrichment analysis and determined by adjusted p-values (< 0.05), FDR values (< 0.25), and normalized enrichment scores (|NES|> 1).

Immunoinfiltration analysis of pivotal biomarkers
We compared the concentration of immune cell accumulation in both carcinomatous and non-carcinomatous tissues and examined the quantity of immune cells in various samples using the CIBERSORT algorithm.CIBERSORT, a deconvolution algorithm, combines the labeled genomes of different immune cell subpopulations to estimate the ratios of 22 immune cells in tissues [33].Nonparametric correlation analysis (Spearman) served to determine the relation of biomarkers and immune infiltrating cells.

miRNA-hub gene network prediction and construction and transcription factor-hub gene network
miRNet is an accessible, web-based tool aimed at assisting in the clarification of microRNA (miRNA) function by combining user data with the existing body of knowledge through web-based functional analysis [34].The miRNet database [35] was used to establish miRNA gene interactions for central genes.Finally, these genes and miRNAs were mapped by the plug-in Cytoscape.TFs for regulatory hub genes were obtained through the miR-Net database, and an adjustment of P-value < 0.05 was regarded as significant.

Identification of potential drugs
The Drug-Gene Interaction Database is a web-based source that provides data on drug-gene interactions and druggable genes from both publishes, databases, and other web-based resources that can be used to evaluate the identification of drugs that interact with these genes [36].We used DGIdb to forecast drugs and molecular comminations that can interact with key genes and to download correlative data.Drug-hub gene interaction networks were mapped using Cytoscape software.

Immunohistochemistry
The expression of eight nuclear genes was inspected using immunohistochemistry.We found that the expression of RRM2 varied widely.Thirty patients with liver neoplasms from the People's Hospital in the Guangxi Zhuang Autonomous Region provided samples.Following paraffin fixation, the samples were divided into serial sections and exposed to incubatie a rabbit RRM2 polyclonal antibody at 4 °C for overnight.The plots were then colored with hematoxylin and eosin (HE).Immunohistochemical staining was performed using the Universal Two-Step Detection Kit (Mouse/Rabbit Enhanced Polymer Detection System) (ZSGB-BIO, PV-9000).The amount of positively stained cells during immunohistochemical labeling was ascertained with antigen content, distribution density, labeling method, and drug sensitivity rate.
We took several representative images using an Olym-busX21 microscope.Each image was then analyzed for general morphometry using ImageJ.Optical density and positive area data were obtained from normal and cancerous tissue by measuring selected stained areas using ImageJ parameters.The higher the optical density was, the more positive the expression and higher the average value.Eventually, statistical methods were used to confirm whether there was an evident distinction in RRM2 expression between the normal and cancer groups.

Statistical analysis
The Wilcoxon rank-sum test was used because normality monitoring made it clear that the specimen failed (p < 0.05).The Wilcoxon rank-sum test represented that the normal controls had markedly reduced scores than the carcinomatous tissue.The diversity (p < 0.001) was statistically significant.Statistical analysis was conducted with the R program, with the following levels of significance: ns, p < 0.05, *p < 0.05, **p < 0.01, and ***p < 0.001.

Weighted gene co-expression networks analysis (WGCNA) and module analysis
We performed WGCNA to investigate the correlation between clinical information and key genes as well as mRNA co-expression network analysis.Genes with obvious expression discrepancies (P < 0.05) were picked as entries for the WGCNA (Fig. 2a, b).In the WGCNA method, β = 5 was the best soft power value for GSE54236 (Fig. 2c), and β = 5 was the best soft power value for GSE177477 (Fig. 2d).28 modules were identified for GSE54236 (Fig. 2e, g), while 8 modules were shortlisted for GSE177477 (Fig. 2f, h).The modules associated with hepatocellular carcinoma and new coronavirus pneumonia were evaluated, and the modules with larger MSs were considered to have greater correlations with disease progression.For HCC, the magenta module exhibited the most significant positive correlation (r = 0.64), and the green module represented the most pronounced negative association in the GSE54236 dataset (r = −0.61).For COVID-19, the magenta module indicates the strongest positive relation (r = 0.72) and the most significant inverse relationship (r = −0.78) is shown in blue in the GSE177477 database.Therefore, these four modules were selected for further analysis.

Functional annotations and analysis
We applied Venn diagram analysis to recognize the DEGs among the four modules, and 223 DEGs were identified (Fig. 3a).These genes may be engaged in COVID-19 and hepatocellular carcinoma mechanisms.GO analysis indicated that the above 223 genes might be associated with DNA replication, regulation of DNA metabolic process, mitotic cell cycle phase transition, and so on.KEGG analysis indicated that these genes might be consistent with cell cycle, vitamin digestion and absorption, nucleotide metabolism, DNA replication, and so on.These results indicate that DNA replication and the cell cycle are jointly concerned with the occurrence and development of these two diseases (Fig. 3b, c).

PPI network construction with module analysis
PPI networks for the shared DEGs were built in Cytoscape, which consisted of 150 vertices and 1232 interaction sets (Fig. 4a).Three tightly connected gene modules including 41 common DEGs and 766 interaction pairs were obtained through Cytoscape's MCODE plugin (Fig. 4b-d).Besides, we did the enrichment analyses of GO and KEGG to ascertain the biological functions and pathways connected with the 41 genes.GO analysis and KEGG results indicated that these genes have relation with DNA replication (Fig. 4e, f ).

Selection and analysis of hub genes
Using the Cytoscape plug-in cytoHubba with seven algorithms (MCC, MNC, Degree, Closeness, Radiality, EPC, and DMNC) to analyze the 223 genes, we figured out the first 20 hub genes (Table 1).After taking the intersection of the Venn diagrams, we found 8 common central genes, including RRM2, TPX2, DTL, CDT1, TYMS, CDCA5, CDC25C, and HJURP (Fig. 5a).Based on the GeneMANIA database, we had an analysis of co-expression network and associated functions.These genes indicated a complicated PPI network with 93.77% co-expression, 5.09% colocalization, 0.71% physical interactions, and 0.43% genetic interactions (Fig. 5b).The full names and related functions are shown in Supplementary Table S1.We identified eight common central genes, namely RRM2, TPX2, DTL, CDT1, TYMS, CDCA5, CDC25C, and HJURP, using the intersection of eight common genes obtained from cytoHubba and 41 central genes obtained from MCODE (Fig. 5c).The GO analysis showed that these eight genes were enclosed mainly in the positive regulation of the mitotic cell cycle, chromosomal region, and chromosome formation (Fig. 5d), while the KEGG results suggested that progesterone-mediated oocyte maturation and metabolism of life activities were important pathways (Fig. 5e).Furthermore, we examined the connections of pathways with genes and the interconnections among the pathways.By the observation, TYMS and RRM2 are not connected to nucleotide metabolism; rather, they have tight relation with pyrimidine metabolism and nucleotide metabolism (Fig. 5f, g).

Identification and validation of the hub genes
We used a boxplot to indicate the expression levels of the eight hub genes in patients with hepatocellular carcinoma.Eight genes were employed differential expression analysis in both carcinomatous and noncarcinomatous tissues.Our studies revealed that all of the gene expression levels in carcinomatous tissues are superior than in non-carcinomatous tissues.Among them, RRM2 and DTL were significantly upregulated (Fig. 6a-h).To check the credibility of these central gene expression levels, we chose an extra dataset of hepatocellular carcinomas and estimated the expression levels of these central genes (Figure S1).The results indicated that hepatocellular carcinoma tissue exhibited significant over-regulation of all central genes as compared non-cancerous tissue.The protein expression levels of the eight genes were examined based on the HPA database.As illustrated in Figure 6i, the TYMS gene had no protein expression.However, superior expression levels of RRM2, TPX2, CDT1, CDCA5, and CDC25C were observed in both carcinomatous tissues, while DTL showed low expression and HJURP represented high expression in both carcinomatous and noncarcinomatous tissues.

Kaplan-Meier plotter identifies eight genes as prognostic markers for HCC survival
Kaplan-Meier plotter was applied to appraise the diagnostic and prognostic efficacy of the eight genes in the clinical context.In principle, the bigger the interval between the two curves (the bigger the bifurcation), the greater the variation in prognosis (endpoint incidence) between the two teams of patients.After evaluating the prognostic value of these eight genes, we found that all of them had P-value less than 0.5.Finally, it was concluded that elevated levels of eight genes (RRM2, TPX2, DTL, CDT1, TYMS, CDCA5, CDC25C, HJURP) in patients with HCC were associated with a poor OS (p < 0.05) (Fig. 7a-h).

Predictive model building and validation
A nomogram is a visually effective representation of the results of risk models and is convenient for predicting outcomes [37].The length of a straight line in a nomogram indicates the impact of different variables and their effect on the outcome.We described the survival period of patients using the K-Mplot database and a closer study of the tables revealed a good prognostic effect for all genes (Fig. 8a) (P-value < 0.05).The ROC plots for the eight genes were analyzed, and the results are presented in the diagram.We further studied the diagram and found that the average under the curves (AUCs) of the eight genes were 0.813, 0.771, 0.798, 0.784, 0.589, 0.851, 0.858, and 0.847.According to the curve position, the whole graph is split into two sections.The AUC is employed to indicate the predictive precision: the AUC is in proportion to the precision of the prognostication.The nearer the cross section of the curve is to the uppermost left-hand angle (the smaller the X and the bigger the Y), the more precise the prediction will be.Based on the data comparison, CDC25C (0.858), RRM2 (0.813), CDCA5 (0.851), and HJURP (0.847) have significant diagnostic value for the diagnosis of HCC (AUC > 0.8 was used as the judgment criterion) (Fig. 8b-i).

GSEA function enrichment analysis
We did a single-gene gene set enrichment analysis(GSEA), a learning-based approach for explaining genome-wide expression profiles [38].GSEA was performed for each gene using the clusterProfilter package, and the top five pathways and the last five pathways that were significantly enriched were considered for representation.The outcomes revealed that eight genes may be engaged in multiple biological processes, containing the cell cycle, DNA replication, and the Fanconi anemia pathway (Fig. 9a-h).

Immune infiltration analysis of hub biomarkers
The enrichment analysis results showed that immunity plays an important role in patients with liver cancer and patients with COVID-19.The CIBERSORT algorithm was used to analyze the abundances of immune cells in different samples.Bar diagrams indicated that there   , 11a).Furthermore, more kinds of immune cells can be detected in hepatocellular carcinoma samples.We found fewer T-cells CD8 and NK cell resting in COVID-19 samples than in control samples (Fig. 10b).Moreover, hepatocellular carcinoma samples had lower T-cell CD8 levels than normal samples, while neutrophils were increased (Fig. 11b).In hepatocellular carcinoma samples, these eight shared genes were strongly negatively associated with plasma cells, and RRM2 was positively associated with neutrophil and NK cell activation, as well as with dendritic cell resting (Fig. 10c-j).In contrast, RRM2 was significantly positively associated with plasma cells and monocytes in COVID-19 samples, while RRM2 was significantly negatively associated with dendritic cell resting (Fig. 11c-j).

miRNA-hub gene network prediction with construction and transcription factor-hub gene network
We aimed to better understand the regulatory roles of hub genes in the pathogenesis of HCC.The miRNANet database was employed to forecast the desired miRNAs of the hub genes.We structured the miRNA gene connection network on the miRNAnet website, which consists of 194 nodes and 367 edges.As shown below, we found that hsa-mir-103a-3p, hsa-mir-107, hsa-mir-129-2-3p, hsa-mir-34a-5p, hsa-mir-147a, hsa-mir-16-5p, and hsa-mir-195-5p interacted with 8 hub genes (Fig. 12a).However, these findings need further validation.Based on the miRNAnet database, we identified nine TFs that may modulate the expression of these genes.The interaction network consisted of two diagnostic genes and nine TFs, including TFDP1, TFCP2, SP1, ESR1, E2F1, CHD8, USF1, YBX1, and TP539.All nine TFs could regulate TYMS, but only E2F1 could regulate RRM2 (Fig. 12b).

Identification of potential drugs that interact with hub genes
The DGIdb database predicted drugs or molecular compounds that may interact with hub genes, and a total of 51 drugs or molecular compounds that may have regulatory relationships with hub genes were screened.The largest number of drugs interacted with TYMS (Fig. 13a), followed by RRM2 (Fig. 13b) and CDC25C (Fig. 13c).These genes did not interacted with TPX2, DTL, CDT1, CDCA5, or HJURP.Therefore, TYMS, RRM2, and CDC25C have significant roles in drug therapy research.There are 13 drugs, including TRIAPINE, CLADRIBINE, and CYT ARA BINE, that interact with RRM2 and may have therapeutic effects by acting on RRM2.

Role of RRM2 in hepatocellular carcinoma and its impact on the tumor microenvironment
Thirty patients with liver cancer from the People's Hospital in the Guangxi Zhuang Autonomous Region provided samples and immunohistochemical screening was 0.312433333, while the average positive area of no-carcinomatous tissue was 232,874.9333m 2 with an average optical density of 0.394316667.Analysis of the t-test results indicated that the significant value was less than 0.001 (p < 0.05), showing that the data from the two groups were very different and that the area of positive coloring in the cancer group was more than that in the normal group.These researches proved the significantly high expression of RRM2 in different databases of hepatic carcinoma tissues and clinical specimens (Fig. 14c).
To confirm the role of RRM2 in the regulation of immunity in LIHC, we analyzed the correlation between RRM2 and immunological features of the tumor microenvironment (TME).

Discussions
Liver injury is a potential complication of COVID-19, and may exacerbate the development of hepatocellular carcinoma due to shared molecular pathways with COVID-19.We generated volcano maps by filtering the gene expression data.This study conducted bioinformatics analyses on two independent datasets, GSE54236 and GSE177477.The liver cancer dataset GSE54236 screened 690 DEGs, while the new coronary pneumonia dataset GSE177477 screened 594 DEGs.Gene co-expression analysis and genomic modular network assays enable the exploration of intergenic relationships [39].Weighted gene co-expression network analysis is an innovative bioinformatics technique used extensively for high-throughput sequencing data analysis.The text discusses the significance of connectivity and genes within nodes, with WGCNA aiming to identify high-order relationships among gene products [40].WGCNA results with higher biological value and reliability are obtained by concentrates on correcting co-expression modules and clinical features [41].The hub genes and modules with biological correlations may function as biomarkers for diagnosis or treatment.Comorbidity mechanisms in a range of diseases have been successfully investigated through the use of the WGCNA method.For example, Zhu et al. used WGCNA to identify ten hub genes involved in the development of Alzheimer's disease and T2DM [42].Bi et al. also employed WGCNA to identify a novel biomarker for distinguishing alcohol-associated HCC from nonalcohol-associated HCC [43].In this study, WGCNA was applied to analyze a total of 3653 genes for COVID-19 and HCC.This study examined the relationships between genes and clinical characteristics by creating co-expression modules with central features.In HCC and COVID-19, a collective of 28 gene modules exhibiting coordinated expression patterns were identified.
This research screened for genes that strongly correlate with HCC and COVID-19.A total of 223 genes commonly found in HCC group and COVID-19 group were ascertain for numerous bioinformatics analyses.GO analysis showed that these genes mainly promoted the regulation of DNA metabolic process, DNA replication, and mitotic cell cycle phase transition.KEGG pathway analysis revealed several enriched pathways, including cell cycle, vitamin digestion and absorption, nucleotide metabolism, and so on.These pathways were further categorized according to the KEGG database.These are primarily associated with cell separation and multiplication.These results propose that patients with COVID-19 may be at an elevated chance of developing hepatocellular carcinoma.According to the MCODE and cytoHubba plug-in in Cytoscape, we ascertained 8 central genes.These genes included RRM2, TPX2, DTL, CDT1, TYMS, CDCA5, CDC25C, and HJURP.These genes were almost upregulated in both HCC and COVID-19 (except for CDT1 and HJURP in GSE54236 and TYMS in GSE177477), indicating their potential importance in the development of HCC and COVID-19.Further validation of the GEPIA database showed that carcinoma survival in HCC patients was associated with 8 central genes.This study indicated that these 8 genes have prognostic and MicroRNAs modulate gene expression by partly or entirely reinforcing the 3' untranslated region (URT) of selected gene miRNAs, triggering miRNA breakdown or interfering with miRNA translation [43,44].This investigation included the construction of a miRNA target gene network and the selection of seven miRNAs (hsa-mir-103a-3p, hsa-mir-107, hsa-mir-129-2-3p, hsamir-34a-5p, hsa-mir-147a, hsa-mir-16-5p, and hsa-mir-195-5p) that interact with the hub genes.These findings indicate that the upregulation of EVA1A by miR-103a-3p potentially acts as a key mediator in inhibiting HCC cell growth and reproduction [45].Hsa-miR-107 regulates the growth of hepatocellular carcinoma cells.The high representation of hsa-miR-107 in tissues affected by hepatocellular carcinoma suggests that aberrantly expressed miR-107 may play the part of a promoter of cancer in this disease [46].MiR-129-3p suppressed the activation of signaling pathways by modulating them and prevented hepatocellular carcinoma cells from proliferating, metastasizing, and infiltrating [47].MiR-34a-5p prevented the growth and development of HCC by targeting VEGFA [48].Researchs have represented that the tissues and cells of liver tumors have down-regulated expression of miR-16-5p.MiR-16-5p suppressed the multiplication, invasion, and metastatic potential of hepatocellular carcinoma cells by specifically binding to and inhibiting IGF1R protein expression [49].Overexpression of hsa-miR-195-5p was found to inhibit the proliferation, invasion, and migration of hepatocellular carcinoma cells by reducing PHF19 expression [50].There is currently limited information available on hsa-mir-147a, and further research is needed.
We predicted that 13 drugs or molecular compounds could be involved in the regulation of hub genes and could be potential drugs for the cure of liver cancer.RRM2 was found to have a positive correlation with most immunomodulators, chemokines, and chemokine receptors in LICH, KIPAN, and KICH.Conversely, it was found to have a negative correlation with most immunomodulators, chemokines, and chemokine receptors in THYM and TGCT.Additionally, RRM2 was found to have a positive correlation with most major histocompatibility complexes (MHCs) in LIHC, KIPAN, and KIRC, as well as in TGCT, NB, and THYM.Analysis using the TISIDB database demonstrated a positive association between the expression of RRM2 and the expression of the chemokines TNFSF4 (ρ = 0.361, p < 0.001), CDC27 (ρ = 0.426, p < 0.001), MICB (ρ = 0.520, p < 0.001), and TAP1 (ρ = 0.355, p < 0.001).Hepatocellular carcinoma tissues had higher levels of RRM2 mRNA than in normal tissue, indicating the potentiality of RRM2 as a biomarker for liver cancer development.
Finally, we acquired clinical samples for immunohistochemical dissection to appraise and verify the accuracy of our findings.Our findings indicate that RRM2 expression is higher in hepatic tumorous tissues compared to healthy ones, which is consistent with our other results.
Due to the recent emergence of neocrown pneumonia, there are limited studies on the association between COVID-19 and liver cancer.Furthermore, a small number of studies have investigated the shared molecular mechanisms between the two through advanced bioinformatics approaches.Liver injury is a known complication of neocoronary pneumonia and may increase the risk of hepatocellular carcinoma.The literature suggests that the 30-day mortality risk of SARS-CoV-2 can be used to inform the prognosis of HCC in the context of SARS-CoV-2 infection.We ascertained shared DEGs and central genes in both cancerous and non-carcinomatous tissue, which could help to further understand the relationship between neocrown pneumonia and hepatocellular carcinoma.However, our study has some limitations.This study is forward-looking and needs to be verified by an outside source.However, most of the data were obtained using a database abstraction approach, without significant support from other studies.Additionally, the new Crowne Plaza disease lacks sufficient data to support it.Future research will center on further validating the roles of the central genes in vitro models.Furthermore, it is crucial to remember that our current study only included the first 8 central genes and lacked detailed information on the molecular mechanisms of both central gene and miRNA regulation in HCC and COVID-19.Additionally, the miRNA-gene interaction networks are based solely on predictions from public databases.Therefore, more research is needed to completely comprehend the roles of central genes and miR-NAs in the development of COVID-19 and HCC.
Similar research methods can be used to analyze the mechanisms of other sequelae caused by COVID-19, advancing a deeper comprehension of the illness and assisting in the creation of all-encompassing treatment and recovery programs to enhance patient survival rates and the living level of survivors.By explaining the mechanisms of many sequelae, attention can be drawn to neocoronavirus pneumonia, increasing awareness of integrated control and preempting the expansion of the epidemic and the virus.In conclusion, our study suggested that miRNA gene regulatory networks may contribute to the pathophysiology of both HCC and COVID-19.We identified some central genes that potentially serve as targets for diagnosis and treatment.However, additional investigation and practical clinical application are required to confirm these results.

Conclusions
In conclusion, COVID-19 vulnerability modules and genes associated with HCC were characterized by coexpression network analysis.This study provides a new idea to study the common molecular mechanism of hepatocellular carcinoma and COVID-19.We identified common DEGs for hepatocellular carcinoma and COVID-19 and performed enrichment and PPI network analyses.Several hub genes, such as RRM2, TPX2, DTL, CDT1, TYMS, CDCA5, CDC25C, and HJURP, have been instrumental in the biological and disease-related mechanisms of hepatocellular carcinoma and COVID-19.Some transcription factors, miRNAs, and drugs that may regulate the hub genes have also been revealed, potentially offering novel therapeutic targets for HCC and COVID-19, but more experimental and clinical practice is needed to verify these findings.RRM2 interacts with TNFSF4, CDC27, MICB and TAPI.Besides, we clinically performed immunohistochemical analysis of RRM2, and the result is consistent with existing research.Therefore, our research provides a novel avenue for exploring therapeutic approaches for liver cancer patients with COVID-19.

Fig. 1
Fig. 1 Volcano maps and heatmap.a Volcano map of GSE177477.b Volcano map of GSE54236.Up-regulated genes are marked in red; down-regulated genes are marked in blue.c Heatmap of DEGs in GSE54236.d Heatmap of DEGs in GSE177477 dataset

Fig. 3
Fig. 3 Venn diagram and module gene enrichment results.a Venn diagram of the crossover of module genes.b Results of enrichment analysis of the GO pathways.c Consequence of enrichment analysis of the KEGG.An adjusted P-value of < 0.05 was regarded as having statistical significance

Fig. 4
Fig. 4 PPI network map and gene module and gene enrichment analysis.a PPI network diagram.b-d Three significant gene clustering modules.e GO enrichment results of the modular genes.f GO enrichment results of the 41 genes

Fig. 5 Fig. 6
Fig. 5 Venn diagram and co-expression network of hub genes.a Venn diagram showing the seven algorithms screened for eight overlapping central genes.b Hub genes and their co-expression genes were analyzed via GeneMANIA.c MOCDE with cytoHubba of overlapping hub genes.d GO enrichment results.e KEGG enrichment results.f Interaction network graph between central genes and pathways.g Interaction network map of pathways

Fig. 10
Fig. 10 Immune infiltration analysis of shared biomarkers in patients with COVID-19.a Histogram of immune cell infiltration.b Box diagram of the proportions of 22 types of immune cells.The diagrams represent the difference in infiltration between the two groups of samples.c-j Correlations between RRM2, TPX2, DTL, CDT1, TYMS, CDCA5, CDC25C, HJURP and infiltrating immune cells

Fig. 11
Fig. 11 Immune infiltration analysis of shared biomarkers in HCC. a Histogram of immune cell infiltration.b Box diagram of the proportions of 22 types of immune cells.The diagrams represent the difference in infiltration between the two groups of samples.c-j Correlations between RRM2, TPX2, DTL, CDT1, TYMS, CDCA5, CDC25C, and HJURP and infiltrating immune cells

Fig. 12
Fig. 12 Integrated miRNA-gene interaction network of eight pivotal genes and the TF regulatory network.a The magenta circles mean eight hub genes.Blue squares mean miRNAs that have relation to hub genes.Abbreviations: miRNA, microRNA.b TF regulatory connectivity map.The green means TFs, and the magenta stands for hub genes

Fig. 14
Fig. 14 Transcription level, immunomodulators, chemokines, and receptors related to RRM2. a Distribution of RRM2 immunological scores in tumor and normal tissues.The ordinate reflects the distribution of the immunological scores in distinct groups, whereas the abscissa indicates the immune cell types.The Wilcoxon test was utilized to contrast the statistical distinctions between two groups, and the Kruskal-Wallis test devoted to affirm the pronounced variations across multiple groups.a Heatmap of immune cell scores.Different hues represent the varied expression distributions in different samples.Significance is denoted by asterisks at levels: *p < 0.05, **p < 0.01, and ***p < 0.001.b Percentages of tumor-infiltrating immune cells in samples.Different colors depict different types of immunological cells.The horizontal axis corresponds to the sample, and the vertical axis represents the quantity of immune cells present in each sample.b Part of the results of the immunohistochemistry experiments.c Mean optical intensity of liver cancer tissue and contiguous tissue in 30 patients with liver cancer.The Cancer Genome Atlas.*P < 0.05, ***P < 0.001.d-g Immunomodulators, chemokines, and receptors associated with RRM2 in liver hepatocellular carcinoma (LIHC).d TNFSF4 e CDC27 f MICB g TAP1

Table 1
Top 20 hub genes rank in cytoHubba