Implication of proliferation gene biomarkers in pulmonary hypertension

Abstract Objective/Background Proliferation is a widely recognized trigger for pulmonary hypertension (PH), a life‐threatening, progressive disorder of pulmonary blood vessels. This study was aimed to identify some proliferation associated genes/targets for better comprehension of PH pathogenesis. Methods Human pulmonary arterial smooth muscle cells (hPASMCs) were cultured in the presence or absence of human recombinant platelet derived growth factor (rhPDGF)‐BB. Cells were collected for metabolomics or transcriptomics study. Gene profiling of lungs of PH rats after hypoxia exposure or of PH patients were retrieved from GEO database. Results 90 metabolites (VIP score >1, fold change >2 or <0.5 and p < .05) and 2701 unique metabolism associated genes (MAGs) were identified in rhPDGF‐BB treated hPASMCs compared to control cells. In addition, 1151 differentially expressed genes (313 upregulated and 838 downregulated) were identified in rhPDGF‐BB treated hPASMCs compared to control cells (fold change >2 or <0.5 and p < .05). 152 differentially expressed MAGs were then determined, out of which 9 hub genes (IL6, CXCL8, CCL2, CXCR4, CCND1, PLAUR, PLAU, HBEGF and F3) were defined as core proliferation associated hub genes in protein proten interaction analysis. In addition, the hub gene‐based LASSO model can predict the occurrence of PH (AUC = 0.88). The expression of CXCR4, as one of the hub genes, was positively correlated to immune cell infiltrates. Conclusion Our findings revealed some key proliferation associated genes in PH, which provide the crucial information concerning complex metabolic reprogramming and inflammatory modulation in response to proliferation signals and might offer therapeutic gains for PH.


| INTRODUC TI ON
Pulmonary hypertension (PH) is a vicious cardio-pulmonary disorder and manifested by progressive increase in pulmonary artery pressure and pulmonary vascular resistance, which ultimately leads to right heart failure and even death. 1 Multiple factors including genetic predisposing genes, 2,3 epigenetic modulations, 4,5 inflammation, 6 altered metabolism 7,8 and environment insults such as hypoxia 9 are reported to cause the remodeling of the pulmonary vasculature, as manifested by overproliferation, anti-apoptosis and high migratory capability of vascular cells. 10 According to WHO classification, PH is categorized into 5 main groups, with pulmonary arterial hypertension (PAH) to be Group 1 PH. 11 As of today, most of therapies are targeted against PAH, including endothelin receptor antagonists, phosphodiesterase type 5 inhibitors, and prostacyclin analogues. 12 These therapies help to improve exercise capacity, hemodynamics and quality of life. However, none of the current treatments are actually curative and long-term prognosis still remains poor.
As proliferation is a critical trigger for PH development, the factors/elements in regulation of proliferation has emerging as a research focus against pulmonary vascular remodeling. For example, low-dose FK506 reverses PH in rats with medial hypertrophy following monocrotaline and neointima formation of Sugen5416/ Hypoxia PH rats (a severe PH model) via rescue of bone morphogenetic protein receptor-2 (BMPR2) signaling dysfunction, 13 which is widely recognized in circulating system and lung tissues of PH patients. 14,15 At the mechanism level, the induction of BMPR2 signaling by FK506 is not only associated with the reversal of endothelial cells (ECs) dysfunction, but also with the inhibition of pulmonary arterial smooth muscle cells (PASMCs) proliferation, as evident by a direct effect of FK506 in inhibiting platelet derived growth factor (PDGF) induced proliferation. PDGF is regarded as the most potent mitogen for PASMCs 16 and the subtype PDGF-BB (dimeric form of PDGF-B) is widely applied as a stimuli for PASMCs functional study in PH. Moreover, a higher PDGF-B and its receptor PDGFR-beta mRNA expression was observed in small pulmonary arteries from patients with idiopathic PAH, and both PDGF-B and PDGFR-beta were localized in PASMCs in small remodeled pulmonary arteries. 17 The elevated PDGF signaling contributes to the pathobiology of PAH 18 by regulation of proliferation and migration of pulmonary vascular SMCs, and its inhibition abolishes pulmonary vascular remodeling in two PH models. 19 Therefore, a better understanding into PDGF-BB mediated gene/metabolite profiling would aid in the discovery of novel target to reverse vascular remodeling. Systemic or local inflammation is another feature of PH. 20,21 The proliferation associated alteration in lung tissues and its impact on immune cell infiltration would also shed some light on the pathogenesis of PH.
Here, we sought to scrutinize the gene/metabolite alterations in PASMCs in response to PDGF-BB, identify some proliferation associated genes/targets by virtue of metabolomics and transcriptomics and unveil the potential link between proliferation and inflammation, which might provide new approaches for the treatment of PH.

| Cell culture and sample collection
Human pulmonary artery smooth muscle cells (hPASMCs) were purchased from the American Type Culture Collection (ATCC; Cat#. PCS-100-023) and cultured at 37°C in the incubator containing 5% For metabolomic study, hPASMCs were seeded into 6-well plates at a density of 0.5 × 10 6 /ml. Twenty-four hours later cells were starved with DMEM/F-12 containing 0.5% fetal bovine serum for another 24 h. Cells were gently washed with PBS followed by addition of fresh blank medium with recombinant human PDGF-BB (rhPDGF-BB) (R& D systems; Cat#. 220-BB) at the concentration of 20 ng/ml or vehicle for 24 h (n = 6/group). Cells were collected, washed, centrifuged at 300 g at 4°C for 10 min and resuspended with PBS to reach a final density of 1 × 10 5 /ml. One milliliter cell suspension per sample was centrifugated at 300 g at 4°C for 5 min.
Supernatants were discarded and cell pellets were stored at −80°C.
Freezing and thawing cycle was avoided before sample processing in case of potential degradation of metabolites.
The expression of CXCR4, as one of the hub genes, was positively correlated to immune cell infiltrates.

Conclusion:
Our findings revealed some key proliferation associated genes in PH, which provide the crucial information concerning complex metabolic reprogramming and inflammatory modulation in response to proliferation signals and might offer therapeutic gains for PH.

K E Y W O R D S
metabolism associated genes, metabolomics, proliferation, pulmonary hypertension, transcriptomics For bulk RNA sequencing, hPASMCs were seeded into 6-well plates and starved as aforementioned. Cells were gently washed with PBS followed by addition of fresh blank medium with recombinant human PDGF-BB (rhPDGF-BB) (R& D systems; Cat#. 220-BB) at the concentration of 20 ng/ml or vehicle BB for 6 h (n = 3/group) before collection. Cells were then washed, centrifuged at 300 g at 4°C for 10 min. Supernatants were discarded and cell pellets were stored at −80°C before further use.

| LC/MS analysis and data process
One milliliter of MeOH:ACN:H 2 O (2:2:1, v/v) solvent mixture was added into the cell samples followed by sonification. Cells were prepared as previously described. 7 Liquid chromatographic separation for processed samples was achieved on a ZORBAX Eclipse Plus C18 column (2.1 × 100 mm, 3.5 μm, Agilent, USA) maintained at 45°C, whereas mass spectrometry was performed on a Nexera X2 system (Shimadzu, Japan) coupled with a Triple TOF 5600 quadrupole-timeof-flight mass spectrometer (AB SCIEX, USA). All the steps for LC/ MS analysis and data preprocessing were depicted in our previous study. 7 We used partial least squares discriminant analysis (PLS-DA) to distinguish the overall difference in metabolic profile between rhPDGF-BB treated-and control hPASMCs. Variables with a variable weight value (Variable Important in Projection, VIP) >1 were

| RNA isolation, RNA sequencing and analysis
Total mRNA was isolated from frozen hPASMCs with Trizol (https://metas cape.org/gp/index.html#/main/step1) and visualized in bar plot with barplot function in R.

| Identification and validation of hub genes
The intersection of MAGs and DEGs in rhPDGF-BB treated-versus vehicle treated cells was regarded as critical proliferation associated genes. Next, protein-protein interaction of those genes were analyzed by STRING (v11.0, https://strin g-db.org). 25 CytoHubba plugin with MCC algorithm or MCODE plugin in Cytoscape (v 3.8.2) 26 were used to find the hub genes of the network identified by STRING.
Hub genes were validated in human lungs from the dataset GSE11 7261 27 with 25 control subjects and 58 PH patients and in rat lungs from the dataset GSE85618 28 including 4 PH rats after hypoxia exposure and 4 rats in normoxic condition. All the datasets were obtained from Gene Expression Omnibus (GEO) database.

| Construction of LASSO regression model
Least Absolute Shrinkage and Selection Operator (LASSO) regression model was constructed with the hub gene panel by glmnet package in R to distinguish PH patients from control subjects. A model index for each sample was generated using the regression coefficients to weight the expression of each hub gene. Samples of dataset GSE11 7261 were randomly assigned to training set (70%) and test set (30%). ROC curves were generated to evaluate the ability of LASSO model to identify PH by ROCR package.

| Immune cell infiltrates and its relationship with CXCR4
To decipher the immune cell heterogeneity of human lung tissues, cell type enrichment analysis from gene expression data of GSE11 7261 for different immune cells was performed based on webtool xCell 29 (https://xcell.ucsf.edu). The relationship of CXCR4 with different cell infiltrates were examined by Pearson correlation analysis.

| Data visualization and statistics
Volcano plot was generated to display the overview of distinguishing metabolites or genes between rhPDGF-BB treated-and vehicle treated hPASMCs. The expression of the distinct metabolites or indicated genes in each individual sample was plotted in heatmap. Hub gene network was visualized by Cytoscape. Correlation of immune cell infiltrates with CXCR4 expression was plotted in scatterplot in R.
Data are presented as the mean ± standard error of the mean (SEM).
Statistical differences between two groups were evaluated with 2-tailed unpaired t test if the samples were normally distributed. Otherwise, Mann-Whitney test was used to detect the difference (GraphPad Prism 8). p < .05 denotes significant differences between two groups.

| Altered metabolic signature in hPASMCs in response to rhPDGF-BB and identification of MAGs
As illustrated in flow chart of Figure Figure 2D. The top 10 metabolite sets were then selected for further analysis in Genecards. According to the results, there were 2701 unique genes identified as MAGs with a relevance score >8.

| Identification of rhPDGF-BB induced proliferation associated hub genes and validation in lungs of human and experimental PH
A total of 152 genes in Figure 3D were then potentially considered as rhPDGF-BB induced proliferation associated MAGs. 411 Protein-protein interactions among the 152 genes were identified and visualized by STRING ( Figure S3). Ten hub genes were identified by cytoHubba plugin ( Figure 4A) and fourteen hub genes were revealed by MCODE plugin (Figure 4B). In particular, there were nine core hub genes (IL6, CXCL8, CCL2, CXCR4, CCND1, PLAUR, PLAU, HBEGF and F3) identified by both methods. The expressions of nine hub genes were then scrutinized in lungs of PH patients from GSE11 7261 dataset ( Figure 4C). CXCR4 encoding C-X-C motif chemokine receptor 4, CCND1 encoding cyclin D1 and HBEGF encoding heparin binding EGF like growth factor (HB-EGF) were significantly increased in PH patients (all p < .05). In addition, Ccl2 encoding C-C motif chemokine ligand 2, Plaur encoding urokinase-type plasminogen activator (u-PA) receptor and Hbegf encoding HB-EGF were elevated in lung tissues of PH rats after hypoxia exposure for 2 weeks compared to that at ambient atmosphere from GSE85618 dataset ( Figure 4D).

| Hub genes-based regression model for PH prediction
In a bid to predict PH, we then sought to construct LASSO regression model based on hub genes. There were 3 genes with

| The correlation of Immune cell infiltration with CXCR4 expression
As inflammation is considered as a critical trigger of PH development, 6 we then investigated the association of different inflammatory cell recruitment in human lung tissues with the expression of the inflammation associated hub gene with non-zero regression coefficient (ie. CXCR4), such that we might get a novel insight into F I G U R E 4 Identification of proliferation associated hub genes and validation in lungs of PH patients and PH rat models. (A) Identification of 10 hub genes by cytoHubba in Cytoscape; each circle represents unique gene and the redder the color is, the higher the MCC score is. (B) Identification of 14 hub genes by MCODE in Cytoscape; each circle represents unique gene. (C) Expression of 9 shared hub genes in lungs of 58 patients with pulmonary hypertension (PH) and 25 control subjects from GSE11 7261 were visualized in box plot. Data represent mean ± SEM. *p < .05, **p < .01, ***p < .001 compared to corresponding control subjects, as analyzed by unpaired t test. (D) The expression of hub genes in lungs of rats under hypoxia for two weeks or in normoxic condition (n = 4 per group). Data represent mean ± SEM. *p < .05, **p < .01, ***p < .001 compared to control rats, as analyzed by unpaired t test the potential role of CXCR4 in the modulation of specific immune cell type. We found that the expression of CXCR4 was positively correlated to immune cell infiltrates as evidenced by a higher ImmuneScore in human lung tissues with higher CXCR4 levels (r = .33, p = .002) ( Figure 6A). There was also a positive correla-  Multi-omics is an emerging strategy that holds promise to discover biomarkers for risk stratification or prognosis and to identify novel therapeutic targets/pathological mechanism in diseased state.

| D ISCUSS I ON
For example, metabolomics is now widely used to discover metabolic perturbations or rapid metabolic alteration in response to PH treatment/surgery. It was reported that the plasma metabolic signature of chronic thromboembolic pulmonary hypertension (CTEPH) in pre-pulmonary endarterectomy (PEA) surgery was distinguishing from that of post-PEA surgery. 30 Metabolites in respond to PEA surgery could serve as suitable noninvasive markers for the evaluation on the therapeutic interventions in the future. In terms of the application of transcriptomics to PH, one of the studies showed that whole blood RNA signature in patients with PAH were distinct from that of control subjects, which was also associated with disease severity and allowed for the identification of patients with poor prognosis. 31 This study provides important markers for risk stratification and prognosis prediction of PH. In addition, the NIH/NHLBI A total of nine proliferation associated hub genes were identified, some of which were cytokine/chemokine related genes.
Although the hub gene IL6 didn't differ PH patients or PH rats from their corresponding controls in our selected datasets, lung-specific IL-6 transgenic mice develop spontaneous PH in normoxia and exhibited exaggerated hypoxia-induced PH. 33 In addition, IL-6 blockade mitigated hypoxia-induced PH and repressed the recruitment of Th17 cells and M2 macrophages in lung tissues after hypoxia exposure. 34 These findings implicate IL6 as a trigger in pathogenesis of PH. CXCL8 (also named IL8) encoding interleukin-8 was reported to be higher in serum of PAH patients and had a negative correlation with cardiac index. Kaplan-Meier analysis showed that levels of interleukin-8 predicted survival in PAH patients, with 5-year survival of PAH patients (interleukin-8 levels of >30 pg/ml) to be 32% compared with 58% for patients with levels ≤30 pg/ml. 35 Nevertheless, the specific role of CXCL8 in the vasculature homeostasis remains elusive. Another inflammatory hub gene, CCL2, was found to be higher in plasma and lung tissues of IPAH patients in previous study by Sanchez et al. 36 although no difference in CCL2 between PH patients and controls were found in dataset GSE11 7261. The divergency might lie in the heterogeneity of PH cohorts among two studies. In their study, PASMCs from PAH patients exhibited stronger migration and proliferation in response to CCL2. This finding suggests a critical role of CCL2 mediated PASMC biological activity in vasculature. In Consistent with the higher CXCR4 expression in lungs of PH patients in our study, CXCR4 was also much higher in lungs of PH rat model. 37 Of note, SMC specific loss of CXCR4 inhibits SMC proliferation and retards the hypoxia-induced PH. 38  PH patients and primarily associated with PH due to left sided heart disease. 43 However, how u-PA and its receptor orchestrate the vasculature in PH development remains elucidated and is worth further investigation. F3 (also named TF) encodes tissue factor, coagulation factor III. Previous study showed that higher TF antigen, and TF mRNA in monocytes were displayed in chronic thromboembolic pulmonary hypertension (CTEPH) patients compared with control subjects. TF was also correlated with inflammatory indicators like CRP, TNFα and CCL2. 44 Similar to the clinical observation, TF mRNA expression had a positive correlation with media hypertrophy (ratio of vessel wall area to total area) and mean pulmonary arterial pressure in rat model of CTEPH. 45 The scrutiny of inflammation-coagulation-thrombosis cycle might open new avenue for the treatment of CTEPH. HB-EGF encoded by HBEGF was identified as a mitogen for SMC decades ago and could be induced by multi-mitogens like PDGF. 46 This is in line with our observation. Besides its role in PASMCs, vascular ECs are able to express HBEGF induced by tumor necrosis factorα. 47 These findings would suggest HBEGF may serve as a crucial proliferative signal amplifier and bridge the crosstalk between EC and SMCs in vessels.
Intriguingly, many of the identified hub genes (i.e., IL6, CXCL8, CXCR4, CCND1, PLAUR, PLAU) are hypoxia-inducible factor-alpha (HIFα, transcription factor modulating adaptive responses to hypoxia) target genes. HIF-1α and HIF-2α are regarded as the main type of HIFα, an important factor triggering PH development. This implicates a common molecular alteration between PDGF-BB induced and hypoxia induced proliferative phenotype in PASMCs. To support this notion, a previous study showed PDGF promoted the metabolic shift toward glycolysis in PASMCs via activation of the PI3K/AKT/mTOR/HIF-1α signaling. 48 This is also in line with another study demonstrating that HIF1α upregulated genes were enriched in glycolysis and NADH regeneration. 49 In addition, cytoskeletal protein paxillin tyrosine phosphorylation were elevated in pulmonary vasculature of hypoxia-induced PH mice, which was abolished by PDGF-BB antagonist (imatinib). In the same study, the increase of paxillin tyrosine phosphorylation in human PASMC was blocked by HIF-1α depletion or by imatinib, suggesting the augmentation/phosphorylation of paxillin in PASMC responded to PDGF-BB could be at least partially regulated by HIF-1α. 50 In addition, HIF-2α-dependent CXCL12 secretion in PHD2 (prolyl hydroxylase-2)-deficient ECs facilitate PASMCs proliferation and HIF-2α-selective inhibitor C76 mitigated the established PH rodent models. 51,52 CXCR4, both HIFα target gene and identified hub gene in our study, is a well-known chemokine receptor for CXCL12. Therefore, it can be assumed that PHD2-deficient ECs could not only induce PASMCs proliferation, but also foster PASMCs migration and even immune cells accumulation in perivascular spaces in response to CXCL12 release from ECs under hypoxia stress. The modulation of these hub genes would be of great potential for therapeutic benefits in PH as they could impact multiple factors driving PH pathogenesis.
There are some limitations in our study. First, only one dataset of PH lung tissues from human (ie. GSE11 7261 including 58 PH patients and 25 control subjects) was included. However, the sample size is much larger compared to most of the datasets of human PH lung tissues. Second, the differentially expressed MAGs in hPASMCs in response to PDGF-BB could not fully represent the genes altered in vivo in PH scenario, and their functional properties in pulmonary vascular remodeling need to be investigated in the future.
In conclusion, we identified a metabolic/gene profile change in hPASMCs in response to rhPDGF-BB. The application of multi-omics allows us to discover hub genes to be responsible for proliferating PASMCs phenotype and to be associated with immune cell infiltrates ( Figure 7). We consider that improved molecular understanding of F I G U R E 7 Schematic overview. This schematic overview (created with BioRender.com) illustrates that proliferation associated hub genes in PASMCs in response to PDGF-BB were explored by virtue of transcriptomics and metabolomics. Those hub genes were demonstrated to be responsible for proliferating PASMCs phenotype and association with immune cell infiltrates into the lung tissues, contributing to the development of pulmonary vascular remodeling the intricate networks associated with proliferation in PH will have major therapeutic implications.

ACK N OWLED G M ENTS
None.

CO N FLI C T O F I NTE R E S T
The authors declare that they have no conflict of interest.

DATA AVA I L A B I L I T Y S TAT E M E N T
The datasets used in the current study are available from the corresponding authors on reasonable request.