Study design and analytical framework

The present study was designed as a retrospective laboratory-based analytical investigation aimed at evaluating AMR and MDR patterns among bacterial isolates recovered from HVS specimens. The study integrated conventional microbiological assessment with advanced statistical and multivariate computational analyses to comprehensively characterize resistance dynamics among clinically important vaginal bacterial pathogens. A total of 220 bacterial isolates were included in the final analytical dataset. The investigation focused primarily on the prevalence of bacterial species, antimicrobial susceptibility patterns, MDR burden, organism-specific resistance profiles, and clustering behavior of resistance phenotypes.

Study population and sample collection

Clinical HVS specimens were collected from female patients presenting with clinical symptoms suggestive of vaginal infection, including abnormal vaginal discharge, irritation, and suspected bacterial vaginitis. Samples were processed consecutively during the study period from January 2024 to December 2024 to minimize selection bias.

Only culture-positive, non-duplicate bacterial isolates with complete microbiological identification and antimicrobial susceptibility testing (AST) records were included. Duplicate isolates obtained from the same patient, incomplete laboratory records, contaminated cultures, and samples without complete antimicrobial susceptibility data were excluded from analysis.

The final dataset consisted of 220 eligible bacterial isolates recovered from HVS specimens. Each isolate represented a single patient episode and was considered an independent observational unit for statistical analysis. For patients showing multiple bacterial growth, the clinically significant predominant isolate was selected according to laboratory interpretation criteria. Only one representative isolate per infectious episode was included to maintain independence of observations.

Because this was a retrospective isolate-based analysis, the study focused on eligible culture-positive bacterial isolates with complete antimicrobial susceptibility records. Complete information regarding total culture-negative specimens was not available from archived laboratory records.

Isolation and identification of bacterial pathogens

All HVS specimens were processed according to standard clinical microbiology laboratory protocols. Initial assessment included direct microscopic examination and Gram staining to evaluate bacterial morphology and inflammatory characteristics.

Samples were inoculated onto blood agar, MacConkey agar, and chocolate agar plates and incubated aerobically at 35–37°C for 18–24 hours. After incubation, bacterial growth was evaluated based on colony morphology, hemolytic characteristics, pigmentation, and Gram staining reactions.

Bacterial identification was performed using conventional biochemical identification methods. S. aureus isolates were confirmed using catalase, coagulase, and mannitol fermentation tests. Gram-negative organisms including E. coli, Klebsiella spp., Citrobacter spp., and Pseudomonas aeruginosa, were identified using oxidase testing, indole production, citrate utilization, urease test, triple sugar iron (TSI) reaction, motility testing, and lactose fermentation characteristics following standard clinical microbiology identification procedures.

Antimicrobial susceptibility testing

AST was performed using the Kirby–Bauer disc diffusion method on Mueller–Hinton agar (Oxoid Ltd., Basingstoke, Hampshire, United Kingdom; prepared according to the manufacturer’s instructions) following standard Clinical and Laboratory Standards Institute (CLSI) interpretive guidelines. Antimicrobial susceptibility interpretation was performed according to CLSI M100 performance standards (2024 edition). Commercially prepared antibiotic discs (Oxoid Ltd., Basingstoke, Hampshire, United Kingdom) with manufacturer-recommended standard concentrations were used for susceptibility testing. Commercially prepared antibiotic discs with standard potencies were used, including penicillin (10 units), erythromycin (15 µg), clarithromycin (15 µg), cefoxitin (30 µg), ciprofloxacin (5 µg), levofloxacin (5 µg), gentamicin (10 µg), amikacin (30 µg), piperacillin (100 µg), ceftazidime (30 µg), meropenem (10 µg), imipenem (10 µg), chloramphenicol (30 µg), vancomycin (30 µg), trimethoprim-sulfamethoxazole (25 µg), and linezolid (30 µg). Following inoculation, Mueller–Hinton agar plates were incubated aerobically at 35–37°C for 18–24 hours using a laboratory incubator (Memmert GmbH + Co. KG, Schwabach, Germany). Quality control of AST was performed using reference bacterial strains, including E. coli ATCC 25922 and S. aureus ATCC 25923, to ensure accuracy and reliability of susceptibility results. Antibiotic susceptibility results were categorized as sensitive, intermediate, or resistant according to standardized zone diameter breakpoints. The antibiotic panel included agents representing multiple clinically relevant antimicrobial classes commonly used for the treatment of bacterial and vaginal infections. The tested antibiotics included penicillin, erythromycin, clarithromycin, cefoxitin, ciprofloxacin, levofloxacin, gentamicin, amikacin, piperacillin, ceftazidime, meropenem, imipenem, chloramphenicol, vancomycin, sulphamethoxazole, and linezolid. Resistance prevalence for each antibiotic was calculated as the proportion of resistant isolates among the total number of tested isolates. Quality control of AST was performed according to CLSI recommendations using reference bacterial strains, including E. coli ATCC 25922 and S. aureus ATCC 25923, to ensure accuracy and reproducibility of susceptibility results.

Definition of multidrug resistance and resistance phenotypes

MDR was defined as acquired non-susceptibility to at least one antimicrobial agent in three or more antimicrobial classes according to the international standardized criteria proposed by Magiorakos et al. [13]. MDR classification was assigned individually to each isolate following evaluation of the complete susceptibility profile. To quantify organism-level resistance burden, an MDR score was calculated for each isolate based on the total number of antibiotics against which resistance was observed.

MRSA status was determined phenotypically using cefoxitin resistance as the surrogate marker for methicillin resistance. Similarly, ESBL production among E. coli isolates was confirmed phenotypically using the combined disc diffusion method based on ceftazidime and ceftazidime-clavulanic acid discs according to CLSI recommendations. An increase of ≥ 5 mm in inhibition zone diameter in the presence of clavulanic acid was interpreted as ESBL-positive. These resistance-associated variables were incorporated into the subsequent statistical and multivariate analyses.

Data processing and management

All microbiological and susceptibility data were compiled into a structured analytical dataset prior to statistical analysis. Data preprocessing involved removal of duplicate entries, correction of formatting inconsistencies, standardization of categorical variables, and harmonization of antimicrobial susceptibility coding. Categorical variables included organism type, MDR status, MRSA status, ESBL status, pregnancy status, and susceptibility outcomes. Antibiotic resistance profiles were transformed into binary numerical matrices for multivariate analyses, clustering procedures, and heatmap generation. The finalized dataset was formatted specifically for statistical computing and visualization within the RStudio environment.

Statistical analysis

All statistical analyses were performed using R statistical software (R Foundation for Statistical Computing, Vienna, Austria) within the RStudio integrated development environment. Descriptive statistical analyses were conducted to summarize bacterial isolate frequencies, antibiotic resistance prevalence, MDR distribution, patient age characteristics, and organism-specific resistance burdens. Continuous variables were summarized using means, medians, and distributional measures, whereas categorical variables were expressed as frequencies and percentages.

Inferential statistical analyses were conducted to evaluate associations between categorical variables. Chi-square tests were applied to assess relationships between organism type and MDR status, age category and resistance status, as well as MRSA and ESBL phenotypes. Statistical significance was considered at a p-value threshold of less than 0.05.

Binary logistic regression analysis was additionally performed to identify potential predictors associated with MDR occurrence. Predictor variables incorporated into the regression model included organism category, patient age, and pregnancy status. Adjusted odds ratios (ORs), standard errors, z-values, and corresponding p-values were calculated to evaluate the relative contribution of each variable to MDR development.

Multivariate and computational resistance analysis

To investigate global resistance structure and inter-isolate similarity patterns, several multivariate analytical approaches were employed. Hierarchical agglomerative clustering was performed using Euclidean distance matrices derived from isolate-level AMR profiles. Ward’s linkage method was applied to identify major resistance clusters among clinical isolates, and dendrogram visualization was subsequently used to interpret clustering relationships and resistance similarity structures.

Heatmap analysis was conducted using binary AMR matrices to visualize organism-specific resistance signatures and co-resistance patterns. Hierarchical clustering of both isolates and antibiotics was integrated into the heatmap framework to facilitate interpretation of resistance clustering behavior across bacterial species.

Principal component analysis (PCA) was further performed to reduce dimensional complexity and identify the major drivers of resistance variability among isolates. Antibiotic susceptibility variables were numerically transformed prior to ordination analysis. PCA loading vectors, explained variance percentages, and isolate clustering distributions were examined to determine the principal contributors to resistance differentiation within the study population.

Data visualization and figure preparation

Publication-quality graphical visualizations were generated to support descriptive, inferential, and multivariate interpretation of the data. Visualization approaches included histograms, horizontal bar charts, MDR distribution plots, organism prevalence graphs, dendrograms, heatmaps, and PCA biplots. Figures were generated using R-based visualization packages and exported in high-resolution formats suitable for thesis preparation and submission to peer-reviewed Q1-indexed journals.

Software and computational packages

Statistical analyses and graphical visualization were performed using R software version 4.4.1 (R Foundation for Statistical Computing, Vienna, Austria) through RStudio version 2026.01.2+418 (Posit Software, PBC, Boston, MA, USA). Data processing and visualization were performed using tidyverse, ggplot2, dplyr, pheatmap, factoextra, cluster, dendextend, reshape2, and scales packages. The software environment and computational packages were used for data organization, statistical analysis, hierarchical clustering, PCA, heatmap construction, and generation of publication-quality graphical outputs.

Ethical considerations

The study utilized anonymized laboratory-derived bacterial isolate data obtained during routine diagnostic procedures. No personally identifiable patient information was accessed during data processing or manuscript preparation. According to institutional policy governing retrospective laboratory-based studies using anonymized secondary microbiological data, formal ethical review and informed consent requirements were waived because no direct patient interaction or intervention was involved. All procedures were conducted in accordance with accepted ethical standards for retrospective microbiological research.

Demographic characteristics of the study population

A total of 220 bacterial isolates recovered from human HVS samples were included in the present investigation. The age distribution of participants ranged from late adolescence to approximately 50 years of age, representing predominantly women of reproductive age. The mean patient age was 32.5 years, while the median age was 33 years, indicating a relatively symmetrical age distribution without substantial skewness (Figure 1). The majority of isolates were recovered from women within the third and fourth decades of life, consistent with the increased susceptibility of reproductive-age women to vaginal bacterial infections reported in previous studies.

Display full size

Figure 1. Age distribution of patients included in the antimicrobial resistance study. The histogram illustrates the frequency distribution of patient ages, while dashed vertical lines indicate the mean and median age values.

Distribution of clinical bacterial isolates

Among the bacterial isolates recovered from HVS samples, S. aureus was the most predominant organism, accounting for 110 isolates (50.0%), followed by E. coli with 85 isolates (38.6%). Lower frequencies were observed for Klebsiella spp., Citrobacter spp., and P. aeruginosa. The predominance of S. aureus and E. coli is in agreement with several previous investigations identifying these organisms as the leading pathogens associated with vaginal infections and abnormal vaginal discharge.

The marked abundance of S. aureus may reflect its ability to colonize the vaginal tract as both a commensal and opportunistic pathogen, whereas the high frequency of E. coli likely reflects fecal contamination and disruption of the normal vaginal microbiota. The recovery of additional Gram-negative organisms further demonstrates the polymicrobial nature of vaginal infections (Table 1 & Figure 2).

Display full size

Figure 2. Frequency distribution of bacterial isolates among clinical samples. Bars represent the absolute number and relative proportion of each bacterial species identified in the study population.

Overall antimicrobial resistance pattern

The overall antimicrobial susceptibility analysis demonstrated extensive resistance against several commonly prescribed antibiotics among the recovered bacterial isolates. The highest resistance frequencies were observed for ceftazidime (89.1%), clarithromycin (83.6%), erythromycin (80.9%), cefoxitin (79.1%), levofloxacin (77.7%), and penicillin (77.3%). Similarly, elevated resistance levels were identified against ciprofloxacin (70.5%), meropenem (55.5%), and gentamicin (50.9%). In contrast, comparatively lower resistance rates were observed for sulphamethoxazole (29.5%), amikacin (20.0%), linezolid (15.5%), imipenem (11.8%), chloramphenicol (11.4%), and vancomycin (10.9%) (Figure 3 & Table 2).

Display full size

Figure 3. Overall antimicrobial resistance pattern among bacterial isolates recovered from high vaginal swab specimens. Bars represent the percentage resistance observed for each tested antibiotic across all clinical isolates.

The overall resistance percentages should be interpreted with caution because not all antimicrobial agents were applicable to every bacterial group. Certain antibiotics, including vancomycin and linezolid, are primarily relevant for Gram-positive organisms, whereas some β-lactam agents are mainly interpreted among Gram-negative bacteria. Therefore, organism-specific resistance patterns provide more clinically meaningful interpretation.

Organism-specific AMR analysis demonstrated substantial variability in resistance patterns among the recovered bacterial pathogens (Table 3). Among Gram-negative isolates, P. aeruginosa exhibited complete resistance to ceftazidime and meropenem (100%), while Citrobacter spp. demonstrated complete resistance to gentamicin (100%). Klebsiella spp. also showed markedly elevated resistance against ceftazidime (90.0%), gentamicin (90.0%), and ciprofloxacin (80.0%). In E. coli, high resistance frequencies were observed for ceftazidime (88.2%), gentamicin (77.6%), and levofloxacin (72.9%). Among Gram-positive isolates, S. aureus demonstrated pronounced resistance against clarithromycin (83.6%), erythromycin (80.9%), levofloxacin (80.9%), cefoxitin (79.1%), and penicillin (77.3%). Conversely, comparatively lower resistance rates were observed for chloramphenicol, imipenem, vancomycin, and linezolid among susceptible organism groups, indicating retained activity of selected reserve antibiotics.

The observed resistance profile demonstrates widespread antimicrobial non-susceptibility among vaginal bacterial pathogens, particularly against β-lactams, macrolides, and fluoroquinolones. Conversely, glycopeptides and selected reserve antibiotics retained comparatively better activity primarily against Gram-positive organisms, particularly S. aureus. Interpretation of antibiotics such as vancomycin and linezolid should therefore be considered organism-specific because these agents are not routinely applicable to Gram-negative pathogens.

Overall multidrug resistance burden

A remarkably high prevalence of MDR was observed among the recovered isolates. Out of the total bacterial population, 209 isolates (95.0%) fulfilled MDR criteria, while only 11 isolates (5.0%) were categorized as non-MDR. The observed MDR prevalence was substantially higher than several previously published reports involving vaginal bacterial pathogens, indicating extensive selective pressure within the study environment.

The MDR score distribution further demonstrated that most isolates were resistant to multiple antibiotic classes simultaneously. The mean MDR score was 5.19, while the median MDR score remained 5, indicating a consistently elevated resistance burden throughout the bacterial population (Table 4 & Figure 4).

Display full size

Figure 4. Distribution of multidrug resistance scores among clinical isolates. Dashed vertical lines indicate mean and median MDR score values.

Organism-specific multidrug resistance burden

The MDR burden varied considerably among bacterial species. S. aureus demonstrated the highest average MDR score, followed by P. aeruginosa and Citrobacter spp. In contrast, E. coli showed a comparatively lower, though still clinically significant, MDR burden.

The elevated MDR score among S. aureus isolates supports the high prevalence of MRSA observed through cefoxitin resistance. Similarly, the increased resistance burden among non-fermenting Gram-negative organisms indicates strong adaptive resistance mechanisms and extensive antimicrobial exposure (Table 3 & Figure 5).

Display full size

Figure 5. Mean multidrug resistance burden among bacterial organisms isolated from HVS samples.

Inferential statistical analysis of MDR associations

Chi-square analysis demonstrated a statistically significant association between bacterial organism type and MDR status (χ² = 12.84, p = 0.012), indicating substantial variation in resistance burden among different bacterial taxa. S. aureus isolates showed significantly higher MDR frequencies compared with several Gram-negative organisms.

Binary logistic regression analysis identified S. aureus isolation as an independent predictor of MDR occurrence (OR = 2.41, 95% CI: 1.18–4.92, p = 0.016). In contrast, patient age and pregnancy status were not statistically significant predictors of MDR burden within the final regression model (p > 0.05).

Hierarchical clustering of antimicrobial resistance profiles

Hierarchical clustering analysis demonstrated substantial heterogeneity among AMR phenotypes of the recovered isolates. Four major resistance clusters were identified based on Euclidean distance metrics, indicating the presence of distinct resistance signatures among vaginal bacterial pathogens.

Certain isolate groups formed compact clusters characterized by short linkage distances, suggesting closely related resistance profiles and potentially shared resistance mechanisms. The clustering architecture additionally demonstrated that MDR isolates were widely distributed across different bacterial taxa, emphasizing the complexity of resistance dissemination within the study population (Figure 6).

Display full size

Figure 6. Hierarchical clustering dendrogram showing relationships among bacterial isolates according to antimicrobial resistance characteristics.

Heatmap visualization of antimicrobial resistance patterns

The AMR heatmap revealed pronounced organism-specific resistance signatures and clear co-resistance patterns among several antibiotics. Strong resistance intensities were consistently associated with cefoxitin, macrolides, ceftazidime, penicillin, and fluoroquinolones, whereas lower resistance frequencies were observed for vancomycin, imipenem, chloramphenicol, and amikacin (Figure 7).

Display full size

Figure 7. Heatmap visualization of antimicrobial resistance profiles among bacterial isolates recovered from HVS samples.

Distinct clustered resistance blocks identified in the heatmap suggest coordinated MDR behavior among subsets of isolates. These patterns likely reflect selective antimicrobial pressure and possible accumulation of multiple resistance-associated phenotypes within vaginal bacterial populations.

Principal component analysis (PCA) of resistance profiles

PCA demonstrated distinct segregation of bacterial organisms according to their AMR characteristics. The first principal component explained 35.5% of total variance, while the second component accounted for 7.1% of overall variability.

S. aureus isolates formed a relatively distinct cluster separated from Gram-negative organisms, reflecting their unique resistance characteristics dominated by β-lactam and macrolide resistance. In contrast, E. coli, Klebsiella spp., Citrobacter spp., and P. aeruginosa demonstrated partial overlap, suggesting similarities in MDR mechanisms among Gram-negative bacteria.

The PCA loading vectors identified cefoxitin, clarithromycin, erythromycin, penicillin, levofloxacin, and ceftazidime as major contributors to variability within the dataset, indicating that these antibiotics were the principal drivers of resistance differentiation among isolates (Figure 8).

Display full size

Figure 8. Principal component analysis (PCA) biplot illustrating clustering of bacterial isolates according to antimicrobial resistance profiles and antibiotic loading vectors.