The Presence of Four Pathogenic Oral Bacterial Species in Six Wild Snake Species from Southern Taiwan: Associated Factors

The oral cavity of snakes serves as a habitat for various microorganisms, some of which may include potential zoonotic pathogens posing risks to hosts and causing wound infections in snakebite victims. Clinical studies on snakebite cases in Taiwan have identified specific pathogens, such as Enterococcus faecalis (Gram-positive), Morganella morganii, Aeromonas hydrophila, and Pseudomonas aeruginosa (Gram-negative). However, the prevalence of these bacteria in the oral cavity of wild snakes remains largely unknown. This study investigated the occurrence of these bacteria in six wild snake species (Naja atra, Bungarus multicinctus, Trimeresurus stejnegeri, Protobothrops mucrosquamatus, Boiga kraepelini, and Elaphe taeniura friesi) from southern Taiwan, along with factors influencing their presence. Oropharyngeal swab samples were collected from a substantial number of wild-caught snakes (n = 1104), followed by DNA extraction, polymerase chain reaction, and gel electrophoresis. The band positions of samples were compared with positive and negative controls to determine the presence of target bacteria in each sample. The overall occurrence rates were 67.4% for E. faecalis, 31.5% for M. morganii, 8.2% for A. hydrophila, and 7.7% for P. aeruginosa. Among snake species, B. kraepelini exhibited dominance in E. faecalis (93.4%), A. hydrophila (17.1%), and P. aeruginosa (14.5%), while male N. atra showed dominance in M. morganii (51.3%). The occurrence of E. faecalis was lowest in winter. The results of multiple logistic regression analyses suggest that factors such as species, sex, temperature, season, and coexisting pathogens may have a significant impact on the occurrence of target bacteria. These findings have implications for wildlife medicine and snakebite management.

Snake oral saliva contains a diverse bacterial community, including harmless environmental strains and potential pathogenic species harmful to humans [10,15,16].Following a snakebite, the venom components may lead to necrosis of affected muscles and skin Microorganisms 2024, 12, 263 2 of 18 tissues, creating an environment conducive to the growth of bacteria carried by snake oral saliva [1,8,9,17,18].Without proper attention and treatment, even bites from non-venomous snakes can result in bacterial invasion, leading to bacterial infections and the development of various chronic diseases, potentially resulting in death in severe cases [1,[19][20][21].Snake antivenom therapy may not prevent the subsequent symptoms caused by bacterial infections.Although antibiotics are commonly used to treat bacterial infections, administering them without understanding the specific infection source may be ineffective and contribute to antibiotic resistance, considering the varying sensitivity of bacterial species to antibiotics [1,10,22].
Targeted approaches offer precision and quantification for known target microbes, while non-targeted approaches provide a more comprehensive and unbiased exploration of the microbial landscape.Researchers should use a combination of both approaches to gain a deeper understanding of microbiological systems.In our investigation, we initially utilized a non-targeted approach to analyze the oral microbiota of seven species of wild venomous snakes in Taiwan, employing high-throughput amplicon sequencing of the full-length 16S rRNA gene [16].Subsequently, we will conduct a targeted approach focusing on specific oral microbes using a molecular genotyping method.The bacterial species associated with secondary infections potentially resulting from snakebites are diverse and include aerobic Gram-positive bacteria, Gram-negative bacteria, and anaerobic bacteria (see Supplementary Table S1).Infections caused by Gram-positive bacteria can manifest with symptoms such as abscess formation, cellulitis, sepsis, meningitis, urinary tract infections, bacterial infections at the bite site, and ocular infections.Enterococcus faecalis, a common pathogenic Gram-positive bacterium found in the gastrointestinal tract, is frequently implicated.Most Gram-negative bacteria are pathogenic, and Morganella morganii is a major contributor to abscess formation and gastrointestinal bacterial infections.It is also frequently detected as a multi-drug-resistant bacterium carried by snakes in their oral cavity.Pseudomonas aeruginosa can lead to soft tissue infections in patients, including necrotizing cellulitis, folliculitis, dermatitis, urethritis, and osteomyelitis, and can induce immunosuppressive reactions in the snakebite site, potentially resulting in nosocomial infections.Infections caused by Aeromonas hydrophila can result in hemorrhagic sepsis and diseases related to diarrhea.Other Gram-negative bacteria can cause various infection symptoms such as gastroenteritis, sepsis, respiratory tract infections, meningitis, diarrhea, fever, and soft tissue infections [8,9,15,[23][24][25][26][27][28].The abundance of the aforementioned four species ranks high among the pathogens identified in our previous study [16].
The bacterial species carried by the oral flora of snakes can impact not only humans but also infect the host snakes themselves [21,26,29].The oral cavity of reptiles may sustain wounds from struggling prey during feeding, tooth loss during hunting, or the ingestion of foreign substances such as artificial waste, strings, or fishing nets.These injuries can lead to swelling and inflammatory symptoms [30][31][32].Inadequate control of environmental temperature and humidity in snake captivity can induce immune suppression in snakes, making them susceptible to the invasion of pathogenic bacteria through the respiratory and gastrointestinal tracts, resulting in systemic infections or even death [33][34][35].Aerobic Gram-negative bacteria are generally considered the most pathogenic bacterial species for reptiles and are the dominant taxa detected in individuals with oral inflammation.Common genera include Pseudomonas, Aeromonas, Proteus, and Escherichia [31,33,36].Among them, Pseudomonas aeruginosa, Providencia rettgeri, and Stenotrophomonas maltophilia are often dominant in the oral cavity of snakes with oral inflammation.Aeromonas hydrophila has historically exhibited the highest occurrence rate and is a significant zoonotic pathogen causing infections and death in snakes, commonly distributed in the oral cavity or respiratory tract of snakes [26,33,35,37,38].To mitigate antibiotic resistance and enhance infection recovery rates, a comprehensive understanding of the carried bacterial species is essential for modifying subsequent treatment approaches and determining the appropriate types of antibiotics to administer [39].
Taiwan is home to six main venomous snake species: Trimeresurus stejnegeri, Protobothrops mucrosquamatus, Naja atra, Bungarus multicinctus, Daboia siamensis, and Deinagkistrodon acutus [46,47].Annually, approximately 800 to 1000 people are bitten by venomous snakes in Taiwan [48,49].Research on these venomous snakes often involves examining the clinical symptoms and wound tissue samples of snakebite patients to identify pathogenic bacterial species.Many previous studies focusing on snakebite patients in Taiwan have primarily used medical records from large hospitals in the northern or central regions for retrospective analysis [23,46,50], with only a few focusing on the southern region [49].In wound samples from snakebite patients in Taiwan, Chen et al. [23] and Huang et al. [46] found that M. morganii was the most prevalent among Gram-negative bacteria, followed by A. hydrophila; Enterococcus strains were the most common among Gram-positive bacteria.Chung et al. [49] analyzed the bacterial composition of six major venomous snake oral samples by the culture method, revealing that the most abundant species among Gram-positive bacteria was E. faecalis, while P. aeruginosa and Proteus vulgaris were the most common among Gram-negative bacteria.Further studies, without cultivating bacteria from oral samples, are necessary to comprehend the occurrence of pathogenic bacteria in the oral cavity of Taiwanese snakes.Additionally, there is still limited exploration of the associations between bacterial species and factors such as snake species, region, season, and the body condition of snakes.Moreover, there is a lack of research on the oral bacterial species in non-venomous snakes in Taiwan, which also harbor potential pathogenic bacteria.
This study aims to investigate the occurrence of pathogenic zoonotic bacterial species in the oral cavities of commonly encountered wild snake species in Kaohsiung, southern Taiwan.The targeted bacterial species include E. faecalis, a Gram-positive bacterium, as well as M. morganii, P. aeruginosa, and A. hydrophila, which are Gram-negative bacteria.Polymerase chain reaction (PCR) technology is employed to detect the presence of these four pathogenic microorganisms.Multiple logistic regression analyses are utilized to examine the association between the occurrence of the four bacterial species and various factors, including snake species, sex, snout-vent length, body mass, season, duration between capture and sampling, temperature and rainfall, altitude, artificial ground area around the capture site, and coexisting bacterial species.The findings of this study can serve as a foundation for antibiotic sensitivity testing of bacterial species and facilitate the development of veterinary medicine for wild snakes in Taiwan.This information can, therefore, be valuable in understanding the treatment approaches for secondary infections in humans suffering from snakebites or bacterial infections, whether in captive or wild snakes, stemming from oral microbes.

Snake Collection
From September 2021 (early autumn) to November 2022 (late autumn), we conducted the collection of six common wild snake species in various administrative districts of Kaohsiung City.The collection comprised four venomous species (N.atra, B. multicinctus, T. stejnegeri, and P. mucrosquamatus), one mildly venomous species (Boiga kraepelini), and one non-venomous species (Elaphe taeniura friesi) (Figure 1).The majority of captured snakes were transported to the animal house at the National Pingtung University of Science and Technology for sample collection and morphological measurements within one to two weeks after capture.Before sampling, continuous access to water was provided, with no administration of food or medication [9].
snakes were transported to the animal house at the National Pingtung University of Science and Technology for sample collection and morphological measurements within one to two weeks after capture.Before sampling, continuous access to water was provided, with no administration of food or medication [9].

Oropharyngeal Sample Collection
Oropharyngeal swab samples were collected from each snake using two cotton swabs [9].Sampling was avoided under the following conditions: excessive blood due to injury in the oral cavity, severe ulceration and purulent discharge from old wounds, inability to ensure the cotton swab's uncontaminated status during sampling, extremely small body size, excessive parasitism in the oral region, poor health condition, or imminent death of the snake [9].A disposable sterile tongue depressor (with a hole in the middle) was used as a mouth opener during sampling (see Figure 1h in [16]) to prevent cross-contamination.After opening the mouth, sterile cotton swabs were employed to sample from the oropharyngeal area.The swabs were rotated to ensure even adherence of saliva to the cotton [9].Following sampling, the sample swabs were placed in a sterilized 1.5 mL microcentrifuge tube and stored at −80 °C.

Oropharyngeal Sample Collection
Oropharyngeal swab samples were collected from each snake using two cotton swabs [9].Sampling was avoided under the following conditions: excessive blood due to injury in the oral cavity, severe ulceration and purulent discharge from old wounds, inability to ensure the cotton swab's uncontaminated status during sampling, extremely small body size, excessive parasitism in the oral region, poor health condition, or imminent death of the snake [9].A disposable sterile tongue depressor (with a hole in the middle) was used as a mouth opener during sampling (see Figure 1h in [16]) to prevent crosscontamination.After opening the mouth, sterile cotton swabs were employed to sample from the oropharyngeal area.The swabs were rotated to ensure even adherence of saliva to the cotton [9].Following sampling, the sample swabs were placed in a sterilized 1.5 mL microcentrifuge tube and stored at −80 • C.

E. faecalis
The forward primer sequence was 5 ′ -GTT TAT GCC GCA TGG CAT AAG AG-3 ′ (located at base position 167-189 of the E. faecalis 16S rRNA, GenBank accession no.NR_113902), and the reverse primer sequence was 5 ′ -CCG TCA GGG GAC GTT CAG-3 ′ (located at base position 459-476 of the E. faecalis 16S rRNA, GenBank accession no.NR_113902).The PCR protocol consisted of an initial denaturation step at 95 • C for 2 min, followed by 26 cycles of denaturation at 95 • C for 30 s, annealing at 60 • C for 1 min, and extension at 72 • C for 1 min.The final extension step was performed at 72 • C for 2 min.The expected product size after gel electrophoresis was 310 base pairs (Figure 2A) [52].
volume of the PCR reaction was 25 µL, with the reagent composition including 15 µL sterile distilled water, 2.5 µL of 10× Taq PCR buffer containing 15 mM MgCl (Cyrusbioscience, New Taipei, Taiwan), 2 µL of 10 µM forward and reverse primers (Purigo Biotech, Taipei, Taiwan) targeting the 16S rRNA gene segment, 0.5 µL of 10 mM dNTPs (ABclonal, Woburn, MA, USA), 2.5 U of Taq DNA Polymerase (Cyrusbioscience, New Taipei, Taiwan), and 2.5 µL of DNA sample.The primer selection and temperature-time configuration of PCR specific to each bacterial species are detailed as follows:

E. faecalis
The forward primer sequence was 5′-GTT TAT GCC GCA TGG CAT AAG AG-3′ (located at base position 167-189 of the E. faecalis 16S rRNA, GenBank accession no.NR_113902), and the reverse primer sequence was 5′-CCG TCA GGG GAC GTT CAG-3′ (located at base position 459-476 of the E. faecalis 16S rRNA, GenBank accession no.NR_113902).The PCR protocol consisted of an initial denaturation step at 95 °C for 2 min, followed by 26 cycles of denaturation at 95 °C for 30 s, annealing at 60 °C for 1 min, and extension at 72 °C for 1 min.The final extension step was performed at 72 °C for 2 min.The expected product size after gel electrophoresis was 310 base pairs (Figure 2A) [52].

M. morganii
The forward primer sequence was 5′-TTT CAG TCG GGA GGA AGG TG-3′ (located at base position 357-376 of the M. morganii 16S rRNA, GenBank accession no.NR_043751), and the reverse primer sequence was 5′-GGG GAT TTC ACA TCT GAC TC-3′ (located at base position 515-534 of the M. morganii 16S rRNA, GenBank accession no.NR_043751).The PCR protocol consisted of an initial denaturation step at 95 °C for 2 min, followed by 30 cycles of denaturation at 94 °C for 1 min, annealing at 65 °C for 45 s, and extension at 72 °C for 2 min.The final extension step was performed at 72 °C for 7 min.The expected product size after gel electrophoresis was 178 base pairs (Figure 2B) [53].

M. morganii
The forward primer sequence was 5 ′ -TTT CAG TCG GGA GGA AGG TG-3 ′ (located at base position 357-376 of the M. morganii 16S rRNA, GenBank accession no.NR_043751), and the reverse primer sequence was 5 ′ -GGG GAT TTC ACA TCT GAC TC-3 ′ (located at base position 515-534 of the M. morganii 16S rRNA, GenBank accession no.NR_043751).The PCR protocol consisted of an initial denaturation step at 95 • C for 2 min, followed by 30 cycles of denaturation at 94 • C for 1 min, annealing at 65 • C for 45 s, and extension at 72 • C for 2 min.The final extension step was performed at 72 • C for 7 min.The expected product size after gel electrophoresis was 178 base pairs (Figure 2B) [53].

P. aeruginosa
The forward primer sequence was 5 ′ -GGG GGA TCT TCG GAC CTC A-3 ′ (located at base position 185-203 of the P. aeruginosa 16S rRNA, GenBank accession no.NR_117678), and the reverse primer sequence was 5 ′ -TCC TTA GAG TGC CCA CCC G-3 ′ (located at base position 1122-1140 of the P. aeruginosa 16S rRNA, GenBank accession no.NR_117678).The PCR protocol included an initial denaturation step at 95 • C for 2 min, followed by 25 cycles of denaturation at 94 • C for 20 s, annealing at 58 • C for 20 s, and extension at 72 • C for 40 s.The final extension step was performed at 72 • C for 1 min.The expected product size after gel electrophoresis was 956 base pairs (Figure 2C) [54].

A. hydrophila
The forward primer sequence was 5 ′ -GAA AGG TTG ATG CCT AAT ACG TA-3 ′ (located at base position 445-467 of the A. hydrophila 16S rRNA, GenBank accession no.NR_074841), and the reverse primer sequence was 5 ′ -CGT GCT GGC AAC AAA GGA CAG-3 ′ (located at base position 1110-1130 of the A. hydrophila 16S rRNA, GenBank accession no.NR_074841).The PCR protocol included an initial denaturation step at 95 • C for 5 min, followed by 30 cycles of denaturation at 94 • C for 2 min, annealing at 56 • C for 2 min, and extension at 72 • C for 2 min.The final extension step was performed at 72 • C for 10 min.The expected product size after gel electrophoresis was 685 base pairs (Figure 2D) [55].

Positive and Negative Controls
Four strains were included in a set of positive and negative control groups during the PCR process.The positive control consisted of a DNA sample from a specific strain of bacteria that had been cultured and purified.After gel electrophoresis and sequencing, the bacterial sequences were verified using the NCBI BLAST (E.faecalis NBRC 100481, M. morganii DSM 14850, P. aeruginosa DSM 50071, and A. hydrophila ATCC 7966).The negative control was sterile distilled water to confirm the absence of contamination [52,55].

Agarose Gel Electrophoresis
The PCR products were separated using a 1.5% agarose gel.For each sample, 4 µL of the PCR product was mixed with 1 µL of 5× prestained loading dye.Electrophoresis was conducted in 1% TAE buffer at 110 volts for 30 min.A 100 bp-prestained DNA ladder was used for size comparison [54].Following electrophoresis, the gel was exposed to UVA light, and images were captured.The DNA ladder and the positive control served as reference points to determine whether bands corresponding to the PCR products of the target bacterial strains appeared at the expected positions and sizes.This method was employed to assess the presence of target bacterial DNA fragments (Figure 2).Additionally, we sequenced the PCR products of select samples to validate the accuracy of species identification by referring to NCBI BLAST.

Statistical Methods
Contingency table analyses, including the X 2 test and Fisher's exact test, were performed using PAST (version 4.12) to examine the significance of associations between variables.In cases where a significant association (p < 0.05) was identified, post-hoc analyses were conducted using the Bonferroni correction to calculate an adjusted p value for each cell in the contingency table.This adjustment aimed to assess the significance of each variable's contribution to deviating from the expected values under the null hypothesis of no association.
Furthermore, multiple logistic regression was employed to investigate the association between bacterial occurrence and morphometric or environmental factors in the sampled snakes [56].The categorical factors include snake species, sex, capture season, and coexistence of other bacterial species; the numerical factors include snout-vent length, scaled mass index, artificial ground area within a 100 m radius of the capture site, average temperature on the capture day, average rainfall for the capture month, elevation of the capture site, and the duration between capture and sampling.The scaled mass index was calculated as follows [57,58]: The b SMA is the slope for the bivariate linear relationship between log 10 (body mass) and log 10 (snout-vent length) following the standardized major axis regression analysis, which was performed using RStudio (version 2022.07.2+576).L 0 represents the mean snout-vent length of all snakes, L i is the snout-vent length of an individual snake, and M i is the body mass of an individual snake.The capture season is differentiated by months, with spring encompassing March to May, summer from June to August, autumn from September to November, and winter from December to February [59].Temperature and rainfall data were sourced from the nearest automatic weather station to the capture site, provided by the Central Weather Bureau's Observation Data Inquire Service (CODiS; https://codis.cwa.gov.tw/,accessed on 15 January, 2023).Land use types at capture sites were obtained from the National Land Surveying and Mapping Center (Taiwan MAP service; https://maps.nlsc.gov.tw/,accessed on 15 January, 2023), including paddy fields, dry fields, orchards, rivers, lakes, sea areas, and forests as non-artificial ground, with all others categorized as artificial ground.Considering the daily movement distances of snakes [60], the artificial ground area within a 100 m radius around the capture point was measured using ArcGIS (version 10.8.1).Multiple logistic regression model analyses were conducted using JMP (version 10.0.0) and IBM SPSS Statistics (version 29.0.0.0), evaluating the impact of all morphometric and environmental factors on the occurrence of target bacterial species and establishing the full model.A reduced model was formed by including only factors with significant impact (p < 0.05).If the AICc (corrected Akaike information criterion) value of the reduced model was smaller than that of the full model, the reduced model was considered the most parsimonious regression model.Examination of the Exp(B) values (i.e., odds ratio) for each factor was conducted, where values greater than 1 indicated a positive association with bacterial occurrence and values less than 1 indicated a negative association.Additionally, the significance of interactions between factors was also tested, and no significant interactions were found.

Comparison of Bacterial
Occurrence by Season, Snake Species, and Sex 3.2.1.E. faecalis E. faecalis was detected in a total of 744 snakes; the overall occurrence rates were 67.4%.There was a significant difference in the occurrence rates of E. faecalis among seasons (X 2 = 55.6, d.f.= 4, p < 0.000001; Figure 3); the lowest occurrence rate occurred in winter 2021 (44.4%), while the highest was observed in autumn 2022 (80.4%).Significant differences were also found in the occurrence rates among snake species (X 2 = 142, d.f.= 5, p < 0.000001; Figure 4).B. kraepelini (93.4%),E. t. friesi (84.5%), and B. multicinctus (83.0%) exhibited significantly higher occurrence rates compared to P. mucrosquamatus (50.0%) and T. stejnegeri (42.6%).There was no significant difference in the overall occurrence rates between all males (69.7%) and all females (65.1%) (Fisher's exact test; p = 0.109; Figure 5); no significant differences were observed in the occurrence rates between males and females within each snake species (Fisher's exact test; all adjusted Ps > 0.05).between all males (69.7%) and all females (65.1%) (Fisher's exact test; p = 0.109; Figure 5); no significant differences were observed in the occurrence rates between males and females within each snake species (Fisher's exact test; all adjusted Ps > 0.05).no significant differences were observed in the occurrence rates between males and females within each snake species (Fisher's exact test; all adjusted Ps > 0.05).

A. hydrophila
A. hydrophila was detected in a total of 90 snakes; the overall occurrence rates were 8.2%.There was no significant difference in the occurrence rates of A. hydrophila among seasons (X 2 = 2.81, d.f.= 4, p = 0.590; Figure 3).However, there was a significant difference in the occurrence rates among snake species (X 2 = 29.9,d.f.= 5, p < 0.00001; Figure 4), with B. kraepelini (17.1%) and B. multicinctus (16.4%) exhibiting significantly higher occurrence rates.The overall occurrence rate in males (8.7%) was not significantly different from that in females (7.6%) (Fisher's exact test; p = 0.583; Figure 5); no significant differences were observed in the occurrence rates between males and females within each snake species (Fisher's exact test; all adjusted Ps > 0.05).

Association Factors of Bacterial Occurrence
The full model of multiple logistic regression incorporated both categorical and numerical factors.Table 1 displays the range of values for the numerical factors.a Refer to the text for explanations.b Artificial ground area within a radius of 100 m from the capture points.c Duration between capture and sampling.

E. faecalis
Factors significantly associated with the occurrence of E. faecalis included snake species, temperature, elevation, season, and coexisting bacterial species (Table 2).Using B. kraepelini as a reference, N. atra, B. multicinctus, T. stejnegeri, and P. mucrosquamatus all exhibited a significantly lower occurrence of E. faecalis (odds ratio < 1).Temperature and elevation showed a significant positive association with bacterial occurrence (odds ratio > 1), but the effect of elevation is actually minor because the B value for elevation is close to 0 (Table 2).In comparison to winter 2021, the occurrence rate in spring 2022 was significantly lower.The presence of the other target bacterial species demonstrated a significant positive association with the occurrence rate of E. faecalis.

M. morganii
Factors significantly associated with the occurrence of M. morganii included snake species, sex, temperature, season, and coexisting bacterial species (Table 3).Using B. kraepelini as a reference, N. atra and T. stejnegeri exhibited a significantly higher occurrence of M. morganii.Compared to males, the occurrence rate in females was significantly lower.Temperature showed a significant negative association with occurrence.The presence of the other target bacterial species demonstrated a significant positive association with the occurrence rate of M. morganii.Table 2. Association statistics for the multiple logistic regression between E. faecalis occurrence and various factors.The table includes the estimated coefficients (B), standard errors (S.E.), Wald statistics (Wald), degrees of freedom (d.f.), probabilities under null hypothesis (P), odds ratios (Exp(B)), and 95% confidence intervals (95% C.I.) of Exp(B).Only significant factors (p < 0.05) are presented.The reference category for subclass factors is indicated within parentheses following the specific factor name.Factors significantly associated with the occurrence of P. aeruginosa included only seasonal and coexisting bacterial species (Table 4).Compared to winter 2021, the occurrence rate in summer 2022 was significantly lower.The presence of the other target bacterial species demonstrated a significant positive association with the occurrence rate of P. aeruginosa.Factors significantly associated with the occurrence of A. hydrophila included only snake species and coexisting bacterial species (Table 5).Using B. kraepelini as a reference, T. stejnegeri, N. atra, and E. t. friesi showed a significantly lower occurrence of A. hydrophila.The presence of the other target bacterial species demonstrated a significant positive association with the occurrence rate of A. hydrophila.

Discussion
The factors contributing to bacterial growth in a snake's oral cavity, aside from venom, may encompass the balance of oral microbiota, the type of prey, environmental factors (humidity, temperature, and substrate composition), the immune system, oral hygiene, and the presence of injury or infection [8][9][10]14,21,[24][25][26]40].This study demonstrates significant associations between the occurrence rates of the four bacterial species and various factors, as summarized in Table 6.Notably, certain factors like sex, snout-vent length, scaled mass index, artificial ground area, rainfall, and the duration between capture and sampling did not exhibit a significant association with bacterial occurrences.In contrast, the presence of other bacterial species in the oral cavity consistently exhibited a positive association with the occurrence rates of the target bacteria, given that sampling of snakes with compromised health conditions or immune activity was avoided in this study.This finding suggests that pathogenic bacterial species may share growth environments with each other, as seen in previous research [61].It also implies that, whether dealing with secondary infections in humans following snake bites or bacterial infections in captive or wild snakes due to oral wounds, multiple pathogenic sources may contribute to the infections.Lin et al. [62] and Yeh et al. [63] indicated that M. morganii and E. faecalis were the most common sources of polymicrobial infections.Both E. faecalis and M. morganii can lead to bacteremia, while A. hydrophila and P. aeruginosa may result in urinary tract infections [23,46].The presence of multiple pathogens in infections implies varying degrees of antibiotic resistance, underscoring the importance of thorough evaluation and testing of antibiotic sensitivity in the future.The associated factors for the occurrence rates of bacterial species are discussed as follows: Table 6.A summary of factors significantly associated with the occurrence of four types of bacteria was determined through multiple logistic regression analyses.The reference category for subclass factors is indicated within parentheses following the specific factor name.Positive associations are denoted by ↑ (0.001 < p < 0.05) or ↑↑↑ (p < 0.001), while negative associations are denoted by ↓ (0.001 < p < 0.05) or ↓↓↓ (p < 0.001).For instance, the occurrence of M. morganii in two snake species (N.atra and T. stejnegeri) is significantly higher than in the reference species (B.kraepelini).Some studies suggest that snake venom possesses antibacterial properties [64][65][66][67][68], hypothesizing that the venom can protect venomous snakes from bacterial infections carried by prey during feeding.However, conflicting research supports the notion that venomous snakes display a greater diversity of oral bacteria compared to non-venomous ones [21,69].In the context of this study, pit vipers (P.mucrosquamatus and T. stejnegeri) exhibit a significantly lower occurrence rate of E. faecalis, a bacterium whose growth has been demonstrated to be inhibited by venom from vipers and cobras [70].In contrast, mildly venomous B. kraepelini and venomous B. multicinctus demonstrate a high occurrence rate of E. faecalis, P. aeruginosa, and A. hydrophila.Furthermore, the occurrence rate of M. morganii in the oral cavity of venomous N. atra and T. stejnegeri surpasses that in other snake species, consistent with prior research on patient wounds [71] and the oral cavities of snakes [16].The reasons behind the heightened bacterial occurrence in the above (mildly) venomous species warrant further investigation.Whether the antibacterial properties of venom are related to only specific venom compositions also remains to be confirmed.

Factor
Our study reveals a positive correlation between the occurrence rate of E. faecalis and environmental temperature, aligning with previous findings that E. faecalis thrives in warmer temperatures [72].Since E. faecalis predominantly inhabits the gastrointestinal tract of animals, the oral microbiota of snakes is influenced by the excretion of prey during feeding [73].The increased activity of snakes due to higher environmental temperatures may lead to an elevated frequency of predation [60], indirectly contributing to the higher occurrence rate of E. faecalis.
M. morganii exhibits the highest occurrence rate among Gram-negative bacteria in this study.This observation aligns with findings from both clinical research on patient wounds and snake oral microbiota studies [23,46,50].Bacteriomic analysis reveals that M. morganii ranks among the 10 most abundant bacterial species in oropharyngeal samples of both N. atra and T. stejnegeri [16].M. morganii, known to inhabit the gastrointestinal tract of animals, is also widely distributed in the natural environment, including soil, water sources, and sewage [74,75].The higher occurrence rate of M. morganii in the oral cavity of male N. atra and B. multicinctus may be attributed to the larger activity range seen in males of some snake species [76,77], increasing the likelihood of contact with pathogens.
P. aeruginosa primarily inhabits human feces and sewage, with less occurrence in natural environments such as soil and streams [78].This study exhibited the lowest occurrence rate among bacterial species, contrary to previous research on Taiwanese snake oral bacteria that used culture-based methods, identifying P. aeruginosa as the Gram-negative species with the highest occurrence rate [49].This discrepancy may stem from the limitation that microbial compositions cultured on agar plates may not entirely represent the original oral microbiota, as indicated by several studies [21,[41][42][43][44][45].Studies have shown that P. aeruginosa outbreaks are more extensive in hospital areas than in non-hospital settings [78].This suggests that the P. aeruginosa identified in clinical studies, particularly in patients with wound infections [23], may predominantly originate from nosocomial infections.Additionally, previous research investigating the seasonal association of P. aeruginosa using human samples demonstrated a peak in its occurrence during cooler months (October to March) [79], which aligns with our findings indicating a higher occurrence rate of P. aeruginosa in winter.
This study has several limitations: (1) The investigation did not analyze the impact of dietary factors on the incidence rates of bacteria in the sampled snake species due to their complex and overlapping feeding habits.Subsequent research could concentrate on snake species with simpler and non-overlapping diets to ascertain whether snake feeding habits significantly influence bacterial occurrence rates.(2) Previous studies [31,33,36] have identified Gram-negative bacteria, particularly P. aeruginosa, P. rettgeri, and S. maltophilia, as dominant strains in the oral cavity of snakes with oral inflammation.However, our study did not sample the oral bacterial community of unhealthy or sick snakes.Future studies should include samples from snakes exhibiting infectious symptoms to explore variations in oral microbiomes under different health conditions.(3) The bacterial strains selected for this study were aerobic or facultatively anaerobic.Future research on the oral bacterial communities of wild snakes should also consider investigating other types of pathogenic bacteria, including anaerobic strains such as B. fragilis and Clostridium spp.

Conclusions
Our research, involving extensive sampling of wild snakes in southern Taiwan, provides insights into the complex factors influencing the presence of oral bacteria in snake populations.Specifically, the prevalence of oral cavity bacteria ranked highest for E. faecalis, followed by M. morganii, A. hydrophila, and P. aeruginosa.Among snake species, B. kraepelini and B. multicinctus dominated in E. faecalis, A. hydrophila, and P. aeruginosa, while male N. atra and T. stejnegeri led in M. morgani.Winter exhibited the lowest occurrence of E. faecalis.Multiple logistic regression analyses underscored species, sex, temperature, elevation, season, and coexisting bacteria as potentially influential factors.The presence of other bacterial species in the oral cavity consistently showed a positive association with the occurrence rates of the target bacteria.While numerous studies have explored bacterial communities in the wounds of patients bitten by vipers and cobras or within the oral cavities of these snakes [8,23,44,46,80], research on other (mildly) venomous [81] snakes as well as non-venomous [21] snakes are notably scarce.To enhance veterinary therapy for such snakes, as well as snakes under different health conditions, additional studies on the bacterial community harbored within these snakes are deemed necessary.Our findings hold significance for both human clinical treatment and wild snake medicine efforts.

Figure 5 .
Figure 5. Sex-based comparison of the occurrence rates of four bacteria (E.faecalis, M. morganii, P. aeruginosa, and A. hydrophila) in the oropharyngeal samples of snakes.

Figure 5 .
Figure 5. Sex-based comparison of the occurrence rates of four bacteria (E.faecalis, M. morganii, P. aeruginosa, and A. hydrophila) in the oropharyngeal samples of snakes.

Table 1 .
Minimum and maximum values of various numerical factors included in the full model of multiple logistic regression.

Table 3 .
Association statistics for the multiple logistic regression between M. morganii occurrence and various factors.Refer to the legend in Table2for an explanation of statistical symbols.

Table 4 .
Association statistics for the multiple logistic regression between P. aeruginosa occurrence and various factors.Refer to the legend in Table2for an explanation of statistical symbols.

Table 5 .
Association statistics for the multiple logistic regression between A. hydrophila occurrence and various factors.Refer to the legend in Table2for an explanation of statistical symbols.