Introduction

Sphingosine kinase 1 (SPHK1) is a lipid kinase that phosphorylates sphingosine to generate the bioactive sphingolipid sphingosine-1 phosphate (S1P). S1P signals intracellularly or can be secreted and activates intracellular signaling via binding several high-affinity G-protein coupled receptor family members, sphingosine-1 phosphate receptor (S1PR) 1–5. The SPHK1/S1P/S1PR pathway is a critical regulator of several physiological processes including Ca2+ homeostasis, cell survival, migration, and inflammation [1]. The S1PR system regulates immune cell trafficking through the circulatory and lymphatic systems as well as modulates vascular barrier function to control the extravasation of immune cells [2, 3]. Under normal conditions, SPHK1/S1P/S1PR is involved in the maintenance of barrier integrity. Conversely, during inflammation, the signaling axis promotes inflammatory signaling by increasing vascular leakiness [4]. Pathological signaling of SPHK1/S1P/S1PR is involved in multiple diseases among which are cancer, diabetes, atherosclerosis, osteoporosis, and pulmonary arterial hypertension (PAH) [5].

In PAH, SPHK1/S1P induces vascular remodeling by promoting the proliferation of pulmonary artery smooth muscle cells (PASMCs) [6]. Hyperproliferation of PASMCs leads to narrowed pulmonary arteries resulting in increased pulmonary vascular resistance which ultimately results in death due to right heart failure. PAH has multiple etiologies and develops as a co-morbidity of diseases that are characterized by excessive or aberrant inflammatory signaling, such as scleroderma, systemic lupus erythematosus, and human immunodeficiency virus (HIV) infection [7,8,9]. Notably, experimental PAH can be induced in animals by exposure to immune-stimulatory triggers such as HIV and interleukin 6 (IL6) [10, 11]. The available pharmaceutical interventions for PAH target mostly vasoconstriction but fail to address the underlying causes of remodeling of pulmonary vascular cells. Hence, there is an urgent need to identify the molecular mechanisms that promote the abnormal proliferation of PASMCs. Therefore, we over-expressed SPHK1 in human PASMCs and performed RNA sequencing (RNAseq) analysis to identify novel pathways that regulate pulmonary vascular cell function.

Materials and Methods

Reagents and Antibodies

Smooth muscle cell growth medium (SMGM2) and supplements were purchased from Lonza (Walkersville, MD). Rabbit anti-signal transducer and activator of transcription 1 (STAT1) and rabbit IgG secondary antibodies were from Cell Signaling Technology (Danvers, MA). Horseradish peroxidase-conjugated β-actin antibody was obtained from Sigma (St. Louis, MO). RIPA buffer, protease inhibitor cocktail, and Super Signal West Femto were from Thermo Fisher Scientific (Waltham, MA).

Cell Culture and Transfection

Human pulmonary artery smooth muscle cells (hPASMCs) were obtained from Lonza. IPAH and control patient cells were obtained from the pulmonary hypertension breakthrough initiative (PHBI) Cell Core (Philadelphia, PA). Cells were cultured in SMGM2 containing growth supplements and 5% FBS (complete medium) at 37 °C in 5% CO2 in a humidified tissue culture incubator. Human SPHK1 (Myc-DDK-tagged) plasmid and pCMV6 empty vector were purchased from Ori Gene (Rockville, MD). For transfection, 300,000 hPASMCs were cultured on 6-well plates followed by transfection of 3.75 μg of plasmid using Xfect transfection reagent (Takara, Mountain View, CA) according to the manufacturer’s protocol. All experiments were performed in low passage cells (≤8). Cells were transfected in a complete medium. The medium was replaced with a complete medium after 4 h, and the cells were collected after 48 h and stored at −80 °C.

RNAseq and Pathway Analysis

Total RNA was purified from transfected hPASMCs over-expressing SPHK1, as well as transfected vector DNA control, by using the RNeasy Plus Mini Kit (RNA >200 nucleotides: Qiagen, Germantown, MD). Three cell replicates were harvested for either the vector control of SPHK1 over-expression. The average RIN for the RNA was 9.4, an indication of good quality preparations. Total RNA was submitted to the Indiana School of Medicine Center for Medical Genomics (Indianapolis, IN) for RNAseq (Illumina) and differentially expressed gene list generation. In brief, the Genomics Core used the KAPA mRNA Hyperprep Kit for library preparation, running on a NovaSeq6000 for 100 base-paired end reads. The Core’s standard RNAseq quantitation work-follow consisted of accessing sequence quality (FastQC and MultiQC (https://www.bioinformatics.babraham.ac.uk/projects/fastqc/ and https://multiqc.info/)), alignment to UCSC hg38 human genome (STAR aligner [12]), mapping quality (NGSUtils, bamutils [13]), read counting (subread/feature Counts [14]), and differential expression (DE) analysis (edgeR [15]). An average of ~29 million reads was uniquely mapped per sample with ~2.3% rRNA contamination. The .bam files have been deposited at GEO (https://www.ncbi.nlm.nih.gov/geo/) under the accession GSE174304.

For gene set enrichment (GSE) analysis and upstream regulator predictions, Enrichr (https://maayanlab.cloud/Enrichr/) and ingenuity pathway analysis (IPA™) were used. Enrichr utilizes the DE gene list (gene symbols only) for its multi-database analyses. Enrichr provides a single-webtool location to query a wide variety of databases for GSE patterns. Section 1A (Tissues) and 1D (transcription factor co-expression) are derived from the ARCHS4 compilation (https://maayanlab.cloud/archs4/). Section 1B (gene ontology (GO) biological process (BP)) is derived from the GO website (http://geneontology.org/). Section 1C is compiled by the Enrichr project with the primary data derived from the ENCODE project (https://www.encodeproject.org/) and the ChIP-X project (https://maayanlab.cloud/Harmonizome/resource/ChIP-X+Enrichment+Analysis).

IPA uses both the gene symbols and the log2 (fold-change) expression ratios. This enables IPA to estimate whether an enriched pathway is activated or inhibited through a z-score calculation. When a pathway generates a z-score > +2.00 its evidence for pathway activation, while z-score < −2.00 indicates pathway inhibition. In addition, z-scores are bias-adjusted based on whether a specific DE gene list has a strong non-random bias to it, e.g., a DE gene list composed of interferon-associated genes. In the data presented here, there appears to be a bias toward immune system-associated genes, so when flagged, bias-corrected z-scores are used. The most recent IPA software was used (version 62089861, building coral, updated 02/17/2021).

Experimental Statistical Analysis

Statistical analysis was performed using Prism 9.0 Software (GraphPad Software, Inc. San Diego, CA). A minimum of three replicates of each experiment was performed and the averaged data are expressed as mean ± SEM. The student’s t-test was used to compare the differences between the means of two groups and P < 0.05 was considered statistically significant.

Results

SPHK1 Over-expression Alters the HPASMC Transcriptome

HPASMCs over-expressing SPHK1 had a total of 473 differentially regulated genes (294 up; 179 down, Supplementary Table 1) that were significantly differentially expressed (FDR, <0.05; Fig. 1A). Among the up-regulated genes was SPHK1 which was increased by 1.5-fold (Fig. 1A, B; Supplementary Table 1). Using traditional statistical GSE categorical analyses (as compiled at https://maayanlab.cloud/Enrichr/ [16]), SPHK1 over-expression modulated expression of genes that regulate multiple cellular processes (Table 1). Enrichr analysis found that the hPASMCs transcriptional pattern after SPHK1 over-expression (GSE adjusted p-values <1.25E−14; Table 1A) was like the DE patterns of Kupffer, respiratory and vascular smooth muscle cells. This supports a tissue-type shift toward both a vascular adherent and immune cell (Kupffer cells) phenotype while maintaining their original vascular smooth muscle characteristics. The top canonical GO BP was related to interferon- and viral infection responses (GSE adjusted p-values <1.20E−15; Table 1B). Querying the transcription factor ChIP-X database in Enrichr (“X” represents a variety of different genome-wide ChIP-type approaches), two of the top predicted transcription factor binding sites for the DE gene set are interferon-regulatory factor 1 (IRF1) and IRF8 (GSE adjusted p-values <1.01E−04; Table 1C). The transcription factor co-expression analysis point to a number of potential immune response regulators (STAT1, STAT2, IRF7, and interferon ι-induced protein 16), which when over-expressed, produces a DE pattern similar to SPHK1 over-expression in hPASMCs (GSE adjusted p-values <1.00E−25; Table 1D). This suggests these transcriptional regulators might be potential regulators of the SPHK1 over-expression gene set signature in hPASMCs.

Fig. 1: RNAseq analysis identified 473 significant differentially expressed genes in hPASMCs after SPHK1 over-expression.
figure 1

A A volcano plot (log2(FC SPHK1/vector) vs. −log10(adjusted p-value)) is shown with the significant genes identified in red triangles. The four genes, SPHK1, IL6, PARP9, and STAT1, discussed in detail are identified as blue triangles and labeled with their gene symbols. B The differential expression found by RNAseq analysis (expressed as log2(FPKM)) for SPHK1, IL6, PARP9, and STAT1 are shown in ascending expression order. The statistical significance of the DE identified by edgeR analysis is indicated by ****(adjusted p-value < 1E−04). For comparison, GAPDH, a typical housekeeping gene, was not significantly different (ns = not significant) and its expression level was 6–7× log2-fold higher (64–128× fold higher) than STAT1

Full size image
Table 1 Gene set enrichment analysis of differentially expressed genes in hPASMCs after SPHK1 over-expression
Full size table

IPA™ summary of the inter-relationship of DE genes and biological functions is shown in Fig. 2. The most striking feature of this figure is the large set of up-regulated immune response genes (especially in the interferon family) and predicted anti-virus infection state of the hPASMCs. In the diseases and functions analysis, the top categories included replication of virus (GSE adjusted p-values = 4.87E−18; predicted to be inhibited, z-score = −3.73) and antiviral response (GSE adjusted p-values = 4.33E−20; predicted to be activated, z-score = +3.23).

Fig. 2: Ingenuity pathways analysis identifies a strong anti-viral phenotype in hPASMCs after SPHK1 over-expression.
figure 2

The summarized global analysis of the DE gene list in hPASMCs after SPHK1 over-expression suggests a strong anti-viral phenotype. The tightly grouped immune response genes in the center (all up-regulated in the DE gene list) are predicted to down-regulate a wide variety of viral infections (DNA and RNA viruses) though there is an edge effect predicting up-regulation of coronavirus pathogenesis. The IPA legend is included for interpreting the figure

Full size image

SPHK1 over-expression leads to an increased cell proliferation gene signature

Our previous studies demonstrated that SPHK1 over-expression leads to increased cell proliferation and mitogen-activated protein kinase (MAPK) signaling (through extracellular signal-regulated kinase (ERKs, MAPK1/3)) [6]. Though these genes were not specifically found to be up-regulated by increased transcription, the DE pattern suggested that there was activation of this pathway (ERK group GSE adjusted p-values = 1.25E−8; predicted to be activated, z-score = +2.3; Supplementary Table 2). Interestingly, IPA gave differential predictions for MAPK3/ERK1(p44) vs. MAPK1/ERK2(p42), with MAPK3 predicted to be an activated upstream regulator (z-score = +2.8) while MAPK1 was predicted to be an inhibited upstream regulator (z-score = −5.3). The predicted upstream regulation pattern observed for G1- > S progression was also supported by IPA predictions (CCND1 GSE adjusted p-values = 1.92E−3; predicted to be activated, z-score = +2.7; cyclin-dependent kinase inhibitor 1A (CDKN1A) GSE adjusted p-values = 1.21E−6; predicted to be inhibited, z-score = −2.2). Several ERK-regulated immediate early genes (FOS, JUN, and EGR1) were predicted to be activated as well (predicted to be activated, z-score > +2.1).

SPHK1 Over-expression Leads to Increased Expression of Genes That Regulate Proinflammatory Signaling Pathways

Analysis indicated that many of the genes that were upregulated by SPHK1 over-expression are involved in the innate and adaptive immune response (Supplementary Table 1). The poly (ADP-ribose) polymerase 9 (PARP9)-STAT1-IL6 genes were up-regulated by SPHK1 over-expression (Fig. 1B, PARP9 (1.5× up, q-value = 2.50E−08); STAT1 (1.6× up, q-value = 2.19E−05); IL6 (1.6× up, q-value = 1.33E−10)). PARP9 promotes activation of STAT1 downstream of interferon signaling. Several cytokines and chemokines were also up-regulated including IL33, CXCL1 (alias GRO1), CXCL8 (alias IL-8), and CXCL6 (Supplementary Table 1). In addition, members of the tumor necrosis factor ligand superfamily (TNFSF), TNFSF10, and TNFSF13B were significantly increased (Supplementary Table 1). To summarize, SPHK1 over-expression leads to increased expression of many genes that can be upregulated by the interferon-inducible immune response.

Predicted Upstream Regulators of Differential Gene Expression Pattern after SPHK1 Over-expression in hPASMCs

The top twenty activated or inhibited IPA predicted upstream regulators are shown in Table 2. The complete IPA upstream regulators output is included as Supplementary Table 2. As shown above, IL6, STAT1, and PARP9 were both up-regulated in hPASMCs after SPHK1 over-expression. In addition, IPA predicted a significant regulatory effect of IL6, STAT1, and PARP9 over-expression predicted from the DE pattern (Figs. 3, 4, and 5 respectively). IL6 and STAT1 can be found in Table 2 showing both a high degree of GSE (the corrected p-values < 1.00E−15) and strongly predicted activation of their pathways (predicted activated z-scores > +4.7). In addition, PARP9 is also a predicted upstream regulator (corrected p-value = 9.70E−10, predicted activated z-score = +3.0, Supplementary Table 2). Thus, IL6, STAT1, and PARP9 are upregulated in hPASMCs over-expressing SPHK1, and the resultant DE gene pattern supports their activating downstream effects as IL6 (cytokine pathway), and PARP9 and STAT1 (transcription regulators).

Table 2 Ingenuity pathways analysis identified potential upstream regulators
Full size table
Fig. 3: Ingenuity pathways analysis identifies A network surrounding IL6 over-expression.
figure 3

IL6 was up-regulated by SPHK1 over-expression in hPASMCs and IPA built a network of interacting genes, which were also significantly up-regulated (orange/red) or down-regulated (blue/green). The observed pattern of significantly DE genes suggests activation of the IL6 network

Full size image
Fig. 4: Ingenuity pathways analysis identifies A network surrounding STAT1 over-expression.
figure 4

STAT1 was up-regulated by SPHK1 over-expression in hPASMCs and IPA built a network of interacting genes, which were also significantly up-regulated (orange/red) or down-regulated (blue/green). The observed pattern of significantly DE genes suggests activation of the STAT1 network

Full size image
Fig. 5: Ingenuity pathways analysis identifies a network surrounding PARP9 over-expression.
figure 5

PARP9 was up-regulated by SPHK1 over-expression in hPASMCs and IPA built a network of interacting genes, which were also significantly up-regulated (orange/red). The observed pattern of significantly DE genes suggests activation of the PARP9 network

Full size image

Global Immune Response Gene Expression Pattern Predicted Downstream Of Interferon-gamma (IFNG) Pathway

IPA network analysis show that IL6, STAT1, and PARP9 map to the interferon-gamma (IFNG) pathway (Fig. 6). Unlike these genes, IFNG itself is not up-regulated by SPHK1 over-expression. However, the pattern of DE genes strongly supports the activation of IFNG signaling in hPASMCs over-expressing SPHK1 (corrected p-value = 3.73E−51; z-score = +7.7). The predicted strong activation is supported by 246 DE genes, from the complete list of 473, that fall within the IFNG regulation network, and by the direction of gene expression of the vast number of these genes following the expected pattern for IFNG activation.

Fig. 6: Ingenuity pathways analysis predicts INFG is an upstream transcriptional regulator.
figure 6

IFNG is predicted to be an important upstream regulator of the observed DE gene pattern in hPASMCs after SPHK1 over-expression, though IFNG itself was not observed to be transcriptionally up-regulated. IPA built a network of genes potentially stemming from IFNG, which were significantly up-regulated (orange/red) or down-regulated (blue/green). The observed pattern of significantly DE genes suggests activation of the extensive IFNG network

Full size image

Increased PARP9 and STAT1 in PASMCs Derived from Patients with PAH

PARP9 and STAT1 protein expression were also increased in PASMCs isolated from idiopathic pulmonary hypertension patients (Fig. 7A, B, P = 0.0007 and P < 0.0001, respectively) suggesting a role for this pathway in the pathogenesis of PAH.

Fig. 7: Idiopathic PAH patients have increased protein levels of PARP9 and STAT1.
figure 7

Western blot showing that PASMCs isolated from patients have increased PARP9 (A) and STAT1 (B) expression. Quantification and normalization to β-actin demonstrate a significant increase in patients (n = 5) when compared to control subjects (n = 7). ***p = .0007 versus control ****P < 0.0001 versus control

Full size image

Discussion

We previously described the role of the SPHK1/S1P/S1PR2 signaling axis in the development of PAH [6]. SPHK1, S1P, and S1PR2 exhibited increased expression in PASMCs isolated from idiopathic PAH patients. S1P mediates intracellular signaling events in smooth muscle cells via ligating S1PR2. Genetic deletion of SPHK1 or pharmacological inhibition of the S1PR2 prevented hypoxia-induced pulmonary hypertension (HPH) in rodent models. Conversely, partial deletion of S1P lyase (Sphl1), which increases the level of S1P, promoted the development of HPH in vivo. Furthermore, SPHK1 over-expression in hPASMCs leads to increase proliferation in vitro. These studies demonstrate that PASMCs has an essential function in promoting the effects of S1P signaling on PAH development.

In the current study, cellular levels of SPHK1 RNA were increased as expected, due to ectopic expression (Fig. 1). Exposure of PASMCs to S1P has been shown to increase the activation of the ERK (MAPK1/3) pathway leading to cell growth [6, 17]. We observed divergent predictions from IPA of the upstream regulator status of MAPK3/ERK1 (p44) and MAPK1/ERK2 (p42) after 48 h of SPHK1 over-expression. ERK1/2 are rapidly and transiently phosphorylated and activated in response to hPASMCs stimulation with S1P. However, the effect of long-term constitutive production of S1P is an area of research that needs deeper exploration. Examining cell cycle regulators showed that SPHK1 over-expression is predicted to promote cell proliferation in G1- > S progression by activation of cyclin D1 (CCDN1) and inhibition of CDKN1A.

A preponderance of the genes that were altered by SPHK1 over-expression was involved in different immune regulatory processes. The SPHK/S1P/S1PR signaling axis has an important role in regulating immune processes in that it is known to regulate both immune cell trafficking and vascular integrity [4, 18]. However, much of the current knowledge in this regard is on endothelial cell inflammatory responses. In resting conditions, the S1P gradient between the circulation (high) and in tissues (low) helps to maintain endothelial barrier protection and to control lymphocyte recirculation through the systemic circulation, and lymphoid vasculature and tissue. During activation of inflammatory signaling, the S1P gradient changes, as does expression of S1PR1 leading to increased vascular permeability as well as increased lymphocyte egress and recruitment to sites of inflammation. This process also involves the signaling of pro-inflammatory cytokines and chemokines [4, 18].

We found that SPHK1 over-expression significantly increased IL6, PARP9, and STAT1 gene expression in hPASMCs (Fig. 1). In macrophage cells, PARP9 is upstream of STAT1 activation as it was demonstrated that its expression is necessary for IFNG induced STAT1 phosphorylation, regulating the expression of INFG induced pro-inflammatory cytokines and chemokines, tumor necrosis factor α (TNFα), interleukin-1β (IL-1β), and chemokine (C–C motif) ligand 2 (CCL2) (alias MCP-1) [19]. In hPASMCs, SPHK1-overexpression increased expression of members of the TNF super ligand family (TNFSF10 and TNFSF13B) and of the interleukin family cytokine, IL33. Expression of chemokines, chemokine (C–X–C motif) ligand 8 (CXCL8, IL-8), CXCL1 (GRO-1), and CXCL6 were also increased. Toll-like receptor 3 (TLR3) expression was also increased and maps to the IL6 signaling network (Fig. 3). Notably, activation of TLR3 leads to IL6 secretion from bronchial epithelial cells [20]. Circulating levels of IL6 are elevated in PAH patients and correlates with worse right ventricular function [21]. Over-expression of IL6 in animals leads to increased muscularization of the pulmonary vasculature and induces spontaneous pulmonary hypertension that is worsened with hypoxia [11, 22]. Furthermore, smooth muscle-specific knock-out of the interleukin 6 receptor (IL6R), provides protection against HPH [23].

In addition to up-regulating mediators of pro-inflammatory signaling, SPHK1 over-expression also enhanced the expression of antiviral genes that support host defense against viral replication. Each of these inflammatory mediators is involved in the IFNG signaling network (Fig. 6). The interferon-induced protein with tetratricopeptide repeat (IFIT) family consists of four members, IFIT1, IFIT2, IFIT3, and IFIT5, that bind to viruses and prevent viral replication by modulating both viral and host cell function [24]. SPHK1 over-expression enhanced the expression of each IFIT member as well as other interferon-regulated antiviral genes. Among the genes that increased were interferon-induced transmembrane protein 1 (IFITM1) which prevents viral entry to host cells, and interferon-induced protein 44-like (IFI44L) which can prevent viral replication without altering viral binding to host cells [25]. All these interferon enhanced genes were associated with up-regulated STAT1 gene expression when pathway analysis was performed (Fig. 4). Pathway analysis also suggests that the effect of up-regulation of these genes on the viral host defense is dependent on the type of virus. For example, replication of Herpesvirus, influenza, and hepatitis C virus was predicted to be inhibited while the Coronavirus pathway was predicted to be activated (Fig. 2). Importantly, anti-SPHK/S1P/S1PR compounds have been suggested as a potential therapeutic approach to treat the SARS-CoV-2 virus (causes COVID-19 disease) as a mechanism to control viral replication [26]. There is a potential link between COVID-19 and pulmonary vascular dysfunction. Patients who succumbed to COVID-19 were found to have thickened pulmonary vascular walls, a major hallmark of PAH [27]. This observation is particularly noteworthy as other viruses, such as HIV, human Herpesvirus-8 (HHV-8) and hepatitis C are also risk factors for developing PAH [9, 28,29,30,31].

Dysregulation of the immune response is clearly associated with the pathobiology of PAH [10, 31,32,33,34,35,36,37,38,39]. The data presented herein suggest that the SPHK1/S1P/S1PR2 signaling axis is important for the regulation of the vascular smooth muscle immune response and highlight the importance of investigating the role of these responses in PAH and other vascular pathologies.