Revealing the action mechanisms of dexamethasone on the birth weight of infant using RNA-sequencing data of trophoblast cells.

Dexamethasone (DEX) could induce low birth weight of infant, and low birth weight has close associations with glucocorticoid levels, insulin resistance, hypertension, and metabolic syndrome in adulthood. This study was designed to reveal the action mechanisms of DEX on the birth weight of infant.Using quantitative real-time polymerase chain reaction (qRT-PCR), trophoblast cells of human placenta were identified and the optimum treatment time of DEX were determined. Trophoblast cells were treated by DEX (DEX group) or ethanol (control group) (each group had 3 samples), and then were performed with RNA-sequencing. Afterward, the differentially expressed genes (DEGs) were identified by R package, and their potential functions were successively enriched using DAVID database and Enrichr method. Followed by protein-protein interaction (PPI) network was constructed using Cytoscape software. Using Enrichr method and TargetScan software, the transcription factors (TFs) and micorRNAs (miRNAs) targeted the DEGs separately were predicted. Based on MsigDB database, gene set enrichment analysis (GSEA) was performed.There were 391 DEGs screened from the DEX group. Upregulated SRR and potassium voltage-gated channel subfamily J member 4 (KCNJ4) and downregulated GALNT1 separately were enriched in PDZ (an acronym of PSD-95, Dlg, and ZO-1) domain binding and Mucin type O-glycan biosynthesis. In the PPI network, CDK2 and CDK4 had higher degrees. TFs ATF2 and E2F4 and miRNA miR-16 were predicted for the DEGs. Moreover, qRT-PCR analysis confirmed that SRR and KCNJ4 were significantly upregulated.These genes might affect the roles of DEX in the birth weight of infant, and might be promising therapeutic targets for reducing the side effects of DEX.

Introduction

Premature delivery is an important factor that affects the morbidity and mortality of perinatals, with a incidence of 5% to 15%.[ Premature infants have immature respiratory system and lack of pulmonary surfactant, thus their lives are seriously threatened.[ For pregnant women in the 24th to 34th week and with the risk of giving birth prematurely, glucocorticoids are recommended for promoting fetal lung maturation.[ However, glucocorticoids used in metaphase and terminal period can significantly reduce the birth weight and head circumference of infant.[ On the other hand, the low birth weight of infant has close association with glucocorticoid levels, insulin resistance, hypertension, and metabolic syndrome in adulthood.[ Therefore, revealing the action mechanisms of glucocorticoids on the birth weight of infant is critical for decreasing the side effects of glucocorticoids.

Glucose is an important nutrient in fetal growth and development, and glucose transporter 1 (GLUT1) and GLUT3 separately play primary and secondary roles in glucose transport in placenta.[ Several studies have found that glucocorticoids can reduce the expression of GLUT in placenta.[ Glucocorticoid dexamethasone (DEX) decreases the expression level of glucocorticoid receptor α (GRα) through reducing the number of GRα+ binucleate cells, which varies with different fetal gender and gestational age.[ DEX is reported to impair glucose transport capacity and induce insulin resistance through reducing insulin receptor substrate-1 (IRS-1), phosphatidylinositol 3-kinase (PI3-K), and protein kinase B (PKB).[ However, the action mechanisms of DEX on the birth weight of infant have not yet been comprehensively reported.

In this study, the trophoblast cells of human placenta treated or untreated by DEX were used for RNA-sequencing. Subsequently, the differentially expressed genes (DEGs) were analyzed and predicted for their potential functions. Furthermore, protein–protein interaction (PPI) network analysis, regulatory network analysis, and gene set enrichment analysis (GSEA) were carried out for deeply exploring the target genes of DEX. In addition, quantitative real-time polymerase chain reaction (qRT-PCR) was applied for confirming the expression changes of key genes. This study might contribute to identify the therapeutic targets for reducing the side effects of DEX.

Methods

Identification of trophoblast cells and determination of the optimum treatment time of DEX

The trophoblast cells of human placenta used in this study were acquired from Cell bank of Chinese Academy of Sciences. Therefore, ethical approval and informed consent were not necessary. The cells were cultured in the Dulbecco's modified Eagle Medium/Ham's F-12 medium (DMEM/F12, GIBCO, Grand Island, NY) containing 10% fetal bovine serum (FBS, GIBCO, Grand Island, NY) and 1% penicillin/streptomycin double-antibody (GIBCO, Grand Island, NY) in a humidified incubator (at 37°C and 5% CO2, Thermo Fisher Scientific Inc., MA). After the cells covered 80% to 90% of a culture dish, they were digested by pancreatin (GIBCO, Grand Island, NY), centrifuged, added with fresh medium, and then transferred to a new culture dish for passage. The trophoblast cells grew well were suspended in 2 mL DMEM/F12 medium (GIBCO, Grand Island, NY) in a 6-well plate and cultured in a humidified incubator (at 37°C and 5% CO2, Thermo Fisher Scientific Inc., MA) for 24 hours. Followed by the cells in DEX group were treated by DEX (final concentration: 100 nM) for 24 hours, and the cells in control group were treated by the same volume of ethanol. Then, total RNA was extracted from the trophoblast cells using RNAiso Plus (Takara, Shiga, Japan) following the manufacturer's protocol, and were determinated using a spectrophotometer (NanoDrop Technologies, Wilmington, Delaware). Subsequently, RNA was reversely transcribed into the first-strand cDNA by a reverse transcription kit (TaKaRa, Tokyo, Japan) and stored for later use.

For polymerase chain reaction (PCR) experiments, the primers of histocompatibility leukocyte antigen-E (HLA-E), GLUT1, GLUT3, and glyceraldehyde-3-phosphate dehydrogenase (GAPDH) were designed using Primer Premier 6.0 software (Premier Software Inc., Cherry Hill, NJ) and then synthesized by Sangon Biotech Co., Ltd (Sangon Biotech Co., Ltd., Shanghai, China) (Table 1). PCR was used to measure the expression level of HLA-E in trophoblast cells, and the amplification products were detected using 1.5% agarose gel electrophoresis. Besides, the expression levels of GLUT1 and GLUT3 in trophoblast cells were evaluated using SYBR green kit (Applied Biosystems, CA). The 20 μL reaction system included 10 μL SYBR Green Mix Buffer (2×), 2 μL cDNA solution (100 ng/μL), 0.4 μL forward primer (10 μM), 0.4 μL reverse primer (10 μM), and 7.2 μL RNase Free ddH2O. The amplifying program was as follows: 95°C for 10 minutes; and 95°C for 15 seconds and 60°C for 60 seconds for 40 cycles; followed by a melting program of 95°C for 15 seconds, 60°C for 1 minute and 95°C for 15 seconds. GAPDH was selected as the reference gene, and each samples had 3 repeats.

Table 1

The primers used for polymerase chain reaction (PCR) experiments.

RNA extraction and RNA-seq library construction

After the trophoblast cells were treated by DEX or the same volume of ethanol for 24 hours (each group had 3 samples), total RNA was extracted using RNAiso Plus (Takara, Shiga, Japan) and were detected using a spectrophotometer (NanoDrop Technologies, Wilmington, DE). After RNA-seq library construction was performed using NEBNext Ultra RNA Library Prep Kit for Illumina (New England Biolabs Inc., Beverly, MA), sequencing data were collected using the Illumina Hiseq 4000 platform (PE150) (Illumina, San Diego, CA). Finally, the sequencing data (accession number: SRP105013) were uploaded into the Sequence Read Archive (SRA) database.

Data preprocessing and DEGs screening

Based on FastQC (http://www.bioinformatics.babraham.ac.uk/projects/fastqc/) method,[ quality control was carried out for the raw data. Based on Tophat (http://tophat.cbcb.umd.edu/)[ and Cufflinks (http://cufflinks.cbcb.umd.edu)[ software, the clean reads were assembled into transcripts. The hg19 reference genome in the University of California Santa Cruz (UCSC, https://genome.ucsc.edu) database and the hg19 RefGenes (a total of 25575 genes) were used as reference genome and reference transcriptome, respectively. Cufflinks software[ was used to assemble and identify novel transcripts for non-RefGenes. Subsequently, sample normalization and calculation of The Fragments Per Kilobase Million (FPKM) for the genes in the 2 groups were performed using Cufflinks software.[ After the FPKM values were performed with log 2 transformation, the t-test method[ in R package (https://www.r-project.org/) was used to identify the DEGs between DEX and control groups. The threshold of P-value < .05 was used for selecting DEGs. Moreover, the hclust method[ in R package was used to performing hierarchical cluster analysis for the DEGs. In addition, principal component analysis (PCA) was conducted using R software.[

Functional and pathway enrichment analysis

Using DAVID database (http://www.david.niaid.nih.gov),[ GO (gene ontology, including 3 items: molecular function (MF), biological process (BP), and cellular component (CC))[ functional and KEGG (Kyoto Encyclopedia of Genes and Genomes)[ pathway enrichment analyses for the DEGs were conducted. The significant terms were selected with the thresholds of gene count ≥ 2 and P-value < .05. To obtain more information, Enrichr method[ (http://amp.pharm.mssm.edu/Enrichr/) was also used for performing GO[ functional, as well as BioCarta,[ KEGG,[ Reactome, NCI-Nature_2016,[ and Humancyc_2016[ pathway enrichment analyses, with P-value < .05 as the threshold.

PPI network analysis

Search Tool for the Retrieval of Interacting Genes (STRING) database, which integrates both known and predicted PPIs of multiple organisms, can be applied for predicting functional interactions of proteins.[ Using STRING database (http://string-db.org/),[ PPIs were predicted for the proteins corresponding to the DEGs. The threshold for PPI prediction was combined score > 0.4. Then, Cytoscape software (http://www.cytoscape.org)[ was used to construct the PPI network. Proteins in the PPI network were defined as network nodes, and the number of interactions involved them was taken as their degrees. In the PPI network, the hub nodes[ were defined as nodes with relative higher degrees.

Regulatory network analysis

Based on the ENCODE_and_ChEA_Consensus_TFs_from_ChIP-X, ENCODE_TF_ChIP-seq_2015 and ChEA_2016 databases in Enrichr method,[ transcription factors (TFs) were predicted for the DEGs, with the cut-off criterion of P-value < .05. Besides, TargetScan software (http://www.targetscan.org/vert_71/)[ was used to predict microRNAs (miRNAs) for the DEGs. The P-value < .05 was set as the threshold. Finally, both TF-target and miRNA-target regulatory networks were visualized by Cytoscape software.[

Mining of deep knowledge

Based on Molecular Signature Database (MsigDB, http://www.broadinstitute.org/gsea/msigdb/index.jsp),[ DEX-induced changes of biological functions, signaling pathways, and targets of TFs and miRNAs were identified by GSEA analysis[ through performing randomized permutation test for whole transcriptome for 1000 times. The threshold was set as adjusted P-value < .05.

Validation of the expression changes of key genes

Using qRT-PCR experiments, the expression levels of serine racemase (SRR) and potassium voltage-gated channel subfamily J member 4 (KCNJ4) were also detected. The primers designed for SRR and KCNJ4 were listed in Table 1. The experiments were performed as described above.

Statistical analysis

The 2−ΔΔCt method[ was utilized for calculating the expression levels of genes. All data were displayed as mean ± standard error of mean (SEM). Data analysis was performed using Graphpad prism (Graphpad Software, San Diego, CA). The threshold for significant difference was set as P < .05.

Results

PCR detected the specific expression of HLA-E, indicating that the cells were trophoblast cells of human placenta. qRT-PCR showed that the expression of both GLUT1 (Fig. 1A) and GLUT3 (Fig. 1B) were downregulated by DEX, especially in trophoblast cells treated by DEX for 24 hours (P < .05). Therefore, the optimum treatment time of DEX was considered as 24 hours and used for the following experiments.

Figure 1

The relative expression of GLUT1 (A) and GLUT3 (B) in trophoblast cells treated by dexamethasone (DEX) for 12 hours or 24 hours, as well as the relative expression of SRR (C) and KCNJ4 (D) in DEX and control groups. DEX = dexamethasone, GLUT1 = glucose transporter 1, KCNJ4 = potassium voltage-gated channel subfamily J member 4.

DEGs analysis

After preprocessing of the raw data, a total of 391 DEGs (174 upregulated and 217 downregulated) were identified in the DEX group in relative to the control group. The number of upregulated genes was less than that of downregulated genes. PCA analysis indicated that the samples could be distinguished obviously (Fig. 2).

Figure 2

Principal component analysis (PCA) shows that dexamethasone and control samples can be separated obviously. PCA = principal component analysis.

Based on DAVID database, the DEGs were conducted with functional and pathway enrichment analyses. Functional enrichment analysis showed that the upregulated genes were enriched in PDZ domain binding (GO_MF; P-value = 2.82E−02; which involved SRR and KCNJ4) and Golgi cis cisterna (GO_CC, P-value = 4.43E−02). However, no significant pathways were enriched for the upregulated genes. For the downregulated genes, the enriched terms mainly included O-glycan processing (GO_BP, P-value = 3.25E−03), cytoplasm (GO_CC, P-value = 5.54E−04), amino acid binding (GO_MF, P-value = 2.27E−02), and Mucin type O-glycan biosynthesis (pathway; P-value = 3.00E−03; which involved polypeptide N-acetylgalactosaminyltransferase, GALNT1) (Table 2 ).

Table 2

The GO terms (A) and pathways (B) enriched for the differentially expressed genes.

Meanwhile, enrichment analysis was also performed using Enrichr method. KEGG enrichment indicated that Glycan biosynthesis pathway was inhibited (Fig. 3A). Besides, NCI-Nature_2016 pathway enrichment showed that epidermal growth factor receptor (EGFR), mitogen-activated protein kinase (MAPK), mechanistic target of rapamycin (mTOR), and CD8 T signal transduction pathways were suppressed after DEX stimulation (Fig. 3B). Moreover, BioCarta pathway enrichment suggested that extracellular signal-regulated kinase (Erk), PI3K and mitogen-activated protein kinase/ERK kinase (Mek) signaling pathways were inhibited following DEX stimulation (Fig. 3C).

Figure 3

The results of KEGG (Kyoto Encyclopedia of Genes and Genomes) (A), NCI-Nature_2016 (B), and BioCarta (C) pathway enrichment analyses based on Enrichr method. KEGG =  Kyoto Encyclopedia of Genes and Genomes.

Network analysis and mining of deep knowledge

For the DEGs, a PPI network involving 200 nodes and 287 interactions was constructed. Importantly, cyclin-dependent kinase 2 (CDK2) and CDK4 had higher degrees in the PPI network. After TFs (such as activating transcription factor 2, ATF2; and E2F transcription factor 4, E2F4) and miRNAs (such as miR-16) were predicted for the DEGs, the TF-target (Fig. 4) and miRNA-target (Fig. 5) regulatory networks separately were constructed.

Figure 4

The transcription factor (TF)-target regulatory network. Red and blue represent upregulated genes and downregulated genes, respectively. Circles and diamonds separately represent genes and TFs. The shade of blue indicates fold change of downregulation. TF = transcription factor.

Figure 5

The microRNA (miRNA)-target regulatory network. Red and blue represent upregulated genes and downregulated genes, respectively. Circles and squares separately represent genes and miRNAs. The size of square indicates −log P-value of target genes.

Using GSEA analysis, DEX-induced changes of biological functions, signaling pathways, and targets of TFs and miRNAs were identified. Especially, the results showed that potassium channel activity was activated (Fig. 6) and ribosome pathway activity was suppressed (Fig. 7) by DEX.

Figure 6

The gene set enrichment analysis (GSEA) result of POTASSIUM_ION_TRANSPORT (A) and the heatmap for the genes enriched in POTASSIUM_ION_TRANSPORT (B). GSEA = gene set enrichment analysis.

Figure 7

The gene set enrichment analysis (GSEA) result of KEGG_RIBOSOME (A) and the heatmap for the genes enriched in KEGG_RIBOSOME (B). GSEA = gene set enrichment analysis.

qRT-PCR analysis to confirm the expression changes of key genes

qRT-PCR showed that the expression of both SRR (Fig. 1C) and KCNJ4 (Fig. 1D) were significantly upregulated (P < .05), which was in accord with the results of DEGs analysis.

Discussion

After the trophoblast cells of human placenta were identified, the optimum treatment time of DEX was determined as 24 hours by qRT-PCR analysis. In this study, a total of 391 DEGs were identified in the DEX group compared with the control group, including 174 upregulated and 217 downregulated genes. Both DAVID database and Enrichr method were used for performing enrichment analysis. The PPI network for the DEGs was consisted of 200 nodes and 287 interactions. In addition, GSEA analysis showed that potassium channel activity was activated and ribosome pathway activity was suppressed by DEX. qRT-PCR analysis confirmed that SRR and KCNJ4 were significantly upregulated in DEX group.

Functional enrichment analysis showed that upregulated SRR and KCNJ4 were enriched in PDZ domain binding. KCNJ2 gene is reported to have expression in human placenta and in all development and differentiation stages of cytotrophoblast cells.[ Lockridge et al[ demonstrated that SRR in peripheral tissues may be critical for glucose homeostasis, and D-serine in β-cells may act as an endogenous islet NMDA receptor (NMDAR) co-agonist. SRR mRNA is detected in placenta, and D-serine transported by amino acid transport system B0 (ATB0) has a higher circulating concentration in the fetus compared with the mother.[ Meanwhile, pathway enrichment analysis showed that downregulated GALNT1 was enriched in Mucin type O-glycan biosynthesis. The family of polypeptide N-acetylgalactosaminyltransferases (GALNTs) consists of 20 members, including GALNT1 to 14 and polypeptide N-acetylgalactosaminyltransferase-like 1 (GALNTL1) to L6.[ Both GALNT1 and GALNT2 are overexpressed in first trimester extravillus trophoblast (EVT) and HTR8/SVneo cells, and the initiating enzyme of O-glycosylation GALNT2 is important for mediating EVT invasion.[ These declared that DEX might affect glucose transport in placenta through regulating GALNT1 enriched in Mucin type O-glycan biosynthesis and SRR and KCNJ4 enriched in PDZ domain binding.

Our results indicate that CDK2 and CDK4 were hub nodes in the PPI network. The simultaneous cytoplasmic mislocalization of CDK2 and p27 functions in regulating cell proliferation, invasion, and migration of human trophoblast cells.[CDK2 expression patterns indicate the decreased cyclin-dependent kinase activity along with placental maturation, and the interactions of CDK2, cyclin E, E2F1, and p27 may mediate aberrant trophoblastic proliferations.[ The knockouts of CDK2 and CDK4 can lead to embryonic lethality through inducing heart defects, and result in retinoblastoma protein (Rb) phosphorylation via regulating E2F-inducible genes.[ The high activity of CDK2 contributes to quickening G1 phase progression and establishing embryonic stem cell (ESC)-specific cell-cycle structure in mouse ESCs.[ Therefore, CDK2 and CDK4 might also be implicated in the roles of DEX in trophoblast cells.

Moreover, TFs ATF2 and E2F4 and miRNA miR-16 were predicted for the DEGs. Expression pattern changes of ATF2, cAMP response-element modulator protein (CREM), and cAMP response-element binding protein (CREB) have correlations with myometrial quiescence in the process of fetal maturation and uterine activation at term.[E2F4 cooperates with pRB in placental development, and is essential for the development of multiple embryonic tissues.[ Zhang et al[ deem that E2F4 interacted with pRB mediates cell cycle exit and plays a cell-intrinsic role during fetal erythropoiesis. miR-16 and miR-21 expression levels are obviously decreased in low birth weight infants and is correlated with fetal growth.[ These suggested that ATF2, E2F4, and miR-16 might also be targets of DEX during fetal maturation.

In conclusion, a total of 391 DEGs were identified in the DEX group. Besides, DEX might affect the birth weight of infant through mediating GALNT1, SRR, KCNJ4, CDK2, CDK4, ATF2, E2F4, and miR-16. However, our results were obtained from bioinformatics analysis and in vitro experiments of trophoblast cells. Therefore, more in-depth experiments should be conducted to further confirm these findings.

Table 2 (Continued)

The GO terms (A) and pathways (B) enriched for the differentially expressed genes.