Ramayo-Caldas et al. Animal Microbiome https://doi.org/10.1186/s42523-020-00038-4 (2020) 2:18 Animal Microbiome RESEARCH ARTICLE Open Access Gut eukaryotic communities in pigs: diversity, composition and host genetics contribution Yuliaxis Ramayo-Caldas1* , Francesc Prenafeta-Boldú2, Laura M. Zingaretti3, Olga Gonzalez-Rodriguez1, Antoni Dalmau4, Raquel Quintanilla1 and Maria Ballester1 Abstract Background: The pig gut microbiome harbors thousands of species of archaea, bacteria, viruses and eukaryotes such as protists and fungi. However, since the majority of published studies have been focused on prokaryotes, little is known about the diversity, host-genetic control, and contributions to host performance of the gut eukaryotic counterparts. Here we report the first study that aims at characterizing the diversity and composition of gut commensal eukaryotes in pigs, exploring their putative control by host genetics, and analyzing their association with piglets body weight. Results: Fungi and protists from the faeces of 514 healthy Duroc pigs of two sexes and two different ages were characterized by 18S and ITS ribosomal RNA gene sequencing. The pig gut mycobiota was dominated by yeasts, with a high prevalence and abundance of Kazachstania spp. Regarding protists, representatives of four genera (Blastocystis, Neobalantidium, Tetratrichomonas and Trichomitus) were predominant in more than the 80% of the pigs. Heritabilities for the diversity and abundance of gut eukaryotic communities were estimated with the subset of 60d aged piglets (N = 390). The heritabilities of α-diversity and of the abundance of fungal and protists genera were low, ranging from 0.15 to 0.28. A genome wide association study reported genetic variants related to the fungal α-diversity and to the abundance of Blastocystis spp. Annotated candidate genes were mainly associated with immunity, gut homeostasis and metabolic processes. Additionally, we explored the association of gut commensal eukaryotes with piglet body weight. Our results pointed to a positive contribution of fungi from the Kazachstania genus, while protists displayed both positive (Blastocystis and Entamoeba) and negative (Trichomitus) associations with piglet body weight. Conclusions: Our results point towards a minor and taxa specific genetic control over the diversity and composition of the pig gut eukaryotic communities. Moreover, we provide evidences of the associations between piglets’ body weight after weaning and members from the gut fungal and protist eukaryote community. Overall, this study highlights the relevance of considering, along with that of bacteria, the contribution of the gut eukaryote communities to better understand host-microbiome association and their role on pig performance, welfare and health. Keywords: Gut eukaryotes, Fungi, Protists, Pig, Host genetics * Correspondence: yuliaxis.ramayo@irta.cat 1 Animal Breeding and Genetics Program, Institute of Agrifood Research and Technology (IRTA), Torre Marimon, 08140 Caldes de Montbui, Spain Full list of author information is available at the end of the article © The Author(s). 2020 Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. Ramayo-Caldas et al. Animal Microbiome (2020) 2:18 Page 2 of 12 Introduction fungal and 227 protist ASVs were identified. The pig gut The gut microbiome harbors thousands of species of ar- mycobiota was dominated by yeast from the Kazachsta- chaea, bacteria, viruses and eukaryotes such as protists nia genus, in particular by the species K. slooffiae and K. and fungi that contribute to host biology. The gut bovina (Fig. 1). Other yeasts associated to the Candida eukaryotic communities of monogastric species show a genus were identified, such as C. glabrata, C. albicans higher interindividual variability, but less abundance, di- and C. (Diutinia) catenulata, which were present in 2– versity and richness than their bacterial counterparts [1]. 7% of the animals (Additional file 1 Table S1). Besides, a The so called gut mycobiome of healthy humans is dom- moderate prevalence (15% of pigs) of Tilletia puccinel- inated by fungal genera like Saccharomyces, Candida, liae was also observed (Additional file 1 Table S1). Penicillim, Aspergillus and Malassezia [2], while protist Along with the ascomycete Saccharomyces arboricola, a genera such as Blastocystis, Entamoeba and Enteromonas number of basidiomycetous yeasts were also detected have been reported in the human gut across different with relatively low incidence (< 2%), such as Sporobolo- geographical location [3–5]. There is currently consider- myces roseus, Trichosporon dohaense, Debaryomyces pro- able interest in understanding the mechanisms through sopidis, Pichia sporocuriosa, Filobasidium globisporum which gut commensal eukaryotic communities contrib- and Vishniacozyma victoriae. Finally, the pig intestinal ute to host homeostasis and health [6]. The role of fungi mycobiota included a third and more diverse group of and protists in the human gut remains poorly under- cosmopolitan fungal species that are generally catego- stood, but some authors suggest that commensal species rized as “soil fungi”. In most cases they have a relatively may induce innate immune response in the host [7], and low occurrence in gut (< 2%) and are primarily associ- they could also have other potential benefits [5, 8, 9]. ated with a general saprotrophic ecophisiology (i.e. As- For example, Blastocystis spp. and non-pathogenic Ent- pergillus spp., Aureobasidium pullulans, Cladosporium amoeba spp. have been associated with a healthy and tenuissimum, Mucor circinelloides and Penicillium highly diverse gut microbiome ecosystem [6, 10, 11], polonicum). while Tritrichomonas musculis modulates the intestinal With regard to protists, the prevalent ASV corre- immune system and increases host protection against sponded to an unidentified species of the subclass Tri- mucosal bacterial infections in mice [12]. chostomatia (superphylum Alveolata in the SAR The few studies about gut eukaryotes conducted in supergroup). Furthermore, the species Neobalantidium pigs (Sus scrofa domesticus) used a limited sample size coli from this class was found in 82% of the animals and mainly focused on isolated members of the (Additional file 1 Table S1). Likewise, with a relatively eukaryote communities or parasite identification under high incidence and average abundance (up to 99 and pathogenic conditions [13–18]. In swine, Kazachstania 29%, respectively), a number of Blastocystis spp. (super- spp. and members of the Saccharomycetaceae family are phylum Stramenopiles in the SAR supergroup). Other predominant in the gut mycobiota [13–16], while gut detected commensal trichomonads included Tetratricho- protist community of healthy pigs is dominated by Blas- monas buttreyi (94% incidence), Hypotrichomonas imi- tocystis spp., Tritrichomonas spp. and Balantidium coli. tans (67% incidence), as well as the Tetratrichomonas Environmental factors such as diet seems to play a key strain PEKPR (10% incidence), and with a very low inci- role in modulating the structure and composition of dence and abundance (< 1%), Tritrichomonas suis. Fi- both eukaryote and prokaryote communities. However, nally, one representative of the ameboid protist species little is known about the genetic control of the gut eu- Entamoeba gingivalis (supergroup Amoebozoa) was also karyotes, since published studies have been focused on detected with a relatively high incidence (77%) but a low the role of host genetics in shaping the gut bacterial abundance (< 1%). communities [19–21]. The main goal of this study is to characterize com- Host and environmental factors modulating the diversity mensal fungi and protists inhabiting the gut of healthy of pig gut microbial eukaryotes pigs, explore the putative host genetic control over di- The modulatory effect of host and environmental factors versity and composition of gut eukaryotes communities, over the diversity and composition of gut eukaryotic and analyze their association with piglets body weight. communities in pigs was evaluated with different ap- proaches. Results from the first PERMANOVA with the Results whole dataset indicated that the combination of farm Composition of the pig gut microbial eukaryotic and animal age represented the most significant effect communities shaping the gut eukaryotic communities (p < 0.0001), The pig gut microbial eukaryote communities from 514 explaining around 44 and 42% of the total variability of healthy pigs were analyzed by sequencing the 18S rRNA fungal and protist communities, respectively. The same gene and the ITS2 region. After quality control, 492 data structure was recovered by the principal Ramayo-Caldas et al. Animal Microbiome (2020) 2:18 Page 3 of 12 Fig. 1 Circular barplot representing the relative abundance of the most abundant fungal Genus (a) and species (b) as well as protists Genus (c) and species (d) coordinates analysis (PCoA) based on the Bray–Curtis (Table 1). Regarding sex effects, gut fungi alpha-diversity distance, which showed two clusters that match with differed between sexes in weaned piglets (p-value = 0.02), sample farms origin (Additional file 2 Figure S1). The di- whereas no differences between castrated males and fe- versity index also revealed important differences be- males mycobiota diversity were observed at 190 d of age tween farms/ages (Kruskal–Wallis test). Samples taken (Table 1). at 190 days of age in the experimental farm showed a Afterwards, we assessed the degree of host genetic con- significantly higher protist (p = 0.013) and fungal (p = trol over the diversity and composition of gut protist and 0.006) alpha-diversity (Additional file 2 Figure S1). In fungal communities by estimating their heritabilities. The contrast, weaned piglets raised in commercial conditions estimated heritability of the α-diversity index was low, showed higher protist beta-diversity. A second PERMA- with posterior mean estimates of ranging from 0.161 to NOVA within farm/age was performed. Results of 0.188 for fungi and protists, respectively (Table 2). The weaned piglets in the commercial farm indicated an im- posterior highest density regions at 95% encompassed portant effect of the batch (p-value< 0.0001) on alpha- values between 0.067 and 0.299, indicating a probably very diversity of both fungal and protist communities, limited but not negligible heritability for these traits. In whereas the pen effect in the second dataset (experimen- the same way, the heritabilities of most representative fun- tal farm) seems to affect exclusively protist diversity gal and protist genera abundance reached low to moderate Table 1 Resume of the permutational multivariate analysis of variance Dataset Community Bacht R2 (P-value) Sex R2 (P-value) Weaned piglets (n = 405) Fungi 0.12 (p < 0.001) 0.01 (p = 0.02) commercial farm Protist 0.07 (p < 0.001) 0.02 (p = 0.60) Finishing pigs (n = 109) Fungi 0.09 (p = 0.35) 0.01 (p = 0.24) experimental farm Protist 0.12 (p = 0.03) 0.01 (p = 0.55) Ramayo-Caldas et al. Animal Microbiome (2020) 2:18 Page 4 of 12 Table 2 Mean (standard deviation) and highest density regions at 95% of marginal posterior distributions for heritabilities of pig gut eukaryotes microbial traits Community Level Trait Heritabiliy (SD) HDR Fungi α-diversity Shannon-index 0.161 (0.053) 0.067–0.266 Genus Kazachstania 0.209 (0.064) 0.095–0.333 Specie Kazachstania_slooffiae 0.216 (0.071) 0.093–0.353 Protist Alpha-diversity Shannon-index 0.188 (0.060) 0.077–0.299 Genus Blastocystis 0.281 (0.101) 0.096–0.483 Trichomitus 0.233 (0.077) 0.085–0.383 Tetratrichomonas 0.159 (0.051) 0.070–0.258 Neobalantidium 0.273 (0.090) 0.108–0.443 Species Trichomitus batrachorum 0.240 (0.075) 0.100–0.382 Tetratrichomonas buttrey 0.158 (0.048) 0.070–0.256 Neobalantidium coli 0.256 (0.086) 0.100–0.425 Blastocystis sp. ATCC 50177 0.245 (0.086) 0.086–0.409 Blastocystis sp. CK86–1 0.203 (0.066) 0.081–0.325 Blastocystis sp. SY94–3 0.205 (0.065) 0.085–0.334 (Table 2). The abundance of the protist Blastocystis genus synonymous, and 4% mapped within non-coding transcript was the most heritable (posterior mean 0.281, SD 0.101), variants (Additional file 3 Table S2). Most of these associ- followed by the abundance of Trichomitus (0.273, SD ated SNPs (164 out of 174) were identified on SSC6 (Table 0.090) and Neobalantidium (0.233, SD 0.077) genera. 3), and were mainly associated with the relative abundance Also, a heritability around 0.209 was estimated for the of two species of Blastocystis genera: CK86–1 and ATCC abundance of the preponderant fungal genus Kazachsta- 50177 (Fig. 2). The remaining significant SNPs located on nia, and lightly larger (0.216) for the abundance of SSC17 (46.77–46.99 Mb interval) and SSC18 (two intervals: Kazachstania slooffiae species. 4,53–4,57 Mb and 25,85–25,88 Mb) were also associated to the relative abundance of members of Blastocystis (Table Identification of host genetic regions linked with the 3). Regarding diversity, the aforementioned 141.91–145.39 composition and diversity of pig gut microbial eukaryotes Mb region of SSC6 resulted also associated with the fungi Results from GWAS revealed few and weak association sig- Shannon diversity-index (Fig. 2). Finally, no significant asso- nals between the host genome and the gut microbial eu- ciations with either the relative abundance of fungi karyotes composition and diversity. We identified a total of (Kazachstania genera or Kazachstania slooffiae) or the di- 174 SNPs as significantly (at chromosome-wide level, versity of protist were observed. FDR < 0.05) related with the abundance of protist commu- nity (Additional file 3 Table S2), located in seven intervals Genes and pathways associated to pig gut microbial distributed across three pig Sus scrofa chromosomes: SSC6, eukaryotes SSC17 and SSC18 (Table 3). The 32.18% of these SNPs A total of 229 protein-coding, 74 long non-coding RNA were intronic variants, 56.89% were intergenic, 5.75% were and two miRNAs (mir138–2 and ssc-mir-186) (Add- located upstream/downstream of genes, 1.15% were exonic itional file 4 Table S3) were annotated within the host Table 3 Resume of the chromosomal intervals identified by the GWAS Chr Start (bp) End (bp) Size (Mb) Associated microbial trait 6 68,999,679 69,466,937 0.467 Blastocystis sp. SY94–3 6 120,437,893 121,498,044 1.060 Blastocystis sp. ATCC 50177; Blastocystis sp. CK86–1 6 122,778,463 132,592,988 9.814 Blastocystis sp. ATCC 50177; Blastocystis sp. CK86–1 6 141,914,558 145,392,721 3.478 Blastocystis sp. ATCC 50177; Blastocystis sp. CK86–1; ITS α-diversity 17 46,764,111 46,988,559 0.224 Blastocystis sp. 18 4,538,602 4,570,818 0.032 Blastocystis sp. ATCC 50177 18 25,850,358 25,885,196 0.035 Blastocystis sp. ATCC 50177 Ramayo-Caldas et al. Animal Microbiome (2020) 2:18 Page 5 of 12 Fig. 2 Significance across the whole genome of the association between SNP markers and the relative abundance of (a) Blastocystis sp. ATCC 50177, (b) Blastocystis sp. CK8166 and (c) Fungi Alpha diversity. The x-axis represents the pig autosomal chromosomes (1–18), and the y-axis the significance level represented as the –log10 (P-value) genomic intervals associated to the gut microbial eu- CDH5, LPAR3, PIK3C3, CCL17, CCL22, CDH5, CKLF karyotes composition and diversity. A detailed explor- and MMP15). Finally, the functional annotation revealed ation of annotated genes revealed that the Interleukin 23 that these genes belong to a variety of physiological pro- Receptor gene (IL23R) holds two intronic SNPs, cesses and gene networks related to Cell-To-Cell Signal- rs343576943 (p-value = 1,8*10− 5) and rs81392455 (p- ing and Interaction, Cellular Assembly and value = 1,8*10− 5) that were associated with fungal α- Organization, Cellular Function and Maintenance, De- diversity. Furthermore, the variants rs81391297, velopmental Disorder, Hereditary Disorder, Cell- rs81391299 and rs81299019, mapped within the 122.8– mediated Immune Response, Humoral Immune Re- 132.6 Mb interval on SSC6 and associated with the rela- sponse, Immune Cell Trafficking, and Protein Synthesis tive abundance of Blastocystis spp., were located on in- (Additional file 5 Table S4). Metabolic pathways most tron eight of the Phosphatidylinositol 3-Kinase Catalytic significantly enriched by the list of candidate genes in- Subunit Type 3 (PIK3C3) gene. On SSC17, the intronic clude Relaxin Signaling, Leptin Signaling in Obesity, variant rs331945396 (SSC17:46832698–46,832,698) in IGF-1 Signaling, PXR/RXR Activation and HIF1α Signal- the Hepatocyte Nuclear Factor 4 Alpha (HNF4A) gene ing pathways. It should be noted the overrepresentation was also associated with the relative abundance of Blas- of immune-related pathways such as Gα12/13 Signaling, tocystis genera (p-value = 9,02*10− 5). Other genes related IL-9 Signaling, IL-23 Signaling, IL-12 Signaling, ERK/ to the immune system were annotated within the afore- MAPK Signaling, Agranulocyte and Granulocyte Adhe- mentioned chromosomal intervals, as for example the sion and Diapedesis (Additional file 6 Table S5). TNF Receptor Superfamily Member 9 (TNFRSF9), the Interleukin 12 Receptor Subunit Beta 2 (IL12RB2) or the Association between pig gut eukaryotes and host Phosphatidylinositol-4,5-Bisphosphate 3-Kinase Catalytic performance Subunit Delta (PIK3CD), as well as other members of Most relevant fungal and protist genus associated with immune-related pathways (e.g. CDH1, CDH11, CDH3, body weight are summarized in Fig. 3. The result Ramayo-Caldas et al. Animal Microbiome (2020) 2:18 Page 6 of 12 suggested association between piglets’ body weight and found in all studied animals across farms, sex and ages, the abundance of Kazachstania fungal genera as well as observation that is in agreement with the ubiquitous de- protist members of Entamoeba, Trichomitus, Tetratri- tection of this yeast from the gut of healthy pigs, and chomonas and Blastocystis. The Shapley Additive Ex- further supports the hypothesis that pigs intestinal tract planation (SHAP) values showed a nonlinear association is the primary ecological niche for K. slooffiae [13]. The between Kazachstania and Blastocystis abundance and closely related K. bovina was detected in the 45% of the piglets body weight (Additional file 7 Figure S2). Fur- studied pigs, and might also play a significant role in the thermore, Entamoeba on one side, and Trichomitus on gut. Several transient fungi that probably arise from the the other, were positive and negatively associated with feeding or the environment were observed but with low body weight (Additional file 7 Figure S2). incidence. These basidiomycetous yeasts have often been isolated from the phyllosphere [26] and are also known Discussion to cause human opportunistic infections in certain cases The analysis of gut microbial eukaryotic communities in [27]. The plant-gut association for certain yeast is more 514 healthy animals have allowed characterizing com- evident in a second group of filamentous fungi that are mensal fungi and protists inhabiting the porcine gut characterized as specialized plant biotrophs. The most tract. In agreement with previous reports in pigs [22–24] representative of those is the basidiomycete T. puccinel- our results evidenced that pigs gut mycobiota is domi- liae that, surprisingly, occurred in 15% of the studied nated by the yeast Kazachstania spp. which corresponds pigs. Species from Tilletia genera are smut fungi that in- to the teleomorphic state of Candida slooffiae and Can- fect various grasses from the Poaceae family and encom- dida bovina [25]. Particularly the species K. slooffiae was pass plant pathogens of economic importance in the Fig. 3 Relevance of fungal and protist genera regarding piglets body weight, according to the Shapley Additive Explanation (SHAP) values. The x- axis represents the SHAP value, the y-axis shows the genera names in order of importance from top to bottom, and the value correspond to the genera mean SHAP value Ramayo-Caldas et al. Animal Microbiome (2020) 2:18 Page 7 of 12 production of cereals and forage grasses. So far, T. pucci- but not negligible heritabilities were obtained for a par- nelliae has only been isolated from weeping alkaligrass ticularly complex trait such as the α-diversity index in (Puccinellia distans), a common ruderal grass in Europe both protist and fungi, which allows hypothesizing that and North America [28]. Plants infected by other Tilletia pigs can have a limited genetic tendency to harbor more species do not pose a toxicity risk for humans but con- or less diverse eukaryotic communities in its gut. These taminated grains might be derived to animal feed due to results should be taken with caution because of the lim- off-flavors [29]. Such association between feed and fun- ited sample size, but as a whole they allow inferring a gal gut might also apply for the plant pathogens Ustilago limited and taxa specific genetic control over diversity hordei and Mycosphaerella tassiana, found in about 3% and composition of gut eukaryotes. Besides, GWAS re- of the animals. Former results are in agreement with vealed genetic variants associated at chromosomal level previous reports on humans and non-human primates with fungal α-diversity (on SSC6) and the abundance of [2, 30–32], and suggest that, differently to pig gut pro- Blastocystis spp. (on SSC6, 17 and 18). Two intronic karyotes [33], the pig gut mycobiome may lack a stable SNPs associated with fungal α-diversity were mapped on core. Consequently, it is expected that a large proportion the IL23R gene. Remarkably, this gene, which is highly of the fungi detected in pig fecal samples may be tran- expressed in Th17 cells, seems to play an important role sients from dietary or environmental origin. Meanwhile, in the proliferation and survival of these cells, which are the prevalence (> 80% of the 514 pig samples) of four critical for host defense against fungal infections [39]. In protist genera (Blastocystis, Neobalantidium, Tetratri- fact, genetic variants in the IL23R gene have been associ- chomonas and Trichomitus) was observed, which is in ated with different susceptibility to fungal infections [40] agreement with previous reports in pigs [17, 34]. Neoba- and IL-23R deficient mouse was susceptible to systemic lantidium is a world-wide parasitic-opportunistic human infection with C. albicans [41]. Furthermore, IL-23R had pathogen that is transmitted through the fecal-oral been reported to be involved in the diversity of ileum route, particularly when there is a close contact with microbiota in humans [42]. In the same chromosome, pigs, which are its asymptomatic reservoir hosts [35]. PIK3C3 gene gathered two SNPs associated to Blastocys- Blastocystis is a ubiquitous protozoan of human and tis abundance. PIK3C3 is required for T cell homeostasis pigs, commonly found in healthy populations [6, 36]. [43–45] and reported to play a relevant role in maintain- However, several species in Blastocystis are causative ing gut homeostasis [46]. Other genes involved in the in- agents of diarrhea in humans through the fecal-oral in- testinal epithelial integrity and gut homeostasis were fection route, being the pig a common reservoir [37]. A annotated on SSC17, such as HNF4A, a relevant regula- number of protists in the Trichomonadida order (super- tor that mediates microbial control of intestinal gene ex- group Excavata) that are common in the digestive tract pression [47]. Also a number of genes related to the of pigs [38] were also found in the present study. immune system were annotated within the aforemen- Colonization with this species has not been associated to tioned chromosomal intervals, as the TNFRSF9 gene that virulence in mammals, but the full pathogenic potential contributes to the development, survival and activation of T. batrachorum has yet to be explored. of T cells [48, 49], the IL12RB2 and PIK3CD genes, and The analysis of host intrinsic and extrinsic factors a plethora of members of the immune-related pathways ascertained the significant modulatory role of the envir- Gα12/13 Signaling, IL-9 Signaling, IL-23 Signaling, IL-12 onmental effects gathered in batch and farm factors (e.g. Signaling, ERK/MAPK Signaling, Agranulocyte and climate, management conditions, diet) on gut eukaryotes Granulocyte Adhesion and Diapedesis. Several studies, composition and diversity. Also important differences most of them performed in cell lines, have tried to shed between the diversity of gut eukaryotic communities in light on Blastocystis–host interactions. These studies weaned and finishing pigs were evidenced, whereas have shown that Blastocystis are able to disrupt intestinal scarce differences due to the sex of the animal were ob- barrier integrity and function, and that they modulate served. We explored the putative host genetic control of the host immune response through the degradation of the diversity and composition of both gut protist and IgA, the suppression of iNOS, and the induction of pro- fungal communities. Despite the relatively low heritabil- inflammatory cytokines such as IL-12 and IL-23, among ity estimates (posterior means from 0.158 to 0.281), they others [50, 51]. Furthermore, Blastocystis subtype- represent the first reported evidence of certain host gen- dependent upregulation of pro-inflammatory cytokines etic determinism of gut eukaryotic profile in pigs. The in macrophages through the activation of MAPK path- abundance of the protist Blastocystis genus was the most ways has been described [52]. Overall, it is noteworthy heritable, with around 28% of their variability explained to highlight that in our study genes related to pathways by the genetic variability of the host, followed by the previously associated with Blastocystis experimental in- abundance of Trichomitus and Neobalantidium genus, fections were identified as being associated with Blasto- and the preponderant fungal genus Kazachstania. Lower cystis relative abundance. Although further genetic and Ramayo-Caldas et al. Animal Microbiome (2020) 2:18 Page 8 of 12 functional analyses are needed to better understand how included after weaning in a typical transition diet. An al- the biological mechanisms of the host could modulate ternative explanation may be also the recognized [6, 10, the diversity and abundance of fungi and protists, and 11, 62] positive effects of Blastocystis and Entamoeba over their relationship with the gut bacterial counterparts. the microbiota diversity and richness, which is turn is as- These findings are in line with recent reports in mice sociated with host gut health. [53], and suggest that polymorphisms located in (or in In summary, our findings highlight the relevance of linkage disequilibrium with) genes related to immune considering the gut eukaryotic communities to better system and gut homeostasis may modulate the diversity understand the porcine gastrointestinal microbiome and and composition of the eukaryote gut communities in its contributions to host performance and health. At that pigs. The evolutionary mechanism beyond the observed point we cannot neglect the contribution of other genera association may be explained by indirect relationships that might be relevant at different ages, under different between mutual favorable selection at host-genome and environmental conditions of for other phenotypic traits. microbial level [54, 55]. Although far to be fully under- We are also aware of some limitations of our study, as stood, some examples of host-genetic microbiota associ- for example the limited sample size or the lack of de- ations have been reported (reviewed in [55]), that similar grees of freedom to estimate all effects, but also of meth- to our results identify candidate genes mainly related to odological constraints such as primers choice or the use host-response against pathogens, sensing microbes, cell of less curated (compared to bacterial) protists and fun- signaling pathways and innate immune system. However, gal databases, that together with the reduced amount of in agreement with previous reports centered in gut bac- public available reference genomes compromise the ac- terial [53, 56], our findings also indicate a limited and curacy of the taxonomic classification. To the best of target specific host genetic control over the composition our knowledge, this study represents the largest effort to of the gut eukaryote communities in pig. characterize the gut fungal and commensal protist com- On another level, we provided evidences of the associ- munities in pigs, but further larger studies including ex- ation between piglets body weight after weaning and perimental validations and alternative meta-sequencing members from the fungal and protists gut eukaryote com- approaches are needed to unveil the role of the host- munity. A previous study with growing piglets, reported associated microbial communities in pigs production no significant effects on growth performance parameters performance, welfare and health. upon oral supplementation with K. slooffiae [57]. How- ever, a recent detailed metabolic study conducted by the Conclusions same authors [58] under laboratory conditions with cul- The diversity and composition of gut commensal tures of K. slooffiae, concluded that this yeast produces eukaryotic communities in healthy pigs at two ages have peptides and short chain fatty acids that might benefit the been characterized. The porcine gut mycobiota is domi- gut health, and provides an additional protein source that nated by yeast, with a high prevalence of Kazachstania, and contains essential amino acids and other useful growth a common set of four protist genera (Blastocystis, Neoba- factors for the animal [58]. Consequently, it seems plaus- lantidium, Tetratrichomonas and Trichomitus) persisted ible that the observed relationships between Kazachstania through the majority of animals. Our results point towards and piglets body weight partially relies on the K. slooffiae a minor and taxa specific genetic control over the diversity positive effect on gut health through the production of es- and composition of the pig gut eukaryotic communities, sential amino acids and short fatty acids that are absorbed but we describe associations with genes functionally related by the host or employed by others microbial members of to the immune system and gut homeostasis that might have the pig gut microbiome. On the other hand, the impact of an effect in modulating the fungi α-diversity and the abun- members of the gut eukaryotic community to growth per- dance of Blastocystis ssp. We provide also evidences of the formance and the digestive enzyme activities have been re- associations between piglets body weight after weaning and ported in shrimp [59]. Furthermore, Blastocystis spp. and members from the gut fungal and protist eukaryote com- non-pathogenic Entamoeba spp. have been associated munity. Overall, our findings highlight the relevance of con- with a healthy and diverse gut microbiome ecosystem [6, sidering, together with that of bacteria, the contribution of 10, 11]. Likewise, mechanisms of horizontal gene transfer gut eukaryotic communities to better understand the host- of genes associated with cellulose degradation and carbo- microbiome association and its role on pig performance, hydrate use from bacterial to protist such as Blastocystis welfare and health. have been documented [60–62]. Therefore, we hypothesize that observed positive relationships between Methods piglets body weight at 60 days and members of Blastocystis Sample collection, DNA extraction and sequencing and Entamoeba may be explained by the capacity of Blas- Fecal samples were collected at two ages from 514 tocystis to degrade carbohydrates, which are routinely healthy Duroc pigs belonging to the same commercial Ramayo-Caldas et al. Animal Microbiome (2020) 2:18 Page 9 of 12 outbred line but allocated in two different farms. A total Shannon index [70], and Beta-diversity was assessed of 405 weaned piglets (204 males and 201 females) dis- using the Whittaker index [71]. To identify environmen- tributed in seven batches were sampled in a commercial tal (farm, batch or pen) and host-covariates (sex, age, farm at 60 ± 8 days of age and mean body weight around body weight) that may modulate the diversity, structure 18.58 kg (SD 3.13), after 4 weeks receiving the same and profile of the eukaryote communities, we run a Per- transition-based diet (Additional file 8 Table S6). The mutational multivariate analysis of variance (PERMA- remaining 109 pigs (50 castrated males and 59 females) NOVA) test using the adonis function from vegan [69]. were raised under intensive standard conditions at IRTA Significance levels were determined after 10,000 permu- experimental farm (IRTA, Monells, Spain), and fecal tations and the multiple comparison tests were per- samples were collected at 190 ± 5 days of age, when they formed using False Discovery Rate (FDR). fed a finishing standard diet. Both groups of pigs were genetically connected. Animal care and experimental Genotype data and estimation of heritability procedures were carried out following national and insti- The putative host genetics determinism of gut tutional guidelines for the Good Experimental Practices eukaryotic profile in pigs was initially assessed by esti- and were approved by the IRTA Ethical Committee. mating the heritability (h2) of both fungi and protists DNA was extracted with the DNeasy PowerSoil Kit alpha-diversity, as well as of their taxa abundance. For (QIAGEN, Hilden, Germany), following manufacturer’s these analyses, we used the samples taken at 60 ± 8 days instructions. Extracted DNA was sent to the University of age in the commercial farm. The Porcine 70 k GGP of Illinois Keck Center for Fluidigm sample preparation Porcine HD Array (Illumina, San Diego, CA) was used and sequencing. The Internal Transcribed Spacer 2 to genotype 390 out of 405 animals sampled in the com- (ITS2) region was amplified using primers ITS3: 5′- mercial farm. Parameters estimation was performed GCATCGATGAAGAACGCAGC-3′ and ITS4: 5′- using the following Bayesian mixed model implemented TCCTCCGCTTATTGATATGC-3′ [63]. Protist-specific with the bglr R package [72]: primers [64], F-566: 5′-CAGCAGCCGCGGTAATTCC- 3′ and R-1200: 5′-CCCGTGTTGAGTCAAATTAAGC- yijk ¼ sex j þ bk þ ui þ eijk 3′ were used to amplify the 18S rRNA gene fragment. Amplicons were paired-end (2 × 250 nt) sequenced on where yijk corresponds to the alpha-diversity (fungi or an Illumina NovaSeq (Illumina, San Diego, CA, USA). protist) or taxa (genera or specie) clr-transformed abun- dance of the th individual of sex j in the kth batch; sexj Bioinformatics and statistical analysis and bk correspond to the systematic effects of jth sex (2 Sequences were analysed with Qiime2 [65], barcode se- levels) and kth batch effect (7 levels), respectively; ui is quences, primers and low-quality reads (Phred scores of the random genetic effect of individual i, distributed as < 30) were removed. The quality control also trimmed case ∼ (0, 2 ) being the genomic relationship sequences based on expected amplicon length and re- matrix calculated using the filtered autosomal SNPs moved chimeras. Afterwards, sequences were processed based on the methodology of Yang et al. [73]; finally, eijk into Amplicon Sequences Variants (ASVs) at 99% of is the random residual term, with a distribution e∼ (0,I identity. ASVs present in less than three samples and σ 2e ). The model was run using a Gibbs sampler with representing less than 0.005% of the total counts were 30,000 iterations and a burn-in of 2000 rounds. Posterior filtered out. Samples with less than 6000 (fungi, n = 21 sample mean and standard deviation of the heritability ð samples) or 10,000 (protists) reads were also excluded. h2 ¼ σ 2u =ðσ 2u þ σ 2e ÞÞ were obtained from the resulted ASVs were classified to the lowest possible taxonomic posterior distributions. level based on SILVA v123 database [66] for 18S rRNA genes, and the UNITE QIIME version (release Genome wide association study (GWAS) 18.11.2018) for fungi [67]. Subsequently, we excluded Quality control was performed with plink [74] to ex- those ASVs not taxonomically classified as protists or clude single nucleotide polymorphisms (SNPs) with fungi. Moreover, we verified the fungi taxonomic assig- minor allele frequencies < 5%, rates of missing genotypes nation following the recommendation of Nilsson [68] by above 10%, as well as SNPs that did not map to the por- a manual examination of the most abundant fungal cine reference genome (Sscrofa11.1 assembly). Then, to ASVs against the International Nucleotide Sequence identify SNPs from the host genome associated with the Database (http://www.insdc.org/). Before the estimation alpha-diversity as well as protists and fungi relative of diversity indexes, samples were rarefied at 6000 abundances, genome-wide association studies (GWAS) (fungi) and 10,000 (protists) reads of depth, to allow an were performed between 42,608 SNPs and the alpha- equal depth. Diversity metrics were estimated with vegan diversity or the centered log ratio (clr) transformed gen- R package [69]. Alpha-diversity was evaluated with the era and species abundance tables. Only the genera and Ramayo-Caldas et al. Animal Microbiome (2020) 2:18 Page 10 of 12 species fully taxonomically classified and present in Supplementary information more than 80% of the samples were included in the ana- Supplementary information accompanies this paper at https://doi.org/10. 1186/s42523-020-00038-4. lysis. For this propose, the genome-wide complex trait analysis (GCTA) software [73] was employed using the Additional file 1 Table S1. Description of the prevalence and relative following model at each SNP: abundance of the pig gut eukaryotes communities. Additional file 2 Figure S1. Sample distribution of fungi (A) and protist (B) communities. Diversity indexes of fungi (C) and protist (D) yijk ¼ sex j þ bk þ ui þ sli al þ eijk communities. Blue color represents the finishing pigs (experimental farm at 190 days) and red the weaned piglets (commercial farm at 60 days). Additional file 3 Table S2. Description of the genetic markers where yijk corresponds to the alpha-diversity (fungi or identified as significantly according to the genome-wide association protist) or taxa (genera or specie) clr-transformed abun- studies. dance of the th individual of sex j in the kth batch; sexj, Additional file 4 Table S3. Description of the genes annotated within the genomic intervals. bk and ui are, respectively, the effects of sex, batch and Additional file 5 Table S4. Functional annotation of the identified infinitesimal genetic effect described in the previous candidate genes. model; sli is the genotype (coded as 0,1,2) for the lth SNP Additional file 6 Table S5. Metabolic pathways significantly enriched of individual i, and al is the allele substitution effect of by the list of candidate genes. SNP l on the analysed trait. A SNP was considered to be Additional file 7 Figure S2. Patterns of association between the piglets significantly associated if the corresponding p-value was body weight at 60-days old and the fungal and protist abundance ac- cording to the Shapley Additive Explanation (SHAP) values. The x-axis rep- lower than 0.05, after Benjamini–Hochberg [75] correc- resent the genera abundance after Cumulative Sum Scaling (CSS) tion for multiple testing at chromosome level. normalization and log transformation. Y-axis represents the SHAP values of the piglets’ body weight. Additional file 8 Table S6. Description of the ingredients and nutrient Gene functional classification and canonical pathway content of the diets used during the experiment. analyses Functional classification and pathway analyses of the an- Abbreviations notated candidate genes were carried out using the In- QIIME: Quantitative insights into microbial ecology; ITS: Internal transcribed spacer; ASV: Amplicon sequence variant; clr: Centered log ratio genuity Pathways Analysis software (IPA; Ingenuity transformation; GWAS: Genome-wide association studies; XGBoost: eXtreme Systems, http://www.ingenuity.com). Significance levels Gradient Boosting; FDR: False discovery rate for enrichment of each canonical pathway in the list of Acknowledgements candidate genes were calculated using Fisher’s exact test, The authors warmly thank all technical staff from Selección Batallé S.A, for and the resulting p-values were corrected for multiple- providing the animal material and their collaboration during the sampling, as test using the Benjamini and Hochberg algorithm [75]. well as to Juan Pablo Sanchez for valuable discussions and comments on the manuscript. The cut-off for considering an enrichment as significant was established at a corrected p-value < 0.05. Authors’ contributions YRC, MB, AD, and RQ designed the study. YRC, MB, OG and RQ performed the sampling. MB and OG carried out DNA extraction. YRC, FPB and LZ Association between pig gut eukaryote communities and analyzed the data. YRC, FPB, LZ, MB and RQ interpreted the results and wrote the manuscript. The authors read and approved the final version of body weight the manuscript. The associations between fungal and protist genera abundance and piglest body weight at 60-days was evalu- Funding YRC was funded by Marie Skłodowska-Curie grant (P-Sphere) agreement No ated using the dataset from the commercial farm (n = 6655919 (EU). LMZ is recipient of a Ph.D. grant from Ministry of Economy 405). In a first step we implemented XGBoost, a fast al- and Science, Spain associated with ‘Centro de Excelencia Severo Ochoa gorithm which incorporates a penalization term in the 2016–2019’ award SEV-2015-0533 to CRAG. MB is recipient of a Ramon y Cajal post-doctoral fellowship (RYC-2013–12573) from the Spanish Ministry of loss function to prevent overfitting (Chen & Guestrin, Economy and Competitiveness. Part of the research presented in this publi- 2016) using the XGBoost R packages. Afterwards, the cation was funded by Grants AGL2016–75432-R, AGL2017–88849-R awarded importance of the variables was evaluated with Shapley by the Spanish Ministry of Economy and Competitiveness. The authors be- long to Consolidated Research Group TERRA (AGAUR, ref. 2017 SGR 1290 Additive Explanation (SHAP) algorithm. SHAP (Lund- and 2017 SGR 1719). berg & Lee, 2017) is an extension of the coalition game from the game theory (Shapley, 1953), which is used in a Availability of data and materials context of multicollinearity. Basically, the model evalu- The raw sequencing data employed in this article has been submitted to the NCBI’s sequence read archive (https://www.ncbi.nlm.nih.gov/sra); BioProject: ates the importance of each variable by measuring a con- PRJNA608629. ditional contribution, i.e., the effect caused in the prediction by including that variable. SHAP values were Ethics approval Animal care and experimental procedures were carried out following obtained and visualized using the ShapXGboost R pack- national and institutional guidelines for the Good Experimental Practices and age [76]. were approved by the IRTA Ethical Committee. Ramayo-Caldas et al. Animal Microbiome (2020) 2:18 Page 11 of 12 Consent for publication 19. Chen C, Huang X, Fang S, Yang H, He M, Zhao Y, et al. Contribution of host Not applicable. genetics to the variation of microbial composition of cecum lumen and feces in pigs. Front Microbiol. 2018;9:2626. 20. Camarinha-Silva A, Maushammer M, Wellmann R, Vital M, Preuss S, Competing interests Bennewitz J. Host genome influence on gut microbial composition and The authors declare no competing interests. microbial prediction of complex traits in pigs. Genetics. 2017;206(3):1637. 21. Crespo-Piazuelo D, Migura-Garcia L, Estellé J, Criado-Mesas L, Revilla M, Author details Castelló A, et al. Association between the pig genome and its gut 1 Animal Breeding and Genetics Program, Institute of Agrifood Research and microbiota composition. Sci Rep. 2019;9(1):8791. Technology (IRTA), Torre Marimon, 08140 Caldes de Montbui, Spain. 2GIRO, 22. Hu J, Nie Y, Chen J, Zhang Y, Wang Z, Fan Q, et al. Gradual changes of gut IRTA, Torre Marimon, 08140 Caldes de Montbui, Spain. 3Centre for Research microbiota in weaned miniature piglets. Front Microbiol. 2016;7:1727. in Agricultural Genomics, CSIC-IRTA-UAB-UB Consortium, Bellaterra, Spain. 23. Urubschurov V, Janczyk P, Pieper R, Souffrant WB. Biological diversity of 4 Animal Welfare Subprogram, IRTA, 17121 Monells, Spain. yeasts in the gastrointestinal tract of weaned piglets kept under different farm conditions. FEMS Yeast Res. 2008;8(8):1349–56. Received: 3 March 2020 Accepted: 29 April 2020 24. Van Uden N, do Carmo Sousa L. Quantitative aspects of the intestinal yeast Flora of swine. Microbiology. 1962;27(1):35–40. 25. Kurtzman CP, Robnett CJ, Ward JM, Brayton C, Gorelick P, Walsh TJ. Multigene phylogenetic analysis of pathogenic candida species in the References Kazachstania (Arxiozyma) telluris complex and description of their 1. Hooks KB, O'Malley MA. Contrasting strategies: human eukaryotic versus ascosporic states as Kazachstania bovina sp. nov., K. heterogenica sp. nov., K. bacterial microbiome research. J Eukaryot Microbiol. 2019;67(2):279–95. pintolopesii sp. nov., and K. slooffiae sp. nov. J Clin Microbiol. 2005;43(1): 2. Nash AK, Auchtung TA, Wong MC, Smith DP, Gesell JR, Ross MC, et al. The 101–11. gut mycobiome of the human microbiome project healthy cohort. 26. Kemler M, Witfeld F, Begerow D, Yurkov A. Phylloplane yeasts in temperate Microbiome. 2017;5(1):153. climates. In: Buzzini P, Lachance M-A, Yurkov A, editors. Yeasts in natural 3. Parfrey LW, Walters WA, Lauber CL, Clemente JC, Berg-Lyons D, Teiling C, ecosystems: diversity. Cham: Springer International Publishing; 2017. p. 171– et al. Communities of microbial eukaryotes in the mammalian gut within 97. the context of environmental eukaryotic diversity. Front Microbiol. 2014;5: 27. Hoog GSD, Guarro J, Gene J, Figueras MJ. Atlas of clinical fungi. 2nd ed; 298. 2000. 4. Hamad I, Raoult D, Bittar F. Repertory of eukaryotes (eukaryome) in the human gastrointestinal tract: taxonomy and detection methods. Parasite 28. Bao X, Carris LM, Huang G, Luo J, Liu Y, Castlebury LA. Tilletia puccinelliae, a Immunol. 2016;38(1):12–36. new species of reticulate-spored bunt fungus infecting Puccinellia distans. 5. Laforest-Lapointe I, Arrieta M-C. Microbial eukaryotes: a missing Link in gut Mycologia. 2010;102(3):613–23. microbiome studies. mSystems. 2018;3(2):e00201–17. 29. Preugschat K, Kersten S, Ettle T, Richter W, Karl H, Breves G, et al. Effects of 6. Chabé M, Lokmer A, Ségurel L. Gut protozoa: friends or foes of the human feeding diets containing increasing proportions of bunt-infected wheat gut microbiota? Trends Parasitol. 2017;33(12):925–34. (Tilletia caries) on performance and health of pigs. Arch Anim Nutr. 2014; 7. Underhill DM, Iliev ID. The mycobiota: interactions between commensal 68(1):55–62. fungi and the host immune system. Nat Rev Immunol. 2014;14(6):405–16. 30. Raimondi S, Amaretti A, Gozzoli C, Simone M, Righini L, Candeliere F, et al. 8. Parfrey LW, Walters WA, Knight R. Microbial eukaryotes in the human Longitudinal survey of fungi in the human gut: ITS profiling, phenotyping, microbiome: ecology, evolution, and future directions. Front Microbiol. and colonization. Front Microbiol. 2019;10:1575. 2011;2:153. 31. Hallen-Adams HE, Kachman SD, Kim J, Legge RM, Martínez I. Fungi 9. Lukeš J, Stensvold CR, Jirků-Pomajbíková K, Wegener Parfrey L. Are human inhabiting the healthy human gastrointestinal tract: a diverse and dynamic intestinal eukaryotes beneficial or commensals? PLoS Pathog. 2015;11(8): community. Fungal Ecol. 2015;15:9–17. e1005039. 32. Mann AE, Mazel F, Lemay MA, Morien E, Billy V, Kowalewski M, et al. 10. Audebert C, Even G, Cian A, Safadi DE, Certad G, Delhaes L, et al. Biodiversity of protists and nematodes in the wild nonhuman primate gut. Colonization with the enteric protozoa Blastocystis is associated with ISME J. 2019;14:609–22. increased diversity of human gut bacterial microbiota. Sci Rep. 2016;6(1): 33. Xiao L, Estellé J, Kiilerich P, Ramayo-Caldas Y, Xia Z, Feng Q, et al. A 25255. reference gene catalogue of the pig gut microbiome. Nat Microbiol. 2016; 11. Tito RY, Chaffron S, Caenepeel C, Lima-Mendez G, Wang J, Vieira-Silva S, 1(12):16161. et al. Population-level analysis of Blastocystis subtype prevalence and 34. Pakandl M. The prevalence of intestinal protozoa in wild and domestic pigs. variation in the human gut microbiota. Gut. 2019;68(7):1180–9. Vet Med (Praha). 1994;39(7):377–80. 12. Chudnovskiy A, Mortha A, Kana V, Kennard A, Ramirez JD, Rahman A, et al. 35. Schuster FL, Ramirez-Avila L. Current world status of <em>Balantidium Host-protozoan interactions protect from mucosal infections through coli</em>. Clin Microbiol Rev. 2008;21(4):626. activation of the Inflammasome. Cell. 2016;167(2):444–456.e414. 36. Wang W, Bielefeldt-Ohmann H, Traub RJ, Cuttell L, Owen H. Location and 13. Summers KL, Frey JF, Ramsay TG, Arfken AM. The piglet mycobiome during pathogenic potential of Blastocystis in the porcine intestine. PLoS One. the weaning transition: a pilot study1. J Anim Sci. 2019;97(7):2889–900. 2014;9(8):e103962. 14. White JK, Nielsen JL, Madsen AM. Microbial species and biodiversity in 37. Rivera WL. Phylogenetic analysis of Blastocystis isolates from animal and settling dust within and between pig farms. Environ Res. 2019;171:558–67. human hosts in the Philippines. Vet Parasitol. 2008;156(3):178–82. 15. Arfken AM, Frey JF, Ramsay TG, Summers KL. Yeasts of burden: exploring 38. Mostegl MM, Richter B, Nedorost N, Lang C, Maderner A, Dinhopl N, et al. the Mycobiome–Bacteriome of the piglet GI tract. Front Microbiol. 2019;10: First evidence of previously undescribed trichomonad species in the 2286. intestine of pigs? Vet Parasitol. 2012;185(2–4):86–90. 16. Shieban F. Studies on intestinal protozoa of domestic pigs in the Teheran 39. Khader SA, Gaffen SL, Kolls JK. Th17 cells at the crossroads of innate and area of Iran**this project is supported in part by funds for the endemic adaptive immunity against infectious diseases at the mucosa. Mucosal diseases research project of plan organization and of the Department of Immunol. 2009;2(5):403–11. Health Science, Teheran University School of Medicine. Br Vet J. 1971;127(3): 40. Carvalho A, Cunha C, Di Ianni M, Pitzurra L, Aloisi T, Falzetti F, et al. iii–v. Prognostic significance of genetic variants in the IL-23/Th17 pathway for 17. Wylezich C, Belka A, Hanke D, Beer M, Blome S, Höper D. Metagenomics for the outcome of T cell-depleted allogeneic stem cell transplantation. Bone broad and improved parasite detection: a proof-of-concept study using Marrow Transplant. 2010;45(11):1645–52. swine faecal samples. Int J Parasitol. 2019;49(10):769–77. 41. Nur S, Sparber F, Lemberg C, Guiducci E, Schweizer TA, Zwicky P, et al. IL-23 18. Matsubayashi M, Suzuta F, Terayama Y, Shimojo K, Yui T, Haritani M, et al. supports host defense against systemic Candida albicans infection by Ultrastructural characteristics and molecular identification of Entamoeba suis ensuring myeloid cell survival. PLoS Pathog. 2020;15(12):e1008115. isolated from pigs with hemorrhagic colitis: implications for pathogenicity. 42. Zakrzewski M, Simms LA, Brown A, Appleyard M, Irwin J, Waddell N, et al. Parasitol Res. 2014;113(8):3023–8. IL23R-protective coding variant promotes beneficial bacteria and diversity in Ramayo-Caldas et al. Animal Microbiome (2020) 2:18 Page 12 of 12 the Ileal microbiome in healthy individuals without inflammatory bowel 65. Bolyen E, Rideout JR, Dillon MR, Bokulich NA, Abnet CC, Al-Ghalith GA, et al. disease. J Crohns Colitis. 2018;13(4):451–61. Reproducible, interactive, scalable and extensible microbiome data science 43. Willinger T, Flavell RA. Canonical autophagy dependent on the class III using QIIME 2. Nat Biotechnol. 2019;37(8):852–7. phosphoinositide-3 kinase Vps34 is required for naive T-cell homeostasis. 66. Quast C, Pruesse E, Yilmaz P, Gerken J, Schweer T, Yarza P, et al. The SILVA Proc Natl Acad Sci U S A. 2012;109(22):8670–5. ribosomal RNA gene database project: improved data processing and web- 44. Parekh VV, Wu L, Boyd KL, Williams JA, Gaddy JA, Olivares-Villagómez D, based tools. Nucleic Acids Res. 2013;41(Database issue):D590–6. et al. Impaired autophagy, defective T cell homeostasis, and a wasting 67. Community U. UNITE QIIME release for fungi. Version 18.11.2018. UNITE syndrome in mice with a T cell-specific deletion of Vps34. J Immunol. 2013; community; 2019. https://doi.org/10.15156/BIO/786334. 190(10):5086–101. 68. Nilsson RH, Anslan S, Bahram M, Wurzbacher C, Baldrian P, Tedersoo L. 45. McLeod IX, Zhou X, Li Q-J, Wang F, He Y-W. The class III kinase Vps34 Mycobiome diversity: high-throughput sequencing and identification of promotes T lymphocyte survival through regulating IL-7Rα surface fungi. Nat Rev Microbiol. 2019;17(2):95–109. expression. J Immunol. 2011;187(10):5051. 69. Jari Oksanen F, Blanchet G, Friendly M, Kindt R, Legendre P, McGlinn D, 46. Zhao S, Xia J, Wu X, Zhang L, Wang P, Wang H, et al. Deficiency in class III et al. Vegan: community ecology package. R package version 2.5–3. 2018. PI3-kinase confers postnatal lethality with IBD-like features in zebrafish. Nat https://www.CRANR-projectorg/package=vegan. Commun. 2018;9(1):2639. 70. Shannon C. A mathematical theory of communication; 1984. 47. Davison JM, Lickwar CR, Song L, Breton G, Crawford GE, Rawls JF. Microbiota 71. Whittaker RH. Evolution and measurement of species diversity. Taxon. 1972; regulate intestinal epithelial gene expression by suppressing the 21(2/3):213–51. transcription factor hepatocyte nuclear factor 4 alpha. Genome Res. 2017; 72. Pérez P, de los Campos G. Genome-wide regression and prediction with the 27(7):1195–206. BGLR statistical package. Genetics. 2014;198(2):483–95. 48. Ward-Kavanagh LK, Lin WW, Šedý JR, Ware CF. The TNF receptor 73. Yang J, Lee SH, Goddard ME, Visscher PM. GCTA: a tool for genome-wide superfamily in co-stimulating and co-inhibitory responses. Immunity. 2016; complex trait analysis. Am J Hum Genet. 2011;88(1):76–82. 44(5):1005–19. 74. Purcell S, Neale B, Todd-Brown K, Thomas L, Ferreira MAR, Bender D, et al. 49. Schwarz H, Tuckwell J, Lotz M. A receptor induced by lymphocyte activation PLINK: a tool set for whole-genome association and population-based (ILA): a new member of the human nerve-growth-factor/tumor-necrosis- linkage analyses. Am J Hum Genet. 2007;81(3):559–75. factor receptor family. Gene. 1993;134(2):295–8. 75. Benjamini Y, Hochberg Y. Controlling the false discovery rate: a practical 50. H. M, K.S T, H. Y. Blastocystis: pathogen or passenger? An evaluation of 101 and powerful approach to multiple testing. J R Stat Soc Ser B Methodol. years of research, vol. 4; 2012. 1995;57(1):289–300. 51. Ajjampur SSR, Tan KSW. Pathogenic mechanisms in Blastocystis spp. — 76. Liu Y. SHAP for XGBoost in R: SHAPforxgboosthttps://github.com/liuyanguu/ interpreting results from in vitro and in vivo studies. Parasitol Int. 2016;65(6, SHAPforxgboost; 2019. Part B):772–9. 52. Lim MX, Png CW, Tay CYB, Teo JDW, Jiao H, Lehming N, et al. Differential Publisher’s Note regulation of proinflammatory cytokine expression by mitogen-activated Springer Nature remains neutral with regard to jurisdictional claims in protein kinases in macrophages in response to intestinal parasite infection. published maps and institutional affiliations. Infect Immun. 2014;82(11):4789–801. 53. Khan AA, Yurkovetskiy L, O’Grady K, Pickard JM, de Pooter R, Antonopoulos DA, et al. Polymorphic immune mechanisms regulate commensal repertoire. Cell Rep. 2019;29(3):541–550.e544. 54. Zilber-Rosenberg I, Rosenberg E. Role of microorganisms in the evolution of animals and plants: the hologenome theory of evolution. FEMS Microbiol Rev. 2008;32(5):723–35. 55. Wang J, Chen L, Zhao N, Xu X, Xu Y, Zhu B. Of genes and microbes: solving the intricacies in host genomes. Protein Cell. 2018;9(5):446–61. 56. Rothschild D, Weissbrod O, Barkan E, Kurilshikov A, Korem T, Zeevi D, et al. Environment dominates over host genetics in shaping human gut microbiota. Nature. 2018;555(7695):210–5. 57. Urubschurov V, Büsing K, Freyer G, Herlemann DPR, Souffrant W-B, Zeyner A. New insights into the role of the porcine intestinal yeast, Kazachstania slooffiae, in intestinal environment of weaned piglets. FEMS Microbiol Ecol. 2016;93(2).. 58. Urubschurov V, Büsing K, Souffrant WB, Schauer N, Zeyner A. Porcine intestinal yeast species, Kazachstania slooffiae, a new potential protein source with favourable amino acid composition for animals. J Anim Physiol Anim Nutr. 2018;102(2):e892–901. 59. Dai W, Yu W, Zhang J, Zhu J, Tao Z, Xiong J. The gut eukaryotic microbiota influences the growth performance among cohabitating shrimp. Appl Microbiol Biotechnol. 2017;101(16):6447–57. 60. del Campo J, Bass D, Keeling PJ. The eukaryome: diversity and role of micro- eukaryotic organisms associated with animal hosts. Funct Ecol. 2019. Early view, Online Version of Record before inclusion in an issue. 61. Todaka N, Inoue T, Saita K, Ohkuma M, Nalepa CA, Lenz M, et al. Phylogenetic analysis of cellulolytic enzyme genes from representative lineages of termites and a related cockroach. PLoS One. 2010;5(1):e8636. 62. Eme L, Gentekaki E, Curtis B, Archibald JM, Roger AJ. Lateral Gene transfer in the adaptation of the anaerobic parasite Blastocystis to the gut. Curr Biol. 2017;27(6):807–20. 63. White TJ, Bruns T, Lee S, Taylor J. PCR protocols amplification and direct sequencing of fungal ribosomal RNA genes for phylogenetics: Elsevier; 1990.. 64. Hadziavdic K, Lekang K, Lanzen A, Jonassen I, Thompson EM, Troedsson C. Characterization of the 18S rRNA Gene for designing universal eukaryote specific primers. PLoS One. 2014;9(2):e87624.