1. Introduction

The importance of the intestinal microbiome for host homeostasis as well as the pathogenesis and treatment of a wide range of inflammatory, immune, chronic, and malignant diseases has become increasingly evident (Honda and Littman, 2016; Rooks and Garrett, 2016). So far, however, most studies have focused on the role of bacteria; another class of commensals, namely the mycobiota, might also contribute substantially, especially conditionally pathogenic fungi such as Candida albicans. A previous study revealed that the human anti-fungal Th17 immunity and pathology rely on cross-reactivity against C. albicans (Bacher et al., 2019), and other studies also showed that the host’s defense against Klebsiella pneumoniae (Kp) is up to Th17 cells and IL-17R signaling (Amezcua Vesely et al., 2019; Nakamoto et al., 2019). Both C. albicans and K. pneumoniae infections are known risk factors for ankylosing spondylitis (AS) (Li et al., 2019). Therefore, we assume that in AS, mycobiota may also synergistically deteriorate the Th17 pathology. Ankylosing spondylitis is the prototype of inflammatory rheumatic diseases grouped under the term spondyloarthritis (SpA) and represents the end phenotype of the axial SpA group (Braun and Sieper, 2007; Sieper and Poddubnyy, 2017). The pathogenesis of AS is not completely understood but likely involves a complex interplay between genetic predisposition involving the human leukocyte antigen, namely HLA-B27, and environmental factors such as mechanical stress and the microbiome (Braun and Sieper, 2007; Sieper and Poddubnyy, 2017). Up to 10% of AS patients are reported to have inflammatory bowel disease (IBD), and 70% show signs indicative of subclinical intestinal inflammation (Simone et al., 2018; Watad et al., 2018; Fragoulis et al., 2019). From genome-wide association studies, it is evident that over 10% of the gene pathways are shared between IBD and AS (Jostins et al., 2012). A recent study showed a clear link between the severity of subclinical gut inflammation and MRI-determined sacroiliac joint involvement (Van Praet et al., 2014). It has long been suggested that intrinsic barrier dysfunction permits non-specific innate immune activation with systemic translocation of bacterial adjuvants (McGonagle et al., 2001). All these factors indicate that an ecological imbalance of microbiota is a dominant risk factor in AS. In AS, the tight junctions between intestinal epithelial cells are prone to increased permeability, leading to a leaky gut (Reid et al., 2011; Koning et al., 2015). Several molecules increased in AS may point toward the idea of a leaky gut, such as lipopolysaccharide (LPS)-binding protein, fatty acid-binding protein, and zonulin (Ciccia et al., 2016; Scher et al., 2016; Ciccia et al., 2017). In recent times, specific alterations have been proposed in the composition of the gut microbiota, rather than the non-specific innate activation theory, and have been associated with a range of immune-mediated disorders (Scher et al., 2016). It has been postulated that carrying HLA-B27 is a predisposing factor for gut dysbiosis, followed by a leaky gut and the subsequent systemic entrance of microbial antigens and adjuvants, which may act as a trigger for enthesitis (Gill et al., 2018). These adjuvants may activate entheseal stromal and resident immune cell populations, leading to the activation of the IL-23/IL-17 axis and secretion of pro-inflammatory cytokines—resulting in enthesitis, osteitis, and local joint inflammation (Costello et al., 2013; Babaie et al., 2018).

In the last 5 years, more than 100 published papers have focused on the interaction between the microbiome and AS (Costello et al., 2013; Xu et al., 2016; Simone et al., 2018; Watad et al., 2018; Huang et al., 2020; Zhou et al., 2020). Bacteroides coprophilus, Parabacteroides distasonis, Eubacterium siraeum, Acidaminococcus fermentans, and Prevotella copri were enriched in AS patients, and pathway analysis revealed increased oxidative phosphorylation, LPS biosynthesis, and glycosaminoglycan degradation (Zhou et al., 2020). Patients with AS have a distinct fecal bacterial signature, which is linked to fecal calprotectin levels, a marker of intestinal inflammation (Klingberg et al., 2019), suggesting a local interplay between intestinal bacteria and gut inflammation in AS. While most of the previous studies were focused on the role of bacteria, few studies have focused on mycobiota, a prototypic and best studied subtype that showed that the gut mycobiota of AS patients was characterized by higher levels of Ascomycota (especially the Dothideomycetes class) and decreased abundance of Basidiomycota, mainly contributed by the death of Agaricales. Compared to healthy controls, AS patients showed decreased ITS2/16S biodiversity ratios and altered bacterial–fungal interkingdom networks (Li et al., 2019). Compared with nonsteroidal anti-inflammatory drugs, treatment with biological agents induced obvious changes in the gut mycobiota of AS patients, and this result was highly associated with disease activity indices, including the AS disease activity index, C-reactive protein, erythrocyte sedimentation rate, and Bath AS disease activity index (Li et al., 2019). In addition, altered mycobiota in AS patients was also found to be associated with the degree of radiographic damage (Li et al., 2019), and in the same study identified a disease-specific interkingdom network alteration in AS, suggesting that fungi, or the interkingdom interactions between bacteria and fungi, may play an essential role in AS development. In our laboratory monkey farm management, we found that some of the Cynomolgus macaques had similar clinical symptoms as AS patients, based on the imageological diagnosis and the clinician’s diagnosis. Therefore, in this study, we aimed to further explore the roles and interactions between bacteria and fungi in AS to prove more prototypic and best-studied subtypes of gut mycobiota and microbiota for diagnosis and treatment.

2. Results2.1. Altered mycobiota in AS monkeys

Fecal samples from 10 AS monkeys with clinical and radiographic assessments (Jia et al.,2022) and 10 healthy controls were collected to systematically characterize thegut mycobiota and microbiota in AS monkeys by 16S rRNA and ITS2 DNAsequencing.. By sequencing ITS2, alpha diversity index difference analysis showed ACE, Chao1, Simpson, and Shannon were significant differences (Table 1, p<0.05), which indicated that the richness and diversity of mycobiota in AS were different from that of the healthy controls. The beta diversity of principal component analysis (PCA) and principal coordinate analysis (PCoA) could successfully distinguish between the two groups (Figures 1A, B). The histogram in Figure 1C showed that the mycobiota structure at the genus level in AS was different from that in the control group, and the heatmap at the genus level (Figure S1A) showed that most of the relative abundance of mycobiota was decreased in the AS, while a few of Kazachstania and Acaulium were increased (Table 2, p<0.05). The result of ANOVA at species level of mycobiota was shown in Table 3. More details of mycobiota structure are presented in Figure S2, which shows that AS and healthy controls could be successfully differentiated at all fungi levels, as the phyla of Rozellomycota, Basidiomycota, Mortierellomycota, Olpidiomycota, Mucoromycota, and Ascomycota (Figure S2A); the class of Saccharomycetes, Agaricomycetes, Archaeorhizomycetes, Sordariomycetes, Orbiliomycetes, Geminibasidiomycetes, Umbelopsidomycetes, Ustilaginomycetes, Rozellomycotina cls Incertae sedis, Malasseziomycetes, and Dothideomycetes (Figure S2B); the order of Helotiales, Thelephorales, Saccharomycetales, Cystofilobasidiales, Polyporales, Filobasidiales, and Dothideales (Figure S2C); the family of Saccharomycetaceae, Mrakiaceae, Nectriaceae, Debaryomycetaceae, Filobasidiaceae, Ophiocordycipitaceae, Thelephoraceae, and Pichiaceae (Figure S2D); the genera of Itersonilia, Kazachstania, Trichothecium, Plectosphaerella, Pestalotiopsis, Apiotrichum, Saccharomyces, Ophiocordyceps, Pichia, Candida, and Naganishia (Figure S2E); the species of Itersonilia perplexans, Kazachstania pintolopesii, K. heterogenica, C. albicans, Nawawia malaysiana, Verticillium dahliae, Curvularia americana, Saccharomyces cerevisiae, Pichia kluyveri, Apiotrichum montevideense, Penicillium herquei, Ganoderma sichuanense, Ophiocordyceps sinensis, and Naganishia albida (Figure S2F). The LEfSe (line discriminant analysis [LDA] effect size) results show that the K. pintolopesii, Saccharomycetaceae, Kazachstania, Saccharomycetales, and Saccharomycetes were specially enriched in the AS group (Figure 1D; more details regarding phylum, class, order, family, genus, and species are given in Figure S2).

Table 1 Alpha diversity index analysis of fecal mycobiota of Cynomolgus macaques (X¯SEM).

Figure 1 Altered mycobiota biodiversity and composition in AS monkeys. (A) Principal Component Analysis (PCA). Note: Dot represents each sample. Different colors represent different groups. The abscissa represents the first principal component, and the percentage represents the contribution value of the first principal component to the sample difference. The ordinate represents the second principal component, and the percentage represents the contribution of the second principal component to the sample difference. (B) Principal-coordinate analysis (PCoA) using the binary Jaccard distance algorithm. The horizontal and vertical coordinates are the two characteristic values that lead to the largest difference between samples, reflecting the major influence degree in the form of a percentage. (C) Global composition of fungal microbiota at the genus levels showing in bar chart. The more details see also in Figures S1, S2. (D) LEfSe (Line Discriminant Analysis (LDA) Effect Size) of AS and health monkeys. Note: The figure shows species whose LDA Score is greater than the set value (the default setting is 4.0). The length of the bar chart represents the size of the impact of different species (i.e., LDA Score). Different colors represent species in different groups.

Table 2 ANOVA (Analysis of Variance) at genus level of mycobiota of Cynomolgus macaques (p <0.05).

Table 3 ANOVA (Analysis of Variance) at species level of mycobiota (p <0.05).

2.2. Altered bacterial microbiota in AS monkeys

We analyzed the bacterial fraction of the microbiota using high-throughput sequencing of the bacterial 16S rRNA gene. Alpha diversity index difference analysis showed only ACE was significant differences (Table 4, p<0.05), which indicated that the species abundance in normal group was higher than that in AS group, but the diversity of bacterial microbiota was no different between two groups. The beta diversity of principal component analysis (PCA) and principal coordinate analysis (PCoA) could successfully distinguish between the two groups (Figures 2A, B). The histogram in Figure 2C showed that the bacterial microbiota structure at the genus level in AS was different from the control group, and the heatmap at the genus level (Figure S3A) showed that Prevotella 2, Blautia, Faecalibacterium, Lachnospiraceae UCG-008, an uncultured bacterium of Erysipelotrichaceae, Ruminococcaceae UCG-003, an uncultured bacterium of Lachnospiraceae, Negativibacillus and Mogibacterium were increased, while the others were decreased (Table 5, p<0.05). The result of ANOVA at species level of bacterial microbiota of Cynomolgus macaques was shown in Table 6. More details with respect to phylum, class, order, family, genus, and species are shown in Figure S4. A LefSe analysis was adopted to identify the bacterial groups that showed considerable and significant differences between the AS and control groups. As shown in Figure 2D, the comparison between the AS and control groups revealed that Prevotalla 2, an uncultured bacterium of Prevotella 2, f Ruminococcaceae, Clostridia, and Clostridiales were specially enriched in the AS group.

Table 4 Alpha diversity index analysis of fecal bacterial microbiota of Cynomolgus macaques (X¯ ± SEM).

Figure 2 Altered bacterial microbiota biodiversity and composition in AS monkeys. (A) Principal Component Analysis (PCA). Note: Dot represents each sample. Different colors represent different groups. The abscissa represents the first principal component, and the percentage represents the contribution value of the first principal component to the sample difference. The ordinate represents the second principal component, and the percentage represents the contribution of the second principal component to the sample difference. (B) Principal-coordinate analysis (PCoA) using the binary Jaccard distance algorithm. The horizontal and vertical coordinates are the two characteristic values that lead to the largest difference between samples, reflecting the major influence degree in the form of a percentage. (C) Global composition of bacterial microbiota at the genus levels showing in bar chart. The more details see also in Figures S1, S2. (D) LEfSe (Line Discriminant Analysis (LDA) Effect Size) of AS and health monkeys. Note: The figure shows species whose LDA Score is greater than the set value (the default setting is 4.0). The length of the bar chart represents the size of the impact of different species (i.e., LDA Score). Different colors represent species in different groups.

Table 5 ANOVA (Analysis of Variance) at genus level of bacterial microbiota of Cynomolgus macaques (p <0.05).

Table 6 ANOVA (Analysis of Variance) at species level of bacterial microbiota of Cynomolgus macaques (p <0.05).

2.3. AS monkeys showed altered bacterial-fungal associations

Fungi and bacteria coexist in the human gut and interact with each other. Although they have been shown to contribute actively to health or disease, few studies have investigated whether the fungal microbiota in AS patients is perturbed. The alpha diversity indexes of Chao1, Ace, Shannon, and Simpson were significantly different (Table 1, p<0.05), which indicated that the richness and diversity of fungi and bacteria in AS were different from the healthy controls. The beta diversity of PCA and PCoA could successfully distinguish between both groups (Figure 1A). In addition to composition differences, we found that the microbiota network (which includes bacteria and fungi) at the genus level in AS monkeys was notably different from that in healthy controls (Figure 3). Specifically, the density of the network in AS monkeys was remarkably higher than that of the healthy individuals (Figure 3, r >0.1, p<0.05, top 50), which suggested an alteration of the entire ecosystem in the guts of the AS monkeys. A higher Spearman’s correlation was found in AS monkeys (104 correlation) compared with healthy controls (90 correlation) (Figure 3). According to the report of Li et al. on AS patients in northeast China, these networks were different (Li et al., 2019), which may be due to the limiting bacterial but not fungal environment of laboratory monkeys or the high temperature, humidity, and climatic conditions in southern China. However, these results need further epidemiological validation. Interestingly, in AS monkeys, we observed a positive correlation between the abundance of Candida and Ruminococcaceae UCG-002, Sarcina, and the p-2534-18B5 gut family group, and a negative correlation with Cladosporium. Blautia showed a negative correlation with Dialister, Itersonilia, Prevotella 7, Saitozyma, and a positive correlation with Agathobacter. The Christensenellaceae R-7 group showed a negative correlation with Alternaria, Cladosporium, Lachnospiraceae UCG-008, Ruminococcaceae NK4A214 group, and uncultured bacterium Ruminococcaceae. Kazachstania showed a negative correlation with Agathobacter, Cladosporium, Prevotella 9, and Candida and a positive correlation with Itersonilia. Lachnospiraceae UCG-008 showed a negative correlation with Itersonilia, the Ruminococcaceae NK4A214 group, and Saitozyma and a positive correlation with the uncultured bacterium Erysipelotrichaceae. Prevotella 7 showed a negative correlation with Clostridium sensu stricto 1 and Pestalotiopsis and a positive correlation with Alternaria, Penicillium, and Saitozyma. Treponema 2 showed a negative correlation with Holdemanella, Prevotella 1, and Saitozyma. Lactobacillus showed a negative correlation with Candida and the Prevotellaceae NK3B31 group and a positive correlation with Penicillium. Megasphaera showed a positive correlation with Dialister and a negative correlation with Prevotella 2 and Ruminococcaceae UCG-003. More details are given in Figure 3 (r >0.1, p<0.05, top 50). Taken together, these results suggest a complex relationship between bacteria and fungi in the gut microbiota. Further studies are needed to better understand the specific communication and interaction between the host and its resident bacteria.

Figure 3 Bacterial–fungal association analysis in the AS and health monkeys. Note: Network graph is a form of correlation analysis. A Spearman rank correlation analysis was conducted according to the abundance and change of each species in each sample, and the data with correlations greater than 0.1 and P<0.05 were screened to construct a correlation network. Based on network graph analysis, the coexistence relationship of species in environmental samples can be obtained, as can the interaction situation and important pattern information of species in the same environment, and the formation mechanism of phenotypic differences between samples can be further explained. (A) Correlation analysis of fecal microorganisms (fungi and bacteria) in healthy monkeys. (B) Correlation analysis of fecal microorganisms (fungi and bacteria) in AS monkeys. Note: The circle represents species, and the size of the circle represents the average abundance of species. The lines represent the correlation between two species, the thickness of the lines represents the strength of the correlation, and the color of the lines: orange represents a positive correlation and green represents a negative correlation.

2.4. K. pintolopesii induces multiple AS immune activation phenotypes in mice

K. pintolopesii was purchased from CBS strains in the Netherlands under platform number Bio-109890. For more information, visit http://www.biobw.org/International-strain/bio-109890.html. Then, K. pintolopesii (7.0 × 108 CFU/ml) was administered to mice by intragastrical for the last 7 days. Our results showed that the serum levels of IL-17A were significantly upregulated (Figure 4A, p<0.05), and the immunohistochemical staining of the colon showed that K. pintolopesii could obviously upregulate the expression of the IL-17A protein (Figure 4E). Further, the expression of SAA1/2, IL-23, and IL-17 were significantly upregulated (Figures 4B–D, p<0.05). The RNA sequencing of colonic tissues showed that about 866 mRNAs should have a change in expression, with 211 mRNAs upregulated and 655 downregulated (fold change [FC] ≥1.5 and FDR<0.05). The KEGG analysis of DEGs showed that the pathways of Cytokine–cytokine receptor interaction, Cell adhesion molecules (CAMs), Natural killer cell-mediated cytotoxicity, Herpes simplex infection, T-cell receptor signaling pathway, Measles, Influenza A, Osteoclast differentiation, Staphylococcus aureus infection, Viral myocarditis, Hematopoietic cell lineage, Autoimmune thyroid disease, Allograft rejection, Antigen processing and presentation, Inflammatory bowel disease (IBD), Intestinal immune network for IgA production, Graft-versus-host disease, and Primary immunodeficiency (Figure 4F) were the top enriched pathways. The GO classification of DEGs showed that the biological processes of immune system, inflammatory response (Figure S5A, p<0.05); cellular components of extracellular exosomes and cytosol (Figure S5B, p<0.05); and molecular functions of antigen binding, immunoglobulin receptor binding, MHC class Ⅱ protein binding, and CD8 receptor binding were enriched (Figure S5C, p<0.05). Furthermore, most pathways were related to infection immunity and inflammation.

Figure 4 K pintolopesii induces multiple AS immune activation phenotypes in mice. (A) K pintolopesii upregulates the levels of IL-17 in serum (p<0.05). (B) K pintolopesii upregulates the expression of SAA1/2 in colon tissues (p<0.05). (C) K pintolopesii upregulates the expression of IL-23 in colon tissues (p<0.05). (D) K pintolopesii upregulates the expression of IL-17 in colon tissues (p<0.05). (E) The expression of IL17A in colon tissue was detected by immunohistochemistry. (F) KEGG analysis of differentially expressed genes showed that the pathways of cytokine–cytokine receptor interaction, cell adhesion molecules (CAMs), natural killer cell-mediated cytotoxicity, Herpes simplex infection, T-cell receptor signaling pathway, measles, influenza A, chemokine signaling pathway, osteoclast differentiation, autoimmune thyroid disease, and inflammatory bowel disease (IBD) were enriched pathways (p<0.05). #p < 0.05, ##p < 0.01, ###p < 0.001.

Deep analysis showed that the mRNAs of Gzma (p = 1.39E−09), Gzmb (p = 1.20E−06), Pik3r5, Map3k14, Ctsw, Casp12, Traf1, Prf1, Pik3cd, Jun, Csf2rb, Csf2rb2, Bcl2a1b, and Bcl2a1d showed changes induced by K. pintolopesii. Granzyme B (GzmB) and Granzyme A (GzmA) are apoptotic proteinase-like serine proteases that can activate caspase-3 and caspase-8 and be transferred to target cells by cytotoxic lymphocytes to trigger cell apoptosis, pyroptosis, and kill virus-infected cells and tumor cells (Jiang et al., 2020; Zhang et al., 2020). This indicated that the necroptosis, pyroptosis, and/or apoptosis delay of neutrophils, lymphocytes, natural killer cells, and/or CD8 cytotoxic T killer cells were triggered by K. pintolopesii, which eventually induces autoimmune diseases or other inflammatory diseases. The changes of Il12rb1, Cd40, Il12rb2, Cxcl16, Ccl9, Ifnar2, Ltb, Csf1r, Il18rap, Il18r1, Il2ra, Bmp2, Il1b, Pdgfc, Csf3r, Cxcl9, Ccr9, Cd27, Il21r, Il4ra, Il2rg, Pdgfd, Il10ra, Crlf2, Ccl5, Ccl2, Inhbb, Il18, Cxcr6, Ccr2, Cxcr3, Cx3cr1, Cxcl2, Xcr1, Il2rb, Csf2rb, Csf2rb2, Ccl28, and Ccr5 mRNAs were enriched into the cytokine–cytokine receptor interaction; and previous studies have demonstrated that the chemokine receptor gene CXCR6 was highly expressed in +CPI colitis T cells and a variety of human cancer cells (Luoma et al., 2020). However, CXCL16, a CXCR6 ligand gene regulated by IFN and TNF signaling pathways, also showed significantly high expression in bone marrow cells and mast cells (Xiao et al., 2015), suggesting that CXCR6 inhibition could be a potential therapeutic target to reduce the metastasis of cancer cells and effectively improve intestinal inflammation. Thus, K. pintolopesii may be a trigger for the AS complication of intestinal inflammation. The changes of Itgb2, Lck, Lcp2, Ptpn6, Cd244a, Cd247, Vav2, Cd48, Gzmb, Fyn, Pik3r5, Ppp3cc, Ifnar2, Zap70, Klrk1, Klrd1, Klrc1, Tyrobp, Lat, Itgal, Nfatc1, Rac2, Vav3, Prf1, Icam1, Pik3cd, Mus_newGene_6834, Fcgr4, Prkcb, Sh3bp2, Fcer1g, and Mus_newGene_6925 were enriched into the Natural killer cell mediated cytotoxicity; and the changes of Itgb2, Cd86, Ltb, Ctla4, Cd28, Il1b, Itgal, Ccl5, Ccl2, H2-Aa, Icam1, H2-DMa, Il18, Jun, H2-Eb1, H2-Ab1, H2-DMb1, Mus_newGene_3650, Mus_newGene_6834, and Mus_newGene_6925 were enriched in rheumatoid arthritis, all of which indicated that K. pintolopesii is a trigger of AS-induced multiple AS immune activation phenotypes. Further details of other pathways are shown in Table S1.

Gene Set Enrichment Analysis (GSEA) using DEXSeq (version 1.12.2, default parameters) showed that positive regulation of GTPase activity, innate immune response (Figure S6A), inflammatory response (Figure S6B), G-protein coupled receptor signaling pathway, positive regulation of gene expression, response to lipopolysaccharide, positive regulation of ERK1 and ERK2 cascade, regulation of cell shape, cellular response to lipopolysaccharide, and immune response were the top 10 enriched biological processes (p<0.0015). The integral component of plasma membrane (Figure S6C), external side of plasma membrane, membrane raft, transcription factor complex (Figure S6D), lamellipodium, cell–cell junction, actin filament, MHC class II protein complex, Golgi medial cisterna, multivesicular body, and immunological synapse were the top 12 enriched cellular components (p<0.0016). GTP binding, GTPase activity, carbohydrate binding (Figure S6E), antigen binding (Figure S6F), integrin binding, receptor activity, protein binding, peptide antigen binding, transmembrane signaling receptor activity, virus receptor activity, cytokine receptor activity, and T-cell receptor binding were the top 12 enriched molecular functions (p<0.0016).

DEU analysis (DEXSeq was used to analyze DEU and performed FDR as the criterion<0.05, Table 7) showed that the exon E020 of Nlrp1b, a key mediator of programmed cell death, was changed (p = 2.07E−07). Overexpression of this gene was reported to induce apoptosis and pyroptosis, and this gene mediates diseases including palmoplantar carcinoma, multiple self-healing, and autoinflammation with arthritis and dyskeratosis, so K. pintolopesii could influence the exon E020 of Nlrp1b to mediate programmed cell death. The E003 (p = 3.55E−17) and E005 (p = 1.29E−15) of Bri3 were different, and a previous study has shown that Bri3 participates in TNFα-induced cell death (Wu et al., 2003). The E001 (p = 5.79E−07) and E004 (p = 7.38E−10) of Efhd2, which were demonstrated to regulate B-cell receptor (BCR)-induced immature and primary B-cell apoptosis (Kroczek et al., 2010) and play a role as a negative regulator of the canonical NF-κB activating branch (Avramidou et al., 2007) and control spontaneous apoptosis through the regulation of BCL2L1 abundance (Avramidou et al., 2007). Furthermore, there were another 23 exons in different genes that also changed (Table 7, FDR<0.05). Taken together, these data show that K. pintolopesii could influence cell pyroptosis and induce inflammation.

Table 7 DEU Analysis at exon levels in K. pintolopesii treated mice, DEXSeq was used to analysis DEU and performed FDR as the criterion<0.05.

2.5. K. pintolopesii induces PANoptosis of macrophage RAW264.7

RNA sequencing of RAW264.7 co-cultured with K. pintolopesii showed that about 2,622 mRNAs changed in expression, with 1,016 mRNAs upregulated and 1,606 downregulated (fold change [FC] ≥1.5 and FDR<0.05, Figure S7A). The differentially expressed genes related to apoptosis were Daxx, Akt2, Map2k1, Bcl2l1, Gadd45b, Nfkbia, Fos, Gadd45g, Ctsl, Tnf, Pmaip1, Ddit3, Pidd1, Xiap, Casp8, Traf1, Ctso, Ctsk, Nfkb1, Pik3r3, Itpr1, Itpr2, Ctsc, Casp3, Birc3, Atm, Map2k2, Mcl1, Atf4, Bcl2, Sptan1, Actg1, Csf2rb2, Ctsf, Bcl2a1b, Bcl2a1d, Bcl2a1a, and Gm49392 (Figure 5A); related to necroptosis were Irf9, Cybb, Sqstm1, Tnfaip3, Hsp90aa1, Dnm1l, Eif2ak2, Tnf, Xiap, Casp1, Casp8, Traf5, Trpm7, Il1a, Tlr3, Birc3, Tyk2, Pygm, Nlrp3, Pla2g4c, Tlr4, H2afv, Sntb2, H2afx, Pla2g4a, Bcl2, and Hist3h2a (Figure 5A); and related to pyroptosis were NLRP3, NLRC4, NAIP6, NAIP5, NAIP2, CASP8, CASP1, and AIM2 (Figure 5A), which suggested that K. pintolopesii may induce PANoptosis in RAW264.7 cells.

Figure 5 Genes related to PANoptosis were screened from differentially expressed genes in the transcriptome of RAW264.7 co-cultured with K. pintolopesii. (A) Heatmap of differentially expressed genes related to PANoptosis.

The KEGG analysis of DEGs showed that the pathways of Herpes simplex virus 1 infection, HIF-1 signaling pathway, NF-kappaB signaling pathway, TNF-signaling pathway, ribosome, apoptosis, legionellosis, small cell lung cancer, p53 signaling pathway, platinum drug, cellular senescence, toxoplasmosis, Epstein–Barr virus infection, human T-cell leukemia virus 1 infection, endocytosis, microRNAs in cancer, FoxO signaling pathway, valine, leucine and isoleucine degradation, rheumatoid arthritis, fluid shear stress, and atherosclerosis were the top enriched pathways (Figures 6A, S7B, C). The KEGG analysis did not enrich the pyroptosis-related pathways because the KEGG data did not include the cell apoptosis pathway yet. The GO classification of DEGs showed that the biological processes of regulation of apoptotic processes, positive regulation of NF-kappaB transcription (Figure S8A, p<0.05), cellular components of nucleus, cytoplasm, and cytosol (Figure S8B, p<0.05), and molecular functions of peptide antigen binding and T-cell receptor binding were enriched (Figure S8C, p<0.05). Furthermore, most pathways were related to the infection immunity and inflammation. There is no doubt that K. pintolopesii can obviously influence the function of macrophages.

Figure 6 The effect of K. pintolopesii on macrophage RAW264.7 (A) KEGG analysis of transcriptome of RAW264.7 treated by K. pintolopesii. (B) CCK8 was used to detect of RAW264.7 treated by the supernatant of K. pintolopesii or K. pintolopesii (C) The content of TNF-α of the supernatant of K. Pintolopesii co-cultured RAW264.7 was determined by Elisa. (D) Western blot was used to detect the GSDMD protein level of RAW264.7 after co-culture with K. pintolopesii. (E) The levels of IL-17A in the serum of mice infected by the products (culture supernatant) secreted by K. pintolopesii through intraperitoneally were determined by Elisa. #p < 0.05, ##p < 0.01, ###p < 0.001.

In addition tests showed that the number of living cells was decreased significantly (p<0.05) after 6 h co-cultured with the supernatant of K. pintolopesii or K. pintolopesii using CCK8, and many floating cells (Figure 6B); and the content of TNF-α of the supernatant of K. Pintolopesii co-cultured RAW264.7 was significantly increased (Figure 6C, p<0.05), also the GSDMD protein levels (Figure 6D), which indicated that K. pintolopesii can cause the pyroptosis. Furthermore, the levels of IL-17A in the serum of mice infected by the products (culture supernatant) secreted by K. pintolopesii intraperitoneally were significantly higher (Figure 6E, p<0.05). This suggested that K. pintolopesii and their secreted products could induce severe inflammatory reactions in vivo and in vitro.

2.6. PANoptosis products promote the proliferation, biofilm formation, and morphological changes of K. pintolopesii

To determine the next reaction of PANoptosis induced by K. pintolopesii, PANoptosis products were added into the culture of K. pintolopesii. Unexpectedly, the colony of K. pintolopesii showed positive correction to the PANoptosis products (Figure 7A), which indicated that the PANoptosis products could promote the proliferation of K. pintolopesii. The morphology of K. pintolopesii were changed in the PANoptosis products and TNF-α (Final concentration of 1,000 pg/ml) could change (Figure 7B), as the same of C. albicans changed from the yeast state to the mycelium state. Those findings suggest that chronic inflammation may cause the transformation of conditioned pathogenic fungi from yeast to mycelium, thereby increasing their pathogenicity. The mechanism of morphological changes of opportunistic pathogens caused by host inflammatory factors is still unclear, and further experiments are needed to elucidate the mechanism. It is suggested that the exotoxin of K. pintolopesii (compared to C. albicans) may be an important substance in activating the IL-17RA pathway.

Figure 7 Effect of PANoptosis products on K. pintolopesii. (A) Plate count. PANoptosis products significantly promotes K. pintolopesii proliferation. (B) PANoptosis products and TNF-α (Final concentration of 1000 pg/mL) could change the morphology of K. pintolopesii and Candida albicans.

3. Discussion

In this study, fecal samples from 10 AS monkeys with clinical and radiographic assessments and 10 healthy controls were collected to systematically characterize the gut mycobiota and microbiota in AS by 16S rRNA gene- and ITS2-based DNA sequencing. Our results showed that the gut mycobiota of K. pintolopesii, I. perplexans, Saccharomycetaceae, Coniochaeta, and Coniochaetaceae were especially enriched in AS. The microbiota of AS monkeys was characterized by increased abundance of the Prevotalla 2 genus, Ruminococcaceae family, Clostridia class, and Clostridiales order. Compared to healthy controls, decreased ITS2/16S biodiversity ratios and altered bacterial-fungal interkingdom networks were observed in AS monkeys. All these findings implied that fungal and bacterial dysbiosis are closely linked to AS and that their communication and interaction play an important role in autoimmune response.

In this paper, the LEfSe analysis showed that the biomarker fungi and bacteria of Saccharomycetaceae, Saccharomycetes, Kazachstania, and K. pintolopesii were rich in AS, while Sordariomycetes, Lactobacillaceae, Lactobacillus, Hypocreales, L. horakii, Trichoderma, Hypocreaceae, T. virens, Agaricomycetes, Russulaceae, and Lactarius were rich in healthy control monkeys (Figure 1C), which indicated that fungal and bacterial dysbiosis are serious occurrences in AS, and the changes of fungi showed more significant and discriminability than bacteria. There is growing research on the effects of fungi on human health, including obesity (Kim et al., 2017; Borges et al., 2018), cancer (Wang et al., 2018; Wu et al., 2019; Zhai et al., 2020), immunity, and inflammation (Leonardi et al., 2018; Li et al., 2019; Richard and Sokol, 2019). Despite the availability of new anti-fungal drugs, infections caused by pathogenic fungi continue to constitute a serious medical problem.

Studies have shown that fungal dysbiosis is correlated with the severity of IBD and, in combination with bacterial dysbiosis, could aggravate DSS-induced colitis. Fungi induce strong ecological changes to the gut microbiome, promote strong local and systemic immunological changes, and affect the immunophenotypes of the intestinal and pulmonary inflammatory responses (van Tilburg et al., 2020). These aspects indicate that fungi can modulate the severity and immunophenotype of the immune response when co-colonized with bacteria. In plant sciences, considerable communication and interaction has been recognized between fungi and bacteria (Li et al., 2019; Jiang et al., 2020; Fitzpatrick et al., 2020), while in humans, there are inadequate data.

A Chinese cohort of 84 AS patients showed that the gut microbiota was perturbed in untreated AS patients with diagnostic potential, and some AS-enriched species might be triggers of autoimmunity by molecular mimicry. More importantly, different inflammatory arthritis share some common microbial signatures (Zhou et al., 2020) such as A. fermentans, B. coprophilus, P. distasonis, E. siraeum, and P. copri as markers specifically enriched in the AS gut, and increased pathways including oxidative phosphorylation, LPS biosynthesis, and glycosaminoglycan degradation.

Genome-wide association studies of this highly genetic disease have implicated specific immune pathways, including the interleukin (IL)-17/IL-23 pathway, control of NF-κB activation, amino acid trimming for major histocompatibility complex (MHC) antigen presentation, and other genes controlling the CD8 and CD4 T cell subsets (Smith, 2015). Bone changes, including systemic bone loss, joint erosion, and calcaneal formation, reflect the combined effects of IL-23 and IL-17 in psoriatic arthritis (PsA) and AS (Babaie et al., 2018; Gravallese and Schett, 2018; Sieper et al., 2019). IL-17A directly or indirectly promotes osteoclast production by producing or inducing the expression of the nuclear factor receptor activator-B ligand (RANKL) (Lin et al., 2014; Ganesan and Rasool, 2017), and IL-23 multiplies the effects on osteoclasts (Ganesan and Rasool, 2017; Shukla et al., 2018). IL-17A has different effects on the maturation of precursor cells into osteoblasts, which depend on the differentiative stage of the precursor cells (Jo et al., 2018). IL-17A blockers have been shown to inhibit bone erosion, delay systemic bone loss in PsA and AS, and form mossy hyperplasia in PsA (Kampylafka et al., 2018; Samarpita et al., 2018). Therefore, in this study, we focus on this pathway to try to identify the pathogenic factors derived from microorganisms.

One suggestion for studying host–parasite interactions of C. albicans in rodents is that C. albicans is not a natural colonizer of mucosal surfaces in these animals and that the equivalent of normal colonizers in rodents is K. pintolopesii (Naglik et al., 2008). C. albicans is an important commensal fungus in the human body, and it is also the main fungal inducer of the Th17 response in humans (Bacher et al., 2019). K. pintolopesii can naturally colonize rodents and monkeys, so it is of great significance to further elucidate the function and pathogenicity of K. pintolopesii in rodents and monkeys. In this study, we found that K. pintolopesii would be an opportunistic and/or pathogenic inducer in AS monkeys. We also found that K. pintolopesii upregulated the levels of IL-23 and IL-17 in mice and could obviously damage the colonic mucosa and its integrity. The KEGG analysis of DEGs of colonic tissues showed that the pathways of cytokine–cytokine receptor interaction, cell adhesion molecules (CAMs), natural killer cell mediated cytotoxicity, T-cell receptor signaling pathway, chemokine signaling pathway, osteoclast differentiation, phagosome, NF-κB signaling pathway, rheumatoid arthritis, TNF signaling pathway, autoimmune thyroid disease, inflammatory bowel disease (IBD), Staphylococcus aureus infection, intestinal immune network for IgA production, antigen processing and presentation, and B-cell receptor signaling pathway (Figure 4F) were the top enriched pathways that were reportedly associated with AS, as determined in the colon of K. pintolopesii-infected mice. Deep analysis showed that Gzma and Gzmb would be a trigger induced by K. pintolopesii for the pyroptosis of neutrophils, lymphocytes, natural killer cells, and/or CD8 cytotoxic T-killer cells. DEU analysis showed that the exon E020 of Nlrp1b, E003 (p = 3.55E−17) and E005 (p = 1.29E−15) of Bri3, E001 (p = 5.79E−07), and E004 (p = 7.38E−10) of Efhd2 would also be the influencing genes on cell pyroptosis and induce inflammation by K. pintolopesii.

>As the first line of defense of immune response, macrophages can maintain the homeostasis of the body by identifying and removing pathogens, killing target cells, antigen presentation, immune regulation, and other functions (Vega-Perez et al., 2021). However, excessive aggregation and activation of macrophages can also lead to tissue damage, so their homeostasis plays an important role in tissue homeostasis (Uderhardt et al., 2019; Locati et al., 2020). In our tests, we found that there were 2,622 DEGs in the RAW264.7 co-cultured with K. pintolopesii (Figure 6 and Figure S6), which indicated that the K. pintolopesii infection induced violent immune responses. Deep analysis showed that genes involved inapoptosis, necroptosis and/or pyroptosis were changed (Figure 6A and Figure S6), and some of them were participate in two or all three pathways (Figure 6A), so we conclude that the K. pintolopesii induces PANoptosis of macrophage RAW264.7. Previous findings indicated that inflammatory disease-induced Th17 cells require the assistance of IL-1β and IL-23. In vitro, to culture the maintain tissue homeostasis Th17 cells requiring the stimulation of TGF-β and IL-6 or IL-21 (Honda and Littman, 2016; Belkaid and Harrison, 2017), and recent reports showed that SAAs could induced pathogenic Th17 cells (Lee et al., 2020).We also found the expressions of IL-17RA, IL-18RAP, IL-23A, and Il-6RA were upregulated (Figure S6A), and the NF-kappaB signaling pathway, TNF-signaling pathway were activated (Figures 6B, C); and the expressions of SAA1, IL-17A, and IL-23 were also upregulated in K. pintolopesii infected mice (Figure 4). All which indicated that the K. pintolopesii infected may through induced the PANoptosis of macrophage, then induces the pathogenic Th17 cells, but the action of basis and mechanism still more further studies.

In summary, our results showed that the gut mycobiota of K. pintolopesii, I. perplexans, Saccharomycetaceae, and Coniochaeta were especially enriched in AS. In addition, the microbiota of AS monkeys was characterized by increased abundance of the Prevotalla 2 genus, Ruminococcaceae family, and Clostridiales order. Compared to healthy controls, decreased ITS2/16S biodiversity ratios and altered bacterial–fungal interkingdom networks were observed in AS monkeys. The characteristic gut mycobiota of K. pintolopesii induces multiple AS phenotypes in mice, such as fungal and bacterial dysbiosis, rheumatoid arthritis, inflammatory bowel disease, and varied immune responses by activation of NK cells, platelets, T cells, leukocytes, and/or B cells. All these findings imply that fungal and bacterial dysbiosis have a close link to AS and that their communication and interaction plays an important role in the autoimmune response. The exotoxin of K. pintolopesii may activate the IL-17RA pathway, which needs further experimental evidence. It is an interesting discovery that macrophage PANoptotic products and host TNF-a can promote the proliferation and morphological transformation of opportunistic pathogens, but more experiments are needed to clarify the specific mechanism.

3.1. Study limitations

Our study has some limitations. First, although we have evaluated the effects of the characteristic gut mycobiota of K. pintolopesii on inducing multiple AS phenotypes in normal mice, the study lacks functional tests on AS model animals. Hence, the mechanism of action of K. pintolopesii on AS is still unknown. Therefore, future studies should aim to use targeted strains in germ-free mice with multi-omics studies to reveal the communication and interaction mechanism of action between these AS characteristic pathogenic organisms for further validation.

4. Materials and methods4.1. Key resources table

4.2. Cell cultures

RAW264.7 cells were obtained from the National Collection of Authenticated Cell Cultures. Cells were cultured in DMEM medium with 10% fetal bovine serum as routine. All cells were grown in a humidified incubator at 37°C with 5% CO2.

4.3. CCK-8 detection

A total of 500 μl of RAW264.7 cell suspension was prepared in a 24-well plate. Plates were pre-incubated in a cell incubator for 12 h with 80 μl of 108 CFU/ml K. pintolopesii (washed with PBS) or 80 μl of culture supernatant of K. pintolopesii. The control group was added with 80 μl of LB medium. Incubate in the incubator for 6 h; replace the fresh cell culture medium and add 50 μl of CCK-8 solution to each well (be careful not to generate bubbles in the well, which will affect the reading of the OD value).

4.4. Animals

Ten-week-old, adult, male C57BL/6 mice obtained from the Center of Laboratory Animal of Guangdong Province (SCXK [Yue] 2018-0002, SYXK [Yue] 2015-0102) were pair-housed in plastic cages in a temperature-controlled (25 ± 2°C) colony room under a 12/12-h light/dark cycle. The animals had ad libitum access to food and water. All experimental protocols were approved by the Center of Laboratory Animals of the Guangdong Institute of Microbiology. All efforts were made to minimize the number of animals used.

4.5. Histopathology and immunostaining

The arthrosis and colon tissues were removed and fixed in 4% paraformaldehyde at pH 7.4 for further pathological experiments. These tissue samples were embedded in paraffin sections after fixation, washing, dehydration, transparency, dipping in wax, and embedding. Four samples from each group were fixed in a 4% paraformaldehyde solution and prepared as paraffin sections. Sections were stained with hematoxylin and eosin (H&E). Immunostaining was performed using paraffin-embedded sections (3-μm-thick) and a two-step peroxidase-conjugated polymer technique (DAKO Envision kit, DAKO, Carpinteria, CA). Slides were observed under a light microscope.

4.6. Microbiome analysis

Total DNA was extracted from 250 to 500 mg of each fecal sample using the QIAamp DNA stool minikit (Qiagen, Germany). The purity and concentration of metagenomic DNA were measured with a NanoDrop 2000 spectrophotometer (Thermo Fisher Scientific, USA).

4.6.1. For 16S

Construction and sequencing: after extracting the total DNA of the sample, primers were obtained according to the design of the conservative region, and the end of the primer was added with a sequencing connector. The microbial 16S rRNA genes were amplified using the forward primer 5′-ACTCCTACGGGAGGCAGCA-3′ and the reverse primer 5′-GGACTACHVGGGTWTCTAAT-3′. PCR amplification was performed, and the products were purified, quantified, and homogenized to form a sequencing library. The established library was firstly inspected by the library, and the qualified library was sequenced by Illumina HiSeq 2500. The original image data files obtained by high-throughput sequencing (Illumina HiSeq and other sequencing platforms) were analyzed and converted into original sequencing Reads by Base Calling, and the results were stored in FASTQ (fq for short) file format, which contained the sequence information of sequencing sequences (Reads) and the corresponding sequencing quality information.

4.6.2. For ITS1

Construction and sequencing: after extracting the total DNA of the sample, primers were obtained according to the design of the conservative region, and the end of the primers was added with a sequencing connector. The microbial ITS genes were amplified using the forward primer ITS1F 5’-CTTGGTCATTTAGAGGAAGTAA-3’ and the reverse primer ITS2 5’-GCTGCGTTCTTCATCGATGC-3’. The rest of the steps are the same as for 16S rDNA sequencing.

4.6.3. Bioinformatics analysis

Sequences of the V3–V4 region of 16S rRNA genes and ITS1 were detected using an Illumina HiSeq 2500 platform (Biomarker Technology Co., Ltd., Beijing, China). OTUs present in 50% or more of the fecal samples were identified as core OTUs. The observed species and Shannon and Simpson indices were calculated with QIIME (version 1.9.1). The abundance and diversity of the OTUs (beta diversity) were examined using principal-coordinate analysis (PCoA) with weighted or unweighted UniFrac analysis in R software. The linear discriminant analysis (LDA) and effect size analysis (LEfSe) were used with the Kruskal–Wallis rank sum test to detect features with significantly different abundances between assigned taxa, and the linear discriminant analysis was performed to estimate the effect size of each feature. The bacterial groups with an LDA score of 4.00 were presented as the most significantly abundant group in the indicated group.

4.7. RNA isolation and sequencing

RNA concentration and purity were measured using NanoDrop 2000 (Thermo Fisher Scientific, Wilmington, DE). RNA integrity was assessed using the RNA Nano 6000 Assay Kit of the Agilent Bioanalyzer 2100 system (Agilent Technologies, CA, USA). A total of 1 μg RNA per sample was used as input material for the RNA sample preparations. Sequencing libraries were generated using the NEBNext UltraTM RNA Library Prep Kit for Illumina (NEB, USA) following the manufacturer’s recommendations, and index codes were added to attribute sequences to each sample. Briefly, mRNA was purified from total RNA using poly-T oligo-attached magnetic beads. Fragmentation was carried out using divalent cations under elevated temperature in the NEBNext First Strand Synthesis Reaction Buffer (5X). First-strand cDNA was synthesized using random hexamer primers and M-MuLV reverse transcriptase. Second-strand cDNA synthesis was subsequently performed using DNA polymerase I and RNase H. Remaining overhangs were converted into blunt ends via exonuclease and polymerase activities. After adenylation of the 3’ ends of DNA fragments, NEBNext adaptor with a hairpin loop structure were ligated to prepare for hybridization. To select cDNA fragments of preferentially 240 bp in length, the library fragments were purified with the AMPure XP system (Beckman Coulter, Beverly, USA). Then 3 μl of USER Enzyme (NEB, USA) was used with size-selected, adaptor-ligated cDNA at 37°C for 15 min, followed by 5 min at 95°C before PCR. Then PCR was performed with Phusion High-Fidelity DNA Polymerase, Universal PCR primers, and Index (X) Primer. Finally, PCR products were purified (using the AMPure XP system), and library quality was assessed on the Agilent Bioanalyzer 2100 system. The clustering of the index-coded samples was performed on a cBot Cluster Generation System using the TruSeq PE Cluster Kit v4-cBot-HS (Illumia) according to the manufacturer’s instructions. After cluster generation, the library preparations were sequenced on an Illumina platform, and paired-end reads were generated.

4.8. Western blotting analysis

Briefly, global brain tissue was dissected from treated mice (purchased from Beijing HFK Bioscience Co., LTD [Certificate No: SCXK (Jing) 2014-0004]) and proteins were extracted with radioimmunoprecipitation assay (RIPA) lysis buffer (Thermo Scientific™ T-PER™ Tissue Protein Extraction Reagent, 78510). The proteins were separated by sodium dodecyl sulfate-polyacrylamide gel electrophoresis and transferred onto polyvinylidene fluoride membranes. After blocking with 5% nonfat dry milk in Tris-buffered saline (20 mM Tris–HCl, 500 mM NaCl, pH 7.4) with 0.2% Tween-20 (Aladdin, T104863), the membranes were probed with antibodies overnight at 4°C, followed by incubation with a horseradish peroxidase-conjugated goat anti-mouse (Servicebio, G2211-1-A) or goat anti-rabbit (Servicebio, G2210-2-A) IgG secondary antibody (1:2,000). Band intensity was quantified using ImageJ software (NIH).

4.9. Gram staining

Put the strain in nutrient broth medium, shake overnight at 37 °C, prepare OD600 = 1 bacterial solution, and dilute with nutrient broth 1:100. After dilution, add it to a 96-well polystyrene-sterilized microplate (180 μl/well) inoculate three wells. The control group was added with 20 μl of sterile nutrient broth, and the supernatant group was added with 20 μl of co-culture supernatant. In the TNF-α group, 20 μl 10,000 pg/ml TNF-α was added to bring the final concentration of TNF to 1,000 pg/ml. All groups were cultured at