Identification of Comorbid Genes of Psoriasis and Coronary Artery Disease Using Ustekinumab with Astragalus

Abstract

Background

Systemic inflammation leading to coronary atherosclerotic heart disease (CAD) increases the risk for Psoriasis (PSO). Despite progress in diagnostic and therapeutic strategies for CAD, it is hard to identify many patients in time, and their lives are threatened.

Purpose

To explore the comorbidity genes of psoriasis (PSO) and coronary artery disease (CAD) under the intervention of Ustekinumab combined with Astragalus, and to provide new ideas for clinical treatment.

Methods

The PSO datasets (GSE13355, GSE14905) and CAD datasets (GSE113079, GSE23561) were obtained from the Gene Expression Omnibus database. Key genes were screened through differential expression analysis, weighted gene co-expression network analysis (WGCNA), gene enrichment analysis, protein–protein interaction (PPI) network construction, and immune infiltration analysis.

Results

A total of 2,412 PSO differential genes and 4,443 CAD differential genes were identified. WGCNA analysis yielded 8,809 PSO-related module genes and 8,389 CAD-related module genes. Cross-analysis identified 186 genes with consistent expression trends. Eventually, aurora kinase B (AURKB), C-X-C motif chemokine ligand 11 (CXCL11), late cornified envelope 3D (LCE3D), tumor necrosis factor ligand superfamily member 10 (TNFSF10), and CD247 were determined as key comorbidity genes, enriched in the proteasome and olfactory conduction pathways.

Conclusion

AURKB, CXCL11, LCE3D, TNFSF10, and CD247 were indicated as comorbid genes for PSO and CAD, and they therefore might be involved in the progression of these diseases.

Introduction

Psoriasis (PSO) is an autoimmune disorder characterized by erythematous, scaly, and indurated plaques, epidermal hyperplasia, and inflammatory infiltration of leukocytes (Duan et al., 2024). Apart from pain, PSO also affects joints and nails, and it usually coexists with systemic diseases such as cardiovascular disease, chronic kidney disease, hyperlipidemia, metabolic syndrome, and so on (Guo et al., 2023; Patel & Kumthekar, 2022). PSO patients also have psychiatric disorders and psychosocial distress, like depression. Intrinsic risk factors for PSO include obesity, diabetes, dyslipidemia, hypertension, and mental stress (Lee et al., 2018). On the other hand, coronary atherosclerotic heart disease (CAD) is characterized by coronary atherosclerotic lesions due to vascular lumen stenosis, hypoxia, or necrosis (Fu et al., 2019; Huang et al., 2019; Kawai et al., 2024). It is the result of atherosclerotic plaque accumulation, which prevents oxygen and nutrients from reaching the heart. Furthermore, inflammation is a key process for CAD development (Figueiredo et al., 2023; Kong et al., 2022). Epidemiologically, age and smoking are known risk factors for CAD (Khera & Kathiresan, 2017; Ng et al., 2020). Furthermore, systemic inflammation leading to CAD increases the risk for PSO (Bacchiega, 2019; Egeberg et al., 2017; Gelfand et al., 2006; Lerman et al., 2017; Shiba et al., 2019; Wu et al., 2015). Despite progress in diagnostic and therapeutic strategies for CAD, it is hard to identify many patients in time, and their lives are threatened (Hüffmeier et al., 2009). Thus, it is necessary to investigate the molecular mechanism for CAD and discover its key genes to develop diagnostic and therapeutic methods.

IL-12 and IL-23 are the main interleukins that can be used as pharmacological targets in treating PSO (McInnes et al., 2013; Ritchlin et al., 2014). Ustekinumab is considered an IL-12/23 inhibitor as it binds to the p40 and blocks the binding to the receptor, thereby reducing the inflammatory response. Risankizumab is known to hinder the process of pro-inflammatory cytokines and chemokine release by binding to the p19 subunit of IL-23, thereby alleviating inflammation and PSO (Baeten et al., 2018). Astragalus, as a traditional Chinese medicine, has been proven to have immunomodulatory and anti-inflammatory effects (Ng et al., 2020; Patel & Kumthekar, 2022). Recent studies have shown that astragalin can alleviate the inflammatory response in the PSO model by inhibiting the nuclear factor-kappa B pathway (Ritchie et al., 2015), but its role in the comorbidity of PSO and CAD remains unclear. Herein, this study investigated the comorbid genes for PSO and CAD through bioinformatics intervention of Ustekinumab plus Astragalus. The findings might provide new insights into the pathological mechanisms of PSO and CAD, as well as clinical treatment.

Materials and Methods

Data Source

PSO and CAD datasets were downloaded from the Gene Expression Omnibus database. The GSE13355 dataset (training set) included a total of 122 samples: 64 normal skin tissue samples (normal control group) and 58 PSO skin tissue samples under Ustekinumab intervention (experimental group).

The GSE14905 dataset (validation set) consisted of 54 samples: 21 normal skin tissue samples (normal control group) and 33 PSO skin tissue samples under Ustekinumab intervention (experimental group).

The GSE113079 dataset (training set) had 141 samples: 48 normal blood samples (normal control group) and 93 CAD blood samples under Ustekinumab intervention (experimental group). The detection was carried out on Agilent-067406 with Human CBC lncRNA + mRNA microarray V4.0.

The GSE23561 dataset (validation set) consisted of 15 samples, including 9 normal blood samples and 6 CAD blood samples treated with Ustekinumab. The Human 50K Exonic Evidence-based Oligonucleotide array was adopted for detection, along with the molecular formulas for Ustekinumab.

Differential Expression Analysis and Enrichment Analysis

The differentially expressed genes (DEGs) between PSO and normal in GSE13355, and CAD samples and control samples in GSE113079 were obtained from DEA through the limma package (Ritchie et al., 2015). The DEA results were presented on a volcano plot drawn by the ggplot2 software package, and heat maps were drawn by the pheatmap software package. The clusterProfiler software package was adopted to explore the functions and enrichment pathways for differential genes to analyze PSO enrichment (Yu et al., 2012).

Identification of Weighted Gene Co-expression Network Analysis and Candidate Genes

Weighted gene co-expression network analysis (WGCNA) was performed on GSE13355 and GSE113079, and the data with maximum relevance were screened from the WGCNA package (Langfelder & Horvath, 2008). The samples were clustered with a soft threshold (β). Module clustering trees were constructed based on a dynamic tree cut. Finally, the correlation between modules and features was calculated to determine key modules. The DEGs in PSO and CAD were intersected, respectively, to identify DEGs with consistent expression trends in the two diseases. The intersecting genes between DEGs and WGCNA were acquired. To explore the interrelationships among crossed genes, the STRING database was used to establish a protein–protein interaction (PPI) network. Next, the ggplot2 software package was used to detect the expression of overlapped genes with reciprocal relationships in the training and validation sets for PSO and CAD, aiming to identify candidate genes.

Key Gene Identification and Gene Set Enrichment Analysis

The candidate genes for PSO and CAD were included in the least absolute shrinkage and selection operator (LASSO) regression analysis, and characteristic genes were obtained through the glmnet software package (Ye et al., 2021). Key genes were identified by overlapping the characteristic genes of the two diseases with a co-expression network constructed using the GeneMANIA database (http://www.genemania.org/). To explore the underlying biological mechanisms, gene set enrichment analysis (GSEA) was performed for each key gene. Then, Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analysis was conducted using clusterProfiler (Yu et al., 2012) and the org.hs.egg.db package with thresholds NOM p < .05, Q < 0.25, |NES| > 1.

Immune Infiltration and Drug Prediction

Subsequently, the immune cell infiltration was compared by t-test to obtain differentially expressed immune cells, and the association between key genes and immune cells was analyzed using the ggcorrplot package. The key genes predicted by microRNA (miRNA) were obtained from the miRDB and miRTarbase databases, and key genes predicted by transcription factors (TF) were collected from the ChEA3 database and used to construct a miRNA-mRNA-TF network. As drugs related to key genes from the Comparative Toxicogenomics Database were detected, the intersection of predicted drugs and common predicted drugs was screened to construct a network of key genes and drugs.

Reverse Transcription Quantitative Polymerase Chain Reaction

The paraffin tissue RNA was extracted by the miRNeasy FFPE Kit (Lot number: 217504; Qiagen) and subjected to reverse transcription quantitative polymerase chain reaction. With glyceraldehyde-3-phosphate dehydrogenase (GAPDH) as an endogenous control, data were analyzed by software and calculated using the 2^−∇∇CT method. The specific primer design is shown in Table 1.

Table 1.

Primer Sequences.

Gene	Forward 5′-3′	Reverse 5′-3′
GAPDH	AGTATGTTGTCACCGCTGG	TAACCTGTCTATACGGAGGGT
vWF	GAGCCACTACGGGATTC	AAGACGAGCTTAGGAGCG
CD55	ACTCTACGGGACGAGC	GAGACGACTTAGGAGA
CD138	GATCGAAGCTAGCAGAGAGTA	GAGCTAGAGCTCGCAGAG

Tissue Staining

Tumor tissues were sliced into sections and stained with 0.1% Sirius Red F3B (Solarbio, China) and 1.3% picric acid or Masson’s trichrome saturated aqueous solution (Solarbio) to assess collagen deposition. Ten non-overlapping regions (n = 7) were scored semi-quantitatively relative to total area. The tissue specimens were dehydrated and fixed before sectioning, and then sections were stained with hematoxylin and eosin. Inflammatory cell markers were analyzed by immunohistochemical (IHC) staining and photographed with a Leica DM750 microscope, an EC3 digital camera, and LAS EZ software (200× magnification). The primary antibodies were CD3 and CD15 (1:5,000). Immunofluorescence (IF) images were acquired with a ZEISS fluorescence microscope using AxioVision Rel 4.8 software (400× magnification). Positive signals from IF images were quantified from at least five random fields using ImageJ 1.52 software.

Statistical Analysis

R language was adopted for all bioinformatics analyses, and limma was used to identify DEGs, Venn to intersect multiple gene sets, and WGCNA to screen modules related to disease scoring. Gene correlation analysis was conducted using the Spearman rank correlation coefficient, with a significance threshold of r > 0.3 and p < .05. The multiple hypothesis test was corrected by the Benjamini–Hochberg method, and a false discovery rate <.05 was considered significant.

Results

Immunomodulatory Effect of Ustekinumab Combined with Astragalus on Psoriasis

To evaluate whether IL-17 inhibition affects the immunomodulatory effects of PSO, we evaluated the condition of immune cell infiltration before and after treatment with Ustekinumab plus Astragalus. The amount of some infiltrating cell subtypes decreased at week 24, but it was clear that the Ustekinumab plus Astragalus treatment failed to resolve the psoriatic synovial inflammation at all. The number of von Willebrand factor (vWF)-positive cells decreased slightly, with no differences in the number of inner CD55+ cells or CD138+ plasma cells (Figure 1).

Figure 1.

The Immunomodulatory Effect of Ustekinumab Combined with Astragalus on Psoriasis. (A) Immunohistochemical (IHC) of Immune Cell Expression. (B). Reverse Transcription Quantitative Polymerase Chain Reaction (RT-qPCR) of Immune Cell Expression Among WF Cells, CD55+ Cells, or CD138+ Plasma Cells. Scale Bar = 200 µm.

Identification of Psoriasis and Coronary Artery Disease Genes

A total of 2,412 genes were screened from the GSE13355 cohort, including 1,258 upregulated genes and 1,154 downregulated genes in PSO (Figure 2A). We observed 4,443 DEGs from the GSE113079 cohort, including 2,348 upregulated and 2,095 downregulated in CAD (Figure 2B). The DEGs for PSO were involved in the regulation of skin development and fatty acid metabolism in graphene oxide biological processes (BP), secretory granule lumen in graphene oxide cellular components (CC), and molecular functions (MF) of graphene oxide, chemokine receptor binding, and peptidase regulatory activity (Figure 2C). Graphene oxide enrichment analysis results showed that DEGs of CAD were mainly enriched in BP, CC, and MF to positively regulate cellular defense, collagen-containing extracellular matrix, cell lysis particles, and cytokine receptor binding, respectively (Figure 2D).

Figure 2.

Identification of Psoriasis (PSO) and Coronary Artery Disease (CAD) Genes. (A) Volcano Plot of Differentially Expressed Genes (DEGs). (B) Heat Map of DEGs. (C) GO Enrichment Analysis of DEGs. (D) Kyoto Encyclopedia of Genes and Genomes (KEGG) Enrichment Analysis of DEGs.

Selection of Nine Genes as Candidate Genes

With no abnormal samples found in the PSO samples in the GSE13355 dataset, the optimal soft threshold (β) was finally set as 6 (Figure 3A). By constructing a module clustering tree, a minimum of 100 genes/module was set, and 10 modules were obtained (Figure 3B). The turquoise module was positively correlated with PSO and was identified as the key module containing 8,809 genes (Figure 3C). As such, we found no abnormal samples in the CAD samples of the GSE113079 dataset, and the optimal threshold (β) was finally determined as 7 (Figure 3D). With the minimum number (n = 100) of genes to cluster into a module, six modules were obtained (Figure 3E).

Figure 3.

Candidate Genes. (A) Selection of Soft Threshold β. (B) GSE13355 Module Clustering Diagram. (C) Grouping Samples into Phenotypes, Heat Map of Relationship Between Gene Modules and Traits. (D) Selection of Soft Threshold β. (E–F) Module Clustering Tree.

Acquisition of Nine Candidate Genes

The intersection of the up- and downregulated DEGs of PSO and CAD was identified, resulting in a total of 186 DEGs with consistent expression trends (Figure 4A). A Venn diagram depicted the intersection of gene fragments with consistent expression trends and key modules for PSO and CAD, and indicated 106 overlapped genes (Figure 4B). To further screen genes, a PPI network of crossed genes was constructed. The network indicated 19 vertices and 30 edges (Figure 4C). Through analysis, we verified 19 genes in PSO and nine in CAD, and thereby considered AURKB, CCL19, CXCL11, CD247, DLGAP5, LCE3D, REN, SPRR2B, and TNFSF10 genes as candidate genes (Figure 4D and 4E).

Figure 4.

Obtaining Candidate Genes. (A) Venn Diagram of Differentially Upregulated Genes. (B) Venn Diagram of Intersected Genes. (C) Protein–Protein Interaction (PPI) Network of Intersected Genes. (D) Verification of Gene Expression in Psoriasis Dataset. (E) Verification of Gene Expression in Coronary Artery Disease (CAD).

Key Genes Related to the Proteasome and Olfactory Transmission Pathway

In order to further identify the key genes, LASSO regression analysis of selected candidate genes was carried out, and indicated six characteristic genes related to PSO (AURKB, CXCL11, LCE3D, SPRR2B, TNFSF10, and CD247) and eight characteristic genes related to CAD (AURKB, CXCL11, LCE3D, CCL19, REN, DLGAP5, TNFSF10, and CD247) (Figure 5A). We then carried out an intersection analysis of featured genes for the two diseases, and the results pointed out five key genes (AURKB, CXCL11, LCE3D, TNFSF10, and CD247) (Figure 5B). GeneMANIA results demonstrated that key gene-related functions were mainly related to cytokine activity (Figure 5C). With regard to PSO, GSEA demonstrated that AURKB was enriched in 75 KEGG pathways, including proteasome and pyrimidine metabolism; CXCL11 in 74 KEGG pathways, including proteasome and toll-like receptor signaling pathways; LCE3D in 76 KEGG pathways, including proteasome and pyrimidine metabolism; TNFSF10 in 73 KEGG pathways, such as proteasome and pyrimidine metabolism; and CD247 in 78 KEGG pathways, including proteasome and cell cycle (Figure 5D). GSEA demonstrated that AURKB in CAD was enriched in 20 KEGG pathways, such as olfactory transduction and spliceosome; CXCL11 in 43 KEGG pathways, such as olfactory transduction and spliceosome; and LCE3D in 61 KEGG pathways, such as olfactory transduction and RNA degradation; TNFSF10 in 8 KEGG pathways, such as ribosome and spliceosome; and CD247 in 30 KEGG pathways, such as olfactory transduction and calcium signaling pathways (Figure 5E). Collectively, the above data demonstrated that these key genes were coenriched in the proteasome and olfactory transduction pathways.

Figure 5.

Key Genes Associated with Proteasome and Olfactory Transmission Pathway. (A) LASSO Regression Analysis. (B) Venn Diagram of Key Genes. (C) Co-expression Network of Key Genes. (D) Kyoto Encyclopedia of Genes and Genomes (KEGG) Enrichment Analysis of Key Genes. (E) Key Gene CD247 KEGG Enrichment Analysis.

Immune Infiltration Analysis of Psoriasis and Coronary Artery Disease

There were 18 types of immune cells that were significantly different in PSO (Figure 6A). Spearman analysis depicted a negative correlation between key genes and B cell memory, T cell regulatory (Tregs), and mast cell resting immune cells, and a positive correlation with CD4 memory activation, T cell follicular helper, and dendritic cell activated immune cells (Figure 6B). The results from immune infiltration analysis also showed a large difference (Figure 6C). More importantly, there was a significant positive correlation between monocyte and immune cell activation (Figure 6D). The predicted miRNAs of key genes were intersected based on the miRDB and miRTarbase databases (Figure 6E). We constructed a miRNA-mRNA-TF network of 93 nodes and 144 edges (Figure 6F).

Figure 6.

Immune Infiltration Analysis of Psoriasis (PSO) and Coronary Artery Disease (CAD). (A) Heat Map of Immune Cell Infiltration. (B) Box Plot of 22 Types of Immune Cell Infiltration in Disease and Healthy Samples. (C) Correlation Plot of Key Genes and Differential Immune Cells. (D) Heat Map of Immune Cell Infiltration. (E) Box Plot of 22 Types of Immune Cell Infiltration in Two Types of Samples. (F) Correlation Diagram of Key Genes and Differential Immune Cells.

Discussion

To identify the potential comorbid genes for PSO and CAD under Ustekinumab intervention, our study evaluated cellular infiltration and molecular pathways to explore how Ustekinumab blockade affects these comorbid genes.

This study demonstrated that cellular infiltration and vascular congestion were significantly increased in patients with comorbid PSO and CAD following Ustekinumab intervention within this mechanistic trial design. Ustekinumab treatment (12 weeks) importantly reduced the number of infiltrating immune cells, which was consistent with a previous study in which sustained synoviocyte infiltration was observed both clinically and sonographically in patients with PSO after effective treatment. Therefore, it can be assumed that identifying comorbid genes for PSO and CAD will benefit future clinical and basic research, and the following results and conclusions are obtained. First, gene differential expression analysis (DEA) was performed on the two disease samples and the healthy samples in the training set. As a result, we selected 2,412 DEGs from the PSO dataset and 4,443 from the CAD dataset. WGCNA was conducted to screen the module genes with the strongest correlation with phenotypes of the two diseases, as the DEGs and module genes in the two diseases were intersected to obtain 106 genes. This study then verified the expression of key genes to obtain nine candidate genes. The enrichment analysis of DEGs identified common pathways among the differential genes. Next, among the nine candidate genes obtained above, the LASSO machine indicated characteristic genes in the two diseases, namely AURKB, CXCL11, LCE3D, TNFSF10, and CD247. Previous studies have also shown that the KEGG pathway and CXCL11 signaling pathway were over-represented in both the total group and ACR20 responders, as CXCL11 was found to be involved in Th17 cell differentiation (Deodhar et al., 2019; Kavanaugh et al., 2014) and considered as a potential therapeutic target. Studies have also noted activation of the CXCL11 signaling pathway in spondyloarthritis (SpA) synovium, and inhibition of CXCL11 is thereby expected to alleviate the inflammatory response of immune cells, skin, and synovial fibroblasts in SpA. Blocking of CXCL11 with rapamycin inhibited arthritis in the SpA rat model, bone remodeling, and osteogenic differentiation of synovial fibroblasts from SpA patients in vitro (Kavanaugh et al., 2016; Nagai et al., 2013). Based on recent data (Roller et al., 2012), we expect that the CXCL11 pathway may also participate in the regulation of innate-like cells and adaptive cells (Th17) or synovial fibroblast cells, though this finding still needs to be confirmed by further experiments (Zhang & Liu, 2002). One limitation of our study is that all these identified indicators have not been verified in patient samples, and their relationship with tumor development and severity remains unclear, which requires further investigations.

Conclusion

This study used DEA, WGCNA, machine learning, enrichment analysis, and drug discovery methods to identify several comorbid genes in CAD and PSO, indicating their role in the development and pathogenesis of CAD and PSO, implying that these identified genes are novel indicators for the diagnosis of and therapeutic targets for the treatment of these diseases.