Complementary Contribution of Fungi and Bacteria to Lignocellulose Digestion in the Food Stored by a Neotropical Higher Termite

Lignocellulose digestion in termites is achieved through the functional synergy between gut symbionts and host enzymes. However, some species have evolved additional associations with nest microorganisms that collaborate in the decomposition of plant biomass. In a previous study, we determined that plant material packed with feces inside the nests of Cornitermes cumulans (Syntermitinae) harbors a distinct microbial assemblage. These food nodules also showed a high hemicellulolytic activity, possibly acting as an external place for complementary lignocellulose digestion. In this study, we used a combination of ITS sequence analysis, metagenomics, and metatranscriptomics to investigate the presence and differential expression of genes coding for carbohydrate-active enzymes (CAZy) in the food nodules and the gut of workers and soldiers. Our results confirm that food nodules express a distinct set of CAZy genes suggesting that stored plant material is initially decomposed by enzymes that target the lignin and complex polysaccharides from fungi and bacteria before the passage through the gut, where it is further targeted by a complementary set of cellulases, xylanases, and esterases produced by the gut microbiota and the termite host. We also showed that the expression of CAZy transcripts associated to endoglucanases and xylanases was higher in the gut of termites than in the food nodules. An additional finding in this study was the presence of fungi in the termite gut that expressed CAZy genes. This study highlights the importance of externalization of digestion by nest microbes and provides new evidence of complementary digestion in the context of higher termite evolution.


INTRODUCTION
Termites (Blattodea: Isoptera) are eusocial cockroaches comprising over 3,000 species widespread in tropical and subtropical regions. These insects play an important role in ecosystems through the decomposition of plant lignocellulose (Brune, 2014;Mikaelyan et al., 2015), having a positive impact in nutrient cycling and soil-water dynamics through the construction of their nests and tunnel networks in soil (Neupane et al., 2015;Siebers et al., 2015;Jouquet et al., 2018).
All termites live in nests containing hundreds to thousands of individuals divided in castes (reproductive individuals, workers, and soldiers) that cooperatively conduct the tasks of the colony. Nests protect termites against the environment and predators and offer a place for feeding, reproduction, and nursery (Bignell and Eggleton, 2000). Higher termites (Termitidae) have evolved remarkable nesting strategies including the construction of epigeal mounds which can reach a prominent size in some species (Korb, 2011). Mounds are usually made of micro-aggregates of soil and feces forming porous structures that control the microclimate (Korb, 2003;Singh et al., 2019) and could function as an external reservoir of microorganisms (Fall et al., 2007;Manjula et al., 2014) or even like a fermenting chamber (Korb, 2011;Schmidt et al., 2014). Some authors have hypothesized that externalization of digestion by nest microorganisms, as seen in the subfamilies Macrotermitinae and Sphaerotermitinae (Garnier-Sillam et al., 1989;Aanen et al., 2002), was the first step in the evolution of Termitidae and might have driven the loss of gut protozoans and the acquisition of specialized lignocellulolytic bacteria resulting in higher efficiency in plant decomposition (Brune, 2014;Chouvenc et al., 2021).
Cornitermes cumulans (Kollar, 1932) (Termitidae: Syntermitinae) is a mound-building termite in pastures and savannas of central South America and in some areas its nest density is exceptionally high compared to other termite species (Coles De Negret and Redford, 1982;Buschini, 2006). The hard epigeal mounds are made externally of soil with an internal core of soft carton made from fecal material and soil particles (Redford, 1984). Although feeding in situ has been reported (Coles De Negret and Redford, 1982), C. cumulans is a harvester termite, cutting and transporting grass and litter into its mound. The plant material is then stored in nodules made with feces and saliva ( Figure 1A) and subsequently consumed by workers (Lima and Costa-Leonardo, 2007). Food nodules showed a high hemicellulolytic activity possibly mediated by bacteria and fungi. It is possible that the pre-processing of plant material by the action of microbes in the food nodules before gut passage acts as a complementation of biomass digestion in C. cumulans. However, only partial taxonomic information for fungal strains was reported because less than 10% of reads were assigned to order level at the time .
In this study, we assessed the taxa and guilds that differentiate fungal communities associated to the food nodules of C. cumulans by using deep taxonomic and guild assignment of datasets generated through internal transcribed spacer (ITS) sequencing. To investigate the potential complementation of lignocellulose digestion among the termite host, the gut symbionts and the microbiota of food nodules, we then performed metagenomic and metatranscriptomic analyses on this tri-partite system to assess: (i) the repertoire and functional contribution of lignocellulolytic enzymes; (ii) how these enzymes are linked taxonomically to the termite holobiont (host, gut, and nest microbiota); and (iii) the differences of enzyme functionality across sterile castes.

ITS Datasets, Sequence Analysis and Guild Assignment
Internal transcribed spacer reads were obtained from raw Illumina datasets of three (3) C. cumulans colonies (CC1, CC2, and CC3) deposited in the European Nucleotide Archive (ENA) with accession no. PRJEB17080 . Libraries were processed using the UPARSE pipeline (Edgar, 2013). Paired end reads were first merged using fastq_mergepairs from USEARCH package version 8.1.1803. Reads with a minimum overlap of 50 bp and a maximum expected error of 0.5 were mapped to fungal ITS2 region using ITSx software (Bengtsson-Palme et al., 2013). Reads were further compared to ITS2 UCHIME reference dataset to filter chimera sequences using chimera UCHIME (Edgar et al., 2011), also implemented in USEARCH package. The filtered reads were then subjected to clustering into OTUs (operational taxonomic units) at 97% of sequence similarity according to UPARSE-OUT algorithm. The OTU table was generated by mapping the reads from each sample back to the OTUs. Taxonomic assignment was performed using sintax command as implemented in USEARCH version 10.0.240 software using RDP Warcup training set v2 (Deshpande et al., 2016) database. Fungal genera with sintax cutoff > 0.7 were considered to discriminate taxa that significantly differentiate fungal communities. Relative abundances were calculated as the number of reads per taxon. OTUs were further classified into ecological guilds using FUNGuild, only considering assignments that were either "highly probable" or "probable" (Nguyen et al., 2016). Fungi were then assigned to saprotrophs (undefined saprotrophs and dung saprotrophs), symbiotrophs (endophytic and lichenized fungi) and plant pathogens. We also used publicly available fungal datasets generated through Illumina sequencing of the ITS2 region, that included samples of: (1) soil from pastures located nearby our sampling sites, which commonly contain C. cumulans mounds (de Oliveira et al., 2020); (2) the carton nest core of the sympatric non-food-storing termite Procornitermes araujoi (Syntermitinae) which builds mounds with outer walls made of thin layer of soil with an internal layer lined with fecal material (Supplementary Figure 1; Moreira et al., 2018); (3) different species of bark beetles (Rassati et al., 2019;Seibold et al., 2019); and (4) nests of a leaf-cutting bee (Rothman et al., 2019).

Diversity and Community Structure Analyses
We used R version 3.4.4 (R Core Team, 2018) to conduct analyses using different software packages. Alpha-diversity Alpha-diversity indices estimated for fungal ITS amplicons. Graphics report the observed operational taxonomic units (OTUs) richness, Chao1 estimator of richness, and Shannon's and Simpson's indices. Boxplots report median, upper and lower quartiles, and maximum and minimum values. Different letters denote significant differences (p-values < 0.05 in GLMM). Samples codes: CC1 = colony 1, CC2 = colony 2, and CC3 = colony 3. F, food nodules; N, carton nest core; W, worker gut; S, soldier gut. estimates were calculated using the function plot_richness in the phyloseq package (McMurdie and Holmes, 2013). Rarefaction was assessed using vegan package (Oksanen et al., 2018). Generalized linear mixed models (GLMM) were performed to check for overall significant differences of α-diversity estimates and the total number of OTUs among samples. The type of substrate was included as fixed effects, while the colony number was included as random effect. Multiple comparisons were conducted by general linear hypotheses using the glhtfunction in the multcomp package (Hothorn et al., 2008). A redundancy analysis (RDA) (999 permutations) of Hellingertransformed abundance data using the vegan package (Oksanen et al., 2018) was used to evaluate whether nest and gut substrates differ in fungal composition at the OTU level. The anova function from the package vegan were used to test the significance of RDA models and axes, respectively. We used Calypso web-server (Version 8.58), an online platform (Zakrzewski et al., 2017), to construct taxonomic heatmaps and the RDA of fungal metanalysis. The transmission of fungi between nest substrates and the gut of workers and soldiers was evaluated by analyzing the number of intersecting OTUs using the UpSetR package (Conway et al., 2017). We used linear discriminant analysis (LDA) effect size (LEfSe) to detect taxa that significantly differentiate fungal communities (Segata et al., 2011). To assess LEFSe, we used Galaxy web application, an online platform for the evaluation of multiple microbial community composition data (Afgan et al., 2016). The LDA score threshold was 4.0. Relative abundances less than 0.01% of the total reads were omitted from further analysis. We then identified variation in taxa using GLMM. To identify OTUs that were significantly different among samples, we used DESeq2 package (Love et al., 2014). Statistical significance was defined as p ≤ 0.05 in both statistical methods. Only OTUs identified by both LEfSe and DESeq2 were discussed. Plots were constructed with the packages ggplot2 (Wickham, 2016) and ggord (Beck, 2017).

DNA Extraction and Metagenomic Sequencing
The whole gut of 200 workers and soldiers and 50 mg of nest substrates (the carton nest core and food nodules) from colony CC1 were extracted in 2 ml tubes containing 1 mL of lysis buffer (500 mM NaCl, 50 mM Tris-HCl, pH 8.0, 50 mM EDTA, and 4% sodium dodecyl sulfate (SDS). Gut samples were obtained by anesthetizing termites on ice for 5 min and dissected with fine forceps. All the samples were kept at −20 • C until DNA extraction and purification using a bead-beating protocol (Yu and Morrison, 2004). DNA integrity was assessed by agarose gel electrophoresis (1.0 w/v) and the quantification was performed using NanoDrop spectrophotometer by measuring the absorbance at 260 nm. For metagenomic sequencing purposes, a library was constructed, using the Nextera library preparation kit (Illumina), according to the manufacturer's instructions. The prepared library was validated and quantified using the Agilent bioanalyzer 2100 system with a 12,000 DNA assay kit (Agilent) and Kapa Biosystems next-generation sequencing library qPCR kit (Kapa Biosystems), respectively. Sequencing was performed at LNBR NGS facility (CNPEM -Campinas, São Paulo, Brazil) using an Illumina HiSeq 2500 platform and applying the paired-end protocol (2 × 100-bp paired ends).

Metatranscriptome Library Preparation and Sequencing
Total RNA (10 µg) was extracted from 200 termite guts and 50 mg of the nest substrates of colony CC1 using Trizol reagent protocol (Invitrogen). Trizol/Chloroform step was performed twice. The total RNA was purified using the RNeasy Plant Mini Kit (Qiagen) following manufacturer's instructions. The quality of RNA was verified using RNAnano chip Bioanalyzer 2100 (Agilent). Good quality RNAs (RIN > 8.0) were submitted to rRNA depletion using the RiboZero rRNA Removal Kit under manufacturer instruction with a slight modification: at depletion step, we used a blend of rRNA removal solution from RiboZero Gold and RiboZero Bacteria kits (1:1) aiming to deplete both prokaryotic and eukaryotic rRNAs. Finally, the depleted RNAs were purified using Ampure XP beads, following the manufacturer's instructions, and kept at −80 • C until library preparation. A total of 50 ng of depleted RNA were used for library preparation using the TruSeq Stranded Total RNA Sample Preparation kit (Illumina), following the manufacturer's instructions. Quality control, quantification, and sequencing were performed as described before.

Metagenome and Metatranscriptome Data Analysis
Metagenome (MG) and metatranscriptome (MT) raw reads were quality checked using FastQC 1 and quality-filtered using Trimmomatic v.0.36 (Bolger et al., 2014). MT reads were also inspected using SortMeRNA to remove rRNA reads, and then both MG and MT reads were taxonomically classified using Kaiju (Menzel et al., 2016). Next, the whole set of MG trimmed 1 https://www.bioinformatics.babraham.ac.uk/projects/fastqc/ reads were de novo co-assembled using IDBA_UD (version 1.1.1) (Peng et al., 2012). Gene prediction and annotation of was performed using Prokka v.1.11 with the meta parameter (Seemann, 2014) and carbohydrate-active enzymes (CAZy) annotation was performed using dbCAN database (Zhang et al., 2018). Taxonomy of predicted protein coding genes were inferred according the easy-taxonomy function of MMseqs2 suite (Steinegger and Söding, 2017) to identify putative bacterial, fungal, and termite-host genes by LCA method. Both, MG and MT reads were mapped to the complete set of predicted genes recovered from the MG de novo assembly using Kallisto v. 0.46.1 (Bray et al., 2016) to estimate genes and transcripts abundance. Normalized expression of protein coding genes was estimated with the Transcripts Per Kilobase Million (TPM) method and used to analyze the abundance of expression of the predicted fungal, bacterial, and host putative CAZy genes. Homology to peptide pattern (Hotpep) was used to infer CAZy function (Busk et al., 2017).

Which Fungal Taxa Are Associated to the Gut and the Nest of Cornitermes cumulans?
We detected 968 fungal OTUs (2,163,723 sequences from the ITS region) associated to the gut and the nests of C. cumulans. Of these, 579 OTUs (1,370,000 sequences, 63.3% of the total) were annotated to the food nodules, 774 OTUs (739,691 sequences, 34.2%) to the carton nest core and 197 OTUs (54,032 sequences, 2.5%) to the gut of workers and soldiers. Deep taxonomic analyses using the Warcup database successfully assigned 88 and 84% of fungal reads at the order and family level, respectively. Overall, nest and gut samples yielded the phyla Ascomycota (99.9% of sequence reads) and Basidiomycota, representing 13 classes, 42 orders, 98 families, and 212 genera (Supplementary Tables 1, 2). Rarefaction curves indicated adequate sampling of fungi for a valid comparison among fungal communities (Supplementary Figure 2). The average number of OTUs was significantly higher in nest samples, the majority classified as undefined saprotrophs (on average 76.5% of sequence reads) and symbiotrophs (endophytes 4.8% and lichenized 2.5%) of the orders Pleosporales, Sordariales, and Xylariales (Ascomycota) using FUNGuild (Figures 1B,C and Supplementary Figure 3).

Fungal Composition Differs Between Gut and Nest Substrates
Observed richness was significantly lower in gut assemblages ( Figure 1D). On the other hand, food nodules and carton core samples of C. cumulans nests did not differ in richness and shared approximately 99% of their sequence reads (293 OTUs), but estimated diversity was significantly lower in food nodules. Constrained ordination showed that community composition of food nodules formed a distinct cluster separated from carton nest core (F = 2.84, df = 1, p = 0.002) and gut samples (F = 6.61, df = 1, p = 0.001) (Figures 2A,B) and (C) Linear discriminant analysis (LDA) combined with effect size measurements (LEfSe, score > 4) and DESeq2 (p ≤ 0.05) analyses to identify the fungal OTUs responsible for the differences between nest substrates. Boxplots report median, upper and lower quartiles, and maximum and minimum values of the relative abundance of taxa identified by both DESeq2 and LEfSe analyses. A p-value < 0.05 was considered significant in GLMM analyses. See Figure 1 for samples codes.

Fungal Transmission Between Nest Substrates and Termites
The gut mycobiota of workers and soldiers was represented by 11 classes, 26 orders, 52 families, 90 genera, and 197 fungal OTUs (Supplementary Table 1) with the prevalence of the genus Phoma (Pleosporales) (Supplementary Tables 2, 3). Diversity and composition of fungal gut assemblages were not influenced by the type of caste (RDA, F = 0.84; df = 1; p = 0.540) ( Figure 2B) and no features with significant differences were found between workers and soldiers. However, we detected a high percentage (between 89 and 99%) of transmission of fungal OTUs between nest substrates and the gut of workers and soldiers. The majority of shared OTUs belongs to Sordariales and Pleosporales. Additionally, workers and soldiers shared between 97.5 and 98.4% of their sequence reads (Supplementary Figure 4).

Termites Nest Exhibits Different Fungal Assemblages
Analyses of the variation in community structure showed that the composition of Ascomycota fungi differed significantly among termite nests (C. cumulans and P. araujoi), soil samples and substrates associated with two phylogenetically distant insects (RDA, F = 15.47; df = 5; p = 0.001) (Supplementary Figure 5A). The nest substrates of termites formed a cluster, clearly separated from other fungal assemblages. The more prominent genera associated to termite nest substrates, like Daldinia (Xylareaceae) and Lecythophora (Coniochaetaceae), were found in lower abundances in other assemblages. However, Thielavia and Chaetomium (Chaetomiaceae) were similarly abundant in both soil and termite nests. We also found that fungal communities associated to bark beetles and bee nests formed separated clusters with their prominent genera present in lower abundances in termite nests (Supplementary Figure 5B).

Bacteria Dominated Microbial Assemblages in Gut and Nest Substrates
The microbial assemblages from colony CC1 were further investigated by metagenomics (MG) and metatranscriptomics (MT) sequencing. Taxonomy analyses of the sequences showed a higher abundance of bacteria over fungi in both nest and gut samples (Figures 3A,B). Actinobacteria and Proteobacteria were dominant in the nest substrates, whereas Firmicutes and Spirochetes were prevalent in the termite gut, confirming our previous results with 16S rRNA sequencing . Ascomycota was present in both nest and gut samples. However, a fungus of genus Neocallimastix (Chytridiomycota: Neocallimastigales) was only found in the gut of workers (F) Detail of the food storage behavior of C. cumulans workers on the nest walls using blue-marked sugar-cane (red arrow). and soldiers. Bacillus (Firmicutes), Treponema (Spirochetes), Streptomyces (Actinobacteria), and Aspergillus (Ascomycota) were the most abundant microbial genera in the samples (Supplementary Tables 6, 7 and Supplementary Figure 6A). We did not obtain RNA from nest core samples; however, the MG showed a similar profile between the carton nest and food nodules ( Figure 3A). The MT sequences from fungi were less abundant (Figure 3B), suggesting a lower metabolic activity from fungi over bacteria. Although our data was derived from a single colony (CC1), the fungal taxa responsible for the expression of CAZy transcripts showed similar abundances among the three colonies used in ITS sequencing, for food nodules (H = 0.74; df = 2; p = 0.689; Kruskal-Wallis test) and the gut of workers and soldiers (H = 1.11; df = 5; p = 0.953) (Supplementary Table 2).

Gut Symbiotic and Host Contribution to Lignocellulose Digestion
The gut microbiota of C. cumulans yielded 88,454 bacterial and fungal 3,370 genes. Of these, 121 were assigned to bacterial and 34 to fungal CAZy. GH5_2 and GH5_4 genes are the highest expressed GH families associated to the gut microbiome. Additionally, the non-catalytic carbohydratebinding module gene (CBM22) is also abundant. Although the functional lignocellulolytic profile overlapped partially with food nodules (Figure 3D), it was linked to different microbial taxa. We also detected a high expression of transcripts associated with endoglucanase (EC 3.2.1.4), xylanase (EC 3.2.1.8), and endomannanase (EC 3.1.2.78) which were linked to Firmicutes, Chytridiomycota and Ascomycota, respectively (Supplementary Figure 6A). On the other hand, we identified 2,262 genes in the termite C. cumulans, but only 21 CAZy gene transcripts. A starch phosphorylase (EC 2.4.1.1), an alcohol-2-dehydrogenase  and the carbohydrate esterases (CE4 and CE10) were predicted in the MG assembly (Supplementary Table 8 and Figure 3D). The total CAZy expression levels associated to cellulose and xylan degradation were higher in the gut of termites, and nearly equally expressed in workers and soldiers ( Figure 3E).

DISCUSSION
Termites rely on a mutualistic association with gut symbionts to decompose the plant biomass; however, the digestive efficiency is achieved by the functional complementation with host endogenous lignocellulases (Tokuda and Watanabe, 2007;Scharf, 2008;Scharf et al., 2011). Moreover, termites also interact with the surrounding soil microbiota and at least two extant lineages (Macrotermitinae and Sphaerotermitinae) engaged in an external association with nest microorganisms, where each component of the termite holobiont has a collaborative role in the decomposition of lignocellulose (Poulsen et al., 2014;Chouvenc et al., 2021).
In this study, we present evidences of complementary digestion by external nest microorganisms in a neotropical grassharvester termite C. cumulans (Syntermitinae) that store plant material inside its nests ( Figure 3F). Our findings revealed that plant biomass (previously cut and packed in the food nodules) is initially decomposed by enzymes that target the lignin and complex polysaccharides as cellulose, xylan, arabinoxylan, and starch. We hypothesized that pre-processing of the plant material reduces the recalcitrance of the lignocellulose before gut passage. Since lignin breakdown by peroxidases needs O 2 (Breznak and Brune, 1994), food nodules are a suitable environment for this process. Peroxidases and dehydrogenases of bacterial origin are known to participate in lignin degradation (Pelmont et al., 1989;Tamboli et al., 2011;Chauhan, 2020). The catalase-peroxidase (EC 1.11.1.21) assigned to Actinobacteria in this study was able to oxidase lignin when incubated with lignocellulosic material (de Gonzalo et al., 2016). The final step in the lignocellulose digestion process involves the passage of the pre-processed plant material in the gut of workers where it is again targeted by xylanases, possibly for overcoming any xylan residuals and making the cellulose more accessible to the endoglucanases produced by the gut microbiota. Although we found that endoglucanase transcripts were expressed in the food nodules, our results showed higher expression of cellulases in the gut of termites.
Other CAZy predicted in the gut were the non-catalytic carbohydrate-binding module CBM22 that promotes the association of xylanases and the carbohydrate esterases of C. cumulans (CE4 and CE10) (Sermsathanaswadi et al., 2017) that could be acting in the breakdown of hemicelluloses (Franco Cairo et al., 2016). The presence of chitinases was predicted in the food nodules and linked to bacteria of the phylum Firmicutes (EC 3.2.1.14) and to an ascomycete fungus (EC 3.2.1.52). Chitinases could be related to fungi cell-wall growth or the activity of mycolytic bacteria strains (Solanki et al., 2012) found in the food nodules. Chitinase transcripts were not obtained in the gut of C. cumulans; however, other higher termites are known to produce chitinases which can improve the digestibility of fungal substrates and decaying wood (Poulsen et al., 2014;Hu et al., 2019) or increase the permeability of peritrophic matrix in the midgut for the absorption of nutrients.
Although there is functional overlap between the food nodules and gut, CAZy were linked to different taxa. It could be explained by the different conditions that make these places more suitable to different kind of microorganisms. In the food nodules, cellulose and hemicellulose digestion were linked to Proteobacteria and Ascomycota, whereas catalase activity was assigned to Actinobacteria. On the other hand, lignocellulases in the gut of C. cumulans were linked to Firmicutes, Ascomycota, and Chytridiomycota. All these microbes dominated the core microbiota and are known to produce cellulases, hemicellulases (Couturier et al., 2016;Grieco et al., 2019;Matsuzawa et al., 2019Matsuzawa et al., , 2020 and laccases (Murphy et al., 2021). In fungus-growing insects that exploit plant-derived resources, the association between fungi and bacteria result in a multipartite functional metabolism of lignocellulose (Barcoto et al., 2020). The ecological relevance of such microbial interaction in the food nodules of C. cumulans needs to be determined.
Because our data was derived from a single colony, we cannot exclude the possibility of intercolonial variation that may affected the expression levels. Indeed, the relative abundances of some fungi orders of food nodules showed variation among the colonies used for ITS sequencing. However, the taxa responsible for CAZy expression showed similar abundances for bacteria  and fungi (this study).
Fungi were less abundant than bacteria in all the samples and their CAZy expression levels associated to cellulose and xylan were also lower. Nevertheless, a noteworthy observation in this study was that the termite gut of C. cumulans is inhabited by fungi that expressed CAZy gene transcripts. One gene transcript, assigned to GH5_4 family and classified as endoglucanase EC 3.2.1.4 was linked to Neocallimastix (Chytridiomycota: Neocallimastigales). This finding suggests this termite possesses a functional mycobiota in its gut that contributes to the lignocellulose digestion. Neocallimastigales comprises anaerobic highly active lignocellulolytic fungi of the gastrointestinal tract of herbivores, which has also been reported in termites (Lee et al., 2015;da Silva et al., 2017).
Most of the fungi associated to the nests of C. cumulans (Syntermitinae) corresponded to saprotrophic Ascomycota, but other guilds within this same fungal phylum are also present in this microenvironment. Like termites, saprotrophic fungi are important plant-litter decomposers in terrestrial ecosystems using lignocellulolytic enzymes (Crowther et al., 2012) and their incidence inside termite nests could be explained by the addition of feces to the internal walls (Brauman, 2000), creating an appropriate environment for the survival and reproduction of these microorganisms (Holt, 1998;Vesala et al., 2019). Lichenized and endophytic fungi were also detected in the nest substrates but in lower abundances. The presence of these fungi has been reported in dead plant material probably consuming organic compounds released by saprotrophic microbes acting on the decomposition of the plant matter (Ottosson et al., 2015), and the mechanisms by which they colonize termite nests are unknown, but they were possibly transported from foraging parties to the mounds by workers. We found a distinct fungal composition of food nodules in comparison to the carton nest core, characterized by the abundance of saprotrophic Xylariales and Sordariales. Members of these fungal orders are known to colonize several substrates, including leaf litter, vertebrate dung, and the interior of termite nests (Okane et al., 2008;Guedegbe et al., 2009). In fungus-growing termites, Xylariales appear to be opportunistic saprotrophs of the stored plant material and fungus combs either due to the failure of microbe suppression mechanisms or after the nests have been abandoned by the termites (Mueller and Gerardo, 2002;Visser et al., 2011).
Approximately 90% of fungal reads were shared between nest fecal material and the gut of C. cumulans. One possible explanation for this is that termites consume the fecal material of their nests, independently of the food storing behavior. Carton feeding has been reported in other termites , but the precise function of this behavior is not completely understood. Fungal colonization of lignocellulosic substrates is known to enrich its nitrogen content (Jacobson et al., 2015) and because most wood-and litter-feeding termites have a diet poor in nitrogen (Tayasu et al., 1994), consumption of nest carton could also help increase nitrogen intake by the termites. An alternative hypothesis is that nest fungi were originally gut symbionts deposited by the workers together with their feces and saliva on the nest walls. Since we found small differences in the number and nature of fungi of gut samples between colonies within each termite species, it is possible that the fungus is part of the natural gut symbionts and could function in their nests as a supplementary reserve of lignocellulolytic microbiota for these termites (Nalepa et al., 2001;Fall et al., 2004). The large number of shared microbes and CAZy transcripts between sterile castes is consistent with the transmission of microorganisms among nestmates. Soldiers, larvae, and reproductive castes cannot feed themselves and depend on workers for nutrition through trophallaxis (Machida et al., 2001).
In wood-feeding termites, fungi are known to play facultative interactions by attracting individuals to the food source or by improving colony development; however, the effect of fungi may also depend on the termite species (Birkemoe et al., 2018). For example, it is known that the wood decayed by the brown rot fungus Gloeophyllum trabeum attracts the termites Reticulitermes virginicus, R. flavipes (Rhinotermitidae), and Nasutitermes columbicus (Termitidae) (Esenther et al., 1961;Smythe et al., 1967) because of the presence of a particular alcohol in this decayed wood, (3Z,6Z,8E)-dodeca-3,6,8-trien-1-ol (Matsumura et al., 1968(Matsumura et al., , 1969, known to be the only compound or one of the compounds of the trail pheromones of these termite species (Sillam-Dussès, 2010). Besides the potential benefits of nest fungi, we cannot discard the possibility of fungal entomopathogens or commensals to be present in the termite nest. Some representatives of Hypocreales found in our samples could be potential entomopathogens (Araújo and Hughes, 2016), but they were found in low frequency because they may be suppressed by the abundant Actinobacteria (Chouvenc et al., 2018) or by termite gland secretions (Rosengaus et al., 2004) or these fungi are just dormant in this environment (Hajek, 2017).
In this work we improved ITS sequence analyses to bring a more comprehensive understanding of the fungal communities associated to the nest, gut, and food nodules stored by C. cumulans and we bring additional metagenomic and metatranscriptomic analyses of the microbial communities associated with this termite species. Taken together with our previous analyses of the 16S rRNA bacterial communities and the differential lignocellulolytic activity in food nodules , we have provided important new evidence suggesting complementary digestion of the plant material inside food storage structures before the passage through the gut of a Neotropical higher termite species. We have shown that fungi are also CAZy contributors in the termite holobiont. This could explain why the survival of this termite increased significantly when workers were fed with food nodules compared to other plant feedstocks (Janei et al., 2020). From the perspective of termite evolution, this study brings new information about the multicomponent nature of the holobiont in higher termites. The complementary digestive strategy of C. cumulans may resemble the one linking the fungus-growing termites and their mutualistic fungus Termitomyces (Poulsen et al., 2014). The association with nest microorganisms is considered as an evolutionary novelty in higher termites that originated millions of years ago in two sister lineages (Macrotermitinae and Sphaerotermitinae) after transition of gut symbionts in the ancestor of Termitidae (Bucek et al., 2019). However, our findings in C. cumulans supports the hypotheses that externalization of the digestion was an ancestral condition of higher termites that evolved through the incorporation of soil lignocellulolytic microbes into the fecal nest material (Chouvenc et al., 2021). Consequently, this association led to the loss of gut flagellates and the acquisition of novel bacterial gut symbionts before the diversification of Termitidae. Thus, we could speculate that ectosymbiosis with nest microorganisms is a plesiomorphic characteristic of higher termites. The externalization of digestion is also present in phylogenetically related species as Syntermes dirus (unpublished results) (Syntermitinae) and Velocitermes heteropterus (Nasutitermitinae) , emphasizing the need to further investigate the ecological relevance of nest microbiota in Termitidae.

DATA AVAILABILITY STATEMENT
The datasets analyzed for this study can be found in the ENA accession no. PRJEB17080, https://www.ebi.ac.uk/ena/browser/ view/PRJEB17080.

AUTHOR CONTRIBUTIONS
AA, GP, AC-L, and FS designed the study. EM, LM, TA, and DP performed the experiments. AA and GP analyzed data. AA, GP, JC, AR, and DS-D contributed to new analyses and wrote the manuscript. All authors contributed to the article and approved the submitted version.