This study was conducted according to the institutional guidelines for animal experiments of Rakuno Gakuen University (approval no. VH16C7). Seventeen primiparous cows (Holstein Friesian, range of age at first calving 22–27 months), received artificial insemination (AI) or embryo transfer (ET), were enrolled from December 2017 to August 2018 at the Rakuno Gakuen University for this study. All cows were diagnosed as either healthy or suffering from endometritis at 35 DPP, according to the methods previously reported (4). Briefly, clinical endometritis was diagnosed with vaginal discharge scoring as follows: 0 = clear or translucent mucus; 1 = mucus containing flecks of white or off-white pus; 2 = exudate containing <50% white or off-white mucopurulent material; and 3 = exudate containing ≥50% purulent material, usually white or yellow, occasionally bloody. The degree of recovery of the uterus was also assessed by rectal examination as follows: 0 = the diameter of largest uterine horn ≤3.5 cm and cervical diameter ≤4.5 cm; 1 = the diameter of largest uterine horn more than 3.5 to <5.5 cm and cervical diameter more than 4.5 to <7.0 cm; 2 = the diameter of largest uterine horn ≥5.5 cm and cervical diameter ≥7.0 cm. The numbers of epithelial endometrial cells and PMNs were counted to assess subclinical endometritis (3). Endometritis was defined either as total gynecological examination scoring >1, or a PMN ratio >10%. Three cows were excluded due to the following reasons: treatment with antibiotics, abnormal deliveries, and urovagina. The criteria for endometritis, scores and sampling dates of the 14 cows (seven healthy cows and seven cows diagnosed with endometritis at 35 DPP) are summarized in Supplementary Figure 1 and Supplementary Table 1. In addition, parturitions were scored for difficulty and obstetric assistance required as follows: 1 = no assistance (n = 6); 2 = slight assistance with only one person (n = 5); 3 = moderate assistance with 2–3 persons (n = 1); 4 = severe dystocia with veterinary treatment (n = 2); 5 = cesarean section (n = 0).

Both uterine and vaginal samples were collected at the same time as follows; Pre (pre-calving, only for vaginal samples), 0 DPP (within 12 h after calving and after expulsion of the placenta), 7, 21, 28, and 35 DPP (Supplementary Table 1).

The perineum and vulva were cleaned with a 70% ethanol swab and wiped with a paper towel. Uterine samples were collected using the cytobrush technique (3). The cytobrush device (Metribrush; Fujihira Kogyo Inc., Tokyo, Japan), double-guarded with sterilized plastic sleeve and tube, was inserted through the vagina and reached the cervix guided by palpation per rectum. Then the plastic sleeve was pulled back and the brush was moved forward and rolled over the endometrium of the uterus body. The cytobrush was retracted into the tube in uterus, then retrieved through the vagina. The tip of the cytobrush was then cut using sterile scissors, placed in 5 ml sterile saline solution and vortexed vigorously. Vaginal samples were obtained by washing methods as previously described (20), as sampling by brushing of the vaginal fornix was not suitable to extract sufficient bacterial DNA. In brief, 50 ml of sterile saline was injected into the vagina, then flushed back 2 or 3 times. Within 1 h, vaginal washings were centrifuged at 8,000 × g for 15 min and pellets were re-suspended in 5 ml sterile saline. After vortexing, 500 μl of uterine and vaginal samples were used for the culture-based method and the remaining samples were mixed in an equal volume of sterile saline containing 40% glycerol. Mixtures were flash frozen in liquid nitrogen, and stored at −80°C until use.

Bacterial DNA was extracted as described previously (21). After thawing, 2 ml of uterine and vaginal samples were mixed with 20 ml of sterile saline solution and centrifuged at 8,000 g for 15 min. Pellets were then suspended in 800 μl 10 mM Tris-HCl and 10 mM EDTA buffer containing lysozyme (Sigma-Aldrich Co., LCC, Missouri, USA) (Final concentration: 15 mg/ml). Mixtures were incubated at 37°C for 1 h. Purified achromopeptidase (Wako Pure Chemical Inc., Osaka, Japan) (Final concentration: 2,000 U/ml) was added and incubated 37°C for 30 min. Finally, proteinase K (Takara-Bio Inc., Shiga, Japan) (Final concentration: 1 mg/ml) with 20% sodium dodecyl sulfate (Sigma-Aldrich Co., LCC) was added and incubated at 55°C for 1 h. The DNA was isolated with phenol:chloroform:isoamyl alcohol (25:24:1, v/v), washed twice with 75% ethanol, and dissolved in 100 μl TE buffer. RNase A (Nippon Gene Co., Ltd., Tokyo, Japan) (Final concentration: 0.1 mg/ml) was added to the mixture, and incubated at 37°C for 1 h. Subsequently, the DNA was purified using a High Pure PCR Template Kit (Roche Inc., Basel, Switzerland) according to the manufacturer's instructions. Elution was performed in 50 μl of TE buffer and then samples were stored at −20°C until further analysis.

Extracted DNA was amplified, targeting the V3–V4 region of the bacterial 16S rRNA gene. The specific universal primer pair 341F (5′-CCTACGGGNGGCWGCAG) and 805R (5′-GACTACHVGGGTATCTAATCC) was used in this study (22). This primer set included the Illumina MiSeq sequencing adapter (forward primer: AATGATACGGCGACCACCGAGATCTACAC; reverse primer: CAAGCAGAAGACGGCATACGAGAT) and a unique barcode sequence that allowed the samples to be pooled for Illumina MiSeq sequencing. Polymerase chain reaction (PCR) was performed using MightyAmp DNA Polymerase (Takara-Bio Inc.) for 35 cycles. The PCR products were quality-checked on a 2% agarose gel electrophoresis and subsequently purified using SPRIselect beads (Beckman-Coulter Inc., California, USA). The amplified DNA was quantified using a ONEdsDNA System (Promega, Madison, WI, USA) and Quantus fluorometer (Promega). The PCR amplicon libraries were prepared by pooling approximately equal amounts of amplified DNA and sequenced on an Illumina MiSeq platform (Illumina, San Diego, USA), using the 2 × 300, v3 600-cycle kit (Illumina).

Illumina Miseq fastq raw reads were analyzed with QIIME 2 platform version: 2020.2 with default scripts (23). Sequences were demultiplexed and processed by the DADA2 program (24). Before merging, reads were trimmed according to quality threshold and adapter length (Forward: 280 bp, Reverse: 220 bp). After quality filtering steps (DADA2), total, average, minimum and maximum number of non-chimeric reads were 1,901,611, 14,086, 7,565, and 17,893 reads, respectively. Amplicon sequence variants (ASVs) obtained after DADA2 [more accurate than traditional operational taxonomic units (OTUs)] were assigned using the qiime feature-classifier classify-sklearn method (25) against the Greengenes database (version 13.8) with 99% similarity. For diversity analysis, the number of sequence reads was rarefied to the minimum sample reads (7,565 reads) using the qiime diversity core-metrics-phylogenetic method (23). The rarefaction curve confirmed that the sub-sampling was sufficient to detect the bacteria species in all samples. The alpha-diversities were calculated using the Simpson and Shannon index. Beta-diversity analysis, represented by principal coordinate analysis (PCoA), was applied to the resulting weighted distance matrices to generate two-dimensional plots. For differential abundance analysis between healthy and endometritis groups, random subsampling of sequence reads was not performed.

Statistical differences in alpha-diversities between healthy and endometritis groups at sampling time points were tested using the Mann-Whitney U-test. The significance of the groups in the community structure was tested using permutational multivariate analysis of variance (PERMANOVA). Relative abundances of genera between healthy and endometritis groups was tested using Mann-Whitney U-test. To identify unique microbial genera for both healthy and endometritis groups, ANOVA-Like Differential Expression 2 (ALDEx2), estimates of the composition of biological features from the number of reads transformed to central log ratio (clr), was performed (26, 27). The correlations between the abundances of genera in uterine and vaginal microbiota were determined using Spearman's rank correlation coefficient. A linear discriminant analysis (LDA) effect size (LEfSe) (28) approach was used to identify bacterial taxonomy that was significantly different between the “before birth” group (Pre) and “after birth” group (0 DPP) for both healthy and endometritic cows.

All samples used for DNA extraction were also cultured under aerobic and anaerobic conditions within 1 h from sampling. Tryptic soy (TS) agar (BD Bacto, New Jersey, USA) and blood liver (BL) agar (Nissui, Tokyo, Japan), which are widely used culture media in the analysis of human fecal microbiota (29) were adopted in this study. Both TS agar and BL agar were supplemented with 5% defibrinated horse blood. Briefly, for aerobic culture, 50 μl of 10-fold serial dilutions of sample suspensions were spread on TS agar and incubated at 37°C for 48 h. The same suspensions were spread on BL agar and incubated at 37°C for 48 h under anaerobic conditions (10% H₂, 10% CO₂, and 80% N₂) using an anaerobic chamber. After incubation, the numbers of colonies with the same morphology on agar were counted. For several different types of colonies, a maximum of three colonies with the same morphology on one plate were picked from each agar and identified using MALDI-TOF MS (Bruker, Bremen, Germany). To arbitrate the results of MALDI-TOF MS and identify bacteria to species level, 16S rRNA gene sequencing was performed by using the universal primer pair 27F (5′-AGAGTTTGATCCTGGCTCAG) and 1492R (5′-GGTTACCTTGTTACGACTT) if bacterial strains obtained an identification scoring of <2.0 (30). The obtained sequences were compared using the basic local alignment search tool (BLAST) with the National Center for Biotechnology Information (NCBI) databases. A species result was adopted only when all three colonies had the same identification result. Mann-Whitney U-test was used to calculate significant differences in the number of bacteria, converted to log CFU/ml of samples, between healthy and endometritis groups.

The alpha-diversities of both uterine and vaginal microbiota of the healthy and endometritis groups were measured using Simpson and Shannon indices (Figure 1). As shown in Figure 1, Simpson and Shannon indices of vaginal microbiota at 35 DPP were found to be significantly higher in the endometritis group than in the healthy group (Mann-Whitney U-test, p < 0.05). As a common trend within the healthy and endometritis groups, vaginal diversity was the highest before calving and tended to decrease as time passed (Mann-Whitney U-test, p < 0.01; Pre vs. 35 DPP), while uterine diversity was lowest at 7 DPP and recovered as time passed (Mann-Whitney U-test, p < 0.05; 7 vs. 0 DPP and 7 vs. 35 DPP).

Principal coordinate analysis based on weighted UniFrac distance showed that the uterine microbiota of the healthy group was clustered and separated from those of the endometritis group at 28 and 35 DPP (PERMANOVA, p < 0.05) (Supplementary Figure 2).

A total of 3,918 ASVs were observed in samples. The ASVs shared between the uterine and vaginal microbiota of all animals at each time point are presented in Figure 2. Only 6.6–17.9 and 9.0–20.1% of all ASVs were shared between uterine and vaginal samples for the healthy group and endometritis group, respectively. In both the healthy and endometritis groups, the number of ASVs was higher in the vagina than in the uterus for 7 DPP. This trend was reversed after 21 DPP, namely the number of ASVs in the uterine samples was higher than that in the vaginal samples. In the healthy group, the number of ASVs at 0 DPP was almost the same in the uterine samples (207 ASVs, 34.4%) and vaginal samples (287 ASVs, 47.7%), whereas the number of ASVs in the uterine samples (414 ASVs, 60.8%) was higher than that in the vaginal samples (173 ASVs, 25.4%) in the endometritis group.

The mean relative abundances of the top 10 bacterial genera in the uterus is displayed by bar chart (Figure 3). The results showed that uterine bacterial communities in the healthy group during the study period were dominated by genus Bifidobacterium (The relative abundance of genus Bifidobacterium at 0, 7, 21, and 35 DPP were 19.5, 11.5, 20.0, 31.4, and 36.4%, respectively) followed by Streptococcus (The relative abundance of genus Streptococcus at 0, 7, 21, and 35 DPP were 17.0, 17.8, 10.3, 13.1, and 19.5%, respectively). In contrast, genus Clostridium was frequently detected and dominated at 0 DPP in the endometritis group (25.9%). The genus Trueperella was detected in both healthy and endometritis cows. With the exception that the genus Ureaplasma was most frequently predominant in the healthy group, vaginal bacterial communities exhibited the same dynamics as uterine in terms of genus Bifidobacterium, Trueperella, and Clostridium (Supplementary Figure 3).

These bacterial dynamics of the uterine and vaginal microbiota were also confirmed by traditional culture-dependent methods. A total of 1,213 colonies were isolated from 14 cows (healthy cows and cows with endometritis), representing 37 bacterial genera. For the three bacterial genera with the highest detection rates, bacteria count data were calculated from the number of colonies diluted 10-fold and converted to log CFU/ml of samples (Supplementary Table 2). As shown in Supplementary Table 2, the number of T. pyogenes was significantly higher in the endometritis group compared to the healthy group at 21 DPP in uterus (N.D. vs. 4.64 ± 0.23 log cfu/ml, p < 0.05) and vagina (2.62 ± 0.76 vs. 5.03 ± 0.27 log cfu/ml, p < 0.05), and at 35 DPP in the vagina (2.62 vs. 4.37 ± 0.73 log cfu/ml p < 0.05).

To verify the difference in genera between the healthy and endometritis groups, differential abundance analysis was performed using ALDEx2, which estimates the composition of biological features in the sample based on the number of reads.

Except for duplications, in total 10 of the bacterial genera were found to be different between the healthy and endometritis groups in uterine samples (Figure 4). The healthy group had higher abundances of genus Enterococcus at 0 DPP, and genus Clostridium at 35 DPP (p < 0.024 and p < 0.045, respectively). Genera differences were most detected at 28 DPP: Prevotella, Bifidobacterium, Vibrio, Streptococcus, Pseudoalteromonas, Peptoniphilus, and Enterococcus were significantly increased in the healthy group with p < 0.0005 to p < 0.049, whereas Phascolarctobacterium and Eubacterium showed a significant increase in the endometritis group (p < 0.044 and p < 0.045, respectively).

For vaginal samples, genus Histophilus and an unclassified genus of the family Mogibacteriaceae were found to be significantly increased in the endometritis group at pre-calving (p < 0.048 and p < 0.047, respectively), whereas unclassified genera of the families of Lachnospiraceae and Ruminococcaseae were significantly increased in the healthy group at 35 DPP (p < 0.029 and p < 0.043, respectively). There were no significant differences in the genera between healthy and endometritis groups at 0–28 DPP for vaginal samples, and 7–21 DPP for uterine samples.

The correlations between the top 20 abundant bacterial genera of the uterine and vaginal microbiota were analyzed mixing both healthy and endometritis groups (Figure 5A). There was a significantly positive correlation between Trueperella and Helcococcus (r = 0.875; p < 0.001) but significantly negative correlations with Corynebacterium, Porphyromonas, Lactobacillus, Streptococcus, unclassified genus of the family Lachnospiraceae, and unclassified genus of the order Clostridiales. The genus Bifidobacterium was shown to have strong positive correlations with the genera Prevotella and Enterococcus (r = 0.964 and 0.954, respectively; p < 0.001), and moderately positive correlations with Lactobacillus and Streptococcus (r = 0.519 and 0.513, respectively; p < 0.01). In the vaginal microbiota, there were statistically significant positive correlations between the genera Trueperella and Helcococcus (r = 0.816; p < 0.001), as well as between Trueperella and Peptoniphilus (r = 0.679; p < 0.001), whereas Trueperella showed moderately negative correlations with unclassified genus of the family Lachnospiraceae and Ruminococcaseae (r = −0.202 and −0.275, respectively; p < 0.01) (Figure 5B). The genus Ureaplasma was the most predominant and negatively correlated with many other bacteria.

For healthy cows, LEfSe analysis resulted in five genera which were significantly discriminative between “before birth” group (Pre) and “after birth” group (0 DPP) (>3 log 10 LDA score, p < 0.05) (Supplementary Figure 4). Genera Helcococcus, Blautia, and Klebsiella were represented in “before birth” group, while genus Sutterella and an unclassified genus of family Alcaligenaceae, were represented in “after birth” groups. By contrast, nine identified bacterial taxa including genera Corynebacterium, Helcococcus, an unclassified genus of the family Veillonellaceae were unique to “before birth” group of endometritic cows (>4 log 10 LDA score, p < 0.05).