The ICC patients and healthy controls were prospectively recruited at the Renmin Hospital of Wuhan University from June 2021 to February 2022. The Helsinki Declaration and the Rules of Good Clinical Practice were followed during the conduct of this study. It was approved by our hospital’s clinical research ethics committee (WDRY2021-KS029) and registered on the Chinese Clinical Trial Registry Platform. A formal informed consent form was signed by each participant. Through the use of medical records and in-person interviews, the demographic information and clinical traits of patients were acquired. The ICC patients were diagnosed by histopathological examination of specimens from surgical resection or percutaneous ultrasound-guided liver needle core biopsy. Mixed-type liver cancer (hepatocellular carcinoma and ICC) was excluded. ICC patients had not received any routine anti-tumor therapy before the collection of feces samples. The individuals who came to our hospital for a medical examination were chosen as the healthy controls. Participants with other malignancies, digestive disorders, autoimmune diseases, and acute or chronic infectious diseases were excluded. To exclude hepatolithiasis, fibropolycystic disease, primary sclerosing cholangitis, viral hepatitis, and parasitic infections from affecting the intestinal fungus of ICC patients, we excluded ICC patients with these factors (although the number of these patients was small). In addition, none of the included subjects reported the use of probiotics, prebiotics, antibiotics, antifungals, or laxatives within the previous 3 months.

A fecal sample was obtained by the participants, transported to the lab within 2 h, and stored at-80°C until DNA was extracted. The CTAB method was used to extract DNA, and 1% agarose gels were used to measure the DNA concentration. The hypervariable ITS2 region of the ITS gene was amplified using the primers [ITS3-2024F (5’-GCATCGATGAAGAACGCAGC-3′) and ITS4-2409R (5′- TCCTCCGCTTATTGATATGC-3′)] with the barcode. Mix the same volume of 1 x loading buffer (containing SYB green) with the PCR products and operate electrophoresis on 2% agarose gels for detection. The PCR products in the mixture were then purified using Thermo Scientific’s GeneJET Gel Extraction Kit. Using the TruSeq^® DNA PCR-Free Sample Preparation Kit (Illumina, United States), sequencing libraries were made according to the manufacturer’s instructions, and index codes were added. To evaluate the library’s quality, Thermo Scientific’s Qubit@ 2.0 Fluorometer was used. Finally, the library was sequenced, and 250 bp paired-end reads were produced on an Illumina Novaseq 6,000 PE250 platform.

The NCBI Sequence Read Archive database contained the raw sequencing information for all samples (PRJNA898483). Paired-end reads were assigned to samples according to their unique barcode and truncated by cutting off the primer sequence and barcode. Then, FLASH (V1.2.7)¹ was utilized to merge paired-end reads (Magoč and Salzberg, 2011). Quality filtering on the raw tags was conducted to obtain the high-quality clean tags (Bokulich et al., 2013) based on the QIIME (V1.9.1)² quality-controlled process (Caporaso et al., 2010). Clean Tags were compared with the Unite Database³ to detect chimera sequences, and then the chimera sequences were removed. Using the Uparse software (Uparse v7.0.1001;⁴ Edgar, 2013) to cluster all Effective Tags, the sequences with ≥97% similarity were assigned to the same OTUs (Operational Taxonomic Units). The sequence with the highest frequency of occurrence in OTUs was chosen as the representative sequence of OTUs. For each representative sequence, the Unite Database (Kõljalg et al., 2013) was used based on the blast algorithm to annotate taxonomic information (Altschul et al., 1990).

A Venn diagram was produced via the R package “VennDiagram” (Version 2.15.3) to show the shared and unique OTUs between the two groups. Using R software (Version 2.15.3), the rank abundance curves were presented. Using QIIME (Version 1.9.1), the alpha diversity (Shannon, Simpson, Chao1, and ACE) and the beta diversity on weighted unifrac distances were determined. The non-metric multidimensional scaling (NMDS) analysis was conducted using the vegan package of R software (Kruskal, 1964), and the principal co-ordinates analysis (PCoA) analysis was completed using the WGCNA, stats, and ggplot2 packages of R software (Minchin, 1987). The R package vegan’s MRPP function was used for multiple response permutation procedure (MRPP) analysis (O'Reilly and Mielke, 1980). The linear discriminant analysis effect size (LEfSe) analysis was carried out using the LEfSe software (Version 1.0) with a default setting of 4 for the linear discriminant analysis (LDA) score screening (Segata et al., 2011). Spearman analysis of the top 20 species in abundance at the genus level with clinical parameters of ICC patients. To further explore the impact of mycobiota community change, functional prediction with FUNGuild annotation tools was deployed (A NHN, B ZS, B STB, 2016).

The categorical variables were compared using Fisher’s exact test. The Wilcoxon rank-sum test was used to compare continuous variables between the two groups. A two-sided p-value <0.05 indicated a significant difference.

After a rigorous pathological diagnosis and exclusion process, 40 stool samples from 23 ICC patients and 17 healthy controls were collected and analyzed using ITS2 rDNA sequencing. Benefiting from a uniform sample collection protocol, all fecal samples were yellow and soft. The clinical characteristics of the participants, such as age, gender, and body mass index (BMI), were similar between the two groups (Table 1). Serum levels of albumin and platelets were significantly reduced in ICC patients compared with healthy controls (Table 1), while the remaining serum parameters were not significantly different (Table 1).

First, we explored the intestinal fungal diversity of the ICC patients and healthy volunteers. The rank abundance curves showed a good richness and uniformity of species in both groups (Figure 1A). The Venn diagram displayed that 995 OUTs were shared between the two groups. However, 2043 and 1,628 OTUs were unique to the control and ICC groups, respectively (Figure 1B). The community richness (Chao1 and ACE) and community diversity (Simpson and Shannon) indicators demonstrated a significantly decreasing tendency of mycobiota diversity in the ICC patients compared to healthy controls (ACE: p < 0.001, Chao1: p = 0.001, Simpson: p < 0.001, Shannon: p < 0.001; Figures 1C–F).

To display mycobiota space between the two groups, beta diversity analysis was conducted through PCoA and NMDS according to the weighted_unifrac matrix. The results indicated that the distribution of fecal fungal communities was significantly segregated clustering between the two populations (Figures 1G,H). Anosim and MRPP analyses further confirmed that inter-group differences were greater than intra-group differences and that there were significant differences in fungal community structure between the two groups (Anosim: p = 0.001, MRPP: p = 0.001; Supplementary Tables S1, S2). Therefore, these results suggest the presence of altered intestinal fungal composition in patients with ICC.

The stacked bar plot presented the distribution of the predominant fungi (top 20 most abundant taxa) at the phylum, genus, and species levels in participants from the two groups (Figures 2A–C). Next, the Wilcoxon rank-sum test was used to explore differences. The results showed that there were 6, 22, 56, 122, 190, and 249 significantly different taxa at the phylum, class, order, family, genus, and species levels between ICC patients and healthy controls, respectively (Supplementary Tables S3–S8). Compared to healthy controls, ICC patients had significantly higher relative abundances of Ascomycota and significantly lower relative abundances of Mucoromycota and Basidiomycota (Figure 2D; Supplementary Table S3). Besides, the 10 most abundant of the above differential taxonomic units at genus and species levels were presented in Figures 2E,F, finding a significant enrichment of Monographella, M. nivalis, Candida, C. albicans, and depletion of Saccharomyces, S. cerevisiae, Pichia, P. mandshurica, Mucor, M. circinelloides, Staphylotrichum, S. coccosporum, Actinomucor, A. elegans, Alternaria, A. alternata, Fusarium, F. oxysporum, Humicola, and H. fuscoatra in the ICC group.

The LEfSe with the criteria of an LDA score ≥ 4.0 was further used to identify differences in the fecal mycobiome between the ICC group and healthy controls. The predominant fungi and differential taxa between the two groups were displayed in the structure of the phylogenetic tree (Figure 3A). As shown in Figure 3B, the intestinal mycobiota of ICC patients were enriched with Ascomycota, Saccharomycetes, Saccharomycetales fam Incertae sedis, Candida, C. albicans, Xylariales, Hyponectriaceae, Monographella, and M. nivalis. For healthy controls, the results revealed significant up-regulation of the abundance of 36 taxa, such as Saccharomyces, S. cerevisiae, Actinomucor, A. elegans, Staphylotrichum, S. coccosporum, Mucor, M. circinelloides, Alternaria, A. alternata, Fusarium, F. oxysporum, Pichia, P. mandshurica, Bipolaris, and Bipolaris sp. These differentially abundant taxa can be considered potential biomarkers.

To further explore the correlation between fungal disorders and tumor progression, we further analyzed the differences in fungal taxa between 12 ICC patients with TNM stage III-IV and 11 ICC patients with TNM stage I-II (8th edition of the UICC/AJCC Cancer staging system). The baseline clinical data of the two groups is shown in Supplementary Table S9. We found that, compared to ICC patients with TNM stages I-II, the fungal flora of ICC patients with TNM stages III-IV exhibited significant segregated clustering, indicating that the gut fungal profile changes significantly with cancer progression (Figures 4A,B). The LEfSe method revealed that some fungal taxonomic units, such as Candida, C. albicans, Candida sp., Dipodascus, D. geotrichum, Ustilaginaceae sp., Clavulinaceae sp., and Bipolaris, were significantly enriched in the ICC patients with stage III-IV (Figures 4C,D); however, Sordariomycetes, Xylariales, Hyponectriaceae, Monographella, Monographella nivalis, Annulohypoxylon, and A. multiforme significantly increased in the ICC patients with stage I-II. Thus, these results suggest that the opportunistic pathogen C. albicans may be involved in the progression of ICC.

We also analyzed the correlation between gut fungal microbes and clinical indicators. We found a significant positive correlation between Dipodascus and carcinoembryonic antigen (CEA), while Saccharomyces and Rhizopus showed a significant negative correlation with alpha-fetoprotein (AFP; Supplementary Figure S2). Wickerhamomyces was found to be significantly negatively correlated with alanine aminotransferase (ALT), and Blumeria and Saccharomycopsis were found to be significantly positively correlated with albumin (Supplementary Figure S2).

Because of the lack of a powerful tool for annotating the function of fungi, we concentrated on the functional guilds of the fungal microorganisms instead, using FUNGuild. As is shown in Figure 5A, the stacked bar plot presented the distribution of the top 10 most abundant trophic modes among participants from the two groups. The most abundant mode in ICC patients was pathotroph, while the most common in healthy control patients was saprotroph (Figure 5A). Then, a Wilcox rank-sum test was implemented on the above modes. The results showed that pathotroph was significantly enriched in ICC patients, while saprotroph, pathotroph-saprotroph, pathotroph-saprotroph-symbiotroph, pathotroph-symbiotroph, saprotroph-symbiotroph, symbiotroph were significantly upregulated in healthy controls (Figure 5B). Thus, our results show that the symbiotic ecological associations of fecal fungi are changed in ICC patients and are dominated by pathological parasitism, which can receive nutrients from and adversely affect host cells.