Gene expression in term placentas is regulated more by spinal or epidural anesthesia than by late-onset preeclampsia or gestational diabetes mellitus

Pre-eclampsia (PE) and gestational diabetes mellitus (GDM) are common complications of pregnancy, but the mechanisms underlying these disorders remain unclear. The aim was to identify the extent of altered gene expression in term placentas from pregnant women with late-onset PE and GDM compared to controls. RNAseq identified few significantly differentially regulated genes in placental biopsies between PE, GDM, or uncomplicated pregnancy (n = 10 each group). Five genes were altered in placentas from PE including 4 non-coding genes and Angiopoietin 2 (ANGPT2). No genes were significantly regulated by GDM. In contrast, many genes were significantly regulated by fetal, maternal and delivery-specific variables, particularly spinal and epidural anesthesia. We selected ANGPT2 and Chemokine (C-X-C motif) ligand 14 (CXCL14) to test with qPCR in a larger set of placentas (n = 475) and found no differences between the groups. However, regression analysis revealed a stronger association between placental ANGPT2 and CXCL14 mRNA expression and fetal, maternal and delivery-specific variables than diagnostic group. To conclude, the gene expression in term placentas are highly affected by fetal, maternal and delivery specific variables. Few regulated genes were found in late-onset PE and GDM placentas, which may suggest that these conditions could be more affected by maternal factors.


Results
Clinical characteristics. Table 1 shows the characteristics of the study population chosen for RNAseq (n = 30) and validation with RT-qPCR (n = 475). In the RNAseq cohorts, as expected, women with GDM had a significantly higher BMI and placental weight and women with PE had a higher DBP at visit 4. In the validation cohort, GDM and PE women had higher BMI, BP, more frequent induction of labor, and shorter duration of labor. In addition, GDM women were older, had higher birth weight and placental weight, had more cesarean deliveries, and were more often nulliparous. Gestational age was shorter in PE women. Women with PE had no indication of poor placental function based on birth weight and the relatively long gestational age underlining the mildness of the phenotype. Blood pressure in the PE patients were higher than 140/90, but we did not have access to exact values after visit 4, and we could therefore not include severity of PE women based on blood pressure.  Smokers a n (%) 1 (10) 3 (30) 1 (10) 73 (20) 29 (30) 3 (16) Primipara b n (%) 3 (  Using World Health Organization (WHO) as opposed to International Association of Diabetes and Pregnancy Study Groups (IADPSG) criteria for the diagnosis of GDM, the differences with controls were somewhat attenuated (Supplemental Table 1). Supplemental Table 2 compares some maternal, fetal and delivery specific variables between women with and without epidural or spinal anesthesia demonstrating that these groups were comparable for most variables. However, a higher umbilical cord venous base excess and a higher incidence of cesarean deliveries were observed in the spinal group compared to women without anesthesia.
Highest expressed genes in term placentas. The genes were ordered by their median expression level.

Effects of maternal, delivery-specific, and fetal variables on placental gene expression by RNAseq.
We evaluated potential maternal, delivery-specific, and fetal confounding variables on placental DEGs with FDR < 0.1. Investigating maternal confounding variables using all 30 samples, 16 DEGs were detected between primiparous compared to multiparous women. Of these, 15 were downregulated and 1 was upregulated in the multiparous women (Fig. 2). Some of these genes are involved in placentation and inflammation, i.e. pentraxin 3 (PTX3), chemokine (C-X3-C motif) receptor 1 (CX3CR1), and TSC22 domain family, member 3 (TSC22D3). Nine DEGs were detected between low and high BMI patients at visit 4. Of these, 8 were upregulated and 1 downregulated in the placenta from women with higher BMI, including increased adipsin (CFD). For BMI visit 1 and age of the mother, we did not find any significantly DEGs, FDR < 0.1 (data not shown).
For the delivery-specific variables (Fig. 3), we found 33 DEGs with epidural anesthesia (22 upregulated and 11 downregulated) and 35 DEGs with spinal analgesia (10 upregulated and 25 downregulated). Women who received general anesthesia (n = 1) and pudendal block (n = 1) were excluded from these analyses. Some DEGs were similar between epidural and spinal anesthesia (e.g., keratins), but there were also specific differences in the gene signature depending on the anesthesia type. In particular, 9 transcripts for heat shock proteins (Hsp) belonging to the Hsp70 and Hsp40 family were abundantly expressed and increased in placentas from women who received epidural anesthesia. In addition, myeloperoxidase (MPO), a marker of neutrophil activation, was enhanced in women who received epidural anesthesia. In women who received spinal anesthesia, genes involved in the regulation of immune responses including several immunoglobins (e.g. IGHG2, IGJ, IGCL2, MARCH1) and chemokines (CXCL9/10) were differentially increased. With regard to differences relating to type of delivery, 1 DEG was downregulated in the women who delivered by cesarean, and 7 DEGs were upregulated in the placentas of women with a longer duration of labor, while no DEGs were associated with induction of labor. For fetal confounding variables ( Fig. 4), 51 DEGs were associated with offspring's sex. Thirty-seven were upregulated and 14 were downregulated when using males as the reference. Most of these were sex-specific chromosomal genes. We found 4 DEGs between the median high and low placental weights. Three were upregulated and 1 was downregulated in those with the highest placental weight. For birth weight, no significant DEGs were observed (Fig. 4). A list of the 70 most DEGs for the confounding variables can be found in Supplemental Table 3.

Effect of PE and GDM on placental gene expression by RNAseq. Only 5 DEGs were identified
between the control and PE placentas with a FDR < 0.1. Of these, there were 3 downregulated non-coding genes, one downregulated protein-coding gene (Angiopoietin 2 [ANGPT2]), and one upregulated non-protein coding gene. Median counts (25 th -75 th percentile) for ANGPT2 in control vs. PE placentas were 663 (224-933) and 239 (155-386), respectively. No DEGs were detected between GDM and control placentas with FDR < 0.1. The 70 most DEGs for the diagnostic groups based on FDR and unadjusted p-value are listed in Supplemental Table 4.

Validation of RNAseq data by RT-qPCR.
To validate our RNAseq data, we performed RT-qPCR. Samples for RT-qPCR validation included not only the 30 samples included in the RNAseq analysis, but all 475 placental samples from the STORK study with high RNA quality. However, since few genes were differentially expressed, in addition to ANGPT2, we chose to measure CXCL14 which was on the top of the lists for both PE and GDM (unadjusted p-value < 0.0006 and <0.001 respectively), and may be known for potentially regulate trophoblast outgrowth at the maternal-fetal interface 19 . As seen in Supplemental Fig. 3, differential (decreased) expression of ANGPT2 in PE placentas was confirmed by qPCR in the same 30 samples used in the RNAseq experiment. However, in the larger sample set, although ANGPT2 expression was decreased in PE placentas compared with controls, the difference did not reach statistical significance (P = 0.74). CXCL14 expression was not increased by Prediction of ANGPT2 and CXCL14 expression in placenta in the large cohort. To evaluate the influence of maternal, delivery-specific, and fetal characteristics on ANGPT2 and CXCL14 expression, we finally performed multivariable regression analysis, including diagnostic group (i.e. GDM or PE) to identify the most important predictors (Table 2). Stepwise linear regression identified gestational age, parity and offspring sex as predictors of ANGPT2 expression. Gestational age followed by birth weight and parity were the variables that predicted the gene expression of CXCL14.

Discussion
We conducted RNAseq analysis to compare gene expression between late-onset PE, GDM and control placentas at term. Using this approach, few genes were differentially expressed between these groups. In contrast, maternal, delivery-specific, and fetal variables were associated with altered expression of multiple genes in the term placenta. In particular, spinal and epidural anesthesia were unexpectedly associated with a significant upregulation of heat shock protection and immunomodulatory genes. When using qPCR in a larger cohort to validate the most highly differentially regulated genes in late-onset PE and GDM placentas identified in the RNAseq experiment (namely, ANGPT2 and CXCL14), we noted that maternal, delivery-specific, and fetal variables were stronger predictors of gene expression than the diagnostic group. These results suggest that maternal, delivery, and fetal variables substantially contribute to placental gene expression and should be considered when designing and evaluating such studies.
Samples chosen for RNAseq were well matched by clinical criteria and the RNA was of acceptable quality with a small percentage of contamination by maternal decidual and nucleated blood cells (~2%). This is similar to data published by Sõber et al. following a single wash with phosphate buffer saline and stabilization in RNA-later solution 14 . Furthermore, the most abundantly expressed genes identified in the term placenta were ones known to be involved in placental function, with good overlap with transcripts identified in previous studies 14 . For the  20 , respectively. Other DEGs associated with offspring sex in the current study were either non-functional genes expressed at low levels or sex specific genes (i.e. TSIX, TXNLGY, DDX3X, TTTY14, TBL1X, TBL1Y and TTTY10) not reported in the prior manuscripts. Totally from the 51 DEGs, 28 were located on the Y-Chromosome while 12 were located on the X-Chromosome. Taken together, the similarities in gene expression, related to abundance and offspring sex, with previous studies as well as the low maternal contamination supports the validity of our placental RNA samples.
A major finding in this study was the significant impact of epidural and spinal anesthesia on placental gene expression. To our knowledge, this has not been previously reported. Since women receiving epidural may be more likely to experience hyperthermia 21 , the increased expression of mRNAs coding for members of the HSP family could reflect the heat shock response, a switch in transcription and translation to preferentially express HSP's that may preserve essential proteins 22 . Other relevant causes of increased HSP expression due to different kinds of stress could be infection, inflammation, hypoxia, injury and also adrenaline, which is a component in epidural anesthesia 23,24 . In particular, we identified a higher expression of the genes encoding HSP70, a HSP frequently studied in placental disease that may also confer harmful effects on binding to toll-like receptors 25 . Of note, increased expression of HSP70 protein has been demonstrated in the placenta of women with PE 26,27 and in relation to delivery-specific variables, such as preterm labor 28 . In contrast, we found no significant association between HSP expression and mode of delivery, induction and duration of labor, or presence of PE. Several studies suggest that underlying maternal inflammation may confer increased risk of fever following epidural anesthesia 29,30 . However, apart from increased MPO mRNA levels in women receiving epidural, possibly reflecting  neutrophil activation during the hyperthermic response 31 , few inflammatory genes were differentially regulated. In contrast, the use of spinal anesthesia was associated with enhanced mRNA expression of several inflammatory genes, including IGJ, CXCL9 and CXCL10, all chemokines specific for T cell signaling through the CXCR3 receptor. Activation of CXCR3 signaling has been implicated in the onset of labor 32 . Raman et al. identified enhanced IGJ expression as the most distinguishable feature of chronic placental inflammation in addition to CXCL9 and other inflammatory components from T and B cells 33 34 and in the spinal cord in experimental animal models in vivo 35 . The use of anesthesia may reflect underlying obstetric conditions and we cannot excluded that the differences in mRNA expression observed for use of epidural and spinal anesthesia may results from longer acting mechanisms and systemic stress responses of the mother.

. Activation of opioid receptors increases the expression of CXCL10 in circulating immune cells in vitro
For maternal confounding variables, the expression of genes involved in inflammation and placentation (i.e., CX3CR1 36 , TSC22D3 37 and PTX3) 38 were associated with parity, indicating that women who have had a prior pregnancy may respond differently to the presence of invading trophoblast compared to nulliparous women. Adipsin expression, which is shown to be abundant in adipose tissue, was increased in the placentas of women with high BMI at term, and placentas from obese women have been shown to secrete higher levels of adipsin 39 .
Surprisingly few genes were differentially expressed in the placentas of women with late-onset PE at term and none in the placentas of women with GDM. We found only five significantly differentially regulated genes in PE placentas based on FDR < 0.1. The result of two recent meta-analyses found few common genes regulated in PE placentas compared to controls 10,12 . Using unsupervised clustering of 7 microarray data sets, Leavey et al. 3 found a large degree of co-clustering of PE and control samples suggesting no distinct gene signature in PE. Direct comparison with previous studies is difficult since these mainly have focused on early-onset PE, a separate entity from late-onset PE as demonstrated in a recent micro-array study 40 . Furthermore, Sõber et al. identified an extensive shift in the placental transcriptome profile in late-onset PE using RNAseq 14 . For one specific DEG, ZDHHC8P1, we found increased expression in late-onset PE patients which were in contrast to their results were it was downregulated in their late-onset PE. However, their patients had a markedly more severe phenotype according to the ACOG 2013 criteria (hypertension with additional symptoms 41 compared to the milder PE phenotype in our patients. This could imply that milder forms of preeclampsia that meet diagnostic criteria may in fact involve no underlying placental pathology and therefore may technically represent gestational hypertension with concomitant kidney dysfunction. It's possible these cases may respond well to standard antihypertensive therapy and not require the termination of pregnancy. As for GDM, previous studies evaluating gene expression have been performed in mixed cohorts (i.e. women with both PE in addition to GDM) 15 , in women on insulin therapy 15,17 , using different diagnostic criteria, or with very small sample numbers 16 , making comparisons difficult. Since our patients were not on insulin, it may be reasonable to assume that they have a milder GDM phenotype compared to other studies 15,17 . Importantly, excluding some differences in BMI and placental weights (Table 1), the major confounding maternal, delivery-specific and fetal variables were similar between our three comparison groups (late-onset PE, GDM and controls). Furthermore, when validating ANGPT2 and CXCL14 gene expression in the larger cohort, maternal, fetal and delivery-specific covariates were by far the strongest determinants of their expression. Since the objective of expression analysis is to identify unique expression patterns for distinct disease entities and not necessarily the effect of confounding conditions, our results highlight the importance of well characterized phenotypes to allow for accurate matching and relevant comparison between groups.
Several limitations of this study should be emphasized. The quality of the RNA could have been better. Five out of 30 RNAseq samples had RIN values below 7, and additional 8 samples had values below 8. However, we judged them to be of acceptable quality as judged by electropherograms. The placenta is a heterogeneous tissue and the biopsy was taken from the placental parenchyma. Other parts and specific cells of the placenta may be more important when investigating genes regulated in PE and GDM. Only one biopsy was taken from each placenta and the collection was performed by different individuals. While the technique for placental biopsy was standardized, we cannot exclude the possibility that biopsies were taken at different sites and that regional differences may have accounted for differences in DEGs 42 . Also, only term placentas were used in our analysis. While this ensures that appropriate control placentas will be available for comparison, it does represent the endpoint of a chronic disease state. It is likely that mRNA expression profiles at an earlier stage will be more reflective of the origin and pathogenesis of the disease and may be more likely to yield a distinctive fingerprint. This is especially true of PE 3 . Also, as we studied late-onset PE, our results cannot be applied to the spectrum of PE. Finally, we have incomplete information about the precise timing from when the placenta was delivered to when the biopsy was taken and the tissue stored at −80 °C.
In conclusion, the current study shows that confounding variables have a greater influence on gene expression in term placenta than do underlying disease states (such as late-onset PE and GDM). We found few genes significantly differentially regulated between late-onset PE, GDM and control placenta, suggesting that maternal factors may be more important than placental (fetal) factors, in the genesis of these pregnancy-specific conditions.

Methods
Study population. The STORK study was a prospective cohort study with a longitudinal design in which 1031 low-risk women of Scandinavian heritage who planned to deliver at Oslo University Hospital, Rikshospitalet, between 2002 and 2008 were followed throughout their pregnancy. Details about the study have been previously published 43 . Exclusion criteria included multiple pregnancy, known pre-gestational diabetes, severe chronic medical conditions (such as lung, cardiac, gastrointestinal or renal diseases), and pregnancies complicated by major fetal malformations or aneuploidy. Patients were routinely monitored for the development of pregnancy complications, such as PE, GDM, or preterm labor. Demographic and clinical data (including use of anesthesia in labor, route and mode of delivery, and neonatal outcome) were abstracted from the medical records. The difference in the formula used at our hospital between the epidural and spinal anesthesia is adrenaline in the epidural anesthesia, in addition to opioids, both fentanyl (or sufentanil in spinal) and bupivacaine in both.
Written informed consent was obtained from all participants. All clinical investigations were conducted according to the principles in the Declaration of Helsinki. The study was approved by the Regional Committee for Medical Research Ethics of Southern Norway in Oslo, Norway.
Clinical Diagnosis of GDM. A 75 g oral glucose tolerance test (OGTT) was performed in the morning after an overnight fast on all women at 30-32 weeks of gestation. Venous blood samples collected into tubes containing EDTA were analyzed at point of care using an Accu-Check Sensor glucometer (Roche Diagnostics GmbH, Mannheim, Germany). Additional venous blood samples were allowed to clot for 30 min. The serum was separated by centrifugation for 10 min at 3000 g and stored at −80 °C. Glucose levels were also measured from frozen serum samples collected at 30-32 weeks using the hexokinase method (Hitachi Modular P800, Roche Diagnostics, Mannheim, Germany) at an accredited clinical chemistry laboratory at Oslo University Hospital Rikshospitalet, as previously reported 44 . GDM was diagnosed on a 75 g OGTT using both the new IADPSG criteria and the old WHO criteria as follows: (1) IADPSG criteria: fasting plasma glucose (FPG) of 5.1-6.9 mmol/L (92-124 mg/dL) and 1 h plasma glucose ≥10.0 mmol/L (≥180 mg/dL) or 2 h plasma glucose 8.5-11.0 mmol/L (153-198 mg/dL); and (2) WHO criteria: 2 h plasma glucose ≥7.8 mmol/L (140 mg/dL) 45 , as previously reported 44 .
Clinical Diagnosis of PE. PE was diagnosed by new-onset blood pressure ≥140/90 mmHg and significant proteinuria (urinary total protein/creatinine ratio >30 or +1 on urine dipstick). All cases were diagnosed after 34 weeks gestation (late-onset PE).
Collection, storage and RNA extraction of placental biopsies. Placental biopsies were collected after vaginal or cesarean delivery. One biopsy was taken per placenta. The STORK study set included 475 placental samples, with good RNA quality, in which 19 were from women with PE, 103 with GDM, and 353 controls (overview, Supplemental, Fig. 1). Blocks of 2-4 cm were taken from the placental parenchyma, briefly washed in phosphate buffer saline, snap frozen in liquid nitrogen, and stored at −80 °C until RNA isolation. Half of the biopsy was homogenized in TRIzol reagent (Invitrogen, Life Technologies) on ice with a tissue grinder (Sigma Aldrich). Total RNA was extracted using TRIzol reagent (Invitrogen, Life Technologies) and purified with RNeasy microkit columns (Qiagen, Netherlands). Purity and concentration of isolated total RNA was measured using Nanodrop ND-1000 Spectrophotometer (Thermo Fisher Scientific Inc., USA) and RNA integrity number (RIN) was estimated using Agilent 2100 Bioanalyzer (Agilent Technologies, USA). RNA samples with a RIN <4 were excluded from subsequent RNAseq and qPCR analysis 46 . RNAseq. Samples with a RIN greater than 4 were preferred for RNAseq analysis (median RIN 8.6; range (25 th , 75 th ), 7.9-9.2). The RIN values and gel/electrophoresis results from the samples used for RNAseq are shown in Supplemental Fig. 2. The samples chosen for RNAseq were matched (based on mothers age and gestational age) between the groups of normal (control), GDM pregnancy and pre-eclamptic pregnancy (10 samples in each group). Sequencing libraries were prepared from 500 ng of total RNA using the TruSeq RNA sample preparation reagents (Illumina, San Diego, California) according to the manufacturer's instructions, with fragmentation for 4 minutes at 94 °C. The libraries were sequenced using 125 bp paired-end sequencing on an Ilumina HiSeq 2000. We recorded an average 22.3 million (range, 20.2-24.6 million) paired reads per sample. Fastq files were generated using bcl2fastq (v1.8.4). Sequence reads were mapped to the reference genome (hg19/GRCh38) using TopHat2 (v2.0.13) and Bowtie2 (v.2.2.3.0). Library sizes and standard deviations for input into TopHat were calculated empirically by aligning 1000000 reads to an index built from human cDNA sequence. Sequence alignment was guided using only previously annotated gene models downloaded from Ensembl (www.ensembl.org; Homo_sapiens.GRCh38.79.gtf). On average, there was 73.6% concordant read pair mapping (range 69.6-76.6), with a mean unique mapping of 94.1%. Raw expression counts were calculated per gene using featureCounts (http://bioinf.wehi.edu.au/featureCounts/) and the same gtf file which was used for the read alignment.