Upregulation of Biomarker Limd1 Was Correlated with Immune Infiltration in Doxorubicin-Related Cardiotoxicity

Zhang, Rui; Hao, Chunshu; Ji, Zhenjun; Qu, Yangyang; Zuo, Wenjie; Yang, Mingming; Zuo, Pengfei; Carvalho, Abdlay; Ma, Genshan; Li, Yongjun

doi:https://doi.org/10.1155/2023/8347759

Mediators of Inflammation

On this page

Abstract Introduction Materials and Methods Results Discussion Conclusions Data Availability Conflicts of Interest Authors’ Contributions Acknowledgments Supplementary Materials References Copyright Related Articles

Research Article | Open Access

Volume 2023 | Article ID 8347759 | https://doi.org/10.1155/2023/8347759

Upregulation of Biomarker Limd1 Was Correlated with Immune Infiltration in Doxorubicin-Related Cardiotoxicity

Rui Zhang,¹Chunshu Hao,¹Zhenjun Ji,¹Yangyang Qu,¹Wenjie Zuo,¹Mingming Yang,¹Pengfei Zuo,¹Abdlay Carvalho,¹Genshan Ma,¹and Yongjun Li¹

Academic Editor: Daniela Caccamo

Received26 Apr 2022

Revised13 Oct 2022

Accepted10 Jan 2023

Published23 Mar 2023

Abstract

Doxorubicin is one of the most common antitumor drugs. However, cardiotoxicity’s side effect limits its clinical applicability. In the present study, Gene Expression Omnibus (GEO) datasets were applied to reanalyze differentially expressed genes (DEGs) and construct weighted correlation network analysis (WGCNA) modules of doxorubicin-induced cardiotoxicity in wild-type mice. Several other bioinformatics analyses were performed to pick out the hub gene, and then the correlation between the hub gene and immune infiltration was evaluated. In total, 120 DEGs were discovered in a mouse model of doxorubicin-induced cardiotoxicity, and PF-04217903, propranolol, azithromycin, etc. were found to be potential drugs against this pathological condition. Among all the DEGs, 14 were further screened out by WGCNA modules, of which Limd1 was upregulated and finally regarded as the hub gene after being validated in other GEO datasets. Limd1 was upregulated in the peripheral blood mononuclear cell (PBMC) of the rat model, and the area under curve (AUC) of the receiver operating characteristic curve (ROC) in diagnosing cardiotoxicity was 0.847. The GSEA and PPI networks revealed a potential immunocyte regulatory role of Limd1 in cardiotoxicity. The proportion of “dendritic cells activated” in the heart was significantly elevated, while “macrophage M1” and “monocytes” declined after in vivo doxorubicin application. Finally, Limd1 expression was significantly positively correlated with “dendritic cells activation’ and negatively correlated with “monocytes” and “macrophages M1’. In summary, our results suggested that limd1 is a valuable biomarker and a potential inflammation regulator in doxorubicin-induced cardiotoxicity.

1. Introduction

Thanks to the advancements in science and technology, numerous treatment schemes have been developed, which have gradually extended the survival times of tumor patients. Chemotherapy drugs including doxorubicin played a vital role in achieving this trend. Doxorubicin, also known as adriamycin, is a member of the anthracycline antibiotics with decades of history [1]. In the antitumor field, doxorubicin is widely applied in clinical settings for the treatment of various types of tumors, from soft tissue and bone sarcomas, cancers of the breast, ovary, bladder, thyroid and lung, acute lymphoblastic leukemia, acute myeloblastic leukemia, Hodgkin lymphoma, etc. [2].

In addition to its unique therapeutic function, this classic agent can lead to a series of side effects in the course of treatment including cognitive impairment [3], hepatotoxicity [4], bone marrow toxicity [5], nephrotoxicity [6], and especially cardiotoxicity [7]. The undesirable effect of cardiotoxicity is especially prominent during the whole treatment period and is usually characterized by irreversible degenerative dilated cardiomyopathy (DCM) and consequently resulting in congestive heart failure (CHF), thus greatly worsening the long-term outcomes of patients [8]. Several mechanisms have been proposed to account for this phenomenon, including the generation of reactive oxygen species (ROS), activation of apoptosis, calcium dysfunction, as well as induction of endothelin-1 (ET-1) and topoisomerase-II [9]. However, understanding this pathological process is still in its infancy, and the specific mechanism remains to be elucidated.

The development of bioinformatics technology has brought great opportunities for fully understanding doxorubicin-induced cardiotoxicity. In the present study, we reanalyzed two datasets GSE23598 and GSE59672 from the GEO database to clarify the differentially expressed genes (DEGs) and the coexpression modules in the doxorubicin-mouse model. Several genes screened out from DEGs and modules were validated in other in vivo model or in vitro model datasets to screen a key gene and evaluate its diagnostic value. According to the analysis results, the hub gene was regarded as inflammation-related, and finally, its correlation with immune infiltration in diseased hearts was calculated.

2. Materials and Methods

2.1. Bioinformatic Datasets

By searching the Gene Expression Omnibus (GEO, https://www.ncbi.nlm.nih.gov/geo/) database, five in vivo experiment datasets GSE23598, GSE59672, GSE81448, GSE97642, and GSE37260 were obtained. In GSE23598, 4 male wild-type (WT) mice aged 6-8 weeks were randomly assigned to the doxorubicin group () and control group () following injection with a single dose of doxorubicin (15 mg/kg, i.p.) or an equal volume of normal saline, respectively. Animals were sacrificed on the 4th day after injection, and hearts were collected for further detection of gene expression profiling by microarray with platform GPL1261. In GSE59672, six male WT mice aged 10-12 weeks were randomly assigned to 2 groups ( in each group) and received similar treatment as the mice from the GSE23598 database, and were sacrificed on the 5th day. The gene expression profiling was performed with the same microarray platform GPL1261. In both GSE81448 and GSE97642, male mice aged 8-10 weeks were injected intraperitoneally with 15 mg/kg doxorubicin ( in each dataset) or an equal amount of PBS ( in each dataset). Mice were euthanized 20 hours after injection, and the expression profiles of the diseased heart were detected both by microarray with platform GPL6887. In GSE37260, rats were treated with doxorubicin () or saline () for 48 hours, and the peripheral blood mononuclear cell (PBMC) was used for detecting transcriptome expression with the bead chip platform GPL6101.

Furthermore, three in vitro experiment datasets GSE42177, GSE154101, and GSE154118 were enrolled. In GSE42177, cardiac myocytes (CMs) were treated with doxorubicin for 0 hr (control group, ) and 8 hrs (doxorubicin group, ). The gene expression profiling was detected by microarray on platform GPL6247. In GSE154101, fibroblasts (FBs) of human pulmonary arterial adventitial were isolated to treat with vehicle (control group, ) and doxorubicin (doxorubicin, ). In GSE154118, smooth muscle cells (SMCs) of the human pulmonary artery were treated with vehicle (control group, ) and doxorubicin (doxorubicin, ). These two datasets detected the expression profile by high-throughput sequencing with GPL18573.

2.2. Identification of Differentially Expressed Genes

Expression profiles of GSE23598, GSE59672, GSE81448, GSE97642, GSE37260, and GSE42177 were downloaded from the GEO database through the “GEOquery” package of R software. The probes corresponding to multiple targets were removed. When multiple probes corresponded to the same molecule, only the probe with the largest signal value was retained. After filtering the data, the combat function of the “SVA” package was used to remove the interbatch difference among GSE23598 and GSE59672 datasets, as well as GSE81448 and GSE97642. Then the standardization situation was viewed through the box chart, and the clustering of samples was checked with the PCA diagram and UMAP diagram. The gene expression profile of GSE37260 was performed with normalized quantiles. Finally, the “limma” package was used to analyze the differentially expressed genes (DEGs) between control and doxorubicin groups in datasets performed with microarray. Besides, DEGs of datasets GSE154101 and GSE154118 were analyzed in the well-known web tool NetworkAnalyst (https://www.networkanalyst.ca/) [10] with the DESeq2 method according to the instruction.

The screening criteria for DEGs are as follows: log₂ fold change (FC) >1 or < -1, and adjusted value (adj. p. val) <0.05. If log₂ FC>1, the expression of DEGs was deemed as upregulated, while < -1 was downregulated.

2.3. Weighted Gene Coexpression Network Analysis

Weighted gene coexpression network analysis (WGCNA) is a systematic biological method for counting the correlation patterns among genes across microarray profiles [11]. The online tool EHBIO platform (http://www.ehbio.com/Cloud_Platform/front/#/) was used to perform WGCNA on the combined dataset (GSE23598 and GSE59672). Modules most related to traits of the WGCNA were enrolled for the following study after the screening.

2.4. Search of Upstream miRNAs

The expression of genes is usually regulated by miRNAs at the posttranscriptional level [12]. Therefore, the upstream miRNAs of the DEGs in GSE23598 and GSE59672 were searched from the miRNET website (https://www.mirnet.ca/) [13]. The miRNA-mRNA networks were subsequently constructed from this web tool.

2.5. Protein-Protein Interaction Network Construction

The “Search Tool for Retrieval of Interacting Genes/Proteins” database (STRING, https://string-db.org/) provides abundant information on protein-protein interaction (PPI) [14], and those DEGs encoded proteins were incorporated into the PPI network with STRING. Then, network visualization was performed using Cytoscape software (https://cytoscape.org/), an open-source software platform [15].

2.6. Potential Therapeutic Drugs Prediction

The Connectivity Map (cMap) database can provide information about the effects of small molecules, genes, or diseases on gene expression signatures [16], and it was applied to search for potential drugs against doxorubicin-induced cardiotoxicity. The online application CLUE (https://clue.io) [17] provides a convenient entrance for the use of the cMAP database. The candidate DEGs were divided into upregulated group and downregulated group and imported into CLUE to calculate the score of all factors including small molecules. When the score is negative, it is considered that small molecules have antagonistic effects. Conversely, when the score is positive, small molecules will be considered to have the effect of exacerbating cardiotoxicity.

2.7. Functional Enrichment Analysis

Gene ontology (GO) analysis mainly includes three aspects: molecular function (MF), biological process (BP), and cellular component (CC) [18]. The Kyoto Encyclopedia of Genes and Genomes (KEGG) database provides biological signaling pathway information enriched by a group of genes [19]. The DAVID database (https://david.ncifcrf.gov/) [20], a comprehensive set of functional annotation tools for understanding the biological meaning of genes, was freely accessed and applied to explore the GO and KEGG results.

2.8. Identifying Biomarker and Hub Gene in Disease Progression

DEGs of the mouse model and the selected WGCNA modules were intersected. The overlapping genes were validated in the merged datasets (GSE81448 and GSE97642), and the three in vitro datasets (GSE42177, GSE154101, and GSE154118) and the rat PBMC dataset (GSE37260) in turn to pick out the hub gene. Gene Set Enrichment Analysis (GSEA) is a knowledge-based approach for analyzing groups of genes that share common biological functions, chromosomal locations, or regulations [21]. The GSEA was analyzed according to groups divided by hub gene with the Sangerbox tools (http://sangerbox.com/Tool), a free online platform. Online database GeneMANIA (http://genemania.org/), a user-friendly web tool for generating hypotheses about gene function, analyzing gene lists, and prioritizing genes for functional assays [22], was used to evaluate the biological function of the hub gene by searching its related genes. GraphPad Prism 8.0 was used to display the expression distribution and perform receiver operating characteristic (ROC) statistical analysis of the hub gene to evaluate its diagnostic efficacy.

2.9. Immune Infiltration Analysis

CIBERSORTx (https://cibersortx.stanford.edu/) [23] is an analytical tool that provides an estimation of the abundances of member cell types in a mixed cell population using a gene expression matrix. The combination of GSE23598 and GSE59672 was applied for analysis in the CIBERSORTx web tool to evaluate the percentages of 22 kinds of infiltrating immunocytes. Statistical analysis was done using Pearson and Wilcoxon test methods. was considered statistically significant. Finally, the correlation between hub gene expression level and scores of immune infiltration was calculated using the Pearson correlation coefficient in GraphPad Prism 8.0.

3. Results

3.1. Identification of DEGs between Control and Doxorubicin Mouse Model

The whole research process is summarized in Figure 1. In accordance with the study design, the batch differences between expression profiles in datasets GSE23598 and GSE59672 regarding wild-type mice were first removed, and the data quality before (Figures 2(a)–2(c)) and after (Figures 2(d)–2(f)) batch removal was assessed and displayed in form of PCA, UMAP, and boxplot. Totally, there were 20813 genes enrolled in the final expression matrix (Table 1). However, only 120 genes were considered to be differentially expressed between 5 control samples and 5 doxorubicin samples according to the criterion: and an adjusted value <0.05 (Figure 3(a)). The expression level of the top 20 upregulated and downregulated genes is shown in the heatmap (Figure 3(b)). Among these 120 DEGs, 54 were upregulated and 66 were downregulated. The number of DEGs under different criteria is represented in Table 1, and detailed information on the 120 DEGs is listed in Supplementary Table 1.

(a)

(b)

(c)

(d)

(e)

(f)

(a)

(b)

(c)

(d)

3.2. Network Analysis of DEGs

MiRNAs are vital regulators of mRNA expression and will form a tight regulatory system with these upregulated or downregulated DEGs. Therefore, the upstream miRNAs were searched with an online tool, and 334 miRNAs were finally revealed to target 109 of the 120 DEGs to form a complex miRNA-mRNA network with 1773 edges (Figure 3(c)). As shown, miR-155-5p, mir-122-5p, and miR-1a-3p have the most downstream target DEGs. Besides, the mRNA-encoded proteins could potentially interact with each other. The PPI network provides information of multiple levels of information regarding the interaction; therefore, the STRING database was used to construct the PPI network, and Cytoscape software was used to visualize the 65 nodes and 205 edges contained in the PPI network (Figure 3(d)).

3.3. Predicting Potential Therapeutic Drugs

The cMAP database provides information about the action of a series of small molecules on the gene expression profiles of several cell lines. By matching the 54 upregulated and 66 downregulated genes to the database, several small molecules including PF-04217903, propranolol, and azithromycin (Table 2) were highlighted as most likely to fight against doxorubicin-induced cardiotoxicity.

3.4. Construction of Weighted Gene Coexpression Network and Identify Key Modules

WGCNA was performed to identify the gene set and construct its connection with phenotype. Hierarchical clustering analysis was performed, and there existed no outlier (Figure 4(a)). Network topology analysis was used to determine candidate power values for relative, balanced scale independence, and mean connectivity in the WGCNA, and finally, the soft-thresholding power was selected as 18 after comprehensive consideration (Figure 4(b)). Totally, there were 8 coexpression modules screened out. The modules containing the most genes were the turquoise and blue modules (Figure 4(c)). The eigengene network indicated that there is good discrimination between different modules (Figure 4(d)).

(a)

(b)

(c)

(d)

A heatmap of module-trait correlations analysis showed that multiple modules were related to doxorubicin-induced cardiotoxicity, while MEturquoise had a maximum positive correlation coefficient and MEblue had a maximum negative correlation coefficient (Figure 5(a)). Finally, genes in the two modules were evaluated for their correlation with the trait (gene significance) and module (module membership), respectively, and displayed in a scatterplot (Figure 5(b)). There was a significant correlation between gene significance scores and module membership scores between the two modules.

(a)

(b)

3.5. Functional Enrichment Analysis of the Modules

To further understand the biological effects and signaling pathway enrichment of MEturquoise and MEblue, GO and KEGG analyses were performed. Genes in MEturquoise were enriched in “cytoplasm,” “nucleotide binding,” “cellular response to insulin stimulus,” etc. of the GO project (Figure 6(a)). “HIF-1 signaling pathway,” “FoxO signaling pathway,” etc. were key biological pathways involved by the MEturquoise (Figure 6(b)). In MEblue, genes were enriched in the “extracellular region,” “heparin-binding,” and “cell adhesion” (Figure 6(c)). Besides, KEGG analysis indicated that the biological pathway enriched in MEblue is different from that in MEturquoise, which included “focal adhesion,” “ECM-receptor interaction,” and “PI3K-Akt signaling pathway,” (Figure 6(d)).

(a)

(b)

(c)

(d)

3.6. Screen out Key Genes to Validate and Evaluate their Expression Distribution

The genes located in the upper right region of the MEturquoise scatterplot and the lower right of the MEblue scatterplot were regarded as potential vital factors. Therefore, genes with gene significance and module membership both scored more than 0.95 and were selected to evaluate their differential expression. As shown in Figure 7(a), 35 genes were selected out from MEturquoise, of which 5 were significantly upregulated, and 9 of 42 genes selected out from MEblue were significantly downregulated in doxorubicin-induced cardiotoxicity (GSE23598 and GSE59672). To assess the role of these 14 hub genes that were key in the progression of cardiotoxicity, we first assessed their expression in the merged dataset of GSE81448 and GSE97642. The PCA, UMAP, and boxplot analyses of these two datasets before and after the combination was displayed in Supplementary Figure 1. Finally, several genes including Limd1 were significantly dysregulated (), but the log₂ FC was all less than 1 (Figure 7(b)).

(a)

(b)

The cardiac tissue mainly consisted of CMs, FBs, SMCs, and immune cells. In order to explore the main effector cell population leading to changes in these 14 key genes, we reanalyzed the datasets (GSE42177, GSE154101, and GSE154118) that evaluate the expression matrix of doxorubicin-treated CMs, FBs, and SMCs, respectively. The analysis results of DEGs were displayed in Supplementary Figure 2 in the form of a volcano plot, and detailed results of the 14 genes including log₂ FC and adj. value were summarized in Table 3. According to Table 3, several genes were upregulated or downregulated, most likely contributed by immune cells. For example, Limd1 was significantly downregulated in CMs induced by doxorubicin but not significantly changed in CFs and SMCs. However, its expression was significantly upregulated in the organizational whole composed of the above elements, which indicates that inflammatory cells may have led to this change trend.

3.7. Limd1 Was a Potential Valuable Circulating Biomarker in the Cardiotoxicity

Due to that the above key genes may be dysregulated in the immune cells, we further evaluated their expression patterns of them in the circulating PBMC of the doxorubicin-related cardiotoxicity model (GSE37260). Among them, the expression of 12 genes was detected, while Evi2a and Mgl2 were nonexistent (Figure 8(a)). Further, Limd1, Spock2, Pcolce, Slamf9, Col15a1, and Nrep were consistent with those in the heart of the mouse model in terms of differential expression trend, while only the dysregulation of Limd1 was significant (). Limd1 was upregulated in the PBMC of the rat treated with doxorubicin (Figure 8(b)), and the area under the ROC curve for Limd1 in diagnosing the cardiotoxicity was 0.8472 (Figure 8(c)), suggesting that Limd1 is a valuable circulating biomarker to judge this diseased state.

(a)

(b)

(c)

3.8. GSEA of Limd1-Associated Gene Set

In order to explicit the biological characteristic involved by Limd1, we ranked genes from the 10 normal and diseased heart samples in the merged dataset by their relative Limd1 expression in the top 50% vs. the bottom 50% for GSEA analysis based on the KEGG database. There were dozens of categories enriched ( but ), and several are functionally closely related to inflammatory cells including the chemokine signaling pathway, cell adhesion molecules pathway, and leukocyte transendothelial migration signaling (Figure 9(a)). These pointed out that the key gene Limd1 was a possible regulator in controlling immunocyte infiltrating tissue.

(a)

(b)

To further explore the potential biological function of Limd1, its interaction protein network was constructed, and the result showed that Limd1 may participate in inflammatory regulation through the NF-kappaB signaling pathway (Figure 9(b)).

3.9. Assessment of Immune Infiltration and the Correlation with Limd1 Expression

The primitive unlog transform microarray matrix regarding two mouse models was merged and removed batch for the subsequent immune infiltration assessment (Figure 10(a)). After being calculated in CIBERSORTx, the composition of immunocytes in the heart tissues of the doxorubicin-induced mouse model was indicated intuitively (Figure 10(b)). The heatmap displayed the correlation between one immunocyte with another (Figure 10(c)), and the results indicated that there exists a positive correlation between “activated dendritic cells” with “regulatory T cells (Tregs)” and “macrophages M2”, and a negative correlation between “activated dendritic cells” with “macrophages M1.”

(a)

(b)

(c)

The percentage data of each immunocyte type between the control and the doxorubicin groups were compared with the Wilcoxon test. As shown in the violin plot, the proportion of “dendritic cells activated” was elevated (), while “macrophage M1’ declined () in the doxorubicin group compared with the control group; the others showed no significance (Figure 11(a)). The correlation between Limd1 expression and the fraction of immunocytes was evaluated. Totally, 3 kinds of immunocytes were found to correlate with the Limd1 level. Among them, “activated dendritic cells” and “macrophage M1” were positively and negatively correlated with Limd1 level, respectively (Figure 11(b)). In addition, although “monocytes” did not change significantly () under doxorubicin, it was negatively correlated with Limd1 level (Figures 11(b) and 11(c)).

(a)

(b)

(c)

4. Discussion

As a widely clinically used drug, doxorubicin has a very powerful antitumor function. However, a growing number of studies support the conception that doxorubicin plays a double-edged sword role during its application. Indeed, the action of doxorubicin on nontargeted tissues may result in cardiotoxicity and eventually worsen to congestive heart failure and death, which will disrupt cancer treatment by needing to control medical dosages of doxorubicin, thus deteriorating the quality of life of patients [24]. Therefore, exploring new strategies against the cardiotoxic effect of doxorubicin would benefit cancer patients to a certain extent.

In this study, we selected two datasets GSE23598 and GSE59672 from the GEO database as the subjects. These two datasets were based on a doxorubicin-induced acute heart injury mouse model performed through intraperitoneal administration of doxorubicin in a dose of 15 mg/kg for 4 or 5 days [25, 26]. After combining and reanalyzing the datasets, we screened out a total of 120 DEGs between the normal myocardial tissue and the toxic myocardial tissue, of which parts were confirmed by previous reports in terms of the expression trend, such as Ctgf [27], Aox1 [28], Alas2 [29], and Ly86 [30]. However, Btg2 was upregulated after doxorubicin treatment in our study, while a previous study found it to be downregulated [31]. On the whole, the differentially expressed genes obtained based on omics were relatively reliable. Furthermore, Cxcl9 was discovered to be downregulated in the heart tissue of the mouse model. Its protein level in plasma was evaluated in breast cancer patients who received doxorubicin, and it was significantly lower in the abnormal decline of the left ventricular ejection fraction (LVEF) group compared with the normal LVEF group [32], indicating that these DEGs are potential circulating biomarkers for clinical application of doxorubicin-related cardiotoxicity in the future. Besides, based on these DEGs, the potential therapeutic drugs against the cardiotoxicity were predicted, and a series of molecules including propranolol were selected. Interestingly, propranolol was reported to have a significant beneficial effect on cardiac injury induced by doxorubicin [33]. Therefore, it is worth expecting whether the other drugs, especially PF-04217903, could relieve this toxicity in the myocardium.

As is well known, miRNAs usually negatively participate in regulating the protein expression of genes through degradation or translational inhibition. We therefore further predicted the upstream miRNAs that target all the 120 DEGs and finally constructed the network containing 334 miRNAs and 109 DEGs targeted by them. In this whole network, miR-155-5p, miR-122-5p, miR-1a-3p, etc. were most closely related to these DEGs. As revealed in the clinical research, circulating miR-122-5p was increased after doxorubicin application in breast cancer patients [34]. In the mouse model, the circulating level of miR-122-5p was higher in cardiotoxicity than in the nontoxicity group after doxorubicin was applied. However, the circulating miR-122-5p was lower in the noncardiotoxicity model compared with the saline-treated model [35]. Besides, the level of circulating miR-122-5p before treatment in patients with myocardial injury after treatment was significantly higher than that in patients without myocardial injury after treatment and could predict adverse cardiac reactions to doxorubicin [34]. Therefore, whether the other miRNAs especially miR-155-5p and miR-1a-3p have clinical application potential is worth exploring.

Genes are usually synergistically involved in disease phenotypes. Accordingly, we extracted the coexpression modules from the whole gene matrix using the WGCNA method and finally determined 8 modules. Among them, the turquoise module (MEturquoise) was most positively correlated with the cardiotoxicity abnormal phenotype, containing 857 genes, while the blue module (MEblue) was most negatively correlated, containing 778 genes. Besides, we performed GO and KEGG analyses, respectively, to fully evaluate the biological roles of the two modules. The biological annotation of GO item demonstrated that those genes in MEturquoise were related to “cytoplasm,” “nucleotide binding,” “cellular response to insulin stimulus,” etc., and genes in MEblue were related to “extracellular region,” “heparin-binding,” and “cell adhesion.” KEGG analysis indicated several biological signaling pathways were significantly enriched in MEturquoise, such as the “FoxO signaling pathway,” which has been validated by a previous study [36]. However, in MEblue, other types of signaling, such as “focal adhesion,” “ECM receptor interaction,” and “PI3K-Akt signaling pathway,” were significantly enriched, and several of them have also been identified involvement in cardiotoxicity [37, 38]. Therefore, these pathways especially those not previously reported, greatly enrich the knowledge regarding the disease and provide potential intervention targets.

The better the correlation between modules and disease traits, the more important function genes have. On these grounds, 35 genes from MEturquoise and 42 genes from MEblue were screened out. The significantly upregulated DEGs in 35 genes from MEturquoise and the significantly downregulated DEGs in 42 genes from MEblue were considered the key ones, and finally, 14 genes were enrolled. To validate the role of these 14 key genes, we first merged and reanalyzed two datasets GSE81448 and GSE97642, which were modeled with the same dose of doxorubicin, but the detection time was advanced to 20 hours after the operation. Limd1, Spock2, Sult1a1, etc. were also discovered to be dysregulated, the same as the previous results. To identify the specific cell populations that lead to changes in the expression of these genes, we reanalyzed another three datasets that detected the CMs, FBs, and SMCs under the intervention of doxorubicin. As the results revealed, dysregulation of several genes including Limd1 and Sult1a1 was not from the above three kinds of cells, and potentially resulted from the inflammatory cells infiltrating into tissue. Therefore, we further evaluate their expression in the PBMC of the rat cardiotoxicity model. In the end, only Limd1 was significantly upregulated, the same as in the heart of the mouse model. Besides, ROC analysis indicated that Limd1 in PBMC had the ability to diagnose this cardiotoxic disease. As is well known, Limd1 (LIM domain-containing protein 1) is a member of the Zyxin proteins and is widely expressed in human tissues. Several previous studies have suggested the suppressor role of Limd1 in tumor diseases [39]. In cardiovascular disease, Limd1 has been reported to be functionally closely related to cardiac fibroblasts [40]. Based on the results, Limd1 was deemed as the hub gene and a valuable biomarker.

Due to that Limd1 was closely related to inflammatory cells in our study, we subsequently determined the biological function of Limd1 by the GSEA tool and protein interactive network. GSEA results of the mouse model suggested that Limd1 may be functionally involved in immunocyte chemotaxis, adhesion, and migration. In the network built by the GeneMANIA database, the Limd1 was depicted to be connected with the NF-kappaB signaling pathway through interacting with Limd1, EGFR, and TRAF6 proteins. As one of the best-understood immune-related pathways, the NF-kappaB signaling pathway is activated by numerous discrete stimuli and participates in regulating activities of the majority of immunocytes including macrophage, dendritic cells, and neutrophils [41], thus exerting a comprehensive role in inflammation.

In fact, the cardiotoxicity induced by doxorubicin is in general characterized by an abnormal inflammatory response, and the NF-kappaB signaling pathway was indeed the mediator [42]. Due to that Limd1 was a highly likely regulator of immune cell infiltrates in the toxic heart according to the above analysis, we therefore finally assessed the immune infiltration in the hearts of mouse model on the 4th day after doxorubicin injection using the CIBERSORTx method. Compared with the hearts in the control group, the doxorubicin-intervened hearts have a significantly higher proportion of “dendritic cells activated” and a lower proportion of “macrophage M1”. Besides, the proportion of monocytes seems to be decreased in the cardiotoxicity model heart, despite not being statistically significant. As a widely studied inflammatory cell type, the proinflammatory M1 macrophage has been found to be increased in doxorubicin-induced cardiomyopathy [43]. The contradiction was due to that the datasets we reanalyzed used a one-time high-dose modeling method, and the literature used intermittent multiple-dosing ways. The one-time modeling indicates that the damage reaches the peak shortly after the beginning, and injury stimulation cannot be maintained at the following. This damage mode was similar to that in other heart diseases, such as acute myocardial infarction (AMI). During the progression of AMI, the response of the monocyte/macrophage system is characterized by the accumulation of proinflammatory monocytes and macrophages over 48 to 72 h, followed by a reparative phase at 4 to 7 days driven by anti-inflammatory macrophages [44]. Despite the M1 macrophage (pro-) was rapidly accumulated and peaking on the fifth day at the damaged site, the ratio of M1 type to the total macrophages (M1 and M2) had decreased significantly earlier than the fifth day [45]. Therefore, our results that the proportion of monocyte/macrophage M1 decreased could be clued from the pattern of myocardial infarction.

Dendritic cells are a kind of bone marrow-derived cells arising from lympho-myeloid hematopoiesis that are responsible for initiating and controlling immune responses, specializing in antigen presentation to drive T-cell priming and differentiation [46]. In research decades ago, the number of interstitial dendritic cells/mm² of a left ventricle section was found to be remarkably increased in animals receiving doxorubicin, compared with saline-treated control [47], which was similar to our analysis. As a matter of fact, dendritic cells are widely involved in the progression of heart diseases. In the mouse model of transverse aortic constriction (TAC)-induced pressure overload, the relative abundance of dendritic cells was found to be increased in the hypertrophic myocardium both at 1 and at 4 weeks post-TAC [48]. Besides, dendritic cells worked as a pathogen to induce autoimmune heart failure under certain conditions [49]. In an acute myocardial infarction mouse model, subcutaneously administration of tolerogenic dendritic cells (tDCs) could reduce infarct size, improve heart function and mouse survival, by timely promoting Treg cell, and activating its mediated macrophage conversion to the reparative M2 to replace the inflammation M1 type [45]. Interestingly, the fraction of activated dendritic cells was positively correlated with that of Tregs and M2 macrophage and negatively correlated with M1 macrophage, as revealed in our results. These suggest a potential pattern that the increasing of dendritic cells may exert its influence on the cardiotoxicity through the subsequent Treg cell and macrophage. Given the biological role of the NF-kappaB signaling pathway in monocyte, macrophage, and dendritic cells, we, therefore, wonder whether Limd1, a potential immunocyte regulator, is associated with the infiltration of these inflammatory cells. Finally, we found that Limd1 expression was negatively correlated with the fraction of monocyte and macrophage and positively correlated with the dendritic cells fraction. Thus, we speculate that Limd1 acts as the regulator of the inflammation system through the NF-kappaB pathway, playing a role in doxorubicin-induced cardiotoxicity.

There are still some limitations to this study. Firstly, the differential expression analyses compared the treatment () and control () groups in the GSE154101 and GSE154118 datasets could lead to some inaccuracy due to the limitations of the sample size. Secondly, the number of sample used for the WGCNA () could also lead to some bias due to the small sample size (because samples is suggested for the WGCNA). Thirdly, the above results were not further validated in our own disease model.

5. Conclusions

In conclusion, Limd1 was upregulated in PBMC and the heart of the doxorubicin-induced cardiotoxicity model and is a valuable circulating biomarker. In particular, we found that Limd1 was functionally correlated with the immune system and significantly correlated with the infiltration of immunocytes including monocytes, M1 macrophage, and dendritic cell. It is worth further exploring to clarify its role in this cardiotoxicity.

Data Availability

The data used to support the findings of this study are available from the first author upon reasonable request.

Conflicts of Interest

The authors declare that there is no conflict of interest regarding this paper.

Authors’ Contributions

Rui Zhang and Chunshu Hao contributed equally to this work.

Acknowledgments

This study was supported by the National Natural Science Foundation of China (grant number 81870213) and Jiangsu Provincial Key Medical Discipline (Laboratory ZDXKA2016023).

Supplementary Materials

Supplementary Table 1. The detailed information of 120 DEGs. Supplementary Figure 1 Remove batch between GSE81448 and GSE9764. Supplementary Figure 2 Volcano plots displayed the upregulated and downregulated genes in doxorubicin-treated CMs, Fibroblasts, and SMCs. (Supplementary Materials)

References

B. Kalyanaraman, “Teaching the basics of the mechanism of doxorubicin-induced cardiotoxicity: have we been barking up the wrong tree?” Redox Biology, vol. 29, pp. 101394–101394, 2020.
View at: Publisher Site | Google Scholar
B. Ongnok, N. Chattipakorn, and S. C. Chattipakorn, “Doxorubicin and cisplatin induced cognitive impairment: the possible mechanisms and interventions,” Experimental Neurology, vol. 324, article 113118, 2020.
View at: Publisher Site | Google Scholar
P. L. Prasanna, K. Renu, and A. Valsala Gopalakrishnan, “New molecular and biochemical insights of doxorubicin-induced hepatotoxicity,” Life Sciences, vol. 250, article 117599, 2020.
View at: Publisher Site | Google Scholar
M. Baxter-Holland and C. R. Dass, “Doxorubicin, mesenchymal stem cell toxicity and antitumour activity: implications for clinical use,” Journal of Pharmacy and Pharmacology, vol. 70, no. 3, pp. 320–327, 2018.
View at: Publisher Site | Google Scholar
A. Li, W. Zhang, L. Zhang et al., “Elucidating the time-dependent changes in the urinary metabolome under doxorubicin-induced nephrotoxicity,” Toxicology Letters, vol. 319, pp. 204–212, 2020.
View at: Publisher Site | Google Scholar
K. N. Timm and D. J. Tyler, “The role of AMPK activation for cardioprotection in doxorubicin-induced cardiotoxicity,” Cardiovascular Drugs and Therapy, vol. 34, no. 2, pp. 255–269, 2020.
View at: Publisher Site | Google Scholar
H. G. Fa, W. G. Chang, X. J. Zhang, D. D. Xiao, and J. X. Wang, “Noncoding RNAs in doxorubicin-induced cardiotoxicity and their potential as biomarkers and therapeutic targets,” Acta Pharmacologica Sinica, vol. 42, no. 4, pp. 499–507, 2021.
View at: Publisher Site | Google Scholar
N. Wenningmann, M. Knapp, A. Ande, T. R. Vaidya, and S. Ait-Oudhia, “Insights into doxorubicin-induced cardiotoxicity: molecular mechanisms, preventive strategies, and early monitoring,” Molecular Pharmacology, vol. 96, no. 2, pp. 219–232, 2019.
View at: Publisher Site | Google Scholar
G. Zhou, O. Soufan, J. Ewald, R. E. W. Hancock, N. Basu, and J. Xia, “NetworkAnalyst 3.0: a visual analytics platform for comprehensive gene expression profiling and meta-analysis,” Nucleic Acids Research, vol. 47, no. W1, pp. W234–w241, 2019.
View at: Publisher Site | Google Scholar
P. Langfelder and S. Horvath, “WGCNA: an R package for weighted correlation network analysis,” BMC Bioinformatics, vol. 9, no. 1, p. 559, 2008.
View at: Publisher Site | Google Scholar
M. Correia de Sousa, M. Gjorgjieva, D. Dolicka, C. Sobolewski, and M. Foti, “Deciphering miRNAs' action through miRNA editing,” International Journal of Molecular Sciences, vol. 20, no. 24, p. 6249, 2019.
View at: Publisher Site | Google Scholar
Y. Fan, K. Siklenka, S. K. Arora, P. Ribeiro, S. Kimmins, and J. Xia, “miRNet - dissecting miRNA-target interactions and functional associations through network-based visual analysis,” Nucleic Acids Research, vol. 44, no. W1, pp. W135–W141, 2016.
View at: Publisher Site | Google Scholar
D. Szklarczyk, A. L. Gable, D. Lyon et al., “STRING v11: protein-protein association networks with increased coverage, supporting functional discovery in genome-wide experimental datasets,” Nucleic Acids Research, vol. 47, no. D1, pp. D607–d613, 2019.
View at: Publisher Site | Google Scholar
P. Shannon, A. Markiel, O. Ozier et al., “Cytoscape: a software environment for integrated models of biomolecular interaction networks,” Genome Research, vol. 13, no. 11, pp. 2498–2504, 2003.
View at: Publisher Site | Google Scholar
J. Lamb, E. D. Crawford, D. Peck et al., “The Connectivity Map: Using Gene-Expression Signatures to Connect Small Molecules, Genes, and Disease,” Science, vol. 313, no. 5795, pp. 1929–1935, 2006.
View at: Publisher Site | Google Scholar
A. Subramanian, R. Narayan, S. M. Corsello et al., “A next generation connectivity map: L1000 platform and the first 1,000,000 profiles,” Cell, vol. 171, no. 6, pp. 1437–1452.e17, 2017.
View at: Publisher Site | Google Scholar
M. Ashburner, C. A. Ball, J. A. Blake et al., “Gene ontology: tool for the unification of biology,” Nature Genetics, vol. 25, no. 1, pp. 25–29, 2000.
View at: Publisher Site | Google Scholar
E. Altermann and T. R. Klaenhammer, “PathwayVoyager: pathway mapping using the Kyoto encyclopedia of genes and genomes (KEGG) database,” BMC Genomics, vol. 6, no. 1, p. 60, 2005.
View at: Publisher Site | Google Scholar
X. Jiao, B. T. Sherman, D. W. Huang et al., “DAVID-WS: a stateful web service to facilitate gene/protein list analysis,” Bioinformatics, vol. 28, no. 13, pp. 1805-1806, 2012.
View at: Publisher Site | Google Scholar
A. Subramanian, P. Tamayo, V. K. Mootha et al., “Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles,” Proceedings of the National Academy of Sciences of the United States of America, vol. 102, no. 43, pp. 15545–15550, 2005.
View at: Publisher Site | Google Scholar
D. Warde-Farley, S. L. Donaldson, O. Comes et al., “The GeneMANIA prediction server: biological network integration for gene prioritization and predicting gene function,” Nucleic Acids Research, vol. 38, suppl_2, pp. W214–W220, 2010.
View at: Publisher Site | Google Scholar
A. M. Newman, C. B. Steen, C. L. Liu et al., “Determining cell type abundance and expression from bulk tissues with digital cytometry,” Nature Biotechnology, vol. 37, no. 7, pp. 773–782, 2019.
View at: Publisher Site | Google Scholar
C. Carvalho, R. X. Santos, S. Cardoso et al., “Doxorubicin: the good, the bad and the ugly effect,” Current Medicinal Chemistry, vol. 16, no. 25, pp. 3267–3285, 2009.
View at: Publisher Site | Google Scholar
Y. Shuai, J. Guo, Y. Dong et al., “Global gene expression profiles of MT knockout and wild-type mice in the condition of doxorubicin-induced cardiomyopathy,” Toxicology Letters, vol. 200, no. 1-2, pp. 77–87, 2011.
View at: Publisher Site | Google Scholar
L. Wang, T. P. Zhang, Y. Zhang et al., “Protection against doxorubicin-induced myocardial dysfunction in mice by cardiac-specific expression of carboxyl terminus of hsp70-interacting protein,” Scientific Reports, vol. 6, no. 1, p. 28399, 2016.
View at: Publisher Site | Google Scholar
C. G. Tocchetti, A. Carpi, C. Coppola et al., “Ranolazine protects from doxorubicin-induced oxidative stress and cardiac dysfunction,” European Journal of Heart Failure, vol. 16, no. 4, pp. 358–366, 2014.
View at: Publisher Site | Google Scholar
M. Schwebe, S. Ameling, E. Hammer et al., “Protective effects of endothelin receptor A and B inhibitors against doxorubicin-induced cardiomyopathy,” Biochemical Pharmacology, vol. 94, no. 2, pp. 109–129, 2015.
View at: Publisher Site | Google Scholar
Z. Wang, J. Gao, H. Teng, and J. Peng, “Effects of Doxorubicin on Heme Biosynthesis and Metabolism in Cardiomyocyte,” Arquivos Brasileiros de Cardiologia, vol. 116, pp. 315–322, 2021.
View at: Publisher Site | Google Scholar
Y. J. Zhang, H. Huang, Y. Liu, B. Kong, and G. Wang, “MD-1 deficiency accelerates myocardial inflammation and apoptosis in doxorubicin-induced cardiotoxicity by activating the TLR4/MAPKs/nuclear factor kappa B (NF-κB) signaling pathway,” Medical Science Monitor, vol. 25, pp. 7898–7907, 2019.
View at: Publisher Site | Google Scholar
Z. Tong, B. Jiang, Y. Wu et al., “MiR-21 protected cardiomyocytes against doxorubicin-induced apoptosis by targeting BTG2,” International Journal of Molecular Sciences, vol. 16, no. 12, pp. 14511–14525, 2015.
View at: Publisher Site | Google Scholar
L. R. Yu, Z. Cao, I. Makhoul et al., “Immune response proteins as predictive biomarkers of doxorubicin-induced cardiotoxicity in breast cancer patients,” Experimental Biology and Medicine, vol. 243, no. 3, pp. 248–255, 2018.
View at: Publisher Site | Google Scholar
H. Kawabata, T. Ryomoto, and K. Ishikawa, “Effect of beta-blocker on metabolism and contraction of doxorubicin-induced cardiotoxicity in the isolated perfused rabbit heart,” Angiology, vol. 51, no. 5, pp. 405–413, 2000.
View at: Google Scholar
S. Gioffré, M. Chiesa, D. M. Cardinale et al., “Circulating MicroRNAs as Potential Predictors of Anthracycline-Induced Troponin Elevation in Breast Cancer Patients: Diverging Effects of Doxorubicin and Epirubicin,” Journal of Clinical Medicine, vol. 9, no. 5, p. 1418, 2020.
View at: Publisher Site | Google Scholar
C. Ruggeri, S. Gioffré, M. Chiesa et al., “A Specific Circulating MicroRNA Cluster Is Associated to Late Differential Cardiac Response to Doxorubicin-Induced Cardiotoxicity In Vivo,” vol. 2018, Article ID 8395651, 2018.
View at: Publisher Site | Google Scholar
A. N. Kavazis, A. J. Smuder, and S. K. Powers, “Effects of short-term endurance exercise training on acute doxorubicin-induced FoxO transcription in cardiac and skeletal muscle,” Journal of Applied Physiology, vol. 117, no. 2014, pp. 223–230, 1985.
View at: Publisher Site | Google Scholar
Z. Cheng, L. A. DiMichele, M. Rojas, C. Vaziri, C. P. Mack, and J. M. Taylor, “Focal adhesion kinase antagonizes doxorubicin cardiotoxicity via p21^Cip1,” Journal of Molecular and Cellular Cardiology, vol. 67, pp. 1–11, 2014.
View at: Publisher Site | Google Scholar
A. Kalantary-Charvadeh, D. Sanajou, M. Hemmati-Dinarvand et al., “Micheliolide Protects Against Doxorubicin-Induced Cardiotoxicity in Mice by Regulating PI3K/Akt/NF-kB Signaling Pathway,” Cardiovascular Toxicology, vol. 19, pp. 297–305, 2019.
View at: Publisher Site | Google Scholar
Z.-Z. Guo, Z.-J. Ma, Y.-Z. He et al., “miR-550a-5p functions as a tumor promoter by targeting LIMD1 in lung adenocarcinoma,” Frontiers in Oncology, vol. 10, pp. 570733–570733, 2020.
View at: Publisher Site | Google Scholar
N. M. Landry, S. G. Rattan, K. L. Filomeno et al., “SKI activates the Hippo pathway via LIMD1 to inhibit cardiac fibroblast activation,” Basic Research in Cardiology, vol. 116, p. 25, 2021.
View at: Publisher Site | Google Scholar
M. G. Dorrington and I. D. C. Fraser, “NF-κB signaling in macrophages: dynamics, crosstalk, and signal integration,” Frontiers in Immunology, vol. 10, p. 705, 2019.
View at: Publisher Site | Google Scholar
F. Yarmohammadi, H. Karbasforooshan, A. W. Hayes, and G. Karimi, “Inflammation suppression in doxorubicin-induced cardiotoxicity: natural compounds as therapeutic options,” Naunyn-Schmiedeberg's Archives of Pharmacology, vol. 394, no. 10, pp. 2003–2011, 2021.
View at: Publisher Site | Google Scholar
T. A. Johnson and D. K. Singla, “PTEN inhibitor VO-OHpic attenuates inflammatory M1 macrophages and cardiac remodeling in doxorubicin-induced cardiomyopathy,” American Journal of Physiology. Heart and Circulatory Physiology, vol. 315, no. 5, pp. H1236–h1249, 2018.
View at: Publisher Site | Google Scholar
I. Andreadou, H. A. Cabrera-Fuentes, Y. Devaux et al., “Immune cells as targets for cardioprotection: new players and novel therapeutic opportunities,” Cardiovascular Research, vol. 115, no. 7, pp. 1117–1130, 2019.
View at: Publisher Site | Google Scholar
E. H. Choo, J. H. Lee, E. H. Park et al., “Infarcted myocardium-primed dendritic cells improve remodeling and cardiac function after myocardial infarction by modulating the regulatory T cell and macrophage polarization,” Circulation, vol. 135, no. 15, pp. 1444–1457, 2017.
View at: Publisher Site | Google Scholar
M. Collin and V. Bigley, “Human dendritic cell subsets: an update,” Immunology, vol. 154, no. 1, pp. 3–20, 2018.
View at: Publisher Site | Google Scholar
J. Zhang, E. H. Herman, and V. J. Ferrans, “Dendritic cells in the hearts of spontaneously hypertensive rats treated with doxorubicin with or without ICRF-187,” The American Journal of Pathology, vol. 142, no. 6, pp. 1916–1926, 1993.
View at: Google Scholar
E. Martini, P. Kunderfranco, C. Peano et al., “Single-cell sequencing of mouse heart immune infiltrate in pressure overload-driven heart failure reveals extent of immune activation,” Circulation, vol. 140, no. 25, pp. 2089–2107, 2019.
View at: Publisher Site | Google Scholar
U. Eriksson, R. Ricci, L. Hunziker et al., “Dendritic cell-induced autoimmune heart failure requires cooperation between adaptive and innate immunity,” Nature Medicine, vol. 9, no. 12, pp. 1484–1490, 2003.
View at: Publisher Site | Google Scholar

Copyright

Copyright © 2023 Rui Zhang et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

PDF Download Citation

Download other formats

Order printed copies

Views

447

Downloads

455

Citations