Cervicovaginal Microbiome Composition Is Associated with Metabolic Profiles in Healthy Pregnancy

Humans have a unique vaginal microbiome compared to other mammals, characterized by low diversity and often dominated by Lactobacillus spp. Dramatic shifts in vaginal microbial communities sometimes contribute to the presence of a polymicrobial overgrowth condition called bacterial vaginosis (BV). However, many healthy women lacking BV symptoms have vaginal microbiomes dominated by microbes associated with BV, resulting in debate about the definition of a healthy vaginal microbiome. Despite substantial evidence that the reproductive health of a woman depends on the vaginal microbiota, future therapies that may improve reproductive health outcomes are stalled due to limited understanding surrounding the ecology of the vaginal microbiome. Here, we use sequencing and metabolomic techniques to show novel associations between vaginal microbes and metabolites during healthy pregnancy. We speculate these associations underlie microbiome dynamics and may contribute to a better understanding of transitions between alternative vaginal microbiome compositions.

many healthy women lacking BV symptoms have vaginal microbiomes dominated by microbes associated with BV, resulting in debate about the definition of a healthy vaginal microbiome. Despite substantial evidence that the reproductive health of a woman depends on the vaginal microbiota, future therapies that may improve reproductive health outcomes are stalled due to limited understanding surrounding the ecology of the vaginal microbiome. Here, we use sequencing and metabolomic techniques to show novel associations between vaginal microbes and metabolites during healthy pregnancy. We speculate these associations underlie microbiome dynamics and may contribute to a better understanding of transitions between alternative vaginal microbiome compositions. KEYWORDS Lactobacillus, longitudinal, metabolome, microbiome, pregnancy, vagina V aginal microbes sustain important physiologies and produce metabolites that can directly or indirectly affect maternal health and infant development during pregnancy. Perturbations to early-life microbiomes and associated metabolic dysfunction have been linked with allergy and autoimmune diseases such as asthma (1)(2)(3)(4). For example, regular prenatal and postnatal farm exposure, i.e., contact with a diversity of microbes during pregnancy and infancy, have been shown to reduce the incidence of chronic health diseases such as asthma and atopy (5). Moreover, recent research has supported the idea of fetal programming, a term describing the process by which the maternal microbiota, as well as maternal antibodies, prepare the infant immune system for the postnatal onslaught of colonizing microbes (6). Others have shown in mice that vaginal dysbiosis, induced by maternal stress, has the potential to negatively affect offspring metabolic profiles (7). Thus, maternal microbes, particularly those of the vaginal tract, are some of the first microbes the offspring will encounter and may be central to the study of early-life microbiome and immune development (8)(9)(10)(11)(12). Indeed, a recent large-scale study of 2,582 women, over 600 of whom were pregnant (a subset of whom were longitudinally sampled), provided evidence for vaginal microbiome restructuring during pregnancy toward a Lactobacillus-dominated community (13). This occurred early in gestation and was associated with a reduced vaginal microbiome metabolic capacity. Postpartum, irrespective of the mode of delivery, the vaginal microbiota resembled that of a gastrointestinal microbiome, likely due to microbial mixing during the birthing process (9), suggesting that both vaginal and gastrointestinal microbial seeding of the neonate occurs.
The human vaginal microbiome maintains low diversity in low-pH conditions and depends on host sugars as carbon sources, with less access to dietary and exogenous nutrients than the gut, skin, or oral cavity. Historically, vaginal microbial communities have been stratified based on hierarchical clustering of the taxa composition (12). Keystone species include Lactobacillus crispatus and Lactobacillus gasseri, which have been associated with maintenance of a simple vaginal microbiome by their production of bacteriostatic and bactericidal compounds (e.g., lactic acid and hydrogen peroxide) and maintenance of a low pH (14)(15)(16) numerically and functionally dominating their respective vaginal communities. A closely related species, Lactobacillus iners, has been associated with health-promoting benefits; however, its genome also encodes the capacity to promote microbiome perturbation by increasing vaginal pH and producing species-specific virulence factors (14,(17)(18)(19). Bacterial vaginosis (BV), the most common gynecological condition in reproductive-age women (20), is characterized by the presence of a more diverse vaginal microbiome and associated with adverse pregnancy outcomes, including preterm birth (21), endometritis (22,23), and spontaneous abortion (24)(25)(26)(27). Recently, vaginal microbial transplants have been successfully implemented as a treatment for intractable BV (28). Despite L. crispatus generally being regarded as a highly beneficial and dominant microbe throughout pregnancy, healthy women from different ethnic groups have markedly different species dominating the vaginal microbiome (15). In fact, many healthy women who lack BV symptoms have vaginal microbiomes dominated by microbes that are associated with BV (29), suggest-ing that taxonomy alone is insufficient to predict health outcomes and that microbial activities, including metabolic productivity, may offer a more contemporary view of microbiome function.
An untargeted, more global assessment of microbiomes and associated metabolites during pregnancy and early life is lacking. To address this gap in knowledge, we collected saliva, urine, and cervical vaginal fluid (CVF) from 18 mothers during each trimester of pregnancy and saliva and urine from offspring through their first year of life. Specifically, we were interested in how maternal CVF microbiome profiles are associated with metabolomic assessments of the same samples. Furthermore, we had the opportunity to examine whether maternal saliva and urine metabolome profiles relate to those of the infant in the first year of life. Here, we present DNA sequencing (amplicon and shotgun) and untargeted metabolomics to characterize microbial and metabolic features of the CVF microbiome throughout pregnancy to determine the composition of the vaginal microbiome from a cohort of healthy Caucasian and Hispanic women, longitudinally sampled throughout a healthy pregnancy. no significant differences in BMI (t test; P ϭ 0.28) or age (t test; P ϭ 0.89) between ethnic groups. Saliva and urine were collected at indicated intervals from each infant up until 1 year of age ( Fig. 1). Saliva, urine, and CVF samples were subjected to metabolomics analysis, whereas only CVF was used for sequence analysis. Sequence analysis included amplicon-based sequencing of the 16S rRNA gene (bacteria and archaea) and ITS2 (fungi) loci and shotgun metagenomic sequencing of the entire microbial community.
Vaginal microbiota support high abundances of Lactobacillus and Bifidobacteriaceae throughout pregnancy. Sixteen individuals (42 total samples) produced sufficient sequence reads for taxonomic assignment using the 16S rRNA gene. Amplicon sequencing stratified cervical samples into those where the most abundant taxon was Lactobacillus spp. (34/42; 81%) or Gardnerella spp. (8/42; 19%) ( Fig. 2A; Fig. S1A). The bacterial taxa in samples with abundant Lactobacillus spp. were significantly less evenly distributed (linear mixed-effects [LME] modeling; P ϭ 0.001; Fig. 2C), with Gardnerella vaginalis being the most abundant in seven of eight samples (88% relative abundance) and a Shuttleworthia taxon being most abundant in one sample (at 23% relative abundance). In samples where Lactobacillus spp. had the highest abundance, a single taxon comprised 50% or more of the sequencing reads (27/34; 79%). L. iners was the most abundant taxon detected in 14 samples from 7 subjects, with a median relative abundance of 79%. In 12 samples from 7 subjects, a Lactobacillus taxon, putatively identified as L. crispatus through metagenomic sequencing ( Fig. S1B and C), had a median relative abundance of 96% and persisted at a relative abundance greater than 90% in subjects 1088, 1120, and 1191. Altogether, the most abundant 3 taxa, L. iners (operational taxonomic unit 1 [OTU_1]), L. crispatus (OTU_2), and G. vaginalis (OTU_3), comprised 66% of the total bacterial sequencing reads.
Whole-genome shotgun sequencing produced, on average, 2.9 million paired-end reads per sample, which decreased to an average of 220,355 paired-end reads per sample following removal of reads that aligned to the human genome. Thirteen individuals (35 total samples) produced sufficient sequence reads for taxonomic assignment, which was concordant with 16S rRNA gene sequence results ( Fig. S1B and C).
Most of the reads classified as L. crispatus or L. iners mapped to a single metagenomic assembled genome, with completeness of 95.7% and 97.1% and redundancy of 0% and 1.4%, respectively.
Twelve samples from seven subjects produced ITS2 sequences ( Fig. 2A); we do not have quantitative data characterizing the abundance of bacterial or fungal biomass. Eleven samples from six subjects contained species of Candida, classified as C. albicans (Fig. S2A), the most abundant fungal taxon in these data. Shotgun metagenomics confirmed these results and allowed for additional identification of reads mapping to  (Fig. S2C). Subject 1,088 was the only participant to deviate from this trend, with a high relative abundance of Aspergillus during the first trimester of pregnancy (Fig. S2A). Alpha diversity indices based on 16S rRNA and ITS2 data, when available, were compared across trimesters. While some subjects exhibited qualitative evidence of compositional shifts in vaginal microbiota with advancing gestation ( Fig. 2A), we did not observe a significant difference in bacterial richness (number of observed OTUs) (LME; P ϭ 0.17), evenness (Pielou's evenness index) (30) (LME; P ϭ 0.46), or phylogenetic diversity (LME; P ϭ 0.21) across trimesters ( Fig. 2D to F).
Highly abundant bacterial taxa were significantly associated with community composition. A nonmetric multidimensional scaling (nMDS) plot of Bray-Curtis dissimilarities showed that vaginal communities clustered by their most abundant bacterium (Fig. 3A). This association between the most abundant bacterial taxa and microbial composition of the sample was significant and explained more than half of the variance using permutational multivariate analysis of variance (PERMANOVA) (R 2 ϭ 56%; P ϭ 0.0001). To account for repeated measures from longitudinal samples from the same individual, we also performed an LME, which required dimensional reduction (LME; R 2 ϭ 69%; P Ͻ 0.0001). Communities with abundant L. crispatus were more similar to each other, sharing more than 90% similarity, in comparison to communities where a different bacterial species was most abundant. While some individuals exhibited a relatively stable microbial community over time, others (6/16) experienced shifts in composition, resulting in a statistically significant change in Bray-Curtis dissimilarity on PC axis 1 between trimesters (Fig. 3B).
The five subjects with Candida detected in at least one of two longitudinally paired samples displayed a significant increase in intersample Bray-Curtis dissimilarity in their bacterial profiles over that interval (e.g., the intervals between trimester 1 to 2, 2 to 3, or 1 to 3), suggesting the presence of Candida may be associated with greater shifts in bacterial composition than those who had no Candida detected (Fig. S2B).  Table S1), with an additional 353 identified compounds. The CVF metabolome as assessed by both mass spectrometry methods did not significantly differ across trimesters (LME; P ϭ 0.6378 for GC-TOF MS; P ϭ 0.3942 for liquid chromatography-tandem mass spectrometry [LC-MS/MS]). This stability was even true for the subset of individuals who exhibited shifts in microbiota composition over trimesters (LME; P ϭ 0.6594 for GC-TOF MS; P ϭ 0.2482 for LC-MS/MS). CVF samples dominated by distinct bacteria exhibited significant differences in metabolic profiles (PERMANOVA; R 2 ϭ 12%; P ϭ 0.0195). A constrained, distance-based ordination plot recapitulated 67% of the community variation observed in the vaginal microbiota (Fig. 4). Superimposed on the ordination plot are GC-TOF MS predictor metabolites, calculated using the distancebased linear models (DISTLM) program in PRIMER-e. Indole-3-lactate (ILA) accounted for 27% of the variation observed in the vaginal microbiota data and was found to be more abundant in vaginal microbiota with abundant L. crispatus (Fig. 4B). Mannitol was also more abundant in samples dominated by L. crispatus (Fig. 4C). In parallel, we found that a pathway for mannitol-1-phosphate production is also more abundant in shotgun metagenomic data sets of CVF samples dominated by L. crispatus (Fig. S3). This linear model identified the top 10 annotated GC-TOF MS metabolites that were associated with variation in the microbial community composition, which are shown in Fig. 4A; these 10 metabolites together might explain almost 57% of the total variation in the microbial community composition. A permutated random forest recapitulated what we found in the DISTLM, identifying mannitol and indole-3-lactate as two of the top variables of importance, specifically for distinguishing microbiomes with high abun- dances of L. crispatus (Fig. S4A). To explore the ability to analyze the metabolome in high throughput, the same sample sets were analyzed using matrix-assisted laser desorption ionization mass spectrometry imaging (MALDI-MSI) (Table S1). Detected ions by MALDI were compared to those identified by GC-MS and LC-MS, and we found that ϳ55% of the metabolites identified had corresponding ions in the MALDI analysis (Table S1).
Metagenomics and functional potential of communities. Distinct functions were associated with each of the vaginal microbial community clusters (PERMANOVA; R 2 ϭ 70%; P ϭ 0.0001; Fig. 5). LEfSe identified several pathways that differed between L. crispatus and G. vaginalis, in particular, an enrichment of ammonia assimilation genes in G. vaginalis (Fig. S4B). Genes involved in mannitol metabolism were enriched in communities where L. crispatus was highly abundant (Fig. S3C and Fig. 4B). Searching the PATRIC database of all sequenced L. crispatus (64 genomes), L. iners (22 genomes), and G. vaginalis (127 genomes) strains revealed annotated genes for mannitol usage and transport for L. crispatus but not for L. iners or G. vaginalis.
Mothers and infants have significantly different saliva and urine metabolomes. Maternal and infant saliva and urine metabolomes were assessed with both GC-TOF MS and LC-MS/MS (lipidomics) in order to study the relationship between maternal and infant metabolomic compositions during early life (see Fig. 1). PERMANOVA analysis of lipidomics data from saliva samples showed the largest difference between mothers and offspring (PERMANOVA; R 2 ϭ 69%; P ϭ 0.001; Fig. S5A). A subset of 50 lipidomics metabolites with high mean abundance, 70% of which were unannotated, showed distinct profiles between mother and offspring salivary metabolomes (Fig. S6A). Likewise, GC-TOF MS salivary metabolomes were also significantly different between mother and offspring, but far less variation was explained (PERMANOVA; R 2 ϭ 12%; P ϭ 0.0001; Fig. S5B). Some metabolites, such as lactulose, were much more abundant in infants and largely absent in mothers (Fig. S6B). Maternal metabolomics profiles  Table S2). The infant metabolome for saliva and urine had little variance attributed to which subject donated the sample, but GC-TOF MS was able to detect a significant change between the infant urine metabolome at birth versus 6 and 12 months of age (PERMANOVA; R 2 ϭ 34%; P ϭ 0.0007; Table S2). Moreover, from lipidomics data, the infant metabolome profile seemed to converge on mothers' metabolomes as they aged, though more samples would be needed to confirm this finding (Fig. S5C). For both saliva and urine, GC-TOF MS-and lipidomicsdetected metabolites were more similar for mother-child pairs than for unrelated individuals (Fig. S7). Mantel tests to determine if intersample relationships were similar between chromatography methods (including both GC-TOF MS and lipidomics) showed a strong correlation between saliva samples and weaker correlations between urine and CVF (Table S3).

DISCUSSION
Exposure to the microbiome in early life is critical for immune and physiological development (1)(2)(3)(4), yet the factors that set this trajectory remain poorly understood. In this study, we followed the vaginal microbiome through the trimesters of pregnancy for 18 women, tracking changes in the bacterial communities with longitudinal samples and capturing their functional potential with metagenomic sequencing and multiple platforms to assess metabolomic profiles. The resolution provided by shotgun metagenomic sequencing allowed us to identify species and characterize the functional gene content of CVF microbiomes. An additional strength of this work is the strict inclusion criteria defining healthy pregnancy (see "Subject information" in Materials and Methods). Moreover, as part of an existing sample cohort, we had the opportunity to measure saliva and urine metabolomes from mothers and children. We aim to establish how the metabolome develops in the first year of life and how maternal-infant saliva and urine metabolomes relate. In our study, most healthy pregnant women exhibited a relatively stable vaginal microbiota throughout the trimesters of pregnancy, dominated by Lactobacillus or, in some cases, more diverse, Bifidobacterium-dominated microbiota. However, a subset of women exhibited compositional shifts in their CVF microbiota as pregnancy progressed, as has also been seen in other larger cohorts (13). We found several strong correlations between particular vaginal communities and metabolites, which may help us understand the physiology underlying distinct vaginal microbiota structures that were evident in our study. Lastly, vaginal microbiota composition predicted which metabolites were present in the CVF samples but not urine or saliva samples from the mothers or infants, suggesting that local microbial metabolism may represent the dominant contributor to the metabolic milieu of the vaginal tract during pregnancy.
Our study supports the results from several other studies that have indicated that the vaginal microbiome is stable during pregnancy (8,31). Specifically, in a longitudinal study that includes 90 women, most retain a microbial community with the same dominant member (in L. crispatus communities, 75% remain stable; in communities with high L. iners abundance, 71% remain stable; and in more diverse communities like those sometimes associated with BV, 58% do not shift) (13). Using metagenomic sequencing to probe microbial community variation, our findings indicate that few bacteria, particularly Lactobacillus species, are highly abundant in the vaginal environment. Indeed, for individuals with vaginal microbiomes numerically dominated by L. crispatus or L. iners, the vast majority of reads mapped to contigs from one strain of L. crispatus or L. iners (Fig. S1B). Of note, the microbiome of some individuals did differ considerably with advancing pregnancy. For instance, the vaginal microbiota of subjects 1180 and 1222 had higher abundances of L. crispatus during the first trimester, but L. iners was more abundant in the remaining trimesters. Brooks et al. (32) demonstrated that shifts in vaginal microbiota structures can be described probabilistically, where shifts from L. crispatus to L. iners are the most likely to occur. This is consistent with the observations made in two individuals from our study; however, due to the small number of samples exhibiting this phenomenon, qualitative assessments were more appropriate than statistical analysis. Of note, vaginal microbiota instability throughout pregnancy was associated with the presence of Candida, a known opportunistic pathogen of the vaginal tract. Since inclusion criteria for this study stipulated no antibiotic treatment, it is unlikely that Candida detection was a result of antibiotic administration. The prevalence of Candida in our cohort is more likely to be reflective of the fact that pregnancy is a known risk factor for candidiasis (33) and to related differences in the vaginal environment, including microbiological colonization. Indeed, L. crispatus has been shown to have anti-Candida activity (34), and 90% (10/11) of samples that were Candida positive came from individuals whose vaginal microbiota were dominated by an organism other than L. crispatus.
A few metabolites were highly indicative of the bacterial community present in each subject and may be useful biomarkers for the type of vaginal microbiota present. The most indicative metabolite was indole-3-lactate (ILA), a tryptophan metabolite whose abundance was correlated with communities having abundant L. crispatus ( Fig. 4B and Fig. S4A). One potential explanation is that L. crispatus produces ILA to competitively exclude the growth of other species (Fig. 6). At physiologically relevant concentrations, ILA has been shown to have antimicrobial properties against both Gram-positive and Gram-negative organisms (35,36). Although the production of lactic acid is generally thought of as a strategy Lactobacillus spp. use to prevent other species from colonizing the vagina, perhaps these organisms also use ILA in a similar or supplementary capacity. Additionally, bacteria-derived ILA (also referred to as indole-lactic acid) has been recently shown to directly move from maternal to infant tissue (6). It has been suggested that indoles may play an important role as a ligand for the human aryl hydrocarbon receptor (AhR), which have diverse functions from immune regulation to metabolism (reviewed in reference 37). Zelante et al. further showed that some lactobacilli produce the related tryptophan catabolite, indole-3-aldehyde (IAld), which provides protection against candidiasis by increasing interleukin-22 (IL-22) production via AhR receptor binding (38). The study also demonstrated that vaginal specific bacteria, such as L. acidophilus, produce IAld in the vaginal environment, which protected against vaginal but not intestinal candidiasis. We measured indole-3-acetate (IAA), the direct precursor to IAld, in our study, but found no difference in its abundance between women dominated by different species of Lactobacillus (data not shown). Because indole-3-lactate can act as a ligand for AhR, we speculate that L. crispatus may regulate the IL-22-AhR response in the vagina, reducing the risk of vaginal candidiasis in the same way IAld does, and potentially activating the AhR response in newborns to prevent early-life candidiasis. Additionally, indole itself may play a role in structuring community composition by selecting for organisms that have adapted to a high abundance of this metabolite, repelling more transient microbes that have not been exposed to higher indole concentrations previously (39).
Whole-genome shotgun metagenomics allowed us to begin to address the functional capacities of these microbiomes. The largest differences between functional capacity appeared to be between communities where L. crispatus or G. vaginalis were the most abundant bacterial taxon. One pathway particularly enriched in Gardnerella communities was the ammonia assimilation cycle (Fig. S4B). Studies have pointed out Gardnerella's preference for ammonia as a nitrogen source (40); moreover, this ability to assimilate ammonia and produce amino acids has been implicated in mutualistic interactions between species of Prevotella, in the context of BV (41). Together, these BV-associated organisms contribute to genital inflammation, which may play a role in the susceptibility of certain diseases, such as HIV (42).
Increased abundance of mannitol when L. crispatus was present is an important and unexpected finding (Fig. 4). Most likely, mannitol contributes to optimizing the tonicity of the vaginal environment, and has recently been considered for this use in developing effective therapeutics for altering the vaginal microbiota (43,44). Even more, mannitol may assist L. crispatus in adhering to the epithelial layer, a strategy the organism may use to competitively inhibit other microbes from colonizing, potentially by drawing out excess water in the mucin layer and altering the mucin structure (45). Irrespective of the biochemistry, these genes, and mannitol in general, represent very specific markers of a community where L. crispatus was most abundant.
Although it is known that homofermentative lactic acid bacteria (LAB) such as L. crispatus (46) convert glucose primarily to lactic acid, it is unclear why mannitol accumulates in this niche. Interestingly, there was no difference in the glucose abundance between the four distinct vaginal communities. Further, metabolomic analysis of our CVF samples failed to capture significant levels of the mannitol precursor fructose; however, previous studies have indicated an appreciable amount of fructose within the cervical mucus of humans (47) and the capability of L. crispatus to utilize fructose as a carbon source (48,49). We speculate that this high extracellular mannitol abundance phenotype may underlie the cell's need to regenerate NADϩ for use in glycolysis. When faced with a limiting amount of pyruvate (or perhaps an upstream glycolytic metabolite) to convert to lactate, homofermentative LAB may be unable to produce sufficient NADϩ to allow glycolysis to continue. To this end, reducing fructose 6-phosphate to mannitol-1-phosphate may be an alternative and vital way L. crispatus regenerates NADϩ for glycolysis (50,51). We did find the gene mannitol-1-phosphate dehydrogenase, responsible for converting fructose 6-phosphate to mannitol-1-phosphate, was highly correlated (R 2 ϭ 0.9) with the relative abundance of L. crispatus. The genes for the conversion of mannitol-1-phosphate to mannitol (presumably an M1P phosphatase) and its subsequent export are currently unknown (51). This may imply that mannitol accumulation is a marker of a cellular switch to NADϩ regeneration by fructose reduction rather than converting pyruvate to lactic acid. Consequently, the decrease in lactic acid production may contribute to community dysbiosis due to a rise in pH (Fig. 6). Future experiments using culturing to elucidate whether these in vivo community data are recapitulated with axenic cultures in vitro are needed.
Furthermore, this study enabled comparison of two metabolomic methods (GC-TOF MS and LC-MS/MS lipidomics) for analysis of the pregnancy and early-life metabolomes. Our data showed that metabolite intensities obtained by GC-TOF MS were more tightly correlated with microbial community composition than those obtained by lipidomics, perhaps indicating that GC-TOF MS is more effective at detecting microbial metabolites than lipidomics, especially during pregnancy. We also show that both the saliva and urine metabolomes of children are more similar to their own mother than to unrelated individuals (Fig. 6). Strikingly, the ability of lipidomics to differentiate mothers from children via saliva was the strongest signal in our metabolomics data (Table S1). The oral microbiome may play a role in this, as there is a well-established community succession in children during early life (reviewed in reference 52), where children begin life with oral microbes that differ from those in adults. Lactulose, detected by GC-TOF MS, was a specific metabolite with increased abundance in infant saliva, which may reflect its use as a treatment for constipation (53) or perhaps even its presence in heated milk (54). Finally, urine metabolomes had a distinct age profile, especially with the lipidomics data. Our data suggest that, over the first year of life, the urine metabolome rapidly converges on the adult metabolome. This is likely a result of the development of the renal system in children (55), along with the development of the gut microbiome and the related metabolites which are processed through the liver and kidneys. Other reasons for this age-related shift include a change in diet and weaning off breastmilk or formula (55). Expectedly, we did not see a strong influence of the vaginal microbiome during pregnancy on the infant saliva and urine metabolome. We suspect that if differences in the vaginal microbiome were to affect the early-life saliva and urine metabolome, those effects would be subtle. The lack of stool samples from the mothers and infants is a limitation of this study, as they may contain a stronger signal of shared metabolomes across mother-infant dyads. Additionally, we explored a high-throughput approach to characterize the metabolomes. By using acoustic deposition in combination with MALDI-MSI, a throughput of ϳ1 s per sample was reached using only 2 l of sample. Of the metabolites identified by GC-MS and LC-MS, ϳ55% had corresponding ions in the MALDI-MSI analysis (Table S1). Future work will focus on confirming these metabolite identifications, but the initial results are promising and indicate that rapid analysis of microbial metabolites using MALDI, an analysis platform routinely used in clinical microbiology laboratories, is feasible (56).
Overall, we provide a broad look at the metabolome during pregnancy and early life, detailing the utility of GC-TOF MS, lipidomics, and MALDI-MSI for saliva, urine, and CVF.
In conclusion, here, we share a high-resolution characterization of the vaginal microbiome, longitudinally sampled throughout healthy pregnancies. We show that, despite the generally accepted view that lactobacilli are indicative of healthy vaginal communities, many women in our healthy pregnancy cohort had non-Lactobacillusdominated communities. The vaginal communities were characterized by a high abundance of one or a few acid-tolerant species, which dictated the physiologic potential and the metabolic profiles of the vaginal microbiome. Many of the metabolites that were specific to these different organisms warrant further investigation, especially considering the recent development of vaginal microbiota transplantation (VMT) as a treatment for BV (28). The metabolites we found to be associated with L. crispatus may be useful as microbiome cultivation approaches are developed to intentionally direct the composition of the vaginal microbiome. For example, indole-3-lactate may support L. crispatus colonization, while mannitol may indicate a shift in metabolism away from fermentation and the production of acid, relaxing the low-pH selection pressure, which normally gives L. crispatus an advantage.

MATERIALS AND METHODS
Subject information. Eighteen women were selected from a larger cohort recruited to address how maternal stress affects child development (57,58) (Table 1). Inclusion criteria for the larger cohort included Ͼ18 years of age, singleton, intrauterine pregnancy, and nondiabetic. Additional inclusion criteria for this study were normal prepregnancy body mass index, vaginal delivery, full-term pregnancy, breastfeeding, and no antibiotics for mother or baby. Generally, these 18 women and their children represented healthy subjects with the most complete sample sets.
Sample collection. Samples were collected at each trimester of pregnancy for women and throughout the first year of life for infants. At each time point, maternal saliva, urine, and cervical vaginal fluid were collected. For infants, urine was collected at birth, 6 months, and 12 months, whereas saliva was sampled at 6 months and 12 months of age. Maternal saliva was collected using a Salivette collection kit, including a small cotton roll contained in a plastic container (Salimetrics, Carlsbad, CA). Mothers were instructed to place the cotton rolls in their mouths until saturated with saliva (approximately 1 to 3 min) and then reseal the swabs in plastic Salivette tubes. Infant saliva was collected using Weck-Cel spears and a swab extraction tube system. Infants were allowed to suck on the spear for 2 min, ensuring saturation. Maternal urine was collected using a sterile collection cup. Infant urine was collected using an adhesive U-Bag attached to the genital region of the infant. A minimum of 2 ml of urine was collected. Cervical vaginal fluid (CVF) was collected by placing three Dacron swabs into the cervix for 10 s to achieve saturation. Each swab was then placed in a plastic vial with 500 l of sterile phosphate-buffered saline (PBS). All samples were initially stored at Ϫ20°C, and then saliva, CVF, and infant urine were subsequently moved to Ϫ80°C storage.
Metabolomics. Prior to processing, samples were thawed from Ϫ80°C storage. Fifty microliters of each sample were subjected to gas chromatography-time of flight mass spectrometry (GC-TOF MS) (59) and liquid chromatography-tandem mass spectrometry (LC-MS/MS, lipidomics). Urine, saliva, and CVF from each time point were sent to the West Coast Metabolomics Center (WCMC) for untargeted metabolomics. GC-TOF MS metabolites were extracted with a mixture of 3:3:2 acetonitrile-isopropyl alcohol-water according to standard operating procedures from the Fiehn Lab at the WCMC (60). LC-MS/MS samples were extracted using a variant of the Matyash method (61). Data were acquired for complex lipids in positive and negative electrospray mode on a Waters CSH column and an Agilent 6530 QTOF mass spectrometer (61). Metabolites were identified by retention time MS/MS matching using MassBank of North America (http://massbank.us) and NIST 17 libraries.
High-throughput metabolomics. All urine, saliva, and CVF were analyzed using matrix-assisted laser desorption ionization mass spectrometry imaging (MALDI-MSI) for high-throughput untargeted metabolomics. Extracted samples in 3:3:2 acetonitrile-isopropyl alcohol-water were diluted 1:2 in water in 384-well plates. Next, an equal volume of MALDI matrix (20 mg/ml of 1:1 2,5-dihydroxybenzoic acid and ␣-cyano-4-hydroxycinnamic acid in 1:3 [vol/vol] H 2 O-MeOH plus 0.2% formic acid) was added. Samples were printed onto a stainless steel blank MALDI plate using an ATS-100 acoustic transfer system (BioSera) with a sample deposition volume of 10 nl. Samples were printed in clusters of four replicates, with the microarray spot pitch (center-to-center distance) set at 900 m. MS-based imaging was performed using an ABI/Sciex 5800 MALDI TOF/TOF mass spectrometer with a laser intensity of 3,500 (arbitrary units) over a mass range of 50 to 3,000 Da. Each position accumulated 20 laser shots. The instrument was controlled using the MALDI-MSI 4800 imaging tool. Surface rasterization was oversampled using a 75-m step size. Average ion intensity for all reported ions was determined using the OpenMSI Arrayed Analysis Toolkit (OMAAT) software package (62).
DNA extraction. Cervical brushes were resuspended in PBS. Two negative extraction controls using sterile PBS were prepared alongside the samples. Aliquots of 100 to 200 l were added to lysing matrix E tubes prealiquoted with 500 of hexadecyltrimethylammonium bromide (CTAB) DNA extraction buffer and incubated at 65°C for 15 min. An equal volume of phenol-chloroform-isoamyl alcohol (25:24:1) was added to each tube, and samples were homogenized in a Fast Prep-24 homogenizer at 5.5 m/s for 30 s. Tubes were centrifuged for 5 min at 16,000 ϫ g, and the aqueous phase was transferred to individual wells of a 2-ml 96-well plate. An additional 500 l of CTAB buffer was added to the lysing matrix E tubes, the previous steps were repeated, and the aqueous phases from paired extractions were combined. An equal volume of chloroform was mixed with each sample, followed by centrifugation at 3,000 ϫ g for 10 min. The aqueous phase (600 l) was transferred to a clean 2-ml 96-well plate, combined with 2 volume equivalents of polyethylene glycol (PEG) and stored overnight at 4°C to precipitate DNA. Plates were centrifuged for 60 min at 3,000 ϫ g. DNA pellets were washed twice with 300 l of 70% ethanol, air-dried for 10 min, and resuspended in 100 l of sterile water. DNA was quantified using the Qubit dsDNA high-sensitivity (HS) assay kit and diluted to 10 ng/l when possible. Although DNA was extracted from CVF, attempts to extract DNA from saliva were unsuccessful, potentially due to the storage swabs trapping the biomaterial.
Amplicon gene sequencing.  [63], and 0.025 U/l Taq DNA polymerase) and 0.4-M barcode-tagged reverse primers (806R) and then amplified in triplicate 25-l reactions, along with a no-template control, for 30 cycles (98°C for 2 min; 98°C for 20 sec, 50°C for 30 sec, and 72°C for 45 sec; repeat steps 2 to 4 29 times; 72°C for 10 min). PCR conditions were identical for ITS2 amplification (primer pair fITS7 [5=-GTGARTCATCGAATCTTTG-3=] and ITS4 [5=-TCCTCCGCTTATTGATATGC-3=]) except for the annealing temperature, which was 52°C. Triplicate reactions were combined and purified using the SequalPrep normalization plate kit (Invitrogen) according to the manufacturer's specifications. Purified amplicons were quantified using the Qubit double-stranded DNA (dsDNA) HS assay kit and pooled at equimolar concentrations. The amplicon library was concentrated using the Agencourt AMPure XP system (Beckman-Coulter), quantified using the Kapa library quantification kit (Kapa Biosystems), and diluted to 2 nM. Equimolar PhiX was added at 40% final volume to the amplicon library; the 16S rRNA amplicon pool was sequenced on the Illumina NextSeq 500 platform on a 153-bp by153-bp sequencing run, and the ITS2 amplicon pool was sequenced on the Illumina MiSeq platform on a 290-bp by 290-bp run.
Shotgun metagenomics sequencing. Sequencing libraries were prepared using the Illumina Nextera kit and methods described in Baym et al. (64). Briefly, DNA from each sample was diluted to 0.5 ng/l and tagmented with the Nextera enzyme (Illumina) for 10 min at 55°C. Following tagmentation, each sample received 1-l forward and 1-l reverse barcodes, which were added via PCR using Phusion DNA polymerase (New England BioLabs). After PCR, the libraries were cleaned of smaller DNA fragments, using AMPure XP magnetic beads (Beckman-Coulter), and pooled by concentration. Libraries were quantified using the Quanti-iT PicoGreen dsDNA kit (Thermo Fisher Scientific), and DNA was run on a gel to check fragment size. These libraries were loaded onto the Illumina Next-Seq 500 at 1.8-pM concentrations and Illumina's midoutput kit for 75-bp paired-end sequencing.

OTU table generation.
Raw sequence data were converted from bcl to fastq format using bcl2fastq v2. 16.0.10. Paired sequencing reads with a minimum overlap of 25 bp were merged using FLASH v1.2.11. Index sequences were extracted from successfully merged reads and demultiplexed in the absence of quality filtering in QIIME (Quantitative Insights into Microbial Ecology, v1.9.1), and reads with more than two expected errors were removed using USEARCH's fastq filter (v7.0.1001). Remaining reads were dereplicated, clustered into operational taxonomic units (OTUs) at 97% sequence identity, filtered to remove chimeric sequences, and mapped back to OTUs using USEARCH v8.0.1623. Taxonomy was assigned with the most current Greengenes database for bacteria (63) (May 2013) and UNITE version 6 for fungi (65). OTUs detected in negative extraction controls (NECs) were considered potential contaminants and filtered by subtracting the maximum NEC read count from all samples; any remaining OTU with a total read count less than 0.001% of the total read count across all samples was removed. Sequencing reads were rarefied to an even depth (28,972 reads for 16S; 91,232 reads for ITS2). To maximize similarity between the raw and rarefied OTU tables, random subsampling was performed at predefined depths 100 times, and the most representative subsampled community was selected based on the minimum Euclidean distance to the other OTU vectors generated for each sample.
16S rRNA gene analysis. Alpha diversity indices and Bray-Curtis dissimilarity matrices were generated in QIIME (66). Linear outcomes were assessed by linear mixed-effects (LME) modeling to adjust for repeated measures using the nlme package (67) in the R environment (68). Variables of P Ͻ 0.05 were considered statistically significant. Data were visualized using Tableau and Adobe Illustrator unless otherwise noted.
Metagenomic analysis. Raw sequences (mean 2,977,881 paired-end reads per sample from 35/38 successfully sequenced samples) were filtered using PRINSEQ v0. 20.4 (69) to filter out sequences that had a mean quality score of 30 or less. Human DNA was next filtered out by aligning the filtered reads to the human genome (hg38) using Bowtie2 v2.2.7 (70) and keeping the reads that failed to align (mean 220,355 paired-end plus 33,744 singleton reads per sample or 10.9% of quality-filtered reads per sample). To analyze functional potential, the reads were run through HUMAnN2 v0.1.9 (71) using default parameters, and differences in pathway abundances were analyzed using LEfSe (72). These reads were also cross-assembled using SPAdes v3.8.2 (73). Each sample was then mapped to this cross-assembly using Bowtie2, samples from the same subject were merged together using Samtools v1.9, and the resulting bam files and the cross-assembly were imported into Anvio4 (74). Taxonomy was assigned to each gene call using Kaiju (75), which subsequently informed a more accurate metagenomics binning of the most abundant microbes present.
Statistical analysis. Unless otherwise noted, statistics were done using the ecological statistics program PRIMER-e v7 (76). Metabolic data were normalized in PRIMER-e by dividing by sum total for each sample. The specific programs used in PRIMER-e were permutational multivariate analysis of variance (PERMANOVA) and distance-based linear models (DISTLM), the former of which calculates the significance and variance explained by a given factor and the latter determines which environmental variables correlate with the biological (microbiome) data. rfPermute (77), an R package for permutated random forests, was also performed to determine which annotated GC metabolites were indicative of microbial composition. PERMANOVA also partitions variance based on each factor, which is done in PRIMER-e by dividing the factor estimate by the sum total estimates of components of variation (ECoV). Traditional R 2 values were also calculated by dividing the sum of squares by total sum of squares. LMEs were carried out as described above; R 2 values for linear mixed models were calculated using the MuMIn package in R (78). Relate tests (analogous to Mantel tests) were used to compare GC-MS and LC-MS data. Bray-Curtis distances were used for all distance-based analyses. To consider repeated measures, linear mixed-effects modeling (nlme package in R) was used to analyze the stability of the microbiome and metabolome through time.
Ethics approval and consent to participate. This study utilized a subset of samples from a larger, longitudinal prospective cohort study designed to analyze the relationship between maternal stress and infant development conducted at the University of California, Irvine (UCI) (57). The University of California's Institutional Review Board approved the protocol, and written, informed consent was obtained from each participant. Research on human subjects was performed in accordance with the Declaration of Helsinki.
Data availability. Sequence data for 16S, ITS2, and shotgun metagenomes were deposited on the National Center for Biotechnology Information (NCBI) sequence read archive (SRA) under the BioProject accession number PRJNA612083. Metabolomics data for all samples can be found in Table S1. R scripts for statistical analysis are published on GitHub at https://github.com/aoliver44/Cervicovaginal-Paper.

SUPPLEMENTAL MATERIAL
Supplemental material is available online only.

ACKNOWLEDGMENTS
The pilot grant awarded to start a UC Center for Pediatric Microbiome Research awarded by the UCI Institute for Clinical and Translational Science (UL1 TR001414) gave rise to this project. A.O. was supported by an NIH T32 training grant (1T32AI14134601A1) from UC Irvine's training program in microbiology and infectious diseases.
We thank Dan Cooper for enthusiasm and strategic advice. We thank the members of the Whiteson Lab, especially Whitney England for insight into metagenomics analysis, and Joann Phan and Tara Gallagher for thoughtful comments during the preparation of the manuscript. Heather Maughan made insightful and clarifying edits that we are grateful for.
We declare that we have no competing financial interests.