Lactobacillus iners dominates the vaginal microbiota of healthy Italian women of reproductive age

41 42 Large sex hormonal fluctua-ons are thought to influence vaginal microbiota, but liNle is known 43 about the impact of small, physiological varia-ons. Here we tracked changes in vaginal microbiota 44 during four key menstrual cycle phases in 61 healthy Italian women from the Women4Health cohort. 45 The microbiota was primarily composed of Lactobacillus, with L. iners being the most abundant. 46 Noteworthy, the high abundance of L. iners contrasts with previous studies in European popula-on s, 47 challenging its proposed pathogenic role, and sugges-ng dis-nct microbio ta profiles within Europe. 48 Individual microbiotas were generally stable, but beta diversity was higher during the follicular 49 phase. Only 11 women exhibited composi-onal shi_s and those mostly occurred between follicular 50 and ovulatory phases. Finally, among the hormones evaluated, 17-beta estradiol had the largest 51 impact on taxa abundance varia-ons. Our study highlights specific features of the Italian popula-on 52 and points to the resilience of the vaginal microbiota to physiological hormonal changes. 53 54 .CC-BY-NC-ND 4.0 International licensemade available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprintthis version posted March 12, 2025. ; https://doi.org/10.1101/2025.03.12.642767doi: bioRxiv preprint 3

55 56 The vaginal microbiota refers to the community of microorganisms that colonize the vagina. Unlike 57 microbiota from other mucosal sites (i.e., the gastrointes-nal tract), it has lower diversity, and it is 58 dominated by a few species, primarily Lactobacillus spp, which play a key role in regula-ng genital 59 health. Producing lac-c acid and an-microbial molecules, these bacteria reduce the risk of 60 coloniza-on by other microorganisms, including poten-al pathogens involved in infec-ons and 61 diseases of the vaginal tract [1]; [2]; [3]; [4]. In general, the prevalence of Lactobacillus species in 62 healthy women increases with puberty and stays stable un-l menopause[5]. 63 The most common species of Lactobacillus are L. iners, L. crispatus, L. gasseri, and L. jensenii. Earliest 64 studies showed that a vaginal microbiota dominated by L. crispatus was associated with a healthy 65 status, whereas L. iners was o_en isolated during infec-ons and therefore a microbiota dominated 66 by it is o_en considered “dysbio-c ”. However, this hypothesis originated by studies with small 67 sample size and a direct pathogenic role for L. iners has not been proved [6]; [7]; [8]. Likewise, a 68 vaginal microbiota dominated by other anaerobic bacteria, such as Gardnerella vaginalis (recently 69 reclassified as Bifidobacterium vaginalis), Prevotella, and Fannyeshea, have been associated with 70 pathological condi-ons of the female genital tract due to the virulence poten-al of some strains [9]; 71 [3]; [10]. 72 The pathogene-c role of L. iners remains under debate since its prevalence can vary substan-ally 73 among different popula-ons [2]; [11]; [12]; [13]; [14]. The vaginal microbiota of African and Hispanic 74 women during the reproduc-ve age is dominated by L. iners , while is much rarer in Northern 75 European women with L. crispatus being the most abundant and prevalent Lactobacillus species. 76 The reasons for this north-to-south gradient in the prevalence of L. crispatus and L. iners are s-ll 77 unknown. It is worth no-ng that these two species can co -exist in the same microbiota, therefore 78 it’s unlikely they are antagonists. 79 The v aginal microbiota has been considered an overall stable ecosystem, although two recent 80 studies conducted on short-term (42 days) and medium -term (18 months) longitudinal follow-up 81 have shown that varia-ons may occur in a small frac-on of women [6]; [15]. These studies evaluated 82 changes in the so-called Community State Types (CST) , a method that classifies the vaginal 83 microbiota into 5 groups depending on the composi-on and the most abundance species [2]; [16]. 84 In the shorter 42-days study, only microbiotas classified as CST-I (dominated by L. crispatus) and CST-85 IV (dominated by anaerobic bacteria) remained stable, while varia-ons were observed in other CST 86 .CC-BY-NC-ND 4.0 International licensemade available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprintthis version posted March 12, 2025. ; https://doi.org/10.1101/2025.03.12.642767doi: bioRxiv preprint 4 types. In the longer 18 months longitudinal study, vaginal microbiota of CST-I and CST-III (dominated 87 by L. iners) remained stable , while CST-IV microbio ta displayed lower temporal stability 88 (approximately 6 months). The underlying causes of these varia-ons, how they occur in popula-on s 89 with different Lactobacillus spp. prevalences, and which species may be responsible for these 90 changes in CSTs are not yet known. Worth men-oning, none of these studies directly characterized 91 sex-hormones, being unable to link composi-onal changes to specific sex hormones fluctua-ons [3]. 92 Here we present vaginal microbiome composi-on and its dynamics in of 61 healthy women from 93 the Italian popula-on cohort Women4Health (W4H), an observa-onal longitudinal study centered 94 on one menstrual cycle set up to explore the interac-ons between sex hormones, vaginal and 95 intes-nal microbiome, and lipid and glucose metabolism in women [17]. We characterized vaginal 96 microbiome composi-on and evaluated dynamics at four -me points pivotal for the menstrual cycle. 97 We also evaluated host factors associated with diversity and iden-fied sex-hormones that are mostly 98 associated with the observed dynamic changes. 99 100

101 102 The Women4Health cohort 103 The Women4Health (W4H) cohort is an ongoing short-term longitudinal study aiming to recruit up 104 to 300 healthy women and to follow them up through a natural menstrual cycle, collec-ng various 105 biological samples and metadata via ques-onnaires [17]. Briefly, women of European origin aged 18 106 to 45 years old, who were not using hormonal contracep-ves and who have not a diagnosis of 107 gynecological, metabolic, or gastrointes-nal diseases, were enrolled in the study at the Ins-tute for 108 Maternal and Child Health – IRCCS Burlo Garofolo in Trieste, Italy. Upon inclusion and signature of 109 an informed consent form, electronic ques-onnaires were used to gather data on family history, 110 clinical and anthropometric measurements, lifestyle and diet . They were then scheduled for 4 111 weekly follow-up appointments to hand in self-collected vaginal swabs, stool samples, saliva (only 112 at first appointment), to undergo a blood withdrawal and to collect an addi-onal ques-onnaire 113 related to their health and lifestyle in days preceding the appointment. Weekly follow-up visits were 114 scheduled according to pivotal phases of the menstrual cycle: the 1st at the follicular phase 115 (between days 6 and 9), the 2nd at the ovulatory phase (between days 12 and 16), the 3rd at the 116 early luteal phase (between days 18 and 22), and the 4th at late luteal phase (between days 24 and 117 .CC-BY-NC-ND 4.0 International licensemade available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprintthis version posted March 12, 2025. ; https://doi.org/10.1101/2025.03.12.642767doi: bioRxiv preprint 5 27). Blood samples were processed within 1 hour a_er collec-on . On serum, we measured 5 sex 118 hormones (beta17-estradiol, luteinizing hormone, follicle s-mula-ng hormone , progesterone and 119 prolac-n). During the first visit, we also measured testosterone and thyroid parameters (thyroid 120 s-mula-ng hormone and thyroxine (T4)) to exclude endocrine disorders. Hormones were measured 121 using 1. 5 ml serum on the ADVIA Centaur CP Immunoassay System instrument from Siemens 122 Healthineers. 123 The study protocol - including informed consent forms - was approved by the Friuli-Venezia Giulia 124 Ethical CommiNee on 09/2021 with prot. N. 0034184/P/GEN/ARCS and modifica-ons amended on 125 09/2023 with prot. N. 0033477/P/GEN/ARCS. 126 Here, we analyzed the vaginal swabs, metadata, and sex -hormone levels in the first 61 enrolled 127 women. Descrip-ve sta-s-cs of the sample are provided in Table S1. 128 129 Vaginal swabs collec=on, processing, and microbial DNA extrac=on 130 Vaginal swabs were self-collected by volunteers before each follow-up visit (the same day or the day 131 before) using the OMNIgene-VAGINAL kits (OMR-130, DNA Genotek), stored at room temperature 132 and handed at IRCCS Burlo Garofolo within 24 hours, where kits were promptly stored, as 133 recommended by the manufacturer, and shipped to the Ins-tute for Gene-c and Biomedical 134 Research - Na-onal Research Council (IRGB -CNR) within 21 days. DNA was extracted from swabs 135 with QIAamp PowerFecal Pro DNA kits (QIAgen) using a custom semiautomated protocol on a 136 QIAcube plarorm. Prior to DNA extrac-on, a bead-bea-ng step on a TissueLyser II (QIAgen) at 7.5 137 Hz for 5 min was performed to improve lysis. The composi-on of vaginal microbiota was established 138 by sequencing amplicons of the V3-V4 and V7-V9 regions of the 16S rRNA gene, along with ITS1 gene 139 of fungi. Library construc-on was performed with a 2-step PCR protocol and dual-barcoding strategy 140 (QIAseq 16S/ITS Region Panel kits, QIAgen). Sequencing was carried out with the MiSeq sequencing 141 system (Illumina) targe-ng an average of 100K 250-bp paired-end reads per sample. Samples were 142 processed in three different batches, each including one QIAseq 16S/ITS Smart Control (QIAgen), a 143 synthe-c DNA construct used both as posi-ve control for library prepara-on and as control for the 144 iden-fica-on of contaminants , and a PCR nega-ve control (dis-lled sterile water). Data obtained 145 were demul-plexed using Qiagen's GeneGlobe so_ware. 146 147 Quality Control (QC) process 148 .CC-BY-NC-ND 4.0 International licensemade available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprintthis version posted March 12, 2025. ; https://doi.org/10.1101/2025.03.12.642767doi: bioRxiv preprint 6 All demul-plexed FastQ underwent quality control screening before bioinforma-c analysis. FastQC 149 v.0.12.1 (see URLs) and Mul&QC v.1.22.2 [18] (see URLs) tools were used to detect outliers based on 150 read quality. The trimmoma&c tool v.0.39 [19] (see URLs) was used for trimming low quality paired-151 end reads. We used the op-ons SLIDINGWINDOW:4:25 to evaluate whether the average quality 152 score was at least 25 per four bases; if not, this part of the read was trimmed. The op-on MILEN:50 153 was used to remove reads shorter than 50 bp. A_er the trimming process, samples were re-analyzed 154 with FastQC and Mul&QC to evaluate their post-process quality. 155 We first processed 185 samples, 4 smart controls and 3 nega-ve controls in two plates and decided 156 to remove i) all samples that had less than 80K reads as input before the trimming process and/or 157 less than 60K reads survived post-trimming (N=10samples), ii) one sample with very high number of 158 ITS reads mapping on Sporobolomyces spp, a known environmental contaminant, and iii) samples 159 whose libraries were prepared together with a nega-ve control that showed contamina-on (>500 160 reads)(N=10). These 21 samples were repeated in a third plate, containing addi-onal 28 samples, 1 161 smart control and 1 nega-ve control. The same quality control filters were applied and only 1 sample 162 was removed from plate 3. Finally, one sample that showed 1.5 million reads and 1.3 million post-163 trimming reads as input, due to an error in library equimolar pooling, was subjected to the 164 subsampling process by randomly selec-ng 187K reads (the average value of reads survived post 165 trimming for Plate 3), through the R package ShortRead v1.62.0 (see URLs). 166 A_er all these QC steps, a total of 212 samples were le_ for downstream analyses with the median 167 number of reads per sample being 138,531 (127,518, 141,450 and 190,072 for plates 1, 2 and 3, 168 respec-vely). 169 170 Characteriza=on of the vaginal microbiota 171 To characterize the bacterial composi-on of the samples, the DADA2 v3.19.0 package [20] (see URLs) 172 was used in the R (v 4.4.1) environment. Using the "filterAndTrim()" func-on , only sequences with a 173 minimum length of 160 bp were selected. With the “removeBimeraDenovo()" func-on , all chimeras 174 were removed. 175 The microorganisms were classified up to the taxonomic level of species with the 176 "assignTaxonomy()" and "addSpecies()" func-on Genome Taxonomy DataBase (GTDB) v.r95 [21] was 177 used to profile a phylogene-cally consistent microbial taxonomy (see URLs). Hypervariable regions 178 V3V4 and V7V9 were processed together. Rare amplicon sequence variants (ASVs) resul-ng from 179 DADA2 were removed; only ASVs that had at least 5 reads in at least 2 samples were maintained for 180 .CC-BY-NC-ND 4.0 International licensemade available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprintthis version posted March 12, 2025. ; https://doi.org/10.1101/2025.03.12.642767doi: bioRxiv preprint 7 subsequent analyses. To increase taxonomic resolu-on for the genus Lactobacillus, we followed a 181 previously described protocol [11]. Using their code and data files provided (see URLs), the genus 182 Lactobacillus was reclassified into different subgenera, allowing the iden-fica-on of more 183 Lactobacillus species. Finally, Lactobacillus species names were replaced by the respec-ve subgenus 184 iden-fied ( Table S2). Next, we selected ambiguous ASVs – i.e. those ASVs that were assigned by the 185 algorithm as Lactobacillus species but belonged to other genera - and uploaded the corresponding 186 sequences on BLAST to iden-fy the most similar species. Since ambiguity remained, we decided to 187 assign these ASVs the label "unclassified" (Table S3). 188 To characterize the fungi composi-on of the samples, the ITS fastQ files were processed with the 189 DADA2 v3.19.0 [20] (see URLs) package on R (v.4.4.1). To this end, “filterAndTrim()” func-on was set 190 as “minLen = 50” and “maxLen = 200”. We profiled samples with the Unite database v.10.0 [22] (see 191 URLs). 192 We derived global metrics of bacterial microbiome composi-on such as alpha and beta diversity 193 using the vegan package v.2.6-4 [23] (see URLs). Alpha diversity was calculated based on taxa counts 194 according to the Shannon -Weiner (H) and Simpson (D ) indices. Beta diversity was calculated 195 according to the Bray-Cur-s dissimilarity index using t axa rela-ve abundances and according to 196 Euclidean distance using taxa rela-ve abundances normalized using the Centered log ra-o (CLR) 197 transforma-on . The Euclidean distance was used in addi-on to Bray-Cur-s because it allowed us to 198 use normalized and corrected data (Table S4). To perform CLR transforma-on we used our own 199 func-on that calculates the transforma-on in parallel among taxa, using the Parallel package v.4.4.1 200 on R (see URLs). 201 In addi-on , we classified our samples according to the Community State Types (CSTs) defini-on of 202 the VALENCIA algorithm [16] (see URLs). Toward this end, it was necessary to manually modify the 203 taxonomy names according to those required by VALENCIA (Table S5). To iden-fy which ASV could 204 be aNributed to Gardenella Vaginalis, given the recent change in nomenclature, we used DADA2 205 with the database SILVA v.138.1 where the nomenclatures of species align to those used in the 206 VALENCIA algorithm. 207 Finally, we used the t-SNE algorithm [24] to visualize the vaginal microbiome composi-on in a two-208 dimensional space and used color-coding labels to describe the observed varia-on in terms of taxa 209 abundance as well as CSTs. For this analysis, we employed the t-SNE algorithm to CLR-transformed 210 data using the R package Rtsne v.0.17 [25] (see URLs). 211 .CC-BY-NC-ND 4.0 International licensemade available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprintthis version posted March 12, 2025. ; https://doi.org/10.1101/2025.03.12.642767doi: bioRxiv preprint 8 Finally, we profiled microbiotas using the same pipeline but subsampling 30K (and 60K) survived 212 reads per sample using the R package ShortRead (see URLs). 213 214 Iden=fica=on of factors associated with beta diversity 215 We evaluated the impact on beta diversity of technical confounders, characteris-cs of the volunteers 216 (as recorded in ques-onnaires) and of sample collec-on. Using the vegan package [23] (see URLs), 217 we first inves-gated whether the first five Bray-Cur-s beta diversity principal coordinate analyses 218 (PCoA) axes were correlated with technical variables, such as microbial DNA concentra-on, library 219 concentra-on, and total reads counts (survived reads). Each PCoA axis was linearly adjusted using 220 only the covariates with which it showed significant correla-on. As the dataset included longitudinal 221 measurements, independence tests were conducted separately for each visit to ensure 222 independence between observa-ons. Namely, we used the Kruskal-Wallis independence test in R 223 through the func-on kruskal.test() of stats package v.4.4.1 [26]; (see URLs). 224 Variables from the ques-onnaires encompassed both categorical and con-nuous data related to 225 demographics (e.g., age), anthropometric measures (e.g., body mass index [BMI], waist-to-hip ra-o 226 [WHR]), lifestyle factors (e.g., smoking, hormonal contracep-ve use), clinical history related to 227 women’s health (e.g., number of pregnancies), and condi-ons surrounding vaginal swab collec-on 228 (e.g., medica-on or supplement use prior to collec-on, sexual intercourse within two days of 229 collec-on, swab -ming, and stool consistency assessed using the Bristol stool scale). Con-nuous 230 variables were categorized to allow applica-on of the independence test. All tested variables are 231 summarized in Table S6. 232 We also run PERMANOVA analyses for the significant associa-ons resul-ng from independence test, 233 using adonis2() func-on from vegan R package [23] (see URLs). 234 235 Sta=s=cal analysis to evaluate vaginal microbiota dynamics 236 To evaluate changes in alpha diversity, beta diversity and taxa rela-ve abundance across the four 237 phases of the menstrual cycle, we used the non -parametric Friedman test. To evaluate changes 238 between two phases (pairwise comparisons) we used the non-parametric Wilcoxon pairwise test for 239 the untransformed data, and the parametric pairwise t-test for the CLR-transformed data. All these 240 compara-ve tests were run using the stat package v.4.4.1 [26] (see URLs). All pairwise comparisons 241 were corrected using the Bonferroni approach. 242 .CC-BY-NC-ND 4.0 International licensemade available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprintthis version posted March 12, 2025. ; https://doi.org/10.1101/2025.03.12.642767doi: bioRxiv preprint 9 To evaluate beta diversity changes (both Bray-Cur-s and Euclidean distance, adjusted for technical 243 covariates) we used two approaches. First, we calculated mean and median beta diversity across all 244 women within each menstrual cycle phase and compared these across the different phases, using 245 the paired Wilcoxon test (Table S7). 246 In a second step, we calculated mean and median beta diversity for each woman among her different 247 collec-ons (up to 4), and compared the distribu-on obtained when comparing her first collec-on to 248 3 collec-ons from different women from the other 3 phases, using again the unpaired Wilcoxon test. 249 All comparisons were corrected using Bonferroni. 250 251 Linear Models 252 We fiNed three groups of linear models to: i) inves-gate whether bacterial abundances vary across 253 the menstrual cycle ; ii) explore the rela-onship between hormones and bacterial abundances 254 varia-ons between two consecu-ve phases and iii) explore the overall rela-onship between 255 hormones and bacteria rela-ve abundances. 256 For all linear models, taxa were filtered using a prevalence > 20% and a mean rela-ve abundance 257 threshold of 0.10 % (Table S8). 258 In the first group of models, we fiNed the following gaussian mixed models: 259 Taxa_adj ~ Visit_number + (1|Woman) (1) 260 where visit_number is a numeric variable coded as 1, 2, 3, 4 according to the four phases of the 261 menstrual cycle, and Taxa_adj are the residuals of the linear models: 262 Taxa_CLR-transformed ~ covariates 263 fiNed using the func-on “LinearRegression()” of sklearn package v.1.5.2 in Python v .3.10.12 [27] 264 (see URLs). The placeholder “covariates” here means that we considered several sets of covariates 265 based on the results of independence tests and technical variables to evaluate if these changes could 266 be confounded by a specific covariate. We opted for stepwise models with increasing number of 267 covariates, rather than running a model with the full set of covariates, given the limited sample size 268 of the study. A total of 10 models were considered (Table S9). 269 From model (1) we used the regression coefficients and the p-values for the variable Visit_number 270 to iden-fy taxa with significant changes . We used permuta-on to derive an empirical p-value and 271 evaluate robustness of significant results (p<0.05). We permuted the column containing taxa 272 abundances and re-run the model 1000 -mes. The propor-on of “fake” p-values less than or equal 273 .CC-BY-NC-ND 4.0 International licensemade available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprintthis version posted March 12, 2025. ; https://doi.org/10.1101/2025.03.12.642767doi: bioRxiv preprint 10 to the original p-value was calculated to provide the empirical p-value. In addi-on, f alse discovery 274 rate (FDR) was instead used to adjust p-values for mul-ple tes-ng (original pvalues, not empirical 275 pvalues), with the func-on “p.adjust()” of stats package and method ‘BH’ [28]. 276 In addi-on, we calculated the interclass correla-on coefficient (ICC) for the model where bacteria 277 were adjusted for the technical covariates (model i) above). This parameter q uan-fies the 278 propor-on of total variance aNributable to between-subject variability, reflec-ng the stability of 279 bacterial abundance across different -me points within subjects. The ICC was calculated using the 280 func-on “icc()” from the package performance v 0.12.4 [28] (see URLs). 281 In the s econd group of models , we again used linear mixed models to inves-gate if changes in 282 bacteria rela-ve abundance between two phases were associated with changes in sex hormones. 283 Firstly, we removed outlier values in the sex hormone variables exceeding 4 standard devia-ons from 284 the mean. Secondly , we calculated the differences between two consecu-ve phases of both sex 285 hormones and bacteria rela-ve abundances and then we normalized these differences using the 286 func-on “ordernorm()” of the R package bestNormalize v 1.9.1 [29](see URLs). 287 The following gaussian mixed models were fiNed: 288 (Hormonet+1 – Hormone t) ~ (Taxa_adj t+1– Taxa_adjt) + (1|Women) (2) 289 where t represents the menstrual cycle phase, and taxa_adj is the CLR-transformed abundance 290 adjusted by technical covariates, calculated as the residuals of the models: 291 Taxa_CLR_transformed ~ Qubit_DNA + Qubit_Library + Total_counts (3) 292 where Qubit_DNA is the DNA concentra-on, Qubit_Library is the library concentra-on, and 293 Total_counts is the total number of survived reads. 294 In the third group of models we inves-gate rela-onship s between taxa rela-ve abundances and sex 295 hormones changes across the en-re menstrual cycle by fing the following: 296 Hormone ~ Visit_number + Taxa_adj + (1|Woman) + Pregnancy_category (4) 297 In both models (2) and (4) we used the same strategies described above for models (1) to derive 298 empirical p-values and to correct p-values for mul-ple tes-ng. To fit all linear mixed models in (1), 299 (2) and (4) we used the lme4 v 1.1.35.5 [30] (see URLs) and the lmerTest v.3.1.3 [28] R packages 300 (see URLs). 301 302 Whole-genome sequencing 303 .CC-BY-NC-ND 4.0 International licensemade available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprintthis version posted March 12, 2025. ; https://doi.org/10.1101/2025.03.12.642767doi: bioRxiv preprint 11 Ten samples with different microbiome composi-on were selected for whole-genome sequencing 304 (WGS). We used the beta diversity matrix to maximize dissimilarity among selected samples, as 305 follows. The first three selected samples were those showing a high rela-ve abundance of L. iners, 306 L. crispatus and L. gasseri, respec-vely, and had the maximum beta diversity between them. Then, 307 we selected seven samples showing maximum distance from the first three samples. Sequencing 308 was performed at Prebiomics srl on the Illumina NovaSeq 6000 sequencing system, and microbiome 309 profiling was carried out using MetaPhlan 4 [31]. An average of 3 million bacterial reads per sample 310 was obtained (Table S10a). Rela-ve abundances from the two data sets (16S and WGS) were 311 compared for taxa with abundances > 0 in at least one of the samples. Spearman's correla-on test 312 was used to check the correla-on s between rela-ve abundances of taxa detected by both 313 approaches and for which at least 3 samples had no-zero value (Table S10b). 314

315 316 Vaginal Microbiome Composi=on 317 A total of 61 women contributed 213 vaginal swabs to this study, of which 212 passed quality control 318 (QC) as described in Methods. Bacterial composi-on was profiled using previously described 319 protocols (see Methods). At the genus level, Lactobacillus was the most common genus, followed 320 by Bifidobacterium and Streptococcus. The two most abundant Lactobacillus spp were L. iners and 321 L. crispatus (Figure S1) and were iden-fied in almost all samples (prevalence > 90%) . When 322 restric-ng to samples with a rela-ve abundance > 60%, L. iners was detected in 92 samples (from 38 323 women) while L. crispatus was detected in 74 samples (from 34 women). In addi-on to Lactobacillus 324 spp, we iden-fied other species, such as Bifidobacterium vaginale and Fannyhessea vaginae (Figure 325 S2 and Table S8). 326 Overall, vaginal microbiota remained stable throughout the menstrual cycle. There were no 327 significant differences in alpha diversity across the phases, as measured by the Shannon (H) and 328 Simpson (D) indices (Figure S3 and Table S11). The average beta diversity was instead higher in the 329 follicular phase compared to ovulatory phase (O) (pwilcoxon = 0.009 and pt-test = 0.007) and early luteal 330 phase (EL) (p wilcoxon = 0.04 and p t-test = 0.03) (Figure 1A and Table S7). Similar results were obtained 331 when using median instead of average (Figure S4 and Table S7), or the overall distribu-on (Figure 332 S5 and Table S12). Nevertheless, these changes were small and did not lead women’s microbiota to 333 .CC-BY-NC-ND 4.0 International licensemade available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprintthis version posted March 12, 2025. ; https://doi.org/10.1101/2025.03.12.642767doi: bioRxiv preprint 12 align to a similar composi-on during a par-cular phase. In fact, the average beta diversity of a 334 woman’s microbiota across all her samples was lower than the average obtained when using random 335 samples from other women (p value << 0.05 for all comparisons), indica-ng that microbiota is largely 336 determined by individual specificity, rather than by menstrual phase (Figure 1B). 337 Analyses of dynamic changes could not be performed for the mycobiome, due to the limited number 338 of samples for which mycobiome could be profiled. Although we sequenced the ribosomal internal 339 transcribed spacer (ITS) region in all 213 samples, only 25 samples had at least 200 reads, cutoff 340 deemed necessary to obtain fungi classifica-on (Figure S6 and Table S13). Twenty-five fungal species 341 were iden-fied, t he most prevalent being Candida albicans , a fungus associated with vaginal 342 infec-ons (candidiasis) but also frequently found as part of the normal vaginal microbiota. The 343 rela-vely low number of reads overall detected in the ITS region may be due to the use of a protocol 344 and kit op-mized for bacterial rather than fungal DNA extrac-on, as reported in previous studies 345 [32]. 346 347 The impact of host phenotypes 348 We then inves-gated the impact on vaginal microbiota of host phenotypes, lifestyle and informa-on 349 on sample collec-on, by correla-ng genus -based and species-based beta diversity PCoA axes with 350 available metadata. Among the 12 variables tested, the variable “having had at least one pregnancy” 351 had the higher impact at both genus and species level, followed by “having used oral hormonal 352 contracep-ves in the past” (Table S14). Other variables that were significant albeit with less impact 353 include age, waist-to-hip ra-o, smoking habits, intercourse in the days preceding the swab, -me of 354 collec-on, and use of medica-on or supplements (Figure S7). Intriguingly, using Permanova test we 355 observed that the microbiota of women who had at least one pregnancy before par-cipa-ng in the 356 study had a lower likelihood to be dominated by L. crispatus (p=0.005), in line with observa-on s 357 from the ISALA study [11], sugges-ng the existence of a poten-al lifelong pregnancy signature 358 (Figure S8). Likewise, using Permanova we observed that the microbiota of older women is less likely 359 dominated by Lactobacillus spp. (p=0.003), in line with the observed decline in their abundance 360 during menopause (Figure S9)[5]; [33]. Similarly, dominance of Lactobacillus spp. is less frequent 361 when swabs are collected a_er feces (p=0.02) (Figure S10), and in those collected in the morning 362 (p=0.003) (Figure S11). 363 364 Frequency of L. iners and L. crispatus along the menstrual cycle 365 .CC-BY-NC-ND 4.0 International licensemade available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprintthis version posted March 12, 2025. ; https://doi.org/10.1101/2025.03.12.642767doi: bioRxiv preprint 13 In contrast to findings from cohort studies on Northern European women, our study found a greater 366 abundance of L. iners than L. crispatus , a microbiome composi-on typically observed instead in 367 Hispanic and African popula-ons [11]; [6]; [12]; [14]; [2] (Table 1). We ruled out the possibility that 368 this result is driven by a poten-al bias from the 16S rRNA gene profiling method, since comparison 369 of the rela-ve abundances of these two Lactobacillus spp in 10 samples showed very high 370 concordance with profiles obtained using shot -gun whole-genome sequencing (WGS) method 371 (Figure S12 and Table S10b). Likewise, we also rule out the possibility that this result is driven by the 372 higher number of reads compared to protocols used in previous studies . In fact, taxa rela-ve 373 abundances were almost iden-cal when using a subsampling approach on our 16S rRNA data to 374 restrict analyses to maximum 30K or 60K reads per sample (Figure S13). 375 In our cohort, L. iners and L. crispatus have a similar prevalence in the follicular phase (98% and 95%, 376 respec-vely), but the average rela-ve abundance was higher for L. iners (average rela-ve abundance 377 40% and 25% for L. iners and L. crispatus, respec-vely). The difference in average rela-ve abundance 378 between these two species was negligible in the other 3 phases ( absolute difference between 379 percentages: 2.03, -1.27, 2.21 in the ovulatory, early and late luteal phases , respec-vely ). 380 Interes-ngly, all microbiotas dominated by L. crispatus in the follicular phase were stable during the 381 menstrual cycle, while a small frac-on ( 20%) of those dominated by L. iners at the follicular phase 382 switched to L. crispatus-dominant - or to other species - in the ovulatory or in the late luteal phase 383 (Figure 2). In addi-on, the microbiota of two women dominated by low prevalent species in the 384 follicular phase became L. crispatus-dominant in the ovulatory phase. Inter-class correla-on (ICC) 385 analyses indicated that these two species explained the largest variance of total variability between 386 women (Table S15). While t hese results suggest a trend for women’s microbiota to switch to L. 387 crispatus during the ovulatory phase, L. iners s-ll remains the most abundant species in the Italian 388 popula-on across all -me points. 389 390 Vaginal microbiota dynamics 391 Given these clear -cut changes observed in main Lactobacillus spp, we sought to thoroughly 392 inves-gate microbiota dynamics along the menstrual cycle. First, we assessed the composi-on and 393 varia-ons along the menstrual cycle using the CST classifica-on (see Methods). Overall, we found 394 that CST-V (L. Jensenii dominated) was the rarest among all samples and phases (2.5%), followed by 395 CST-II (L. gasseri dominated) (10.2%) and CST-IV (defined by moderate to high abundance of different 396 anaerobic species, such as Prevotella spp and G. vaginalis, and low abundance of Lactobacillus spp 397 .CC-BY-NC-ND 4.0 International licensemade available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprintthis version posted March 12, 2025. ; https://doi.org/10.1101/2025.03.12.642767doi: bioRxiv preprint 14 (9.2%) (Figure S14A). In line with rela-ve abundances, the CST-III (L. iners dominated) was the most 398 common (42.9%) and with higher prevalence in the follicular phase, while the second most common 399 CST was CST-I (L. crispatus dominated) (Figure S14B), whose prevalence increased during ovulatory 400 phase (Figure S14B). Most women maintained a stable CST throughout the menstrual cycle: only 11 401 out of 61 women changed CST at any -me point (18 if we consider changes in subCST) (Figure 3A). 402 Interes-ngly, 72% (N=8) of these changes occurred from follicular to ovulatory phases (Figure 3B), 403 and of which 6 were changes from CST-III or CST-IV (N=2) to CST-I (N=4). The remaining changes 404 occurred from early to luteal phase, with 3 women’s microbiota switching to CST-IV. 405 None of the women’s microbiotas classified as CST-I exhibited CST shi_s, indica-ng the high stability 406 of communi-es dominated by L. crispatus. Of the women who undergo CST shi_s, only one changed 407 CST types more than once ( Figure S15). We were unable to detect associa-on with CST/subCST 408 changes and presence of fungi, given the limited amount of informa-on ( only 4 out of 11 women 409 exhibi- ng shits in CST had fungi informa-on at the relevant -me points). Likewise, we did not detect 410 associa-on with CST/subCST changes with having had sex intercourse in the two days prior sample 411 collec-on (p>0.05). 412 Secondly, we examined the longitudinal changes of taxa rela-ve abundances using linear mixed 413 models and controlling for poten-al confounding variables, as described in Methods. We found that 414 the rela-ve abundance of a few bacteria species varies significantly across phases, with the majority 415 remaining stable when different covariates were adjusted. For the genus Lactobacillus, L. iners, L. 416 gasseri and L. jensenii showed significant varia-on across mostly all correc-ons applied, while L. 417 crispatus was not significant in any of the models (Figure 4). Likewise, for the genus Streptococcus, 418 two taxa, namely S. unclassified and S. agalac&ae varied significantly across phases regardless of 419 covariates used. Finally, another bacterial species, Dialister B micraerophilus showed significant 420 varia-on in 3 of the 10 models (Figure 4). 421 422 Longitudinal changes of bacteria abundance associate with sex hormone levels 423 We inves-gated w hether changes in Lactobacillus and Streptococcus species between two 424 consecu-ve phases were associated with corresponding changes in sex hormone levels. Changes in 425 L. iners between ovulatory and early luteal phases were associated with changes in luteinizing 426 hormone (LH) levels (p=0.04), and those between early and late luteal phases in L. jensenii with 427 prolac-n levels (p=0.04) (Table S16 and Figure 5A). These results were however only nominally 428 significant, therefore larger samples size and/or addi-onal cohorts are needed to be confirmed. 429 .CC-BY-NC-ND 4.0 International licensemade available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprintthis version posted March 12, 2025. ; https://doi.org/10.1101/2025.03.12.642767doi: bioRxiv preprint 15 Furthermore, we iden-fied several bacteria whose longitudinal fluctua-ons along the en-re 430 menstrual cycle were significantly associated with levels of sex hormones – for a total of 14 431 associa-on s (Figure 5B and Figure S16 and Table S17). Most of the associa-on s were observed with 432 beta-17 estradiol (BES17). Some were in common between BES17 and LH, including Prevotella genus 433 and species P. &monensis_A, and genera Finegoldia and Anaerococcus genera. Only one associa-on 434 was detected with prolac-n (PRL) , namely with abundance of Bifidobacterium vaginale. All these 435 associa-ons were nega-ve , indica-ng nega-ve correla-ons between sex hormones levels and the 436 rela-ve abundance of these bacteria ( Figure 5B and Figure S16). When considering permuted p-437 values to evaluate robustness of results to outliers, the number of significant associa-ons decreased 438 to 9 out of 14. Notably, 4 of these 9 associa-ons, all linked to BES17, also meet a stricter significance 439 threshold applied a_er correc-ng for mul-ple tes-ng, sugges-ng that these are likely genuine. 440 441

442 We carried out the first short-term longitudinal study on the vaginal microbiota of non-pregnant, 443 healthy women from an Italian popula-on. In the W4H cohort we observed a higher abundance of 444 L. iners compared to L. crispatus, in contrast from earlier studies on cohorts from Northern Europe 445 [6]; [11] (Table 1). Of note, a recent study on young French women [15] reported again a higher 446 prevalence of L. crispatus-dominated CST (CST-I) over the L. iners-dominated CST (CST-III), although 447 the difference in frequency between the two CSTs was less pronounced (40.5 vs 38.1% , while our 448 study has obtained 35% for CST-I and 43% for CST-III) (Figure S14) compared to observa-ons from 449 cohorts based in Belgium (ISALA study) , Danmark and Sweden (MiMens study) [11]; [6]. In th is 450 cohort, 70% of women were enrolled in Bordeaux (south of France). Hence, considered these results 451 and the higher prevalence of L. iners in African and Hispanic cohorts [12]; [2]; [13], our study of the 452 Italian popula-on provides evidence for a north-to-south gradient of increased L. iners abundance 453 even within Europe. 454 The reasons for these geographical differences in microbiota composi-on remain unclear but results 455 from our study and those from the ISALA study offer insights for specula-on. For example, in our 456 study we found that women who previously had children were less likely to have a microbio ta 457 dominated by L. crispatus. Likewise, in the ISALA study having children was nega-vely correlated 458 with L. crispatus and posi- vely with L. iners; this variable was the strongest host factor associated - 459 with opposite direc-on - with the abundance of these two species [11]. Geographic differences may 460 thus be related to the delivery mode , or to p ost-pregnancy health care prac-ces . Another 461 .CC-BY-NC-ND 4.0 International licensemade available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprintthis version posted March 12, 2025. ; https://doi.org/10.1101/2025.03.12.642767doi: bioRxiv preprint 16 contribu-ng factor might be related to different cultural and hygiene prac-ces, such as the soap type 462 used, the frequency of washing, the use of in-mate wipes and even the type of pads, tampons or 463 other periods supplies [34]; [4]. While we do not have this informa-on in our cohort, results from 464 the ISALA study indicate that the use of menstrual pads is nega-vely correlated with L. crispatus 465 abundance, while the use menstrual cups was posi-vely correlated [11]. Intriguingly, the ISALA study 466 also reported a posi-ve correla-on with L. Crispatus abundance with use of hormonal contracep-ve 467

(combina-on pill, vaginal ring or patch). Previously, also Tuddenham S. [35] showed that 468 the propor-on of CST I ( L. crispatus -dominated) was higher among users of oral hormonal 469 contracep-ve than the propor-on of CST I among non-users, albeit this result was significant only in 470 White and not African American. No associa-on was instead seen in these two studies betwe en 471 hormonal contracep-ves use and L.iners. 472 All our volunteers had a natural menstrual cycle and did not use hormonal contracep-ves for at least 473 the past 30 days. This, together with the overall limited use of hormonal contracep-ve in the Italian 474 popula-on ( 19.1% of female s during reproduc-ve age versus ~30% in Northern and Western 475 European countries, according to WHO [36] could explain the observed low abundance of L. 476 Crispatus. Finally, differences in abundance of the two Lactobacillus species may be aNributed to 477 varia-ons in the host -immune system, in the amount and quality of vaginal secre-on, in gene-c 478 variants of gene encoding receptors expressed on the surface of the vaginal epithelial cells, and 479 other gene-cally determined factors specific to each host. 480 The very high prevalence and abundance of L. iners in our cohort of healthy women is striking and 481 challenges the idea that L. iners dominated microbiotas are 'unhealthy' or ‘dysbio-c’. In fact, W4H 482 volunteers have reported no current symptoma-c vaginal infec-ons, so we can rule out the 483 abundance of L. iners in our samples due to either vaginal pathologies or infec-ons. We speculate 484 that the presence of different strains could be the most plausible explana-on of geographic and 485 clinical differences. Unfortunately, with the 16S sequencing approach we (and others) are unable to 486 iden-fy which strains of L. iners are present in the studied vaginal microbiota. However, given its 487 high gene-c variability, is reasonable to assume that certain strains of L. iners may promote healthy 488 condi-ons, while others may increase risk of infec-on , thus explaining its presence as species 489 abundance in both healthy and dysbio-c vaginal condi-ons [37]; [8]. In addi-on, a woman can be 490 colonized by mul-ple coexis-ng strains with different gene-c characteris-cs and metabolic 491 func-ons , enabling the survival of this species in a variety of vaginal microenvironments [38]. A deep 492 .CC-BY-NC-ND 4.0 International licensemade available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprintthis version posted March 12, 2025. ; https://doi.org/10.1101/2025.03.12.642767doi: bioRxiv preprint 17 gene-c study on L. iners strains is needed to beNer understand the ecology of this bacterial species 493 and their differen-al impact on infec-on s and other gynecological diseases. 494 Our results on dynamic changes along menstrual cycle phases agree with previous studies repor-ng 495 that vaginal microbiota is generally stable over -me and support current knowledge on the strong 496 stability of microbiomes dominated by L. crispatus [3]; [15]; [6]. Furthermore, we were able to 497 pinpoint to species responsible for the dynamic changes occurring in a minority of women, namely 498 L. iners, L. gasseri, L. jensenii and Streptococcus spp. The decreased stability conferred by these 499 species could be explained by the produc-on of two dis-nct isoforms of lac-c acid by these vaginal 500 microorganisms. Unlike L. Crispatus, L. iners, Bifidobacterium spp. and Streptoccosus spp. produce 501 the L- enan-omer of lac-c acid that allows easier coloniza-on by other microbes than D- enan-omer 502 [39]; [4]. Furthermore, L. iners possesses addi-onal genes compared to L. crispatus, encoding stress-503 tolerance proteins, iron–sulfur proteins, exogenous L-cysteine transport, and exhibits superior 504 metabolic adapta-on to the changing carbohydrate sources in the vaginal environment [40]; [1]. 505 Intriguingly, we observed that even if only a few women’s microbiotas change along the menstrual 506 cycle, most changes occur from follicular to ovulatory phase, with most of these aligning to an L. 507 crispatus profile. Beta diversity was also lower during the ovulatory and early luteal phases 508 compared to follicular phase, which could be due to a protec-ve mechanism to support a poten-ally 509 fer-lized ovum. It is possible that the fluctua-on of sex hormones plays a key role in vaginal 510 microbiota changes. Our analyses showed that BES17 has the major impact on longitudinal varia-on 511 in taxa abundances, although none of the hormone-taxa associa-on involved Lactobacillus species 512 and thus they cannot fully explain changes we observed during ovulatory phase. Nonetheless, these 513 associa-on s with BES17 are fascina-ng . In fact, we detected significant nega-ve rela-onships with 514 Prevotella spp ., Finegoldia and Corynebacterium genera, taxa that are o_en found at high 515 abundances in women a_er menopause with menopause-associated symptoms, such as vulvo -516 vaginal atrophy and genito-urinary symptoms [41]. The decline in estrogen levels during menopause 517 may play a prominent role in driving microbiota change, as reduced estrogen leads to a decrease in 518 glycogen by vaginal cells. This causes an elevated vaginal pH, crea-ng condi-ons that favor the 519 coloniza-on of mul-ple microbial species adapted to a less acidic environment. 520 Our study has several strengths and novel-es, including the high depth 16S rRNA sequencing, the 521 enrollment of naturally menstrua-ng women, the direct measurement of sex hormones, and the 522 inves-ga-on of a south European healthy popula-on. While our results have not a direct clinical 523 .CC-BY-NC-ND 4.0 International licensemade available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprintthis version posted March 12, 2025. ; https://doi.org/10.1101/2025.03.12.642767doi: bioRxiv preprint 18 implica-on, knowledge of vaginal microbiota dynamics under physiological sex hormones varia-ons 524 in all popula-ons is an essen-al route for inves-ga-ng the poten-al role of microbes also in 525 diagnosis and preven-on of diseases condi-ons such infer-lity and endometriosis. 526 We also recognize the limita-on s of the study. First, we did not collect samples during menses, a 527 stage where important changes in bacterial composi-on may occur. Moreover, we have focused on 528 a single menstrual cycle and therefore we were unable to assess if observed changes were cyclic, i.e. 529 the microbiota of that minority of women experiencing shi_ returns to the original profile a_er 530 menstrua-on. Furthermore, we have no measurement of vaginal pH, a crucial factor in the 531 regula-on of the microbiota. Finally, while we have employed a very deep 16S sequencing approach, 532 we are limited to the rela-ve abundance of bacteria and do not have informa-on on bacterial 533 func-on or strains, features that metagenomic sequencing could provide. 534 In conclusion, our study highlights the importance to extend inves-ga-ons on vaginal microbiota in 535 healthy women from diverse popula-ons, also from the same ethnicity, and the need to assess 536 stability and dynamics with sex hormone changes, including those occurring during a natural 537 menstrual cycle. 538 539 DATA AVAILABILITY STATEMENT 540 All our code used to analyze the data and instruc-ons are available at: hNps://github.com/Sanna-s-541 LAB/Women4Health 542 The human-depleted 16S rRNA and ITS sequences are deposited at the Sequence Read Archive (SRA) 543 repository under accession number PRJNA1222832 [which will be released upon acceptance of the 544 manuscript]. Par-cipant-level personally iden-fiable data cannot be shared in respect to 545 par-cipants' informed consent and are protected under the Italian Personal Data Protec-on Code 546 and European Regula-on 2016/679 of the European Parliament and of the Council (GDPR) that 547 prohibit distribu-on even in pseudo-anonymized form. 548 All data necessary to support the conclusions drawn in this study are available in the Supplementary 549 Tables. 550 551 URLs 552 Code for the reclassifica=on of Lactobacillus sp.: 553 hNps://github.com/LebeerLab/Ci-zen-science-map-of-the-vaginal-microbiome 554 .CC-BY-NC-ND 4.0 International licensemade available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprintthis version posted March 12, 2025. ; https://doi.org/10.1101/2025.03.12.642767doi: bioRxiv preprint 19 Code and files for VALENCIA algorithm 555 hNps://github.com/ravel-lab/VALENCIA 556 557 Tools 558 FastQC: hgps://www.bioinforma&cs.babraham.ac.uk/projects/fastqc/ 559 Trimmoma&c: hgps://github.com/&mflutre/trimmoma&c?tab=readme-ov-file 560 Mul&QC: hgps://github.com/Mul&QC/Mul&QC 561 R packages: 562 vegan: hgps://github.com/vegandevs/vegan/ 563 lme4: 10.32614/CRAN.package.lme4 564 lmerTest: 10.32614/CRAN.package.lmerTest 565 Rtsne: 10.32614/CRAN.package.Rtsne 566 Parallel: hgps://www.rdocumenta&on.org/packages/parallel/versions/3.6.2 567 STAT: hgps://doi.org/10.32614/CRAN.package.STAT 568 DADA2: hgps://www.bioconductor.org/packages/release/bioc/html/dada2.html 569 ShortRead: hgps://kasperdanielhansen.github.io/genbioconductor/html/ShortRead.html 570 Performance: 10.32614/CRAN.package.performance 571 sklearn: hgps://jmlr.csail.mit.edu/papers/v12/pedregosa11a.html 572 573 Databases: 574 GTDB: hgps://gtdb.ecogenomic.org/stats/r95 575 UNITE DB: hgps://doi.org/10.15156/BIO/2959330 576 577 AUTHOR’S CONTRIBUTION 578 Conceptualiza-on : F.D.S., G.G., S.S. 579 Wri-ng - original dra_ Prepara-on: E.V., F.C., A.M., S.S. 580 Wri-ng – review and edi-ng: F.B., M.L.F. 581 Project supervision: S.S. 582 Funding acquisi-on: G.G., S.S. 583 Sta-s-cal and Bioinforma-c analys es: E.V., F.C., S.S. 584 Cri-cal support for bioinforma-c analyses: D.V.Z., V. LF., R.G., J.S., N.K., A.Z., S.S. 585 Samples and data collec-on: S.L., S.C., G.V.B., A.K., F.D.S., D.M, G.G. 586 .CC-BY-NC-ND 4.0 International licensemade available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprintthis version posted March 12, 2025. ; https://doi.org/10.1101/2025.03.12.642767doi: bioRxiv preprint 20 Vaginal samples processing: A.M., F.Cr., S.I., F.B., M.L.F. 587 All authors read and revised the manuscript and approved its final version. 588 589 ACKNOWLEDGMENTS 590 We thank all the volunteers who par-cipated in the study for their -me and commitment. We thank 591 the physicians of the University of Trieste for their contribu-on to the recruitment phase, including 592 but not limited to, Francesco Cracco, Roberta Maria Gen-le, Elena Stefani and Michele 593 Stracquadaini, the Directors of I.R.C.C.S. Burlo-Garofolo and of IRGB-CNR for logis-c support, the 594 Agenzia Regionale Sardegna Ricerche for providing access to laboratories, technologists Davide 595 Murrau, Marco Masala and Michele Marongiu from the IRGB-CNR for IT support. Finally, we thank 596 Dr. Alessandra Meloni and Dr. Ferdinando Coghe from the Azienda Ospedaliero Universitaria di 597 Cagliari and Prof. Stefano Guerriero from the University of Cagliari for having joined the 598 Women4Health project and ini-a- ng a second recrui-ng center in Cagliari star-ng in October 2024. 599 600 FUNDING STATEMENT 601 This study was co-funded by the Italian Ministry of Health through the contribu-on given to the 602 Ins-tute for Maternal and Child Health IRCCS Burlo Garofolo, Trieste - Italy (SD 02/21 to G.G.), by the 603 European Union (ERC Stg 2022 to S.S., acronym SEMICYCLE, GA n.101075624), by the Next 604 Genera-on EU, in the context of the Na-onal Recovery and Resilience Plan, Investment PE8 – Project 605 Age-It: “Ageing Well in an Ageing Society, ” [DM 1557 11.10.2022 to S.S.] and by NutrAGE grant (CNR 606 Project FOE-2021 DBA.AD005.225) . Views and opinions expressed are , however, those of the 607 author(s) only and do not necessarily reflect those of the European Union or the European Research 608 Council. Neither the European Union nor the gran-ng authority can be held responsible for them. 609 In addi-on, S.S. also received a PRIN2022 grant by the Next Genera-on EU funds (DSB.PN004.021 610 2022PMZKEC_LS2_PRIN2022 SANNA, CUP B53D23008300006). 611 612 USE OF AI STATEMENT 613 The authors declare that they have used genera-ve ar-ficial intelligence, specifically Quillbot, 614 ChatGPT and Reverso online plarorms only to check spelling and grammar. The authors declare that 615 ChatGPT was also used to improve R scripts for handling data formang and genera-on of figures. 616 617 .CC-BY-NC-ND 4.0 International licensemade available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprintthis version posted March 12, 2025. ; https://doi.org/10.1101/2025.03.12.642767doi: bioRxiv preprint 21 Table 1. Prevalence and rela=ve abundance of L. iners and L. crispatus in our study and other European cohorts. 618 The table reports the sta-s-cs from previously published studies on European cohorts and those from our study, for L. iners and L. crispatus. The 619 column Time point indicates weather sta-s-cs refer to samples collected in a specific menstrual cycle phase or at any -me during menstrual cycle 620 (or menopause). 621 Frequency (%) Prevalence (%) Average Rela4ve Abundance (%) Study Country Sequencing

N reads (K) N samples Age range Time point CST-I CST-III L.iners L.crispatus L.iners L.crispatus ISALA Belgium 16S (V4) 25 3000 18-73 any 43 28 72 90 24 38 MiMens Denmark, Sweden 16S (V3V4) WGS -- 49 20-28 Follicular phase -- -- 63 76 12 59 Ovulatory phase -- -- 66 83 14 67 i-Predict France 16S (V3V4) 34 241 18-25 any 41 38 -- -- 38 45 W4H Italy 16S (V3V4, V7V9) 155 61 18-45 Follicular phase 22 31 98 95 40 25 Ovulatory phase 31 27 93 93 37 35 622 .CC-BY-NC-ND 4.0 International licensemade available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprintthis version posted March 12, 2025. ; https://doi.org/10.1101/2025.03.12.642767doi: bioRxiv preprint 22 Figures 623 Figure 1. Average beta diversity across menstrual cycle phases. 624 In (A) Average Beta diversity across the four menstrual cycle phases: Follicular phase (F), Ovulatory 625 phase (O), Early Luteal phase (EL), and Late Luteal phase (LL). The black line inside the boxplot 626 represents the median, with the lowest and highest values within the 1.5 interquar-le range 627 represented by the whiskers. The black lines between violins indicate significant differences (*p< 628 0.01). P-values were obtained using paired T-test. In (B) Average beta diversity within the same 629 woman across phases (Self) and average of beta diversity among different women (Random 1 to 630 Random 10). In both panels, beta diversity was calculated using Euclidean distances on transformed 631 taxa adjusted for technical covariates (see Material and methods). 632 633 A) 634 635 636 B) 637 638 639 .CC-BY-NC-ND 4.0 International licensemade available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprintthis version posted March 12, 2025. ; https://doi.org/10.1101/2025.03.12.642767doi: bioRxiv preprint 23 640 641 Figure 2: Distribu=on of the most abundant taxa across menstrual cycle phases. 642 Rela-ve abundance of the 10 most abundant species in the four different phases. In (A) samples 643 are sorted according to the abundance of L. iners, L. crispatus, L. gasseri, respec-vely. Samples in 644 (B), (C), and (D) are ordered as in (A). The grey bars show missing informa-on for those women in 645 the specific week. In the x-axis -tle we indicate for how many samples microbiome profile was 646 available. 647 648 A) 649 650 B) 651 652 653 654 655 656 657 .CC-BY-NC-ND 4.0 International licensemade available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprintthis version posted March 12, 2025. ; https://doi.org/10.1101/2025.03.12.642767doi: bioRxiv preprint 24 C) 658 659 D) 660 661 662 .CC-BY-NC-ND 4.0 International licensemade available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprintthis version posted March 12, 2025. ; https://doi.org/10.1101/2025.03.12.642767doi: bioRxiv preprint 25 Figure 3: Community State Type (CST) across menstrual cycle phases 663 In (A) t-SNE representa-on of the vaginal microbiome composi-on, with samples colored according 664 to their CST classifica-on. Arrows highlight instances where a woman changes CST between 665 consecu-ve phases, with the nota-on start_phase → final_phase displayed above each arrow. In 666 the boNom-le_ corner, a pie chart shows the propor-ons of phases where CST changes occur. The 667 legend provides the color coding for each CST. In (B) Alluvial plot showing Community State Type 668 (subCST) changes between phases. The four columns represent the Follicular (F), Ovulatory (O), Early 669 Luteal (EL), and Late Luteal (LL) phases, respec-vely. In legend of the plot: CH: a change occurred in 670 at least one phase, complete data available (15/61); CH + NA: a change occurred in at least one 671 phase, missing samples in some phases (3/61); NO CH: No changes occurred, complete data 672 available (26/61); NO CH + NA: No changes occurred, but missing samples in some phases (17/61). 673 674 A) 675 676 B) 677 678 679 .CC-BY-NC-ND 4.0 International licensemade available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprintthis version posted March 12, 2025. ; https://doi.org/10.1101/2025.03.12.642767doi: bioRxiv preprint 26 Figure 4: Longitudinal changes in taxa rela=ve abundance across phases 680 This figure shows the results of a linear model inves-ga-ng whether bacterial abundances change 681 across menstrual cycle phases. In each panel, the model coefficients for each taxon are shows with 682 points represen-ng the beta coefficient of the variable -me and the bars the 95% confidence interval 683 of the es-mate. The actual p -values are shown above each point, and the permuted p-values are 684 displayed in blue below each point for actual p -values <0.05, or NA otherwise. Colour coding 685 indicates the significance of the associa-ons according to original p-values. 686 Each panel corresponds to the same model using different correc-on for the taxa, with the addi-on 687 of an increasing number of covariates. In par-cular: model 1 is the basic model with CLR transformed 688 bacteria, model 2 adds to previous model the technical covariates, model 3 adds Age, model 4 adds 689 sex_less_than_two_days, model 5 adds Pregnancy_category, model 6 adds Pill_use, model 7 adds 690 Swab_a_er_feces, model 8 adds Bristol_stool_scale, model adds Swab_morning_or_not and model 691 10 adds WHR_ranges. 692 693 694 .CC-BY-NC-ND 4.0 International licensemade available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprintthis version posted March 12, 2025. ; https://doi.org/10.1101/2025.03.12.642767doi: bioRxiv preprint 27 Figure 5: Rela=onship between sex hormones and bacteria 695 These forest plots display the significant associa-ons between sex hormones and bacterial taxa 696 when considering two consecu-ve phases (A) or in the overall menustral cycle (B). The x -axis 697 represents the model coefficients, indica-ng the strength and direc-on of each associa-on, while 698 the y-axis lists bacterial taxa. Each line illustrates the 95% confidence interval for the coefficient, with 699 colors represen-ng the sex hormone associated with each taxon. The p -value for each model is 700 shown above the corresponding point, and the permuted p -value (in blue) is shown below. 701 Associa-ons that remain significant a_er FDR correc-on, with a threshold of 0.2, are highlighted 702 with an asterisk on the y -axis. In the legend, the color coding for the hormones is as follows: 703 BES17=beta-estradiol, FSH=follicle-s-mula-ng hormone, LH=luteinizing hormone, PRL=prolac-n, 704 and PROG=progesterone. 705 706 A) 707 708 709 B) 710 711 .CC-BY-NC-ND 4.0 International licensemade available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprintthis version posted March 12, 2025. ; https://doi.org/10.1101/2025.03.12.642767doi: bioRxiv preprint 28 712 713 714 715 .CC-BY-NC-ND 4.0 International licensemade available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprintthis version posted March 12, 2025. ; https://doi.org/10.1101/2025.03.12.642767doi: bioRxiv preprint 29 716

717 [1] M. T. France, H. Mendes-Soares, e L. J. Forney, «Genomic Comparisons of Lactobacillus 718 crispatus and Lactobacillus iners Reveal Potential Ecological Drivers of Community 719 Composition in the Vagina», Appl. Environ. Microbiol. , vol. 82, fasc. 24, Art. fasc. 24, dic. 720 2016, doi: 10.1128/AEM.02385-16. 721 [2] J. Ravel et al., «Vaginal microbiome of reproductive-age women», Proc. Natl. Acad. Sci., 722 vol. 108, fasc. supplement_1, pp. 4680–4687, mar. 2011, doi: 723 10.1073/pnas.1002611107. 724 [3] P. G a j e r et al., «Temporal Dynamics of the Human Vaginal Microbiota», Sci. Transl. Med. , 725 vol. 4, fasc. 132, pp. 132ra52-132ra52, mag. 2012, doi: 10.1126/scitranslmed.3003605. 726 [4] J. B. Holm et al., «Integrating compositional and functional content to describe vaginal 727 microbiomes in health and disease», Microbiome , vol. 11, fasc. 1, p. 259, nov. 2023, doi: 728 10.1186/s40168-023-01692-x. 729 [5] P. Łaniewski e M. M. Herbst-Kralovetz, «Connecting microbiome and menopause for 730 healthy ageing», Nat. Microbiol. , vol. 7, fasc. 3, pp. 354–358, mar. 2022, doi: 731 10.1038/s41564-022-01071-6. 732 [6] L. W. Hugerth et al., «Defining Vaginal Community Dynamics: daily microbiome 733 transitions, the role of menstruation, bacteriophages, and bacterial genes», Microbiome , 734 vol. 12, fasc. 1, Art. fasc. 1, ago. 2024, doi: 10.1186/s40168-024-01870-5. 735 [7] J. Ravel et al., «Daily temporal dynamics of vaginal microbiota before, during and after 736 episodes of bacterial vaginosis», Microbiome , vol. 1, fasc. 1, Art. fasc. 1, dic. 2013, doi: 737 10.1186/2049-2618-1-29. 738 [8] M. I. Petrova, G. Reid, M. Vaneechoutte, e S. Lebeer, «Lactobacillus iners: Friend or Foe?», 739 Trends Microbiol. , vol. 25, fasc. 3, pp. 182–191, mar. 2017, doi: 740 10.1016/j.tim.2016.11.007. 741 [9] S. M. Bloom et al., «Cysteine dependence of Lactobacillus iners is a potential 742 therapeutic target for vaginal microbiota modulation», Nat. Microbiol. , vol. 7, fasc. 3, Art. 743 fasc. 3, mar. 2022, doi: 10.1038/s41564-022-01070-7. 744 [10] C. A. Muzny, P . Łaniewski, J. R. Schwebke, e M. M. Herbst-Kralovetz, «Host–vaginal 745 microbiota interactions in the pathogenesis of bacterial vaginosis», Curr. Opin. Infect. 746 Dis., vol. 33, fasc. 1, p. 59, feb. 2020, doi: 10.1097/QCO.0000000000000620. 747 [11] S. Lebeer et al., «A citizen-science-enabled catalogue of the vaginal microbiome and 748 associated factors», Nat. Microbiol. , vol. 8, fasc. 11, Art. fasc. 11, nov. 2023, doi: 749 10.1038/s41564-023-01500-0. 750 [12] M. G. Serrano et al., «Racioethnic diversity in the dynamics of the vaginal microbiome 751 during pregnancy», Nat. Med. , vol. 25, fasc. 6, Art. fasc. 6, giu. 2019, doi: 752 10.1038/s41591-019-0465-8. 753 [13] X. Wei et al., «Vaginal microbiomes show ethnic evolutionary dynamics and positive 754 selection of Lactobacillus adhesins driven by a long-term niche-specific process», Cell 755 Rep., vol. 43, fasc. 4, p. 114078, apr. 2024, doi: 10.1016/j.celrep.2024.114078. 756 [14] O. S. E. Roachford, A. T. Alleyne, e K. E. Nelson, «Insights into the vaginal microbiome in a 757 diverse group of women of African, Asian and European ancestries», PeerJ, vol. 10, p. 758 e14449, nov. 2022, doi: 10.7717/peerj.14449. 759 [15] J. Tamarelle et al., «Vaginal microbiota stability over 18 months in young student women 760 in France», Eur. J. Clin. Microbiol. Infect. Dis. , set. 2024, doi: 10.1007/s10096-024-04943-761 3. 762 .CC-BY-NC-ND 4.0 International licensemade available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprintthis version posted March 12, 2025. ; https://doi.org/10.1101/2025.03.12.642767doi: bioRxiv preprint 30 [16] M. T. France et al., «VALENCIA: a nearest centroid classification method for vaginal 763 microbial communities based on composition», Microbiome , vol. 8, fasc. 1, Art. fasc. 1, 764 nov. 2020, doi: 10.1186/s40168-020-00934-6. 765 [17] F . Busonero et al., «The Women4Health cohort: a unique cohort to study women-specific 766 mechanisms of cardio-metabolic regulation», Eur. Heart J. Open, vol. 4, fasc. 2, Art. fasc. 767 2, mar. 2024, doi: 10.1093/ehjopen/oeae012. 768 [18] P . Ewels, M. Magnusson, S. Lundin, e M. Käller, «MultiQC: summarize analysis results for 769 multiple tools and samples in a single report», Bioinforma. Oxf. Engl., vol. 32, fasc. 19, 770 Art. fasc. 19, ott. 2016, doi: 10.1093/bioinformatics/btw354. 771 [19] A. M. Bolger, M. Lohse, e B. Usadel, «Trimmomatic: a flexible trimmer for Illumina 772 sequence data», Bioinforma. Oxf. Engl., vol. 30, fasc. 15, Art. fasc. 15, ago. 2014, doi: 773 10.1093/bioinformatics/btu170. 774 [20] B. J. Callahan, P . J. McMurdie, M. J. Rosen, A. W. Han, A. J. A. Johnson, e S. P . Holmes, 775 «DADA2: High resolution sample inference from Illumina amplicon data», Nat. Methods , 776 vol. 13, fasc. 7, Art. fasc. 7, lug. 2016, doi: 10.1038/nmeth.3869. 777 [21] D. H. Parks, M. Chuvochina, P .-A. Chaumeil, C. Rinke, A. J. Mussig, e P . Hugenholtz, «A 778 complete domain-to-species taxonomy for Bacteria and Archaea», Nat. Biotechnol., vol. 779 38, fasc. 9, Art. fasc. 9, set. 2020, doi: 10.1038/s41587-020-0501-8. 780 [22] K. Abarenkov et al., «The UNITE database for molecular identification and taxonomic 781 communication of fungi and other eukaryotes: sequences, taxa and classifications 782 reconsidered», Nucleic Acids Res., vol. 52, fasc. D1, Art. fasc. D1, gen. 2024, doi: 783 10.1093/nar/gkad1039. 784 [23] J. Oksanen et al., «Vegan: Community Ecology Package. R package version 2.0-2», gen. 785 2012. 786 [24] L. van der Maaten e G. Hinton, «Visualizing Data using t-SNE», J. Mach. Learn. Res. , vol. 9, 787 fasc. 86, pp. 2579–2605, 2008. 788 [25] J. Krijthe e L. van der M. (Author of original C. code), Rtsne: T-Distributed Stochastic 789 Neighbor Embedding using a Barnes-Hut Implementation. (7 dicembre 2023). 790 Consultato: 16 dicembre 2024. [Online]. Disponibile su: https://cran.r-791 project.org/web/packages/Rtsne/index.html 792 [26] K. Bolar, STAT: Interactive Document for Working with Basic Statistical Analysis. (1 aprile 793 2019). Consultato: 16 dicembre 2024. [Online]. Disponibile su: https://cran.r-794 project.org/web/packages/STAT/index.html 795 [27] F. Pe d r e g o s a et al., «Scikit-learn: Machine Learning in Python», J Mach Learn Res , vol. 12, 796 fasc. null, pp. 2825–2830, nov. 2011. 797 [28] A. Kuznetsova, P . B. Brockhop, R. H. B. Christensen, e S. P . Jensen, lmerTest: Tests in 798 Linear Mixed EMects Models . (23 ottobre 2020). Consultato: 16 dicembre 2024. [Online]. 799 Disponibile su: https://cran.r-project.org/web/packages/lmerTest/index.html 800 [29] R. A. Peterson, «Finding Optimal Normalizing Transformations via bestNormalize», R J., 801 vol. 13, fasc. 1, pp. 294–313, giu. 2021. 802 [30] D. Bates et al., lme4: Linear Mixed -EMects Models using «Eigen» and S4. (11 gennaio 803 2025). Consultato: 14 gennaio 2025. [Online]. Disponibile su: https://cran.r-804 project.org/web/packages/lme4/index.html 805 [31] A. Blanco-Míguez et al., «Extending and improving metagenomic taxonomic profiling with 806 uncharacterized species using MetaPhlAn 4», Nat. Biotechnol., vol. 41, fasc. 11, Art. fasc. 807 11, nov. 2023, doi: 10.1038/s41587-023-01688-w. 808 .CC-BY-NC-ND 4.0 International licensemade available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprintthis version posted March 12, 2025. ; https://doi.org/10.1101/2025.03.12.642767doi: bioRxiv preprint 31 [32] S. Ahannach et al., «Microbial enrichment and storage for metagenomics of vaginal, skin, 809 and saliva samples», iScience, vol. 24, fasc. 11, p. 103306, nov. 2021, doi: 810 10.1016/j.isci.2021.103306. 811 [33] A. L. Muhleisen e M. M. Herbst-Kralovetz, «Menopause and the vaginal microbiome», 812 Maturitas , vol. 91, pp. 42–50, set. 2016, doi: 10.1016/j.maturitas.2016.05.015. 813 [34] A. M. Holdcroft, D. J. Ireland, e M. S. Payne, «The Vaginal Microbiome in Health and 814 Disease-What Role Do Common Intimate Hygiene Practices Play?», Microorganisms , vol. 815 11, fasc. 2, p. 298, gen. 2023, doi: 10.3390/microorganisms11020298. 816 [35] S. Tuddenham et al., «Lactobacillus-dominance and rapid stabilization of vaginal 817 microbiota in combined oral contraceptive pill users examined through a longitudinal 818 cohort study with frequent vaginal sampling over two years», eBioMedicine , vol. 87, p. 819 104407, gen. 2023, doi: 10.1016/j.ebiom.2022.104407. 820 [36] United Nations, Contraceptive Use by Method 2019: Data Booklet . United Nations, 2019. 821 doi: 10.18356/1bd58a10-en. 822 [37] W. Kwak et al., «Complete Genome of Lactobacillus iners KY Using Flongle Provides 823 Insight Into the Genetic Background of Optimal Adaption to Vaginal Econiche», Front. 824 Microbiol. , vol. 11, mag. 2020, doi: 10.3389/fmicb.2020.01048. 825 [38] J. B. Holm, K. A. Carter, J. Ravel, e R. M. Brotman, «Lactobacillus iners and Genital 826 Health: Molecular Clues to an Enigmatic Vaginal Species», Curr. Infect. Dis. Rep., vol. 25, 827 fasc. 4, pp. 67–75, apr. 2023, doi: 10.1007/s11908-023-00798-5. 828 [39] S. S. Witkin, H. Mendes-Soares, I. M. Linhares, A. Jayaram, W. J. Ledger, e L. J. Forney, 829 «Influence of Vaginal Bacteria and d- and l-Lactic Acid Isomers on Vaginal Extracellular 830 Matrix Metalloproteinase Inducer: Implications for Protection against Upper Genital Tract 831 Infections», mBio, vol. 4, fasc. 4, p. 10.1128/mbio.00460-13, ago. 2013, doi: 832 10.1128/mbio.00460-13. 833 [40] J. M. Macklaim, G. B. Gloor, K. C. Anukam, S. Cribby, e G. Reid, «At the crossroads of 834 vaginal health and disease, the genome sequence of Lactobacillus iners AB-1», Proc. 835 Natl. Acad. Sci., vol. 108, fasc. supplement_1, pp. 4688–4695, mar. 2011, doi: 836 10.1073/pnas.1000086107. 837 [41] M. N. Anahtar et al., «Cervicovaginal Bacteria Are a Major Modulator of Host 838 Inflammatory Responses in the Female Genital Tract», Immunity, vol. 42, fasc. 5, pp. 839 965–976, mag. 2015, doi: 10.1016/j.immuni.2015.04.019. 840 841 842 .CC-BY-NC-ND 4.0 International licensemade available under a (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is The copyright holder for this preprintthis version posted March 12, 2025. ; https://doi.org/10.1101/2025.03.12.642767doi: bioRxiv preprint