Simultaneous Analysis of Bacterial and Fungal Communities in Oral Samples from Intubated Patients in Intensive Care Unit

Intubated patients in intensive care units (ICUs) too frequently contract ventilator-associated pneumonia or Candida infections. Oropharyngeal microbes are believed to play an important etiologic role. This study was undertaken to determine whether next-generation sequencing (NGS) can be used to simultaneously analyze bacterial and fungal communities. Buccal samples were collected from intubated ICU patients. Primers targeting the V1-V2 region of bacterial 16S rRNA and the internal transcribed spacer 2 (ITS2) region of fungal 18S rRNA were used. V1-V2, ITS2, or mixed V1-V2/ITS2 primers were used to prepare an NGS library. Bacterial and fungal relative abundances were comparable for V1-V2, ITS2, or mixed V1-V2/ITS2 primers, respectively. A standard microbial community was used to adjust the relative abundances to theoretical abundance, and NGS and RT-PCR-adjusted relative abundances showed a high correlation. Using mixed V1-V2/ITS2 primers, bacterial and fungal abundances were simultaneously determined. The constructed microbiome network revealed novel interkingdom and intrakingdom interactions, and the simultaneous detection of bacterial and fungal communities using mixed V1-V2/ITS2 primers enabled analysis across two kingdoms. This study provides a novel approach to simultaneously determining bacterial and fungal communities using mixed V1-V2/ITS2 primers.


Introduction
Oral bacterial and fungal communities play important roles in various nosocomial infections such as ventilator-associated pneumonia (VAP) [1], hospital-acquired pneumonia (HAP) [2], oral mucositis [3], and systemic infections [4]. In the case of VAP, the microaspiration of oropharyngeal microbes is believed to play a role in the etiology of ventilator-associated pneumonia (VAP) [2], and the prevalence of oral colonization by VAP-associated pathogens such as Pseudomonas aeruginosa, Enterobacteriacea, Acinetobacter, S. aureus, and Streptococcus spp. are markedly elevated in VAP patients [1]. Moreover, fungal load increases rapidly in intubated patients after admission to intensive care units (ICUs) [5] and can account for a substantial proportion of healthcare-associated infections in intubated patients [6,7]. Proper oral hygiene practices are recommended to reduce the risk of HAIs and other hospital-related conditions caused by oral bacteria and fungus [8].
The monitoring of oral bacterial and fungus abundance changes has several merits. First, the data obtained can be used to predict and diagnose infections [9]. For example, changes in the oral microbiome composition can indicate the presence of infections, such as VAP or bloodstream infections [10]. Second, it provides an understanding of the impact of interventions as a comparison of the oral microbiome before and after interventions, such as antibiotics or oral practice, allows clinicians to understand how these interventions impact the microbiome and potentially enable treatment plans to be adjusted accordingly [11,12]. Thus, the efficient monitoring of oral bacterial and fungal changes in ICU patients and other hospital conditions is important.
Next-generation sequencing (NGS) produces large numbers of sequences at unprecedented speeds. Amplicon analysis is commonly used for microbiome studies and has been used in large projects, including the Human Microbiome Project [13]. The 16S rRNA gene is highly conserved among bacteria and is commonly used to differentiate bacteria at taxonomic levels [14], whereas the internal transcribed spacer (ITS) of the nuclear ribosomal DNA is commonly used in fungal studies. The ITS region is located between the 18S and 28S rRNA genes and comprised of the ITS1 and ITS2 regions, the latter of which is located between the 5.8S and 28S rRNA genes and has been suggested to be useful for clinical studies [15].
In this study, we sought to determine whether bacterial and fungal communities can be simultaneously determined using primers targeting the V1-V2 region for bacteria and the ITS2 region for fungi. Libraries were prepared using individual primers (V1-V2 or ITS2) and combined primers (V1-V2/ITS2 mix) for NGS. The results obtained were compared for alpha diversity and relative abundance at the species level. Furthermore, NGS results were compared with real-time PCR results to confirm the relative abundance measured by NGS.

Patients
The study inclusion criteria applied were: (1) an ICU patient, (2) >18 years of age, (3) endotracheal tube use with a holding method, and (4) consent for participation from the patient's family. The intubated patients enrolled had a heterogeneous set of underlying diseases (Table S1). Oral samples were collected at multiple time points to monitor microbiome and mycobiome shifts. Broad-spectrum antibiotics were administered to treat underlying diseases. This study was approved by Institutional Review Board (DAUHIRB-19-072). The ICU patients involved were mostly sedated or had poor levels of consciousness, and thus, family members were considered legal representatives with power of attorney as stipulated by the International Council for Harmonization [16]. During the recruitment period, the study purpose, voluntary nature of participation, confidentiality of information, and procedures used in the study were explained to family members, and each family provided informed consent.

Oral Swab Sample Preparation
Buccal samples were collected from 4 patients during ICU treatment (day 0-day 12) using the Levine technique. A sterile cotton swab was rotated over an area of approximately 1 cm 2 with sufficient pressure [17] and then placed in a conical 15 mL ID-labeled tube and stored at −80 • C.

Real-Time PCR
To quantify bacterial and fungal abundance, real-time PCR was applied using speciesspecific primers for VAP-associated pathogens, including Acinetobacter baumannii, Pseudomonas aeruginosa, S. aureus, S. epidermidis, and C. albicans (Table S2) [3]. The reaction conditions used for DNA amplification were 95 • C for 10 min followed by 40 cycles of 95 • C for 30 s and 60 • C for 30 s. The analyses were conducted using an RT-PCR system (Thermo Fisher Scientific, Waltham, MA, USA). Total bacteria abundances were determined using universal 16S rRNA, and total fungal abundances using 18S rRNA. The relative abundances of targets were calculated using the ∆∆Ct method using 16S rRNA for normalization [3].

Bioinformatic Analysis, Statistical Analysis, and Visualization
Basic microbiome analysis was performed using QIIME2 (version 2020.6) [22] and its associated plugins. DADA2 (Divisive Amplicon Denoising Algorithm 2) was applied to denoise, join the paired reads, and remove the chimera. Choa1 indices, Faith's phylogenetic diversity (PD) indices, and Shannon's indices were used to measure alpha diversities. The Mann-Whitney U test was used to determine the significance of alpha diversities. To assign taxonomy to unique representative sequences, we used pretrained Naive Bayes classifier with the eHOMD [23] and UNITE [21] databases. To visualize the internal interactions and further measurement of the microbial community, Sparse Correlations for Compositional data (SparCC) [24] were used to calculate the correlation coefficient with corresponding p-value between each two species. The network was then visualized by Cytoscape [25], with the nodes denoting the species and connections representing the existence of correlation.

Data Availability
The raw sequencing data were deposited at NCBI GenBank under BioProject ID PRJNA773666 (BioSample SAMN22503998-SAMN22504007).

Taxonomy Comparison
After preprocessing, relative bacterial community abundances obtained using V1-V2 and mixed V1-V2/ITS2 primers were compared. To facilitate the direct comparisons, taxa assigned to the bacteria or fungi kingdom in V1-V2/ITS2 mix data were selected, and relative abundances were calculated for the bacteria and fungi, respectively. Thus, the relative abundance of bacteria from V1-V2/ITS2 was compared with V1-V2 data, and fungal abundance from V1-V2/ITS2 was compared with ITS2 data. When relative abundances were plotted for individual samples depending on the bacterial or fungal kingdom, most taxa showed similar abundances ( Figure 1A,B). However, some discrepancies were observed between the relative abundances of Acinetobacter baumannii, Mycoplasma spp., Pseudomonas aeruginosa, and Streptococcus spp. in Pt [04] samples. Overall, a significant correlation was observed between V1-V2 and mixed V1-V2/ITS2 primers ( Figure 1C).

Taxonomy Comparison
After preprocessing, relative bacterial community abundances obtained using V1-V2 and mixed V1-V2/ITS2 primers were compared. To facilitate the direct comparisons, taxa assigned to the bacteria or fungi kingdom in V1-V2/ITS2 mix data were selected, and relative abundances were calculated for the bacteria and fungi, respectively. Thus, the relative abundance of bacteria from V1-V2/ITS2 was compared with V1-V2 data, and fungal abundance from V1-V2/ITS2 was compared with ITS2 data. When relative abundances were plotted for individual samples depending on the bacterial or fungal kingdom, most taxa showed similar abundances ( Figure 1A,B). However, some discrepancies were observed between the relative abundances of Acinetobacter baumannii, Mycoplasma spp., Pseudomonas aeruginosa, and Streptococcus spp. in Pt [04] samples. Overall, a significant correlation was observed between V1-V2 and mixed V1-V2/ITS2 primers ( Figure 1C).
Next, the relative abundances of fungal communities, as determined using ITS2 primers and mixed V1-V2/ITS2 primers, were compared. When the relative abundances were plotted for individual samples for specific primers, most of the taxa showed similar abundances (Figure 2A,B). The correlation between abundances determined using ITS2 and mixed V1-V2/ITS2 primers was significant ( Figure 2C). Taken together, the results obtained using the V1-V2/ITS2 primers, V1-V2, and ITS2 primers were comparable, indicating that the V1-V2/ITS2 primers simultaneously produced a relative abundance similar to those obtained using individual single primers.  Relative abundances of bacterial species depending on primers: (A) V1-V2, (B) V1-V2/ITS2 primers, and (C) the correlation between abundances determined using V1-V2 and V1-V2/ITS2 primers.
Next, the relative abundances of fungal communities, as determined using ITS2 primers and mixed V1-V2/ITS2 primers, were compared. When the relative abundances were plotted for individual samples for specific primers, most of the taxa showed similar abundances (Figure 2A,B). The correlation between abundances determined using ITS2 and mixed V1-V2/ITS2 primers was significant ( Figure 2C). Taken together, the results obtained using the V1-V2/ITS2 primers, V1-V2, and ITS2 primers were comparable, indicating that the V1-V2/ITS2 primers simultaneously produced a relative abundance similar to those obtained using individual single primers.

Microbial DNA Standard
To test how well V1-V2/ITS2 primers performed in terms of predicting mixed bacterial and fungal abundances, we used a microbial community DNA standard containing a mixture of genomic DNA isolated from pure cultures of eight bacterial and two fungal strains. When relative abundances were plotted versus raw read counts, determined relative abundances of the two fungal strains almost reached 50%, which compared with a theoretical abundance of less than 5%. Thus, read counts of taxa assigned to the fungi kingdom were reduced by a factor of 20. For example, if the initial read count for the fungi kingdom was 20,000 and the read count for the bacteria kingdom was 100,000, the taxa assigned to the fungi kingdom was reduced to 1000. After adjustment, predicted relative abundances closely matched theoretical abundance ( Figure 3A).
Using adjusted relative abundances, we examined the correlation between V1-V2/ITS2 targeting next-generation sequencing (NGS) and real-time PCR (RT-PCR). Bacterial ( Figure 3B) and C. albicans ( Figure 3C) abundances showed significant correlations indicating the adjustment made optimally predicted abundances in clinical samples.

Microbial DNA Standard
To test how well V1-V2/ITS2 primers performed in terms of predicting mixed bacterial and fungal abundances, we used a microbial community DNA standard containing a mixture of genomic DNA isolated from pure cultures of eight bacterial and two fungal strains. When relative abundances were plotted versus raw read counts, determined relative abundances of the two fungal strains almost reached 50%, which compared with a theoretical abundance of less than 5%. Thus, read counts of taxa assigned to the fungi kingdom were reduced by a factor of 20. For example, if the initial read count for the fungi kingdom was 20,000 and the read count for the bacteria kingdom was 100,000, the taxa assigned to the fungi kingdom was reduced to 1000. After adjustment, predicted relative abundances closely matched theoretical abundance ( Figure 3A).

Microbial DNA Standard
To test how well V1-V2/ITS2 primers performed in terms of predicting mixed bacterial and fungal abundances, we used a microbial community DNA standard containing a mixture of genomic DNA isolated from pure cultures of eight bacterial and two fungal strains. When relative abundances were plotted versus raw read counts, determined relative abundances of the two fungal strains almost reached 50%, which compared with a theoretical abundance of less than 5%. Thus, read counts of taxa assigned to the fungi kingdom were reduced by a factor of 20. For example, if the initial read count for the fungi kingdom was 20,000 and the read count for the bacteria kingdom was 100,000, the taxa assigned to the fungi kingdom was reduced to 1000. After adjustment, predicted relative abundances closely matched theoretical abundance ( Figure 3A).
Using adjusted relative abundances, we examined the correlation between V1-V2/ITS2 targeting next-generation sequencing (NGS) and real-time PCR (RT-PCR). Bacterial ( Figure 3B) and C. albicans ( Figure 3C) abundances showed significant correlations indicating the adjustment made optimally predicted abundances in clinical samples.  Using adjusted relative abundances, we examined the correlation between V1-V2/ITS2 targeting next-generation sequencing (NGS) and real-time PCR (RT-PCR). Bacterial ( Figure 3B) and C. albicans ( Figure 3C) abundances showed significant correlations indicating the adjustment made optimally predicted abundances in clinical samples.

Alpha Diversity Depending on Primers
Chao1, Faith PD, and Shannon indices were used to evaluate alpha diversities depending on the primers. Similar Chao1 and Shannon index was observed between V1-V2 and bacterial reads from mixed V1-V2/ITS2 primers and between ITS2 and fungal reads from mixed V1-V2/ITS2 primers ( Figure 4A,B). However, V1-V2/ITS2 primers exhibited significantly higher richness than ITS2 or V1-V2 primers based on Faith PD results, which incorporate phylogenetic differences ( Figure 4C).

Alpha Diversity Depending on Primers
Chao1, Faith PD, and Shannon indices were used to evaluate alpha diversities depending on the primers. Similar Chao1 and Shannon index was observed between V1-V2 and bacterial reads from mixed V1-V2/ITS2 primers and between ITS2 and fungal reads from mixed V1-V2/ITS2 primers ( Figure 4A,B). However, V1-V2/ITS2 primers exhibited significantly higher richness than ITS2 or V1-V2 primers based on Faith PD results, which incorporate phylogenetic differences ( Figure 4C).

Discussion
The human oral microbiome is a most diverse ecosystem, comprising numerous c munities [26]. In this ecosystem, more than 600 bacterial [27] and 100 fungal species compose the microbiome. Microbes residing in the human oral cavity have been w

Discussion
The human oral microbiome is a most diverse ecosystem, comprising numerous communities [26]. In this ecosystem, more than 600 bacterial [27] and 100 fungal species [28] compose the microbiome. Microbes residing in the human oral cavity have been wellreflected in human health and disease [29]. However, the oral fungal microbiome, also known as mycobiome, has not received as much scrutiny despite its ecological and clinical importance [30], as evidenced by its relative lack of attention [28]. As a result, potential interactions between the two kingdoms in the oral microbiome have yet to be thoroughly investigated. To gain a complete understanding of how the human oral microbiota affects health and disease, it is essential to take a comprehensive approach that accounts for both intra-and cross-kingdom interactions within the oral microbiota [31].
The use of high-throughput sequencing for 16S/ITS targeted amplicon high-throughput sequencing and metagenomics has allowed for a closer examination of the bacterial aspect of the human microbiome to better comprehend its function in health and disease [32,33]. The simultaneous detection of bacteria and fungi communities has several advantages over sequencing individual kingdoms. Simultaneous detection can provide relative abundance in the two kingdoms and enable the identification of intra-and interkingdom interactions. In this study, we sought to determine whether the abundances of bacterial and fungal communities can be simultaneously determined using primers targeting the bacterial V1-V2 region and the fungal ITS2 region.
Total read counts obtained using mixed V1-V2/ITS2 primers were significantly higher than those obtained using V1-V2 or ITS2 primers (Table). Since all samples were run simultaneously and equal concentrations of libraries were applied to chips, it is not clear why the V1-V2/ITS2 primer produced significantly higher read counts. Furthermore, the input to nonchimeric read count percentage was significantly higher for the ITS2 primer than the V1-V2 or V1-V2/ITS2 primers. Since all the samples were worked up using the same PCR conditions, PCR cycles are probably not responsible for these differences. On the other hand, primer specificities or target DNA abundances within samples might explain these differences [34].
After preprocessing, the relative abundances of bacterial and fungal communities were compared to determine whether mixing primers influenced results (Figures 1 and 2). When relative abundances obtained using the different primers were plotted for individual samples, most taxa showed similar abundances. However, some discrepancies in the relative abundances of bacteria were observed at the species level in some samples. Nonetheless, when the relative abundances of assigned taxa were compared, high correlations were observed for bacterial and fungal communities. Taken together, mixed V1-V2/ITS2, V1-V2, and ITS2 primers produced comparable results, suggesting that the use of V1-V2/ITS2 primers simultaneously produced relative abundances similar to those obtained using individual single primers.
Next, a microbial DNA standard was used to test the ability of the V1-V2/ITS2 primer mix to simultaneously predict the abundance of bacterial and fungal species (Figure 3). When we plotted relative abundances versus raw read counts, large differences were found between theoretical and predicted abundances. We suggest several possible reasons for this discrepancy. First, the copy numbers of 16S and 18S rRNA for bacteria and fungi are markedly different. For bacteria, the average 16S rRNA copy number is 4.2 [35], while the 18S rRNA copy number for fungi ranges from 38 to 91 [36], and for C. albicans, the copy number of 18S rRNA is 55 [37], which suggests the use of 18S rRNA exaggerated the abundance of fungi. Second, ITS2 primers produced more nonchimeric reads, which would also increase the relative abundance of fungi. Thus, adjustment was required for the fungal strains. When NGS and RT-PCR abundances were plotted using adjusted relative abundance, bacterial and C. albicans abundances showed significant correlations, which indicated the adjustment optimally predicted bacterial and fungal abundances in clinical samples.
Finally, we plotted bacterial and fungal relative abundances obtained using mixed V1-V2/ITS2 primers. When alpha diversity was determined, a similar Chao1 and Shannon index was observed between V1-V2 and bacterial reads from mixed V1-V2/ITS2 primers and between ITS2 and fungal reads from mixed V1-V2/ITS2 primers ( Figure 4A,B). These results support that the mixed V1-V2/ITS2 primer produced results comparable to those produced by individual V1-V2 and ITS2 primers, respectively. Since fungal diversity was significantly lower than bacterial diversity, no difference was observed between V1-V2 and V1-V2/ITS2 primers, suggesting that fungal diversity was completely masked by bacterial diversity. However, V1-V2/ITS2 primers exhibited significantly greater richness than the ITS2 or V1-V2 primers based on Faith PD results, which incorporate phylogenetic differences ( Figure 4C). Thus, Faith's PD seems to save fungal diversity information.
In many studies, the correlations between the abundances of fungal and bacterial species have been analyzed following individual by community analysis [38][39][40]. Furthermore, since information about the relative abundances of bacteria and fungi is limited, it is difficult to merge reported relative abundances directly in a single plot. When we simultaneously plotted bacterial and fungal relative abundances obtained using V1-V2/ITS2 primers for individual samples, the relative abundances of C. albicans and Clavispora lusitaniae were determinable within bacterial populations ( Figure 5A). One of the advantages of using V1-V2/ITS2 primers is that results enable shifts in fungal relative abundances to be better understood than those obtained using individual primers. For example, in a patient [01], an increase in C. tropicalis abundance was observed on days 8 and 9 samples (Figure 2). However, when V1-V2/ITS2 primers were used, the overall relative abundance of fungus was found to decrease on days 8 and 9, indicating the increases in C. tropicalis abundance observed on these days were the result of C. tropicalis surviving during a decrease in total fungal abundance rather an increase in C. tropicalis abundance ( Figure 5). Thus, the simultaneous determination of bacterial and fungal community abundances provides more comprehensive information regarding microbiome changes compared to single analysis.
A microbiome network was constructed to visualize bacterial and fungal interactions and revealed strong inter-and intrakingdom interactions ( Figure 5B). Interestingly, C. albicans interacted positively with Mycoplasma salivarium, Fusobacterium nucleatum subsp. Vincentii, and Desulfomicrobium orale bacteria. There are many studies on interactions between Candida and oral bacteria, such as C. albicans and S. mutans [41,42], which have influenced their synergistic relationship to oral health. It is a very useful study model that could serve as a preventive or therapeutic target for oral diseases. However, the abundance of S. mutans is relatively low in the oral samples of healthy elderly subjects [43], and S. mutans is rarely found in intubated patients [44]. Thus, other oral bacterial species should also be considered in studies on interactions between bacteria and C. albicans. In the present study, C. lustitaniae formed a dense network with several oral bacteria; for example, its abundance was positively correlated with those of Treponema socranskii, Staphylococcus aureus, Saccharibacteria (TM7) spp., and Streptococcus anginosus and negatively correlated with those of Rothia mucilaginosa, Haemophilus parainfluenzae, Lautropia mirabilis, Streptococcus sanguinis, and Neisseria oralis. In a previous study, Rothia, Haemophilus, Lautropia, and Neisseria were more frequently found in healthy individuals, whereas Treponema, Staphylococcus, and Saccharibacteria (TM7) were more frequent under disease conditions, such as periodontitis [45]. The biological relevance of these correlations has not been determined, and thus, studies are required to determine the nature of the relationships between fungal and bacterial species, as they may reveal interactions that enable modulation of the relative abundances of targeted species. This study is limited by a relatively small sample size caused by the requirement for family consent for the collection of oral samples from ICU patients. Additionally, only elderly patients were included, which could have introduced bias due to the limited number of participants. Therefore, we recommend a larger, well-designed study be undertaken to validate our findings in the future.
In the bacterial and fungal microbiome analysis, data regarding technical issues related to primer development are still lacking. In this study, we demonstrated that analyzing bacterial and fungal communities simultaneously in human oral samples can reduce analysis costs and enable the identification of infectious causes across two kingdoms. By performing network analysis, we suggested that bacterial-fungal interactions are present in the oral microbiome. Moreover, we identified C. lustitaniae and C. albicans, among a few other taxa, as potential keystone taxa of the oral microbiome of ICU patients. Overall, our study has facilitated an understanding of the determining factors and cross-kingdom interactions of the human oral microbiome. These findings will lead to the development of more effective approaches for diagnosing complex clinical infections.