Skip to main content

CEBPD modulates the airway smooth muscle transcriptomic response to glucocorticoids

Abstract

Background

CCAAT/Enhancer Binding Protein D (CEBPD), a pleiotropic glucocorticoid-responsive transcription factor, modulates inflammatory responses. Of relevance to asthma, expression of CEBPD in airway smooth muscle (ASM) increases with glucocorticoid exposure. We sought to characterize CEBPD-mediated transcriptomic responses to glucocorticoid exposure in ASM by measuring changes observed after knockdown of CEBPD and its impact on asthma-related ASM function.

Methods

Primary ASM cells derived from four donors were transfected with CEBPD or non-targeting (NT) siRNA and exposed to vehicle control, budesonide (100 nM, 18 h), TNFα (10 ng/ml, 18 h), or both budesonide and TNFα. Subsequently, RNA-Seq was used to measure gene expression levels, and pairwise differential expression results were obtained for exposures versus vehicle and knockdown versus control conditions. Weighted gene co-expression analysis was performed to identify groups of genes with similar expression patterns across the various experimental conditions (i.e., CEBPD knockdown status, exposures).

Results

CEBPD knockdown altered expression of 3037 genes under at least one exposure (q-value < 0.05). Co-expression analysis identified sets of 197, 152 and 290 genes that were correlated with CEBPD knockdown status, TNFα exposure status, and both, respectively. JAK-STAT signaling pathway genes, including IL6R and SOCS3, were among those influenced by both TNFα and CEBPD knockdown. Immunoblot assays revealed that budesonide-induced IL-6R protein expression and augmented IL-6-induced STAT3 phosphorylation levels were attenuated by CEBPD knockdown in ASM.

Conclusions

CEBPD modulates glucocorticoid responses in ASM, in part via modulation of IL-6 receptor signaling.

Background

Asthma is a chronic inflammatory respiratory disease characterized by variable airflow limitation and airway hyperresponsiveness to specific environmental stimuli that affects over 22 million Americans and incurs an annual cost of $81.9 billion in the U.S. [1]. Treatment of asthma according to established guidelines includes use of inhaled glucocorticoids to control symptoms in patients with persistent asthma, and “bursts” or long-term use of oral formulations to treat exacerbations or severe forms of asthma, respectively [2]. Chronic use of glucocorticoids elicits considerable adverse effects and may alter tissue sensitivity [3]. Studies to better understand glucocorticoid responses have thus been undertaken to identify mechanisms of resistance and improve personalized treatment strategies [3].

Airway smooth muscle (ASM) is a prominent asthma-related cell type that is directly involved in airway remodeling and airway narrowing [4, 5]. In addition to reducing inflammation, glucocorticoids reduce asthma symptoms by modulating other ASM-dependent processes, including impaired bronchodilation [6], airway hyperresponsiveness [7], and increased ASM contractility [8]. In cells, glucocorticoids exert some of their effects via direct modulation of gene transcription through glucocorticoid receptor (GR) binding to DNA at glucocorticoid response elements (GREs) [3]. Some targets of glucocorticoids include TNFα-inducible pro-inflammatory genes whose expression is modulated by nuclear factor κB (NF-κB) and interferons [9, 10].

CCAAT/Enhancer Binding Proteins (C/EBPs) are a family of six transcription factors that regulate immune responses, as well as cell growth, arrest and differentiation [11]. One of these proteins, CCAAT/enhancer binding protein δ (CEBPD), has been linked to various conditions with altered inflammatory responses [12], including cancers [13], lipopolysaccharide-induced acute lung injury [14,15,16], pulmonary Aspergillus fumigatus conidia infection [17], atherosclerosis [18], and Alzheimer’s disease [19]. According to gene expression microarray and RNA-Seq studies, CEBPD expression increases with glucocorticoid exposure in ASM [20, 21]. Additionally, exposure to the glucocorticoid dexamethasone increases GR occupancy near CEBPD in A549 cells, suggesting that CEBPD is a primary glucocorticoid-responsive GR target [22]. Activation of CEBPD by inflammatory factors, including interleukin-6 (IL-6) and tumor necrosis factor-α (TNFα), has also been observed in a variety of tissues, indicating that the modulation of inflammation by CEBPD involves complex tissue-specific signaling pathways that may have opposing outcomes depending on cellular context [12, 13, 23]. Inflammatory cytokines such as IL-6, induce the binding of C/EBPs to promoters of acute phase genes to control their transcription [23], and this binding can be inhibited by steroids [24, 25] via the activation of GR and direct interaction between GR and C/EBPs [26].

Previously, we identified hundreds of ASM glucocorticoid-responsive genes, most of which were consistently differentially expressed in cells derived from asthma donors versus donors without asthma [20]. The greatest difference in fold change based on asthma status that we observed among these genes was for CEBPD, suggesting that it may contribute to differences in glucocorticoid responses in people with asthma via complex interactions with signaling pathways involving pro-inflammatory cytokines (e.g., TNFα), which are also differentially expressed in asthma. Here, we sought to characterize the effects of CEBPD knockdown on the ASM transcriptomic response to glucocorticoid and TNFα exposures, as well as its impact on related ASM function.

Methods

Detailed methods are provided in the Additional file 1.

ASM RNA-Seq library construction, sequencing and data analysis

Total RNA was extracted from ASM cells derived from four non-asthma donors that were transfected with CEBPD or NT siRNA, and exposed to vehicle control, the glucocorticoid budesonide (BUD) (100 nM), TNFα (10 ng/ml), or BUD + TNFα for 18 h. Stranded RNA-Seq libraries were prepared and sequenced on an Illumina HiSeq 2500 instrument. RNA-Seq data are available in the Gene Expression Omnibus (GEO) under accession GSE146017.

The RAVED pipeline was used to analyze RNA-Seq data [27]. Differential expression analysis was performed for ten pairwise comparisons: CEBPD siRNA versus NT siRNA under the four exposure conditions (control, BUD, TNFα, BUD + TNFα); TNFα versus control in cells transfected with NT siRNA and CEBPD siRNA; and BUD versus control and BUD + TNFα versus TNFα in cells transfected with NT siRNA or CEBPD siRNA. Genes with Benjamini–Hochberg adjusted p-values (i.e., q-values) < 0.05 were considered significant. Results of individual gene’s expression across samples were visualized as boxplots, where the line in the center represents the median value, the box spans the inter-quartile range, and the whiskers show the minimum and maximum (without outliers) of the normalized read counts.

Weighted gene co-expression network analysis

Weighted gene co-expression network analysis (WGCNA) was performed using the WGCNA R package [28] to identify groups of genes with similar expression patterns. Correlations between the resultant groups of genes and 11 phenotype variables (based on transfection, exposure, and donor status) were obtained.

Ontological category enrichment analysis

Overall gene set enrichment analysis (GSEA) was performed using the fgsea R package [29]. For select WGCNA co-expression groups, enrichment analysis was performed using modified Fisher’s exact tests [30]. Ontological categories with q-values < 0.05 were considered significant.

ChIP-Seq data analysis

CEBPD-binding sites were identified using the brocade pipeline [31] applied to ChIP-Seq data from GEO entry GSE32465, which measured CEBPD binding in HepG2 and K562 cell lines [32].

Immunoblots

ASM cells derived from 6 non-asthma donors were transfected with CEBPD or NT siRNA and exposed to DMSO (control), IL-6 (10 ng/ml, 30 min), BUD (100 nM, 24 h), or BUD + IL-6. Immunoblot signals were developed for CEBPD, IL-6Rα, α Tubulin, STAT3, and phosphorylated STAT3 (pSTAT3) from protein samples and changes in band intensities were assessed with paired Student’s t-tests. The ratios of signals were visualized as barplots of height equivalent to the mean across donors and error bars representing standard errors (SEs) across replicates.

ASM traction microscopy

Primary human ASM cells from non-asthma donors were transfected with CEBPD or NT siRNA. Fourier transfer traction microscopy (FTTM) [33] was used to measure traction forces in cells at baseline or exposed to either the contractile agonist histamine (1 μM) or the β2-agonist isoproterenol (1 μM) for 5 min. Traction forces were normalized to the baseline traction and visualized as barplots of height equivalent to the mean across five donors and error bars representing SEs across replicates.

Results

RNA-Seq data met quality control (QC) considerations

Thirty-two ASM samples corresponding to four non-smoking donors without chronic disease, four exposure conditions and either CEBPD or NT siRNA transfection status were prepared. Expression of CEBPD in knockdown cells compared to their donor-paired NT siRNA-transfected cells was decreased an average of 67% according to RT-qPCR measurements for all but one sample that showed no change in CEBPD expression levels and thus, was excluded from RNA-Seq library preparation, resulting in 31 samples selected for RNA-Seq (Additional file 1: Fig. E1). These RT-qPCR measurements also showed that CXCL8 expression increased more than tenfold with TNFα exposure, demonstrating an expected pro-inflammatory response, and the increased CXCL8 expression was blunted by BUD, demonstrating an expected decrease in inflammatory response with glucocorticoid exposure. The RNA-Seq data for the 31 samples sequenced was deemed of high quality (Additional file 1: Fig. E2A, B, Table E1), and all samples were included in differential expression analyses. Normalized CEBPD read counts decreased by 70% in the knockdown cells compared to NT siRNA transfected cells (Additional file 1: Fig. E2C). In NT siRNA cells, genes differentially expressed in response to BUD exposure included well-known glucocorticoid-responsive genes (e.g., FKBP5, TSC22D3, GLUL, PER1, CRISPLD2) [20, 21]. Genes differentially expressed in response to TNFα exposure included well-known pro-inflammatory cytokines (e.g., IL6, CXCL8) [9, 10].

Fig. 1
figure 1

Gene co-expression groups associated with various exposures. A Heatmap of correlations between eigengenes and 13 experimental conditions in three gene co-expression groups. B Heatmaps of gene expression represented by log2(normalized count + 1) for all 31 samples in three gene co-expression groups. Samples in these group can be distinguished based on CEBPD siRNA status, TNFα exposure status and both. BUD budesonide, NT non-targeting

Overall ASM transcriptomic changes in response to CEBPD knockdown

Comparison of CEBPD siRNA versus NT siRNA samples found 1,617, 1,459, 1,330 and 1,985 differentially expressed genes within control, TNFα, BUD and BUD + TNFα exposure conditions, respectively (Table 1). A total of 3037 genes were differentially expressed under at least one of these exposure conditions, and 588 of the genes were in common across the four exposure groups. The five top-ranked genes influenced by CEBPD knockdown according to lowest q-value for each of the four exposure conditions, representing 13 unique genes, included TNFRSF10D, a TNF receptor gene with an inhibitory role in apoptosis [34], whose expression was increased by CEBPD knockdown, and TXNIP, an NF-κB inhibitor gene [35] whose expression was decreased by CEBPD knockdown (Table 2). GSEA found that 23 ontological categories were significantly enriched in at least one exposure group with CEBPD knockdown (Additional file 1: Table E2, Fig. E3). Seven of these categories were shared across all exposures: peptide chain elongation, ribosome, 3' UTR mediated translational regulation, influenza viral RNA transcription and replication, nonsense mediated decay enhanced by the exon junction, SRP-dependent cotranslational protein targeting to membrane, and influenza life cycle.

Table 1 Number of significantly differentially expressed genes in various comparisons
Table 2 Top differentially expressed genes in CEBPD siRNA versus NT siRNA cells under exposures of control, BUD, TNFα, or BUD + TNFα

Influence of CEBPD knockdown on ASM transcriptomic response to TNFα exposure

When comparing TNFα versus control exposure, there were 2315 and 1953 differentially expressed genes in NT siRNA and CEBPD siRNA cells, respectively, 1515 of which were in common (Table 1). The log2 fold changes corresponding to the differentially expressed genes were broadly similar in NT siRNA and CEBPD siRNA cells (Additional file 1: Fig. E4A). The five top-ranked genes according to lowest q-value whose expression was altered by TNFα in each of the two transfection status conditions, representing seven unique genes, are shown in Table 3. Although these seven genes, which included the cytokines IL32 and IL6, were significantly differentially expressed regardless of CEBPD knockdown status, genes such as IER3 and ICAM1 had reduced levels of expression with CEBPD knockdown and were significantly differentially expressed when comparing CEBPD knockdown versus NT siRNA status within TNFα-exposed cells (q-value = 7.73 × 10–5 for IER3 and 5.47 × 10–4 for ICAM1). GSEA found that the ontological categories overrepresented by genes in response to TNFα exposure were the same regardless of knockdown status (Additional file 1: Table E3, Fig. E5). In contrast, the ontological categories obtained for the CEBPD siRNA versus NT siRNA comparison among TNFα-exposed cells found that the categories smooth muscle contraction and nitric oxide stimulates guanylate cyclase were affected by CEBPD knockdown (Additional file 1: Table E2, Fig. E3). Individual genes that drove these differences in ontological category overrepresentation included ITGA1 and MYL9 for smooth muscle contraction (CEBPD siRNA versus NT siRNA q-value = 2.24 × 10–17 and 3.95 × 10–8, respectively) and GUCY1B3 and MRVI1 for nitric oxide stimulates guanylate cyclase (CEBPD siRNA versus NT siRNA q-value = 6.92 × 10–3 and 1.08 × 10–14, respectively) (Additional file 1: Fig. E6).

Table 3 Top differentially expressed genes in TNFα versus control in NT siRNA or CEBPD siRNA cells

Influence of CEBPD knockdown on ASM transcriptomic response to budesonide exposure

When comparing BUD versus control exposure, there were 470 and 421 differentially expressed genes among NT siRNA and CEBPD siRNA samples, respectively, 276 of which were in common (Table 1). When comparing BUD + TNFα versus TNFα exposure, there were 535 and 474 differentially expressed genes in NT siRNA and CEBPD siRNA, respectively, 264 of which overlapped. Table 4 lists the five top-ranked genes according to q-value for BUD versus control or BUD + TNFα versus TNFα in either transfection condition, yielding 15 unique genes, which include the well-known glucocorticoid-responsive genes GLUL and DUSP1 [20, 21]. Some top-ranked genes, such as IL1B and PTGS2, had expression levels that were highly induced by TNFα exposure (TNFα versus control q-values < 10–10) and therefore, had greater observed differences in expression with glucocorticoid exposure in the BUD + TNFα co-stimulation than the BUD condition. Overall, however, the log2 fold changes of the differentially expressed genes in BUD versus control and BUD + TNFα versus TNFα were broadly similar by transfection status (Additional file 1: Fig. E4B, C).

Table 4 Top differentially expressed genes in BUD versus control and BUD + TNFα versus TNFα in NT siRNA or CEBPD siRNA cells

We found that 23 ontological categories were significantly enriched in at least one of the four comparisons involving BUD (i.e., BUD versus control in: (1) NT siRNA, and (2) CEBPD siRNA samples; and BUD + TNFα versus TNFα in: (3) NT siRNA, and (4) CEBPD siRNA samples; Additional file 1: Table E4 and Fig. E7). Six of these categories were shared across the four comparisons, including smooth muscle contraction, suggesting that contraction-related gene expression changes were influenced by BUD regardless of transfection or TNFα co-stimulation status. In contrast, some categories were enriched only under conditions that also involved CEBPD knockdown and/or TNFα exposure. For example, the cytokine-cytokine receptor interaction category was not enriched in the BUD versus control within NT siRNA condition, but it was enriched in the BUD versus control within CEBPD knockdown and/or TNFα exposure conditions due to the differential expression of genes such as IL6, IL1A, and IL1B observed under the latter conditions. An example of a category that was enriched only with co-stimulation of TNFα was extracellular matrix organization, which was driven by the collagen-related genes COL12A1, COL7A1, COL5A3 and COL13A1 that were differentially expressed with TNFα + BUD exposure but not BUD alone.

Identification of gene co-expression groups and their association with CEBPD knockdown and TNFα exposure status

We selected the 1,365 genes that were differentially expressed with CEBPD knockdown in any exposure condition for WGCNA analysis. Soft-thresholding power (β) of 18 was chosen to generate an unsigned weighted co-expression network (Additional file 1: Fig. E8). Of eight groups of co-expressed genes identified, three that were significantly correlated with CEBPD knockdown or exposure status, but not donor status, were considered further (Additional file 1: Fig. E9). Correlation coefficients showed that Group 1 (composed of 197 genes) was correlated with CEBPD knockdown status only, Group 2 (composed of 152 genes) was correlated with TNFα exposure status, and Group 3 (composed of 290 genes) was correlated with both CEBPD knockdown and TNFα exposure status (Fig. 1A). To a lesser extent, Groups 2 and 3 were correlated with BUD exposure (p < 0.05). Hierarchical clustering using gene expression levels (i.e., log2(normalized counts + 1)) of the genes in each co-expression group clustered the 31 samples according to their transfection and/or TNFα exposure status, consistent with phenotypes that they were associated with (Fig. 1B). Group 3 subjects were not perfectly clustered, which may be due to the slight correlation between their eigengenes and donor status. In terms of ontological categories overrepresented by the eigengenes, Group 1 contained genes involved in regulation of actin cytoskeleton; Group 2 contained genes involved in interferon, cytokine-cytokine, and JAK-STAT signaling; and Group 3 contained genes involved in translation, influenza viral RNA transcription and replication, and JAK-STAT signaling (Additional file 1: Table E5).

JAK-STAT signaling pathway genes co-expressed in response to both CEBPD knockdown and TNFα exposure

The Group 2 (TNFα exposure-associated) JAK-STAT signaling pathway genes that were overrepresented included IL10RB, IL13RA2, IL15RA, IL7R, LEP, STAT1, STAT4, while the Group 3 (CEBPD knockdown- and TNFα exposure-associated) JAK-STAT signaling pathway genes included CCND3, IL24, IL6R, LIF, PIM1, SOCS1, SOCS2, SOCS3, SPRY4. We focused further on the JAK-STAT signaling pathway genes in Group 3 because its co-expression patterns were associated with a combined phenotype that most relates to the influence of CEBPD on asthma glucocorticoid responses. RNA-Seq results for IL6R, SOCS3, SOCS1 and SOCS2—genes known to participate in IL-6 signaling pathways—across all conditions showed that their expression levels differed with TNFα exposure versus control within the NT siRNA samples differently than within the CEBPD siRNA samples, consistent with the Group 3 phenotype (Fig. 2; Additional file 1: Table E6). Analysis of a ChIP-Seq dataset involving HepG2 and K562 cells found putative CEBPD-binding sites near the transcription start sites (TSS) of each of these four genes, suggesting that CEBPD can directly modulate their transcription (Additional file 1: Fig. E10).

Fig. 2
figure 2

CEBPD knockdown influenced ASM expression of genes in the IL-6 receptor signaling pathway. A Expression levels of IL6R, SOCS3, SOCS1, and SOCS2 in response to control, BUD, TNFα, and BUD + TNFα exposures in NT siRNA and CEBPD siRNA transfected cells obtained with RNA-Seq data show that CEBPD knockdown altered their response to BUD, TNFα and/or BUD + TNFα exposure. Boxplots show the median value at the center line, the box spans the inter-quartile range, and the whiskers span the minimum and maximum (without outliers) of normalized read counts (N = 3–4 donors per condition). Individual read count values are displayed as points. B IL-6 receptor signaling pathway diagram indicating known relationships among IL-6R, SOCS and STAT3 proteins. BUD budesonide, NT non-targeting

CEBPD knockdown blunted IL6-induced IL-6R signaling in ASM

Given that IL-6R had the strongest effect among the IL-6 signaling genes with CEBPD knockdown according to RNA-Seq data, we sought to determine the role of CEBPD knockdown on IL-6-induced IL-6R signaling pathways by measuring changes of IL-6R protein levels and downstream phosphorylation of STAT3 (Fig. 2B). The CEBPD siRNA-transfected cells had decreased levels of CEBPD under control exposure relative to NT siRNA samples, a difference that was starker with the BUD exposure induction of CEPBD, suggesting that the knockdown effectively reduced CEBPD protein levels (Additional file 1: Fig. E11). The increase in CEBPD that was elicited with BUD exposure in NT siRNA-transfected cells (p < 0.05) was substantially diminished with CEBPD knockdown (Fig. 3A). In NT siRNA-transfected cells, IL-6R protein levels were significantly increased with BUD versus control exposure (p < 0.001), a change that was abrogated with CEBPD knockdown (Fig. 3B). The expected IL-6R-mediated induction of pSTAT3 by IL-6 was observed in both NT and CEBPD siRNA transfected cells (p-value < 0.005), and while BUD exposure further augmented pSTAT3 levels in NT siRNA transfected cells (p-value < 0.05), this effect was reduced with CEBPD knockdown (Fig. 3C).

Fig. 3
figure 3

CEBPD knockdown blunted IL-6 receptor protein expression and signaling in ASM. A Representative immunoblot lanes for tubulin and CEBPD expression. Immunoblot band density quantification results for CEBPD (normalized to tubulin under NT siRNA control or CEBPD siRNA control within groups) show that its levels were induced by BUD exposure among NT siRNA-transfected samples, an effect that was substantially reduced with CEBPD siRNA transfection. B Representative immunoblot lanes for tubulin and IL-6R expression. Immunoblot band density quantification results for IL-6R (normalized to tubulin under NT siRNA control or CEBPD siRNA control within groups) were significantly increased with BUD versus control exposure among NT siRNA-transfected cells, an effect that not was significant for BUD versus control exposures with CEBPD knockdown. C Representative immunoblot results examining pSTAT3 normalized to total STAT3 expression under NT siRNA control or CEBPD siRNA control within groups. pSTAT3 was induced with IL-6 exposure and further augmented with addition of BUD exposure in NT siRNA transfected cells. The augmented pSTAT3 induced by BUD + IL-6 versus IL-6 was no longer significant with CEBPD knockdown. Barplots of signal ratios are of height equivalent to the mean across donors, and the error bars represent standard errors (SEs) across replicates with N = 6 donors per condition. Intact immunoblot bands are provided in Additional file 1: Fig. E11. BUD budesonide, NT non-targeting, pSTAT3 phosphorylated STAT3

CEBPD influenced baseline ASM contractility

CEBPD siRNA transfection efficiency in ASM cells was confirmed to be 58% for this experiment (Fig. 4A). At baseline, the average cell traction force was significantly higher (p-value < 0.05) in CEBPD knockdown cells (224.4 ± 31.6 Pa, mean ± standard error hereafter) compared to cells transfected with NT siRNA (176.1 ± 30.3 Pa) (Fig. 4B). ASM traction was increased by histamine (contractile agonist) and decreased by isoproterenol (β2 agonist) relative to baseline levels, but these responses did not significantly differ with CEBPD knockdown (Fig. 4C, D).

Fig. 4
figure 4

CEBPD influenced baseline ASM contractility. A CEBPD expression reduction in ASM was confirmed using qPCR 48 h after siRNA transfection. B CEBPD knockdown in ASM cells led to a statistically significant higher traction at baseline (T0). CEBPD knockdown in ASM cells did not alter C the contractile response to histamine (His) or D the relaxation response to isoproterenol (Iso) based on the relative traction values. Barplots are of height equivalent to the mean across five donors, and error bars represent standard errors (SEs) across replicates with N = 5 donors per condition. *p-value < 0.05. NT non-targeting

Discussion

We and others have observed a large number of glucocorticoid-responsive genes in ASM cells [20, 21, 36, 37]. Our most recent study found that CEBPD had the largest difference in glucocorticoid-induced expression changes in ASM from non-asthma donors versus fatal asthma donors (log2 fold change of 1.43 versus 0.48 with budesonide exposure) [20], suggesting that differing levels of CEBPD expression may influence glucocorticoid responsiveness in people with asthma. Our present study characterized the ASM transcriptomic responses to glucocorticoid and TNFα exposures in the context of CEBPD knockdown, which involved performing 10 pairwise differential expression comparisons. Due to the large number of differentially expressed genes observed among these many comparisons, we performed a weighted gene co-expression analysis to identify the groups of genes specifically changed in response to CEBPD knockdown in the context of other relevant exposures, thereby facilitating the identification of individual genes and pathways for validation studies of our main trait of interest.

Pairwise differential expression results revealed many changes in ASM with CEBPD knockdown, a large proportion of which were specific to exposure conditions. Our pathway-level results of these exposures recapitulated known pathways, including that TNFα is involved in innate immunity and interferon signaling, and that glucocorticoids alter cytokine-cytokine receptor signaling and smooth muscle contraction. The large number of differentially expressed genes observed in CEBPD siRNA versus NT siRNA transfected samples under glucocorticoid and/or TNFα exposures included expected findings, such as changes in genes related to the ontological categories interferon signaling and downstream signaling events of B cell receptor, as well as novel findings of relevance to asthma, such as the alteration of genes involved in smooth muscle contraction and nitric oxide stimulates guanylate cyclase. Although many genes from among these comparisons are of interest to better understand asthma, we proceeded with WGCNA to focus on groups of genes with similar changes across exposures and transfection status.

WGCNA identifies gene co-expression groups based on their topological similarity across samples, and it is able to identify relationships of these co-expression groups with multiple phenotypes under consideration [28]. We tailored WGCNA to our study goals by (1) including differentially expressed genes from the CEBPD siRNA versus NT siRNA across the four exposure comparisons, and (2) constructing networks that included connections regardless of the direction of expression changed by CEBPD knockdown (i.e., we used an unsigned correlation network). Although use of a selected set of genes biases the identification of gene co-expression groups, in this case, it allowed us to identify three salient gene co-expression groups with expression patterns corresponding to CEBPD knockdown and/or TNFα exposure status. We verified that individual genes within the groups had results consistent with their grouping: CEBPD was among the genes in Group 1 (CEBPD knockdown-associated) and many cytokine-related genes were among the genes in Group 2 (TNFα exposure-associated). Interestingly, TNFα exposure and CEBPD knockdown resulted in greater transcriptomic effects than budesonide: (1) there were substantially more differentially expressed genes identified in the pairwise TNFα versus control or CEBPD siRNA versus NT siRNA, than BUD versus control conditions, and (2) none of the gene co-expression groups were highly correlated with budesonide exposure status. Therefore, Group 3 was deemed most relevant to our question of understanding the impact of CEBPD on asthma-related gene expression changes. The Group 3 (CEBPD knockdown- and TNFα exposure-associated) JAK-STAT pathway genes included some whose expression was decreased with CEBPD knockdown (e.g., IL6R, SOCS3, SOCS1) and some whose expression changed in the opposite direction (e.g., SOCS2), demonstrating that WGCNA was helpful to identify sets of genes that changed under specific conditions, regardless of the direction of this change. Comparison of JAK-STAT pathway genes in Group 3 versus Group 2 was helpful to identify the TNFα-modulated ones that were also changed by CEBPD knockdown, which led us to select the IL-6R pathway for further study.

Consistent with RNA-Seq results showing that IL-6R transcript levels were substantially reduced with CEBPD knockdown in the BUD versus control exposure comparison, immunoblot results showed that CEBPD knockdown resulted in decreased IL-6R protein expression when comparing BUD versus control exposures. Extension of RNA-Seq results to the protein level also revealed that IL-6 receptor signaling vis-à-vis IL-6-induced pSTAT3 expression remained intact, although the fold-change of pSTAT3 induced with BUD + IL-6 versus IL-6 alone was only statistically significantly different among NT siRNA-transfected cells, suggesting an overall reduction of IL-6R signaling with CEBPD knockdown. Future studies are needed to investigate more detailed mechanisms whereby altered CEBPD expression and its posttranslational modifications affect IL-6R signaling in ASM to influence glucocorticoid responses in asthma.

Smooth muscle contraction pathway genes were enriched in (1) genes differentially expressed with CEBPD knockdown under the condition of TNFα exposure and (2) genes differentially expressed with budesonide exposure regardless of CEBPD knockdown status. The traction microscopy results support a potential modest effect of CEBPD knockdown on ASM contractile force at baseline, however, CEBPD had little effect on ASM excitation–contraction coupling. Together, these results suggest that CEBPD alone is not likely to directly regulate ASM contractility in response to glucocorticoid exposure.

Several limitations of our study are worth noting. First, we did not determine whether CEBPD modulated IL6 signaling via membrane-bound IL-6R or trans-signaling of its soluble form. Because prior studies found that membrane-bound IL-6R was not present in ASM, while its soluble form was [38], it is likely that CEBPD influences the IL-6 pathway via trans-signaling. Of note, a specific asthma phenotype has been proposed to correspond to IL-6 trans-signaling, as patients with increased IL-6 trans-signaling had more exacerbations, eosinophilia, and submucosal T cells and macrophages [39], and a coding genetic polymorphism in the IL-6R gene that promotes trans-signaling has been linked to lung function differences in people with severe asthma [40]. Second, additional experiments are necessary to determine whether some of the transcriptomic effects observed may have resulted from direct protein–protein interactions among CEBPD, NF-κB and GR. Third, our statistically significant findings for the effect of CEBPD on IL-6R and pSTAT3 had modest effect sizes, which may be due to relatively long exposure times and resulting compensation by other C/EBP family members. Additional experiments are necessary to determine the time courses of CEBPD effects, as well as the concomitant role of CEBPA, CEBPB, and related proteins on IL-6 signaling.

In summary, we found that CEBPD knockdown resulted in many ASM transcriptomic changes in response to glucocorticoid and TNFα exposures. Among these, CEBPD knockdown influenced expression of several TNFα-induced JAK-STAT pathway genes, including the IL-6 receptor. Further mechanistic insights regarding these CEBPD-mediated ASM transcriptomic changes may lead to an improved understanding of glucocorticoid responses in patients with asthma.

Availability of data and materials

RNA-Seq data are available in the Gene Expression Omnibus (GEO) (https://www.ncbi.nlm.nih.gov/geo/) under accession GSE146017.

Abbreviations

ASM:

Airway smooth muscle

BUD:

Budesonide

CEBPD:

CCAAT/enhancer binding protein D

ELISA:

Enzyme-linked immunosorbent assay

FTTM:

Fourier transfer traction microscopy

GSEA:

Gene set enrichment analysis

GR:

Glucocorticoid receptor

GRE:

Glucocorticoid response element

IHC:

Immunohistochemistry

NF-κB:

Nuclear factor κB

NES:

Normalized enrichment score

NT:

Non-targeting

TSS:

Transcription start site

WGCNA:

Weighted gene co-expression network analysis

References

  1. Nurmagambetov T, Kuwahara R, Garbe P. The Economic Burden of Asthma in the United States, 2008–2013. Ann Am Thorac Soc. 2018;15:348–56.

    Article  PubMed  Google Scholar 

  2. Chung KF, Wenzel SE, Brozek JL, Bush A, Castro M, Sterk PJ, et al. International ERS/ATS guidelines on definition, evaluation and treatment of severe asthma. Eur Respir J. 2014;43:343–73.

    Article  CAS  PubMed  Google Scholar 

  3. Kan M, Himes BE. Insights into glucocorticoid responses derived from omics studies. Pharmacol Therapeut. 2020;218:107674.

    Article  Google Scholar 

  4. Prakash YS. Airway smooth muscle in airway reactivity and remodeling: what have we learned? Am J Physiol Lung Cell Mol Physiol. 2013;305:L912-933.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  5. Prakash YS, Halayko AJ, Gosens R, Panettieri RA, Camoretti-Mercado B, Penn RB, et al. An Official American Thoracic Society Research Statement: current challenges facing research and therapeutic advances in airway remodeling. Am J Respir Crit Care Med. 2017;195:e4-19.

    Article  CAS  PubMed  Google Scholar 

  6. Slats AM, Sont JK, van Klink RHCJ, Bel EHD, Sterk PJ. Improvement in bronchodilation following deep inspiration after a course of high-dose oral prednisone in asthma. Chest. 2006;130:58–65.

    Article  CAS  PubMed  Google Scholar 

  7. Trifilieff A, El-Hashim A, Bertrand C. Time course of inflammatory and remodeling events in a murine model of asthma: effect of steroid treatment. Am J Physiol Lung Cell Mol Physiol. 2000;279:L1120-1128.

    Article  CAS  PubMed  Google Scholar 

  8. Lakser OJ, Dowell ML, Hoyte FL, Chen B, Lavoie TL, Ferreira C, et al. Steroids augment relengthening of contracted airway smooth muscle: potential additional mechanism of benefit in asthma. Eur Respir J. 2008;32:1224–30.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  9. Kalliolias GD, Ivashkiv LB. TNF biology, pathogenic mechanisms and emerging therapeutic strategies. Nat Rev Rheumatol. 2016;12:49–62.

    Article  CAS  PubMed  Google Scholar 

  10. Bradley J. TNF-mediated inflammatory disease. J Pathol. 2008;214:149–60.

    Article  CAS  PubMed  Google Scholar 

  11. Lekstrom-Himes J, Xanthopoulos KG. Biological role of the CCAAT/enhancer-binding protein family of transcription factors. J Biol Chem. 1998;273:28545–8.

    Article  CAS  PubMed  Google Scholar 

  12. Ko C-Y, Chang W-C, Wang J-M. Biological roles of CCAAT/Enhancer-binding protein delta during inflammation. J Biomed Sci. 2015;22:6.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  13. Balamurugan K, Sterneck E. The many faces of C/EBPδ and their relevance for inflammation and cancer. Int J Biol Sci. 2013;9:917–33.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  14. Yan C, Johnson PF, Tang H, Ye Y, Wu M, Gao H. CCAAT/enhancer-binding protein delta is a critical mediator of lipopolysaccharide-induced acute lung injury. Am J Pathol. 2013;182:420–30.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  15. Yan C, Ward PA, Wang X, Gao H. Myeloid depletion of SOCS3 enhances LPS-induced acute lung injury through CCAAT/enhancer binding protein delta pathway. FASEB J. 2013;27:2967–76.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  16. Do-Umehara HC, Chen C, Urich D, Zhou L, Qiu J, Jang S, et al. Suppression of inflammation and acute lung injury by Miz1 via repression of C/EBP-delta. Nat Immunol. 2013;14:461–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  17. Liu L, Zhong J, Chen B, Wang W, Xi H, Su X. CCAAT/enhancer binding protein (C/EBP) delta promotes the expression of PTX3 and macrophage phagocytosis during A. fumigatus infection. J Leukoc Biol. 2021;111:1225.

    Article  PubMed  Google Scholar 

  18. Lai H-Y, Hsu L-W, Tsai H-H, Lo Y-C, Yang S-H, Liu P-Y, et al. CCAAT/enhancer-binding protein delta promotes intracellular lipid accumulation in M1 macrophages of vascular lesions. Cardiovasc Res. 2017;113:1376–88.

    Article  CAS  PubMed  Google Scholar 

  19. Ko CY, Chang LH, Lee YC, Sterneck E, Cheng CP, Chen SH, et al. CCAAT/enhancer binding protein delta (CEBPD) elevating PTX3 expression inhibits macrophage-mediated phagocytosis of dying neuron cells. Neurobiol Aging. 2012;33(422):e11-25.

    Google Scholar 

  20. Kan M, Koziol-White C, Shumyatcher M, Johnson M, Jester W, Panettieri RA, et al. Airway smooth muscle-specific transcriptomic signatures of glucocorticoid exposure. Am J Respir Cell Mol Biol. 2019;61:110–20.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  21. Himes BE, Jiang X, Wagner P, Hu R, Wang Q, Klanderman B, et al. RNA-Seq transcriptome profiling identifies CRISPLD2 as a glucocorticoid responsive gene that modulates cytokine function in airway smooth muscle cells. PLoS ONE. 2014;9: e99625.

    Article  PubMed  PubMed Central  Google Scholar 

  22. Reddy TE, Pauli F, Sprouse RO, Neff NF, Newberry KM, Garabedian MJ, et al. Genomic determination of the glucocorticoid response reveals unexpected mechanisms of gene regulation. Genome Res. 2009;19:2163–71.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  23. Poli V. The role of C/EBP isoforms in the control of inflammatory and native immunity functions. J Biol Chem. 1998;273:29279–82.

    Article  CAS  PubMed  Google Scholar 

  24. Akira S. IL-6-regulated transcription factors. Int J Biochem Cell Biol. 1997;29:1401–18.

    Article  CAS  PubMed  Google Scholar 

  25. Nishio Y, Isshiki H, Kishimoto T, Akira S. A nuclear factor for interleukin-6 expression (NF-IL6) and the glucocorticoid receptor synergistically activate transcription of the rat alpha 1-acid glycoprotein gene via direct protein-protein interaction. Mol Cell Biol. 1993;13:1854–62.

    CAS  PubMed  PubMed Central  Google Scholar 

  26. Mitani T, Takaya T, Harada N, Katayama S, Yamaji R, Nakamura S, et al. Theophylline suppresses interleukin-6 expression by inhibiting glucocorticoid receptor signaling in pre-adipocytes. Arch Biochem Biophys. 2018;646:98–106.

    Article  CAS  PubMed  Google Scholar 

  27. Kan M, Shumyatcher M, Diwadkar A, Soliman G, Himes BE. Integration of transcriptomic data identifies global and cell-specific asthma-related gene expression signatures. AMIA Annu Symp Proc. 2018;2018:1338–47.

    PubMed  PubMed Central  Google Scholar 

  28. Langfelder P, Zhang B, Horvath S. Defining clusters from a hierarchical cluster tree: the Dynamic Tree Cut package for R. Bioinformatics. 2008;24:719–20.

    Article  CAS  PubMed  Google Scholar 

  29. Sergushichev A. An algorithm for fast preranked gene set enrichment analysis using cumulative statistic calculation. bioRxiv. 2016.

  30. Hosack DA, Dennis G, Sherman BT, Lane HC, Lempicki RA. Identifying biological themes within lists of genes with EASE. Genome Biol. 2003;4:R70.

    Article  PubMed  PubMed Central  Google Scholar 

  31. Diwadkar AR, Kan M, Himes BE. Facilitating analysis of publicly available ChIP-Seq data for integrative studies. AMIA Annu Symp Proc. 2019;2019:371–9.

    PubMed  Google Scholar 

  32. Gertz J, Savic D, Varley KE, Partridge EC, Safi A, Jain P, et al. Distinct properties of cell-type-specific and shared transcription factor binding sites. Mol Cell. 2013;52:25–36.

    Article  CAS  PubMed  Google Scholar 

  33. Butler JP, Tolić-Nørrelykke IM, Fabry B, Fredberg JJ. Traction fields, moments, and strain energy that cells exert on their surroundings. Am J Physiol Cell Physiol. 2002;282:C595-605.

    Article  CAS  PubMed  Google Scholar 

  34. Marsters SA, Sheridan JP, Pitti RM, Huang A, Skubatch M, Baldwin D, et al. A novel receptor for Apo2L/TRAIL contains a truncated death domain. Curr Biol. 1997;7:1003–6.

    Article  CAS  PubMed  Google Scholar 

  35. Sherbet GV. Tumour Suppression and Inhibition of Angiogenesis by TXNIP (Thioredoxin Interacting Protein). Molecular Approach to Cancer Management [Internet]. Elsevier; 2017 [cited 2020 May 22]. p. 101–3. Available from: https://linkinghub.elsevier.com/retrieve/pii/B9780128128961000106.

  36. Masuno K, Haldar SM, Jeyaraj D, Mailloux CM, Huang X, Panettieri RA, et al. Expression profiling identifies Klf15 as a glucocorticoid target that regulates airway hyperresponsiveness. Am J Respir Cell Mol Biol. 2011;45:642–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  37. Misior AM, Deshpande DA, Loza MJ, Pascual RM, Hipp JD, Penn RB. Glucocorticoid- and protein kinase A-dependent transcriptome regulation in airway smooth muscle. Am J Respir Cell Mol Biol. 2009;41:24–39.

    Article  CAS  PubMed  Google Scholar 

  38. Ammit AJ, Moir LM, Oliver BG, Hughes JM, Alkhouri H, Ge Q, et al. Effect of IL-6 trans -signaling on the pro-remodeling phenotype of airway smooth muscle. Am J Physiol Lung Cell Mol Physiol. 2007;292:L199-206.

    Article  CAS  PubMed  Google Scholar 

  39. Jevnikar Z, Östling J, Ax E, Calvén J, Thörn K, Israelsson E, et al. Epithelial IL-6 trans-signaling defines a new asthma phenotype with increased airway inflammation. J Allergy Clin Immunol. 2019;143:577–90.

    Article  CAS  PubMed  Google Scholar 

  40. Hawkins GA, Robinson MB, Hastie AT, Li X, Li H, Moore WC, et al. The IL6R variation Asp(358)Ala is a potential modifier of lung function in subjects with asthma. J Allergy Clin Immunol. 2012;130:510-515.e1.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

Download references

Acknowledgements

Not applicable.

Funding

Funding was provided by National Institutes of Health (NIH) R01 HL133433, R01 HL141992, R01 HL139496, and the Center of Excellence in Environmental Toxicology (P30 ES013508). The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.

Author information

Authors and Affiliations

Authors

Contributions

Designed the study: MK, CK, RAP, QL and BEH. Conducted experiments: MS, XJ, VP, GC, EG, CK, WJ and BL. Analyzed data: MK, ARD, CK and BEH. Wrote the manuscript with input from other co-authors: MK, CK and BEH. Edited and approved the final draft: all authors. All authors read and approved the final manuscript.

Corresponding author

Correspondence to Blanca E. Himes.

Ethics declarations

Ethics approval and consent to participate

The lung tissue from which airway smooth muscle was derived was provided by the National Disease Research Interchange and the International Institute for the Advancement of Medicine and is obtained commercially from deceased anonymous donors. The University of Pennsylvania Committee on Studies Involving Human Beings, the Rutgers Biomedical and Health Sciences Institutional Review Board, and the Harvard Institutional Review Board have determined that use of these cells does not constitute Human Subjects research since all donor tissue was harvested anonymously and de-identified.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no conflicts of interest.

Additional information

Publisher's Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary Information

Additional file 1

: Table E1. RNA-Seq Quality Control Metrics. Table E2. Gene set enrichment analysis results corresponding to the CEBPD siRNA versus NT siRNA comparisons. Table E3. Gene set enrichment analysis results corresponding to the TNFα versus control comparisons. Table E4. Gene set enrichment analysis results corresponding to budesonide exposure. Table E5. Ontological categories enriched within gene co-expression groups. Table E6. RNA-Seq differential expression results for IL6R, SOCS3, SOCS1, and SOCS2 across the 10 comparisons made. Figure E1. Sample quality control prior to RNA-Seq via RT-qPCR of CEBPD and CXCL8. Figure E2. RNA-Seq data quality control. Figure E3. Significant gene set enrichment analysis categories corresponding to the CEBPD siRNA versus NT siRNA comparisons. Figure E4. Overall RNA-Seq results for each exposure were generally consistent with CEBPD knockdown. Figure E5. Significant gene set enrichment analysis categories corresponding to the TNFα versus control comparisons. Figure E6. TNFα-responsive genes whose expression changed with CEBPD knockdown selected from two significantly changed ontological categories. Figure E7. Significant gene set enrichment analysis categories corresponding to budesonide-responsive genes. Figure E8. Selection of soft-thresholding power (β) for weighted gene co-expression network analysis. Figure E9. Correlations between gene co-expression groups and phenotypes. Figure E10. CEBPD-binding sites near select IL-6 signaling pathway genes. Figure E11. Full representative immunoblots and CEBPD/Tubulin expression levels showing adequacy of CEBPD knockdown.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Kan, M., Sun, M., Jiang, X. et al. CEBPD modulates the airway smooth muscle transcriptomic response to glucocorticoids. Respir Res 23, 193 (2022). https://doi.org/10.1186/s12931-022-02119-1

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s12931-022-02119-1

Keywords