Association between changes in salivary microbiota with glucose levels during pregnancy – findings from a pilot study

Background : Microbial shifts that correspond to host variations during pregnancy are vital in health maintenance. Significant changes have been reported in the oral microbiota of pregnant women when compared with nonpregnant women, but little is known about the dynamic shifts in oral microbiota during the pregnancy course. Methods : In this study, changes in salivary microbiota in 81 healthy pregnant women throughout the early stage (G1: 9-14 weeks), middle stage (G2: 21-28 weeks), and late stage (G3: 31-38 weeks) were investigated with 16S rRNA sequencing techniques. Correlations between salivary microbiota and maternal characteristics, including fasting blood glucose (FBG) levels, were evaluated. Results : Alpha diversity indexes were stable throughout pregnancy, but significant changes were found in beta diversity measured by weighted and unweighted UniFrac distances. Fourteen dominant trimester-specific taxa were identified using the LEfSe method, including Bacteroidetes in G1, Proteobacteria in G2 and Firmicutes in G3 at the phylum level. Tax4Fun prediction analysis revealed significant changes in Genetic Information Processing, Environmental Information Processing, Unclassified and Human Diseases in G2 and in Metabolism in G3 when compared to G1. Significant correlations were found between FBG levels and microbial composition, and these correlations were independent of gestational diabetes mellitus (GDM) status. Conclusion : Within the limitations of this study, the dynamic changes in salivary microbiota during pregnancy were characterized, and beyond pregnancy, FBG was also involved in shaping the salivary microbiota.

Pregnancy is a complex biological process involving profound hormonal, immunological and metabolic changes. It is reasonable for the oral microbiota to coordinate to respond to these changes [16]. Relatively few studies have investigated the oral microbiota relevant to pregnancy. Both the quantity and composition of the oral bacteria change significantly in the pregnant state when compared to the nonpregnant state, and more periodontal pathogens, including Porphyromonas gingivalis (Pg), have been detected [17,18]. Four studies investigated the changes in several microbiota, including that of the oral cavity, during pregnancy [18][19][20][21]. Stable microbial maintenance in the saliva and teeth/gum was reported by Daniel B.
DiGiulio et al. throughout pregnancy based on the 16S rRNA amplification technique [19]. In a second study from the same group, a subset of 10 samples from the same dataset were further analyzed using a Metagenomic method. This method found relative stability over time within subjects but high interindividual variability in the saliva microbiota [21]. However, in another small-sample-size study by Wenzhen Lin, differences in supragingival microbial structure between the trimesters were observed [18]. Changes in the microbial community in a compound sample mixed from tongue, hard palate and the gum line between trimesters were also identified by L. Dunlop et al.'s study [20]. In consideration of the limited number of studies and the inconsistent conclusions, shifts in oral microbiota composition during pregnancy need to be further explored.
Metabolic changes are recognized as one of the main factors affecting microbial features during pregnancy [16]. In many ways, the metabolic changes associated with pregnancy are similar to those that occur in metabolic syndrome, including weight gain, elevated fasting blood glucose levels, insulin resistance and low-grade inflammation [9]. Obesity is related to the diversity and structure of the salivary microbiome [22]. The oral microbiota could coordinate with glucose metabolism during pregnancy, considering the varying degrees of associations between gestational diabetes mellitus (GDM) and the saliva microbiota [23][24][25].
Hence, the present study was performed to investigate the changes in the oral microbiota of healthy Chinese women during pregnancy and to analyze the possible relationship between the oral microbiota and general maternal characteristics as well as metabolic parameters, including FBG levels.

Study Population
This study was based on a whole-pregnancy follow-up cohort (The gastrointestinal microbial cohort study during pregnancy). Healthy participants who met the inclusion and exclusion criteria were recruited from the clinic of the Obstetrics Department of Peking Union Medical College Hospital at their first visit during pregnancy. They were all aged from 20 to 45 years, had their first visit before the 13 +6 week of gestation, and committed to regular visits. Subjects who had digestive system diseases; metabolic diseases including hypertension, diabetes and dyslipidemia; or immune system diseases and tumors and subjects who had taken antibiotics or probiotics within the preceding 3 months were excluded. Participants who were diagnosed with eclampsia or preeclampsia during pregnancy and those who failed to provide salivary samples from all three trimesters were also excluded. Finally, 81 pregnant women were included in this study. All of them were free of complications including hypothyroidism, acute fatty liver and intrahepatic cholestasis.
Each participant completed a questionnaire at the beginning of the study, and their age, height, prepregnant weight, prepregnant body mass index (preBMI), gravidity, parity and number of fetuses (N.Fetus) were recorded. Clinical information, including gestational weight gain (GWG), levels of fasting blood glucose and disorders of pregnancy such as gestational diabetes mellitus (GDM), was obtained from their medical records. A GDM diagnosis was made according to the criteria of the World Health Organization [26]. This study was conducted with the informed consent of all participants and was approved by the Ethics Committee of Peking Union Medical College Hospital (JS-1535).

Saliva collection
Saliva samples were collected from each participant three times in the early (G1: 9-14 weeks gestation), middle (G2: 24-28 weeks gestation) and late (G3: 32-38 weeks gestation) stages of pregnancy, respectively. Sampling was performed at the Obstetrics Clinic at 8-10 am. The participants were asked not to brush their teeth on the sampling day. After rinsing the mouth with water, approximately 1.5 mL of unstimulated saliva was collected into a 5 mL sterile Eppendorf tube, transported on ice to the laboratory, and then kept frozen at −80⁰C until use.

DNA Extraction, Amplification, and Sequencing
Total genome DNA from the samples was extracted using the CTAB method as previously described [27]. The DNA concentration and purity were monitored on 1% agarose gels. 16S rRNA V4 hypervariable regions were amplified used specific primers (16S V4: 515F-806R) with barcodes. All PCRs were carried out with a Phusion® High-Fidelity PCR Master Mix (New England Biolabs). After quantification by 2% agarose gel electrophoresis, PCR samples with a bright main band of 300 bp were retrieved with a Qiagen Gel Extraction Kit (Qiagen, Germany) for further experiments. Sequencing libraries were generated using a TruSeq® DNA PCR-Free Sample Preparation Kit (Illumina, USA) following the manufacturer's recommendations, and index codes were added. The library quality was assessed on the Qubit@ 2.0 Fluorometer (Thermo Scientific) and Agilent Bioanalyzer 2100 system. Then, the library was sequenced on an IlluminaHiSeq2500 platform, and 250 bp paired-end reads were generated.

Data analysis
Raw sequencing data were processed according to the QIIME (V1.9.1, http://qiime.org/index.html) quality controlled process. Chimera sequences were then detected and removed with the help of the UCHIME algorithm (UCHIME Algorithm, http://www.drive5.com/usearch/manual/uchime_algo.html) [28]. All of the remaining sequences were clustered into operational taxonomical units (OTUs) by Uparse software (Uparse v7.0.1001， http://drive5.com/uparse/) [29]. Sequences with ≥97% similarity were assigned to the same OTUs. Species annotations were performed against the Silva132 database (http://www.arb-silva.de/). Multiple sequence alignment was conducted using the MUSCLE software (Version 3.8.31) in order to study the phylogenetic relationships of the different OTUs and the differences of the dominant species in the different samples (groups) [30]. The OTU abundance information was normalized using a standard of sequence number corresponding to the sample with the least sequences. Subsequent analyses of alpha diversity and beta diversity were all performed basing on this output normalized data.
Alpha diversity was applied in analyzing the complexity of species diversity for a sample through 4 indexes, including Observed-species, Chao1, Shannon, and Simpson. All these indexes in our samples were calculated with QIIME (Version 1.9.1) and displayed with R software (Version 2.15.3). Changes among the three trimesters were analyzed using two-way ANOVA by SPSS 13.0 software. Alpha diversity indexes were set as the dependent variables. The subject and trimester were treated as the fixed factors and were set as the main effect. Beta diversity analysis was used to evaluate the differences between samples in terms of species complexity. Beta diversity on both weighted and unweighted Unifrac distances was calculated by QIIME software (Version 1.9.1). PCoAs were performed using the ade4 and vegan package and visualized using ggplot2 in R software. For each grouping variable, 95% confidence ellipses were calculated. Differences in microbial communities during pregnancy were assessed by analysis of molecular variance (AMOVA) based on weighted or unweighted UniFrac distances using the AMOVA function [31]. A linear discriminant analysis (LDA) effect size (LEfSe, http://huttenhower.sph.harvard.edu/galaxy) method was used to find biologically relevant features, which highlights biological consistency, statistical significance, and effect relevance [32]. The threshold of the logarithmic LDA score for discriminative features was >4.0. Spearman's rank correlation was utilized to analyze the relationship between salivary microbiota and clinical characteristics. The statistical calculations were performed and visualized using R software (Version 2.15.3). Partial correlation coefficients were calculated using the package "ppcor" [33]. Canonical correspondence analysis (CCA) [34] was performed using the package "vegan" to evaluate the association between maternal parameters and microbial distributions.
Clinical factors screened by variance inflation factor analysis (VIF) [35] were used.
The clinical factors were fitted to the ordination plots using the 'envfit' function of the Vegan package in R [36]. Microbial function prediction was performed using package Tax4Fun in R software against the SILVA database. Differences between groups were detected by Student's t-test. values <0.05 were considered statistically significant.

Basic clinical characteristics
A total of 81 pregnant women were included in this study. Their basic clinical characteristics are shown in Table 1 Changes in Alpha diversity indexes including the observed_species, Chao1, Shannon and Simpson were evaluated during pregnancy. As shown in Table 2, although a decreased tendency in observed species index and Chao1index was found, none of the four indexes showed a significant variation with gestational progress. Beta diversity analysis was performed and demonstrated by PCoA. Unweighted and weighted UniFrac distances were used. AMOVA analysis was first carried out and displayed significant differences between every two gestational trimesters on both unweighted and weighted conditions (Table S1). From the PCoA results, although the G1, G2 and G3 samples were not clearly segregated, a trend of separation was shown by unweighted and weighted UniFrac distances ( Figure 1).
For more details to explain the microbial differences between trimesters, dominant trimester-specific taxa were scanned using the LEfSe method. Taxa with LDA score >4 were displayed in an LDA bar graph, and a circular cladogram was generated to show the differentially abundant taxa. Compared to the taxa in the other two trimesters, four taxa in G1, 3 taxa in G2 and 7 taxa in G3 were significantly enriched. Respectively, they were Bacteroidetes in G1, Proteobacteria in G2 and  (Table S2).
To further investigate the effect of FBG on the salivary microbiota, Spearman analysis of the association between FBG and the top 10 phyla was performed. and Methylococcaceae respectively). Besides, significant correlations were also found between Neisseriaceae and FBG (partial correlation analysis, r=-0.142, p=0.029) (Table S3A). At the genus levels, 31 genera by spearman analysis and 28 genera by partial correlation analysis were identified to be significantly correlated to FBG (Table S3B).

Influence of GDM on the salivary microbiota during pregnancy
FBG levels were not only related to the pregnancy stage but were also influenced by GDM. The salivary microbiota of GDM+ and GDM-women were compared. Among the 81 pregnant women, 14 of them developed GDM in their second trimester. None of them reported a history of probiotics or received any insulin treatment. Both the GDM+ group and the GDM-group displayed a significant increase in FBG levels when gestation entered the third stage. FBG levels in the GDM+ group were significantly higher than they were in the GDM-group at every counterpart stage (Mann-Whitney U-test, p<0.05) (Table S4).
Alpha diversity measured by the Chao1 and Shannon indexes was compared between the GDM+ and GDM-groups (Table S5). No significant differences between these two groups for the Chao1 and Shannon indexes in any trimester were found, nor were any significant changes during pregnancy observed in either group (Mann-Whitney U-test, p>0.05). AMOVA was performed to compare the unweighted and weighted UniFrac distances between the GDM+ and GDM-groups and found no significance in the three trimesters (Table S6). Significant changes during pregnancy were also found in healthy women (Table S7).

Discussion
This longitudinal study with 81 healthy pregnant women was performed to investigate the changes in salivary microbiota during pregnancy and to identify the maternal factors associated with those variations. Significant changes were found in community structures and functional pathways among the three trimesters. A total of 14 trimester-specific dominant taxa were identified and might contribute to the changes in the pathways. Among the maternal characteristics, including age, preBMI, GWG, gravidity, parity, and N.Fetus, a significant correlation between FBG and microbial composition was identified. The effect of FBG on salivary microbiota was independent of GDM.
During pregnancy, dramatic weight gains and hormonal, metabolic, and immune changes occur. Microbial changes would be reasonable and are considered vital for a healthy pregnancy [16]. Only a few previous studies have focused on oral microbial alterations during pregnancy. When measured by alpha diversity, almost all of these studies reported stable maintenance of the oral microbiota [18][19][20]. In our study, no DiGiulio et al. [19]. According to our results, a separation, although not completely clear, of the saliva microbial communities was seen among the three gestational trimesters by PcoA results. This verified the findings of the supragingival microbiota by Wenzhen Lin et al., which were based on 11 Chinese pregnant women and displayed segregation of microbial communities during pregnancy [18]. Significant differences between every two trimesters were confirmed by AMOVA, whether unweighted and weighted UniFrac distances were used. This result was consistent with the findings obtained from a compound sample mixed from tongue, hard palate and gum line [20].
Tax4Fun prediction results demonstrated that a clear variation occurred in the second trimester. Metabolism of other amino acids, Cellular processes and signaling and Neurodegenerative diseases were significantly increased, while Carbohydrate metabolism, Nucleotide metabolism, Replication and repair and Endocrine and metabolic diseases were significantly decreased. When gestation progressed to the third trimester, pathway abundance changed to levels comparable to the first stage.
Energy metabolism was dominant in G1 and Signaling molecules and interaction was most abundant in G3. Four abundance dominant taxa in G1, 3 in G2 and 7 in G3 were identified by the LEfSe method. The shared variation tendency during pregnancy along with the significant correlations between those changed taxa and functional pathways indicated a possible contribution of these taxa to the functional prediction changes. However, predicted by 16s rRNA sequencing analysis, these findings were informative, further analysis such as shot-gun metagenomic sequencing would be required in order to confirm these observations.
But controversy still exists with regards to the FBG changes throughout pregnancy.
According to our result, a significant decrease was found in G2 when compared to G1, this was in accordance with most of the studies [38]. When gestation came to G3, FBG was significantly increased and peaked among the three trimesters. This was a relatively rare report. Only small but not significant increases of FBG in the third trimester than the second trimester were reported previously [39,40]. In addition to gestational weeks, FBG had a significant effect on salivary microbiota distribution.
When measured by glycoalbumin, glucose metabolism was also reported to be related to gut microbial diversity, although only a minor aberrance was found in the gut microbiota during pregnancy [41]. To further investigate the effect of FBG on saliva microbiota, correlation analysis between FBG and the top 10 phyla was performed.
Significant associations were found between Proteobacteria and Firmicutes and FBG levels, and both phyla belonged to the core oral microbiome [1] and were differently expressed during pregnancy, with Proteobacteria at its highest abundance in G2 and Firmicutes at its highest abundance in G3. An effect of FBG levels on changes of salivary microbiota during pregnancy was indicated. The increase in FBG with pregnancy was significantly associated with a decrease in Proteobacteria and an increase in Firmicutes. After controlling for the factors of maternal age, preBMI, GWG, gravidity, parity and N.Fetus, the association between FBG and Firmicutes was not significant, which means that their association might be affected by those confounders. Proteobacteria is involved in inflammatory response and insulin resistance [42,43], and a close relationship between Proteobacteria and glucose metabolism was displayed. Pasteurellaceae, Rhodocyclaceae, Legionellaceae, Methylococcaceae and Neisseriaceae of Proteobacteria phylum were next identified to be responsible for this correlation.
As FBG was also affected by GMD status, we analyzed the influence of GDM on the changes in the microbiota. The variation in salivary microbiota during pregnancy was not related to GDM, and no significant differences in microbial alpha nor beta diversity were found between the GDM+ and GDM-groups. This was in accordance with the results of Mie K. W. Crusell et al. In their study, the link between GDM and salivary microbiota was evaluated. Contrary to the significant changes from late pregnancy to postpartum, no significant variations were found between GDM patients and healthy controls, whether in their third trimester or postpartum, when measured by alpha and beta diversity indexes [24]. However, in another study by Yajuan Xu et al., a significant decrease in salivary microbiota in alpha diversity was found in GDM patients when compared with healthy controls [23]. Given the limitation of the small number of GDM patients in our study, the association between GDM and oral microbiota needs to be further explored.
The present study was subject to several limitations. Although dietary habits may affect oral microbiota, the diets of the participants were not evaluated in the study.
Oral health, especially gingival health, is recognized to be closely related to pregnancy. Gingivitis is more common in pregnancy, and its severity can be aggravated during the progression of pregnancy [6]. Another important weakness in our study is the fact that oral health conditions during pregnancy were not examined or monitored, which could potentially influence our observations. None of the mothers, however, reported oral pain or discomfort during the study, and whether changes in salivary microbiota correlated with oral health status needs to be further investigated.

Conclusions:
In conclusion, significant changes were found in beta diversity of salivary microbiota during pregnancy. The most abundant expression was disclosed for Bacteroidetes in CCA: canonical correspondence analysis; VIF: variance inflation factor analysis.

Consent for publication
Not applicable

Availability of data and materials
The sequence data of the 243 salivary samples has been submitted to the SRA database and the accession number was PRJNA663243.

Ethics approval and consent to participate
This study was conducted with the informed consent of all participants and was approved by the Ethics Committee of Peking Union Medical College Hospital (JS-1535).

Competing interests
The authors declare that they have no competing interests.

Funding
This study was supported by the CAMS initiative for Innovative Medicine XZ and PW wrote the paper. All authors contributed to manuscript revision, read, and approved the submitted version  Values are presented as the mean±SD G1, G2, and G3 stand for early, middle and late pregnancy, respectively.
Changes in alpha diversity indexes during pregnancy were evaluated by two-way ANOVA. Alpha diversity indexes were set as the dependent variables. The subject and trimester were treated as the fixed factors and were set as the main effect. respectively. The differentially expressed function prediction identified by paied-T tests between G1 and G2 (B) , between G2 and G3 (C) and between G1 and G3 (D) . Bold numbers indicate a statistical significance (p < 0.05).