Critically ill patients with COVID-19 show lung fungal dysbiosis with reduced microbial diversity in patients colonized with Candida spp.

BACKGROUND: The COVID-19 pandemic has intensified interest in how the infection affects the lung microbiome of critically ill patients and how it contributes to acute respiratory distress syndrome (ARDS). We aimed to characterize the lower respiratory tract mycobiome of critically ill patients with COVID-19 in comparison to patients without COVID-19.
METHODS: We performed an internal transcribed spacer 2 (ITS2) profiling with the Illumina MiSeq platform on 26 respiratory specimens from patients with COVID-19 as well as from 26 patients with non-COVID-19 pneumonia.
RESULTS: Patients with COVID-19 were more likely to be colonized with Candida spp. ARDS was associated with lung dysbiosis characterized by a shift to Candida species colonization and a decrease of fungal diversity. We also observed higher bacterial phylogenetic distance among taxa in colonized patients with COVID-19. In patients with COVID-19 not colonized with Candida spp., ITS2 amplicon sequencing revealed an increase of Ascomycota unassigned spp. and 1 Aspergillus spp.-positive specimen. In addition, we found that corticosteroid therapy was frequently associated with positive Galactomannan cell wall component of Aspergillus spp. among patients with COVID-19.
CONCLUSION: Our study underpins that ARDS in patients with COVID-19 is associated with lung dysbiosis and that an increased density of Ascomycota unassigned spp. is present in patients not colonized with Candida spp.

Introduction

Coronaviruses are important human and animal pathogens. At the end of 2019, a novel coronavirus was identified as the causative agent of a pneumonia outbreak in Wuhan, China, that subsequently spread worldwide in a global pandemic. Globally, there have been more than 240 million confirmed cases of COVID-19, including nearly 5 million deaths (https://www.who.int/). Critically ill patients with COVID-19 may develop acute respiratory distress syndrome (ARDS), requiring admission in intensive care unit (ICU) and mechanical ventilation, which predisposes them to bacterial and fungal superinfections (Bassetti et al., 2020; Lansbury et al., 2020; Chong et al., 2021). Candida is one of the most frequently isolated pathogens in ICUs, affecting between 6% and 10% of patients (Zhang et al., 2020; Koehler et al., 2019). The estimated mortality rate attributed to invasive candidiasis ranges from 19% to 40%; this mortality is even higher among ICU patients, approaching 70% (Kullberg and Arendrup, 2015; Kullberg et al., 2015; Marra et al., 2011). Recently, it has been reported that fungi are more frequently detected among patients with SARS-CoV-2, with Candida albicans being the most frequently isolated yeast (Calderaro et al., 2021). Furthermore, the wide use of antibiotics and corticosteroids along with the damage exerted by SARS-CoV-2 may allow commensal yeast to invade internal organs, causing deeply invasive infections (Arastehfar et al., 2020; Talento and Hoenigl, 2020; Talento et al., 2020; Posteraro et al., 2020).

Despite the low number of studies on lung mycobiome, growing evidence indicates that the fungal microbiota is altered in critically ill patients (Krause et al., 2016). Fungi found in the human respiratory tract are predominantly from the Dikarya subkingdom, which is composed of the phyla Ascomycota and Basidiomycota. In healthy individuals, the fungal burden is generally low, and the mycobiome appears to be largely composed of environmental fungi or fungi disseminating from the oral cavity (Dickson and Huffnagle, 2015). By contrast, more stable fungal communities can colonize the lung when its physiology is altered. As an example, in most patients with cystic fibrosis, the fungal burden is increased, whereas alpha diversity is reduced and correlates with disease severity (Iliev and Leonardi, 2017). Despite the speculated importance of lung dysbiosis in the genesis of both ventilator-associated pneumonia (VAP) and ARDS, few studies have examined the lung mycobiome in these patient populations (Krause et al., 2016; Dickson et al., 2015). In the lungs of patients with COVID-19, fungal colonization/infection represents a major concern although the clinical significance is debated (Peng et al., 2021). Indeed, Candida albicans has been reported as the most frequently isolated yeast from the lung, whereas COVID-19–associated pulmonary aspergillosis (CAPA) has been reported in a few centers (Bartoletti et al., 2020; Permpalung et al., 2021; Borman et al., 2020). Despite this fact, no data concerning the lung mycobiome in patients with COVID-19 with ARDS have been reported. The goal of this study is to analyze the composition of lung mycobiome in mechanically ventilated patients with COVID-19 with ARDS.

Methods

Population

Twenty-six patients with COVID-19 and 26 patients without COVID-19 were enrolled in this study. All the patients were recovered at IRCCS Sant'Orsola Malpighi University Hospital, Bologna, Italy, from March to April 2020. The study was conducted in accordance with the Declaration of Helsinki. Samples were coded, and analysis was performed with an anonymized database. Informed consent for study participation was obtained from each patient. The study was approved by the local institutional review board (Comitato Etico AVEC; approval number n. 283/2020/Oss/AOUBo). Clinical and demographic characteristics are presented in Table 1 .

Table 1

Demographic and clinical characteristics of patients with and without COVID-19.

Variable	N	Patients without COVID-19 (%, IQR)(n = 26)	Patients with COVID-19 (%, IQR)(n = 26)	p valuea
Age, median (IQR)	52	64 (52–71)	68 (65–69)	0.39
Male, n (%)	52	15/26 (57.7)	19/26 (73.1)	0.24
Time for symptoms onset to BAL collection, median (IQR)	46	16.27 (10–21)	17.91 (12–24)	0.51
Time hospitalization to BAL collection, median (IQR)	45	7.78 (2–11)	12.09 (7–15)	0.99
Mechanical ventilation, n (%)	52	6/26 (23.1)	26/26 (100)	≤0.01
ICU admission, n (%)	46	3/26 (12)	26/26 (100)	≤0.01
ICU death,b n (%)	52	1/26 (3.8)	8/26 (30.8)	0.01
SOFA score (media, IQR)c	26	Not available	3 (2–4)
Immunomodulation, n (%)
Toclizumab	25	1/25 (4.0)	14/25 (56.0)	<0.01
Interferon gamma	25	0/25 (0)	5/25 (20.0)	0.05
Corticosteroids	49	21/24 (87.5)	15/25 (60.0)	0.05
Antiviral, n (%)	50	22/25 (88)	11/25 (44.0)	0.002
Antibiotics, n (%)	50	15/25 (60)	17/25 (68.0)	0.77
Voriconazole, n (%)	50	1/25 (4)	3/25 (12.0)	0.61
Candida spp. (culture) colonization in the respiratory tract, n (%)	50	0/26 (0)	6/24 (25.0)	0.01
Candida spp. (culture) colonization,d n (%)	52	4/26 (15.4)	16/26 (61.5)	0.001
Candidemia (culture), n (%)	52	0/26 (0)	2/26 (7.7)	0.35
Aspergillus fumigatus (culture), n (%)	50	0/26 (0)	2/24 (8.3)	0.13
Galactomannan positive, n (%)	31	3/13 (23.1)	8/18 (44.4))	0.27
Aspergillus PCR positive, n (%)	52	1/26 (3.8)	6/26 (7.7)	0.61
Aspergillus NGS positive, n (%)	52	0/26 (0)	1/26 (3.8)	0.26
Bacterial infections in the respiratory tract (culture), n (%)	52	7/26 (27)	23/26 (88)	0.001

Chi-square of Fisher exact text; Mann-Whitney test for continuous data.

Among patients with COVID-19, 3 of 6 patients (50%) with Candida colonization died versus 5 of 18 patients (27.7%) without colonization (p = 0.317)

SOFA score for colonized versus noncolonized patients was not statistically significant (p = 0.537).

Candida spp. colonization in the respiratory tract and other sites (oral, nasal, respiratory, rectal, genital sites).

BAL, bronchoalveolar lavage; ICU, intensive care unit; IQR, interquartile range; NGS, next-generation sequencing; PCR, polymerase chain reaction; SOFA, Sequential Organ Failure Assessment.

Bronchoalveolar lavage collection and microbiological analysis

The presence of SARS-CoV-2 was detected by the reverse transcriptase polymerase chain reaction (RT-PCR) assay. Briefly, detection of SARS-CoV-2 was performed by real-time RT-PCR following the World Health Organization and/or Centers for Disease Control and Prevention protocols in a QuantStudio S5 Real-time PCR system (Thermo Fisher). All patients with COVID-19 were managed in a dedicated COVID-19 ICU and underwent mechanical ventilation (Table 1). The control group included 26 patients admitted to the hospital, presenting with clinical and radiological findings of pneumonia and with PCR performed on nasopharyngeal swab or bronchoalveolar lavage (BAL) negative for SARS-CoV-2. Radiological findings of the control group were consistent with interstitial pneumonia in 9 patients (35%), or patchy multifocal infiltrates in 8 cases (30%), or other findings in the remaining ones (35%). In our COVID-19 population, 8 of 26 cases of CAPA were classified as probable, according to the recent European Confederation of Medical Mycology and the International Society for Human and Animal Mycology (ECMM/ISHAM) consensus criteria (Koehler et al., 2020). Direct microscopy or biopsy, which are the microbiological criteria to classify CAPA as proven, were not performed in our study. Of the 26 patients with COVID-19, 6 were mechanically ventilated and 3 were admitted in the ICU (Table 1).

BALs were cultured for the isolation of fungal and bacterial pathogens at 32°C for 5 days for filamentous fungi and at 37°C for 2 days for yeasts and bacterial pathogens. In particular, selective media such as CHROMagar Candida and Sabouraud-dextrose (SAB) chloramphenicol agar (Vacutest Kima, Italy) as well as nonselective media such as SAB were used for in vitro fungal growth. Selective media such as Agar Herellea, Agar salt-mannitol, agar chocolate haemophilus and nonselective media such as sheep blood agar were used for bacterial growth. The level of galactomannans (GMs), which are fungal antigens, was measured with a sandwich enzyme-linked immunosorbent assay (ELISA; Platelia Aspergillus; Bio-Rad Laboratories) in BAL specimens following manufacturer's instructions. DNA extraction from 1-mL volume of respiratory material (diluted 1:1 with dithiothreitol) and real-time PCR for Aspergillus spp. targeting rDNA 18S (ELITe MGB kit, Elitech Group, Italy) were performed on ELITe InGenius automated platform. DNA was eluted in 100-µL volume, and the DNA copy number was expressed as copies/mL. Negative and positive controls were used in each run as well as standard 10-fold dilutions of Aspergillus DNA.

Next-generation sequencing

The lower respiratory mycobiome was characterized both in critically ill patients with COVID-19 and in noninfected patients. The lower respiratory bacterial microbiome of these samples has been previously analyzed by Gaibani et al. (Gaibani et al., 2021). Total microbial DNA was extracted from BAL samples using the QIAamp 96 PowerFecal QIAcube HT Kit on the QIAcube HT instrument (QIAGEN, Hilden, Germany) following the manufacturer's instructions. A bead‐beating step with Lysing Matrix E (MP Biomedicals) was performed on a FastPrep24 bead-beater (MP Biomedicals, Irvine, CA) at 6.0 movements per second for 40 seconds before total DNA extraction. Negative controls were PCR-grade water that underwent library preparation steps and next-generation sequencing (NGS). DNA was quantified using the Qubit 4 Fluorometer (Fisher Scientific). Internal transcribed spacer 2 (ITS2) was amplified using the primer set ITS3: 5′-GCATCGATGAAGAACGCAGC-3′ and ITS4: 5′-TCCTCCGCTTATTGATATGC-3′ (White et al., 1990). PCR products were purified with a magnetic bead–based clean‐up system (Agencourt AMPure XP; Beckman Coulter, Brea, CA). Indexed libraries were prepared by limited‐cycle PCR using Nextera technology and further cleaned up with AMPure XP magnetic beads (Beckman Coulter). Libraries were pooled at equimolar concentrations (4 nM), denatured, and diluted to 5 pM before loading onto the MiSeq flow cell. Sequencing on Illumina MiSeq platform was performed by using a 2 × 250 bp paired-end protocol, according to the manufacturer's instructions (Illumina, San Diego, CA). Sequencing reads were deposited in the National Center for Biotechnology Information Sequence Read Archive (NCBI SRA; BioProject ID PRJNA742164)

Data analysis

Paired-end sequenced reads of samples were analyzed, combining PANDAseq2 and QIIME2 version 2018.6 (Bolyen et al., 2018). The Divisive Amplicon Denoising Algorithm 2 (DADA2) (Hall and Beiko, 2018) plug-in was used to remove noise and chimeras and to generate amplicon sequence variants (ASVs). Quality filtering and clustering were performed using VSEARCH (Rognes et al., 2016). High‐quality reads were classified taxonomically using the UNITE reference database version 7.2 (UNITE Community, 2017: UNITE QIIME release. Version 01.12.2017. UNITE Community. https://doi.org/10.15156/BIO/587481). Samples that had less than 1000 reads after Illumina MiSeq sequencing were discarded. The bacterial abundance data were imported into R (version 3.6.1) on Rstudio v1.1.456, where all statistical analyses were performed using R package phyloseq (McMurdie and Holmes, 2013). Using the decontam R package at 1% and 5% stringency (Davis et al., 2018) on our negative controls (ultrapure water samples that underwent the whole library preparation processes), we assessed the absence of detectable contaminant ASVs. We obtained 258 fungal taxa after quality filtering and the removal of unidentified fungal phyla. Sample sequencing reads were not rarefied to avoid introducing unwanted bias because samples reached the total ASV number asymptote at around 800 reads even when having higher sequencing read depths (Supplementary Figure S1). The differences in alpha diversity were evaluated based on the data distribution of metrics, using analysis of variance (ANOVA) and Tukey honestly significant difference tests for normally distributed data or Wilcoxon–Mann-Whitney with Holm-Bonferroni correction method for non-normally distributed data. To check that the nonapplication of rarefaction would have not led us to inexact conclusions, alpha diversity analysis was repeated after rarefying to the lowest sequencing depth, and results matched the ones obtained without rarefaction (Supplementary Figure S2). To compare microbial composition between samples, beta diversity was measured by calculating the Bray-Curtis distance matrix. Principal coordinates analysis (PCoA) was applied on the distance matrices to generate bidimensional plots in R. Dispersion of the PCoA clusters was compared using the betadisper function in R vegan package (Anderson and Walsh, 2013) . The permutational ANOVA (PERMANOVA) test (Anderson and Walsh, 2013), calculated using the function adonis in the vegan package (Oksanen et al., 2014), was performed to determine whether there was a significant separation between different sample groups. The plots were graphed using ggplot2 R packages (Wickham, 2016). Dissimilarity percentage (SIMPER) analysis function (https://github.com/asteinberger9/seq_scripts) based on R packages vegan and dplyr was used to determine the contribution of individual taxa driving the average dissimilarities between groups. A p value <0.05 after false discovery rate (FDR) correction was considered as statistically significant. Linear discriminant analysis (LDA) effect size (LEfSE) algorithm (Segata et al., 2011), a tool which is hosted on the Galaxy web application at https://huttenhower.sph.harvard.edu/galaxy/, was also used to discover bacterial biomarkers associated with patients with COVID-19. The differences in abundance were regarded as significant when the logarithmic LDA score was higher than 2.

Statistical analysis

Given that the microbial diversity of the mycobiome has not been well described in patients with COVID-19 and there is no established change in diversity genotypes that has been previously reported to have clinical effects in this population, an empirical sample size of 26 patients per group was planned based on empirical assessments from previous lung microbiome studies. For descriptive analysis, categorical variables are presented as counts and percentages, continuous variables as means and standard deviations if normally distributed, or as medians and interquartile ranges (IQRs) if non-normally distributed. For group comparisons, the Student t test and Mann-Whitney U test were used for normally distributed quantitative variables or skewed distributions, respectively. The Fisher exact test was used for categorical variables. Shapiro-Wilk and Kolmogorov-Smirnov tests, as well as visual inspection, were used to test data for normality assumptions. Given the number of study subjects per group (n = 26), we opted not to perform multivariable analysis to control for latency age and immunosuppressive treatments as inclusion of more than a couple of variables would likely lead to model overfitting and inaccurate odds estimates with wide confidence intervals (CIs) and poor generalizability of findings.

Results

Clinical characteristics

This study was performed in Italy during the COVID-19 major pandemic wave of March–April 2020 when ICUs were overwhelmed with patients with COVID-19, and few contemporary patients without COVID-19 were available for analysis in health care systems. Thus, in that period, precise patient enrollment based on clinical and demographic features was difficult to plan, especially for control population. Limitations of this study include variation in age, sex, and treatment of the patients, and its strength relies on having the patient samples from a single hospital for more uniform data collection.

Clinical and demographic characteristics of the patients without COVID-19 and of the ones with COVID-19 are presented in Table 1. The median age was 64 (IQR 52–78) for the patients without COVID-19 and 68 (IQR 65–69) for the ones with COVID-19 (p = 0.05). Time from symptom onset to BAL collection as well as from hospitalization to BAL collection was not significantly higher for the patients with COVID-19 than those without COVID-19 (p = 0.52, p = 0.59, respectively). On the contrary, significant differences among the 2 groups of patients occurred for mechanical ventilation and ICU admission. Indeed, only 6 patients without COVID-19 were mechanically ventilated (p ≤ 0.01) (Table 1). Rates of antibacterial and antifungal use were similar between the 2 groups. However, a higher percentage of patients with COVID-19 received immunomodulating therapies including tocilizumab (71% vs 10%, p < 0.01), interferon gamma (20% vs 0%, p = 0.05), and corticosteroids (87.5% vs 60.0%, p = 0.05). Mortality in patients with COVID-19 was significantly higher than in those without COVID-19 (p = 0.01). Moreover, among patients with COVID-19, mortality rate was not significantly different for colonized versus noncolonized patients (p = 0.317).

As presented in Table 1, among the 26 patients with COVID-19, 6 BAL specimens were positive for Candida spp. (3 C albicans, 2 C glabrata, and 1 C tropicalis) and 2 BALs for Aspergillus fumigatus by culture. Patients with COVID-19 were more likely to have positive BAL GMs (8/18, 44%) with an index ≥ 1.00 than those without COVID-19 (3/23, 13%). Of those without COVID-19, none was positive for Candida and Aspergillus spp. by culture; however, 3 patients were positive for BAL GMs. Patients with COVID-19 were more likely to be colonized with Candida spp. (3/23, 13%; p = 0.010). Interestingly, the Candida spp. colonization in different body sites (eg, pharyngeal, rectal, and genital tracts), including the respiratory tract, was higher for patients with COVID-19 than for those without COVID-19 (p = 0.001). Furthermore, 2 patients with COVID-19 developed candidemia, but those without COVID-19 did not develop invasive infections (p = 0.35). Bacterial superinfections were significantly higher in patients with COVID-19 than in those without COVID-19 (p = 0.001); a detailed description of bacterial pathogens is reported in Supplementary Table 1 (Table S1).

Candida colonization is accompanied by decreased alpha diversity of fungi in patients with COVID-19

Comparisons on fungal microbial composition using ITS2 amplicon sequencing in BAL samples from patients with COVID-19 and those without COVID-19 were performed, but no significant differences were found between these groups. We then focused on individuals with COVID-19, in whom we detected microbial composition differences between the patients who were colonized by Candida spp. versus those who were not colonized. The 6 individuals with COVID-19 who were colonized by Candida spp. were sex-, age-, and treatment-matched with the rest of the patients with COVID-19 (Table 2 ). The same analysis could not be applied on the patients without COVID-19 because they were not colonized by Candida spp.

Table 2

Demographic characteristics of patients with COVID-19

Variable	Candida colonization (n = 6)	No Candida colonization (n = 18)	p valuea
Average age (± SD)	67.3% (± 5.1)	65.2% (± 9)	0.55
Sex (Female), n (%)	2 (33.3)	5 (27.8)	1
Antibiotics, n (%)	3 (50)	5 (27.8)	0.37
Steroids, n (%)	3 (50)	11 (61.1)	0.67
Immunomodulants, n (%)	3 (50)	15 (83.3)	0.14
COVID-19 infection, n (%)	6 (100)	18 (100)	1

Fisher exact test; Mann-Whitney test for continuous data.

ITS2 amplicon sequencing analysis revealed a lower fungal diversity (richness and evenness of distribution of fungal taxa: observed species index p value FDR-corrected = 0.0224; Shannon index p value FDR-corrected = 0.0004; inverse Simpson index p value FDR-corrected = 0.006) in patients with COVID-19 with Candida spp. colonization with respect to uncolonized ones (Figure 1 ). A widespread ecological simplification of the fungal microbiome in the BAL samples from patients who were colonized by Candida spp. can be appreciated when looking at the fungal relative abundance (%) distribution in each specimen (Supplementary Figure S3). Moreover, patients with COVID-19 who were colonized by Candida spp. showed a different clusterization of their fungal pulmonary microbiota (PERMANOVA based on Bray-Curtis dissimilarities = 0.002) (Figure 2 ) and a significantly higher relative abundance of Candida spp. (SIMPER, p value FDR-corrected = 0.008) with respect to noncolonized patients. In contrast, patients who were not colonized exhibited increased density of Ascomycota unclassified spp. (LDA LEfSE > 2.00) and a negative association with Candida spp. (Figure 3 ). Patients with COVID-19 colonized by Candida spp. showed a higher bacterial phylogenetic distance among taxa (p value FDR-corrected = 0.049) (Figure 4 ). However, no bacterial genera were specifically associated with Candida colonization.

Figure 1

Alpha diversity indexes of fungal microbiomes. Violin plots showing the comparison of alpha diversity measures between patients who were positive for SARS-CoV-2 (n = 24) with or without Candida spp. colonization (n = 6 and n = 18, respectively). Median, first, third quartile, p values with false discovery rate correction, and outliers are shown.

Figure 2

Principal coordinate analysis (PCoA) on Bray-Curtis distance metric at the OTU level calculated on patients with COVID-19 (n = 24) with Candida spp. colonization (n = 6, blue dots) and without Candida spp. colonization (n = 18, red dots). Each sample is represented by a dot. Axis 1 explained 21.1% of the variation, and axis 2 explained 12.4% of the variation observed.

Figure 3

Linear discriminant analysis (LDA) effect size (LEfSE) analysis. The plot was generated using the online Galaxy web platform tools at https://huttenhower.sph.harvard.edu/galaxy/. The length of the bar column represents the LDA score. The figure shows the microbial taxa with significant differences between the patients with COVID-19 colonized by Candida spp. (green bar) and not colonized by this microorganism (red bar) (LDA score > 2). NEG, negative; POS, positive.

Figure 4

Alpha diversity phylogenetic diversity (PD) whole tree index of bacterial microbiomes. Violin plots showing the comparison between patients who were positive for SARS-CoV-2 (n = 22) with or without Candida spp. colonization (n = 6 and n = 16, red box and blue box, respectively) calculated on their bacterial microbiome. Median, first, third quartile, p value with false discovery rate correction, and outliers are shown. Two bronchoalveolar lavage samples did not have sufficient 16S rRNA sequencing reads and could not be analyzed.

Increase of Ascomycota unclassified spp. in patients with COVID-19 uncolonized by Candida spp

Among patients with COVID-19, the LDA LEfSE analysis found a significant (LDA > 2) increase in Ascomycota unclassified spp. in patients uncolonized by Candida spp. Moreover, Candida spp. was not identified by culture in patients without COVID-19. Among the 31 patients tested for GM in the BAL fluid (18 patients with COVID-19 and 13 patients without COVID-19), Aspergillus spp. was detected by NGS in 1 patient with COVID-19 (Supplementary Figure S3) (sample #7), but in none of the patients without COVID-19 (Figure 5 ). Furthermore, the NGS extracts were further analyzed by a target-specific PCR assay able to detect different species of Aspergillus genus. In this way, we confirmed the positivity of case #7 (Figure 5) and we found 1 positive case among patients without COVID-19 (#1).

Figure. 5

Alluvial graph representing specimens positive for Aspergillus spp. as determined by next-generation sequencing (NGS) and polymerase chain reaction (PCR) both in patients with and without COVID-19.

BAL, bronchoalveolar lavage; NGS BAL, next-generation sequencing on NGS dedicated extracted specimens; PCR BAL, polymerase chain reaction on NGS dedicated extracted specimens; follow-up BAL PCR, polymerase chain reaction on PCR dedicated extracted specimens collected from 0 to 6 days after those for NGS.

Considering that a main concern in patients with COVID-19 with ARDS was pulmonary aspergillosis (CAPA), we further analyzed BAL specimens collected during routine diagnostic procedures for the presence of Aspergillus spp. by PCR. In total, we found 6 positive cases, including case #7 (Figure 5). These specimens were collected from 0 to 6 days after those used for NGS analysis. Interestingly, among the 6 patients with COVID-19, 2 were colonized by C glabrata (Figure 5).

Discussion

The main finding of this study was that COVID-19 infection was associated with lung dysbiosis characterized by a shift to Candida species colonization and a decrease in fungal diversity. We also observed higher bacterial phylogenetic distance among taxa in patients who were colonized.

To date, few studies have evaluated the lung mycobiota using high-throughput sequencing (Soltani et al., 2021). In healthy individuals, the fungal burden is generally low and lung mycobiota appears to be largely composed of environmental fungi such as Cladosporium, Aspergillus, Penicillium and yeasts belonging to the 2 main phyla Ascomycota and Basidiomycota, whereas more stable fungal communities colonize the lung when its physiology is altered (Krüger et al., 2019). As an example, alpha diversity is reduced and fungal burden is increased in patients with cystic fibrosis (Linnane et al., 2021). Once established, the dysbiotic fungal communities seem to persist even in the presence of antibiotic or immunosuppressant therapy (Iliev et al., 2017).

Similar disruptions in lung homeostasis have been reported in patients with bacterial VAP (Fernández-Barat et al., 2020). Krause et al. observed a decrease of fungal diversity in a ventilated patient with pneumonia and an increase of Candida spp., representing 75% of fungal species (Krause et al., 2016). Furthermore, it has been recently reported that fungi are more frequently detected among mechanically ventilated patients with SARS-CoV-2 than in nonventilated ones, with C albicans being the most frequently isolated microorganism (Calderaro et al., 2021; Palabiyikoğlu et al., 2001). In our population of mechanically ventilated patients with COVID-19, we did observe a decrease in fungal diversity and an increase in lung Candida spp. colonization. Interestingly, patients with COVID-19 showed a significantly higher Candida spp. colonization in body sites different from the respiratory tract than in the control population. Furthermore, candidemia was present in 2 patients with COVID-19 and was absent in the control population. The mortality rate in patients with COVID-19 was higher than in those without COVID-19; however, among patients with COVID-19, no significant difference was present in patients colonized by Candida spp. versus those who were not. Similarly, no significant difference was present in the Sequential Organ Failure Assessment score of patients who were colonized versus those who were not. Accordingly, Terraneo et al. (Terraneo et al., 2016) reported that Candida spp. airway colonization in patients with ICU-associated pneumonia does not influence outcomes in these patients. On the contrary, we did not find Candida spp. colonization in the COVID-19 control population, which is probably due to the low rate of mechanically ventilated patients in our case series (23.1%). A possible contribution of Candida colonization to the pathogenesis of bacterial pneumonia has been hypothesized on the basis of observational clinical data (Fromentin et al., 2021), although it is unclear if there is a direct effect or if colonization is a signal of impaired host defenses. Direct bacteria-fungi interactions may play a role in the pathogenesis of VAP, although the precise mechanisms remain to be defined. For example, overgrowth of Candida spp. in the respiratory tract may impair host bacterial phagocytosis (Oliver et al., 2019). In our COVID-19 population, broad-spectrum antibiotic therapy was more commonly used in patients with Candida colonization and positive Aspergillus GMs antigen, but the CI of the odds ratio estimate is relatively large because of missing data. However, in the analysis by Krause et al. (Krause et al., 2016), antibiotic therapy was not reported to influence mycobiota composition. Similarly, we did not find a specific association with Candida colonization of the respiratory tract and concurrent antibiotic therapy in patients with COVID-19. However, patients with Candida colonization did demonstrate higher bacterial phylogenetic distance among taxa, which implies a higher ecological differentiation among bacterial metabolic functions. In our study, bacterial superinfections, as determined by culture, were significantly higher in patients with COVID-19 versus those without COVID-19, while no significant difference was present for the administration of broad-spectrum antibacterial therapy. Finally, we cannot exclude that bacterial superinfections could have influenced the microbial niche and consequently the comparison of those with COVID-19 and those without COVID-19 in terms of Candida colonization. No significant difference was found in patients with COVID-19 versus those without COVID-19 in the time from symptom onset to BAL collection, as well as from hospitalization to BAL collections. So, the patients with and without COVID-19 were homogeneous in their clinical characteristics except for mechanical ventilation and ICU admission. We cannot exclude that this might have had some impact on our findings when comparing those with and without COVID-19.

CAPA must be considered a serious and potentially life-threatening complication in patients with severe COVID-19 receiving immunosuppressive treatment (Machado et al., 2021). In fact, many studies reported an increase in CAPA in intubated patients with COVID-19–related ARDS (Bartoletti et al., 2020; Dewi et al., 2021). In this study, we found that corticosteroid therapy was frequently associated with subsequent positive GM tests. ITS2 amplicon sequencing revealed an increase of Ascomycota unclassified spp. in patients with COVID-19 who were not colonized with Candida spp. Nevertheless, only 1 patient was positive for Aspergillus spp. Considering Aspergillus spp. PCR results, on PCR dedicated DNA extracts collected from 1 to 6 days later than those for NGS, we detected an additional 6 positive cases in patients with COVID-19 that were not positive by NGS.

In conclusion, our study demonstrates that lung fungal dysbiosis is more severe in critically ill patients with COVID-19. In particular, Candida spp. colonization was accompanied by a decreased diversity (richness and evenness of fungal taxa) of fungi in the respiratory tract overall, while in patients not colonized by Candida, the mycobiome was characterized by a higher density of unclassified fungi from the Ascomycota phylum. CAPA is a major concern for patients with COVID-19, and by NGS, we found 1 case positive for Aspergillus spp. among patients with COVID-19 with ARDS. By analyzing BAL collected later than BAL that underwent ITS2 amplicon sequencing, PCR assays confirmed this case together with 6 other positive cases. In this regard, it is possible that a fluctuation of Aspergillus DNA copy number during CAPA evolution as well as the effect of voriconazole therapy in enhancing the presence of free DNA, and the presence of Candida colonization in 2 PCR-positive cases may have influenced the different rate of positivity between PCR and NGS methods. Our study is unique from previous lung microbiome investigations in that we analyzed the fungal community of the lower respiratory tract of critically ill patients with COVID-19 and we compared these patients with contemporary ones hospitalized at the same hospital with non–COVID-19 pneumonia. Limitations of this study reflect the reality of performing studies during the major pandemic wave of March–April 2020 in Italy, when ICUs were overwhelmed with patients with COVID-19 and few contemporary patients without COVID-19 were available for analysis. Hence, enrollment was suboptimal in terms of the number of patients analyzed and variation in clinical and demographic settings. Further studies are needed to confirm our results on this interesting topic.

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest. No animal studies are presented in this manuscript.

Funding source

The study was conducted independently of study sponsors. Ricerca Fondamentale Orientata (RFO) 2020 of Dr. Cricca Monica, Alma Mater Studiorum, University of Bologna.

Author contributions

EV contributed to the study design; carried out the experiments, data analysis, and writing; discussed the results; and commented on the manuscript. PG contributed to the study design, discussed the results, and commented on the manuscript. AC contributed to the study design, discussed the results, and commented on the manuscript. AL carried out the experiment, discussed the results, and commented on the manuscript. MB carried out data collection, discussed the results, and commented on the manuscript. PV discussed the results and commented on the manuscript. TL discussed the results and commented on the manuscript. SA discussed the results and commented on the manuscript. RL carried out data analysis, discussed the results, and commented on the manuscript. MC contributed to data collection, study design, data analysis, and writing; discussed the results; and commented on the manuscript.

Ethical approval statement

The study was approved by the local institutional review board (Comitato Etico AVEC; approval number n. 283/2020/Oss/AOUBo). The participants provided written informed consent to participate in this study. No potentially identifiable human images or data are presented in this study.