molBV reveals immune landscape of bacterial vaginosis and predicts human papillomavirus infection natural history

Abstract

Bacterial vaginosis (BV) is a highly prevalent condition that is associated with adverse health outcomes. It has been proposed that BV's role as a pathogenic condition is mediated via bacteria-induced inflammation. However, the complex interplay between vaginal microbes and host immune factors has yet to be clearly elucidated. Here, we develop molBV, a 16 S rRNA gene amplicon-based classification pipeline that generates a molecular score and diagnoses BV with the same accuracy as the current gold standard method (i.e., Nugent score). Using 3 confirmatory cohorts we show that molBV is independent of the 16 S rRNA region and generalizable across populations. We use the score in a cohort without clinical BV states, but with measures of HPV infection history and immune markers, to reveal that BV-associated increases in the IL-1β/IP-10 cytokine ratio directly predicts clearance of incident high-risk HPV infection (HR = 1.86, 95% CI: 1.19-2.9). Furthermore, we identify an alternate inflammatory BV signature characterized by elevated TNF-α/MIP-1β ratio that is prospectively associated with progression of incident infections to CIN2 + (OR = 2.81, 95% CI: 1.62-5.42). Thus, BV is a heterogeneous condition that activates different arms of the immune response, which in turn are independent risk factors for HR-HPV clearance and progression. Clinical Trial registration number: The CVT trial has been registered under: NCT00128661.

Conflict of interest statement

The authors declare no competing interests.

© 2022. The Author(s).

Figures

Fig. 1. Microbial community features of bacterial vaginosis.

Panel A shows a heatmap of the 15 most prevalent bacterial species that are indicated to the right of the heatmap. Each column represents a participant. Hierarchical clustering separates samples into two primary clades: one dominated by Lactobacillus and one with polymicrobialism. There is a significant tendency of the BV-positive cases to be found in the polymicrobial clade and BV-negative in the one dominated by Lactobacillus based on either the Amsel or Nugent diagnosis (see the “Methods” section) (p < 0.001 for both). Panel B shows the alpha diversity differences between BV diagnosed by either Amsel or Nugent criteria and the microbial communities based on the Chao1, Fisher, and Shannon diversity indices (all p < 0.001), as indicated at the right of the panels. Panels C and D show beta diversity analyses using PCoA and the Jensen–Shannon diversity index for the Amsel BV (panel C) (R2 = 0.25, p < 0.001) and Nugent BV diagnosis (panel D) (R2 = 0.59, p < 0.001). Panel E shows the top 20 microbial markers (based on W-stat) for detecting BV using the “clean” BV status sample set (Amsel+/Nugent+ vs. Amsel−/Nugent−). The y-axis represents the ANCOM W-stat, while the x-axis represents the mean relative abundance difference between BV+ and BV− cases for each bacterial taxon. The size of the circles represents relative abundance. Source data are provided as a Source Data file.

Fig. 2. Cytokines associated with molBV categorical states in the CVT dataset.

Panel A forest plot showing the OR and 95% confidence intervals computed using a linear model between the cytokine quartiles and the three ordinal states of molBV (i.e., BV-negative, BV-intermediate, and BV-positive) derived using 431 biologically independent CVM samples. Only cytokines with an adjusted q-value < 0.001 are presented. ORs in panel A represent the odds of moving to the immediate next ordinal BV state. Panel B shows the correlation network between all cytokines with a Pearson correlation >0.6 using all baseline samples. Panel C is a volcano plot showing the ORs on the x-axis when comparing molecular BV-negative vs. BV-positive and the –log(q-value) on the y-axis. The ratios that had a q-value < 3.77*10-44 (i.e., −log(q-value)>100) are indicated and labeled in red. Panel D shows the pairwise Pearson correlation of the highly significant ratios presented in panel C (colored red). Ratios that had an OR < 1.0 were inverted for symmetry of data presentation. Clusters from the strongest molBV-associated cytokine ratios appear to fall into two primary groups; ones that include IL-1β and those with TNF-α. Panel E presents a box and violin plot for the log(IL-1β/IP-10) ratio for BV-negative (colored in blue), BV-intermediate (colored in yellow), and BV-positive (colored in red) samples. Panel F shows the bacterial species identified by ANCOM predicting IL-1β/IP-10 inflammation (i.e., above the median) in women within the BV-negative group shown in the small figure below with a red border. Comparison is between samples above the median vs. below the median. Panel G shows the ORs and 95% confidence intervals for the top 10 ratio combinations of the 32-cytokines (based on adjusted q-value) when comparing BV-negative to BV-positive women that had IL-1β/IP-10 levels below the cohort median (n = 171 biologically independent samples). Source data are provided as a Source Data file.

Fig. 3. Molecular BV, inflammation, and HR-HPV clearance.

Panel A shows the Kaplan–Meier curves for HR-HPV clearance colored by sustained BV status (i.e., having a molBV value above (red) or below (blue) the cohort median for both measured visits) with the unadjusted p-value presented in the bottom left corner of the plot and number of women at risk at each time point presented in the accompanying table directly below the plot. Panel B shows the Kaplan–Meier curves comparing women with either sustained high (red) or low (blue) IL-1β/IP-10 inflammation markers across the two analyzed visits. Sustained high and low IL-1β/IP-10 refers to women that had IL-1β/IP-10 levels above or below the cohort median for both measured visits, respectively, with the unadjusted p-value presented in the bottom left corner of the plot and number of women at risk at each time point presented in the accompanying table directly below the plot. Time is shown in days at the bottom of the figures and tables. Source data are provided as a Source Data file.