The faecal metabolome of black howler monkeys (Alouatta pigra) varies in response to seasonal dietary changes.

Mammals rely on the metabolic functions of their gut microbiota to meet their energetic needs and digest potentially toxic components in their diet. The gut microbiome plastically responds to shifts in host diet and may buffer variation in energy and nutrient availability. However, it is unclear how seasonal differences in the gut microbiome influence microbial metabolism and nutrients available to hosts. In this study, we examine seasonal variation in the gut metabolome of black howler monkeys (Alouatta pigra) to determine whether those variations are associated with differences in gut microbiome composition and nutrient intake, and if plasticity in the gut microbiome buffers shortfalls in energy or nutrient intake. We integrated data on the metabolome of 81 faecal samples from 16 individuals collected across three distinct seasons with gut microbiome, nutrient intake and plant metabolite consumption data from the same period. Faecal metabolite profiles differed significantly between seasons and were strongly associated with changes in plant metabolite consumption. However, microbial community composition and faecal metabolite composition were not strongly associated. Additionally, the connectivity and stability of faecal metabolome networks varied seasonally, with network connectivity being highest during the dry, fruit-dominated season when black howler monkey diets were calorically and nutritionally constrained. Network stability was highest during the dry, leaf-dominated season when most nutrients were being consumed at intermediate rates. Our results suggest that the gut microbiome buffers seasonal variation in dietary intake, and that the buffering effect is most limited when host diet becomes calorically or nutritionally restricted.

INTRODUCTION

Mammals that consume diets rich in cellulose and other plant fibres depend on gut microbes to extract short chain fatty acids (SCFAs) and additional energy sources from food items (Flint et al., 2008; Flint & Bayer, 2008; Lambert, 1998; Mackie, 2002). In mammals relying heavily on plant foods, the gut microbiome plays an essential role in breaking down phenolics and other plant toxins (Dearing & Kohl, 2017; Greene et al., 2020; Kohl et al., 2014, 2016). Additionally, dietary shifts, environmental changes in food availability, rainfall and temperature, and habitat differences, such as forest structure and anthropogenic changes, appear to more strongly influence the gut microbiome composition of animals that ingest a high proportion of indigestible plant fibre than animals that ingest a diet composed principally of fruit, flowers and/or invertebrates (Frankel et al., 2019; Greene et al., 2019). Moreover, there is evidence of functional and compositional convergence in the gut microbiome of animals with high‐fibre diets, both within and across taxonomic groups (Amato et al., 2019; Hale et al., 2018; Ley et al., 2008; Muegge et al., 2011).

Gut microbial composition can rapidly change in response to dietary variation (e.g., high fat to low fat consumption) (David et al., 2014; Turnbaugh et al., 2009). Several studies have identified seasonal shifts in gut microbial composition, including in wild mice (Apodemus sylvaticus) (Maurice et al., 2015), wild Bale monkeys (Chlorocebus djamdjamensis) (Trosvik, Rueness, et al., 2018), captive and wild giant pandas (Ailuropoda melanoleuca) (Wu et al., 2017; Xue et al., 2015), wild flying squirrels (Pteromys volans orii) (Liu et al., 2019), wild sage grouse (Centrocercus urophasianus) (Drovetski et al., 2019), wild Verreaux's sifakas (Propithecus verreauxi) (Springer et al., 2017), wild tench (Tinca tinca) (Dulski et al., 2020) and wild geladas (Theropithecus gelada) (Baniel et al., 2021; Trosvik, Muinck, et al., 2018). Seasonal shifts in the gut microbiome of these species are linked to changes in the specific foods consumed, food availability or macronutrient composition of the diet (Kartzinel et al., 2019; Orkin et al., 2019; Ren et al., 2017). For example, in African great apes (Gorilla gorilla gorilla and Pan troglodytes troglodytes), the gut microbiome contains more fibre‐degrading taxa and higher cellulose degradation functional potential when individuals are consuming a leaf‐heavy diet (Hicks et al., 2018). In Chinese alligators (Alligator sinensis) and hibernating ground squirrels (Ictidomys tridecemlineatus), mucin‐degrading gut microbial taxa increase in abundance during periods of fasting, when less energy is available from dietary sources (Carey et al., 2012; Tang et al., 2019). Gut microbial community stability and resilience improve host health by reducing long periods of gut microbial dysbiosis (Allaway et al., 2020; Sommer et al., 2017). However, the ability of the microbiome to plastically respond to short‐ and long‐term changes in environmental conditions could also be important for host health and fitness. This plastic response in gut microbiome composition and function influences host phenotype, resulting in variable host physiology and behaviour in different environmental circumstances (Davidson et al., 2018; Moeller & Sanders, 2020; Stearns, 1989; West‐Eberhard, 2003).

Evidence is beginning to emerge that diet‐related gut microbial changes may buffer energy and nutrient availability. In humans, we see marked, rapid shifts in gut microbiome composition and SCFA production immediately after an ultramarathon (Grosicki et al., 2019) or in response to increases of other forms of exercise (Estaki et al., 2016; Keohane et al., 2019), that are consistent with an increase in gut microbial efficiency beneficial for host health. Gut microbial community composition shifts also have been observed in pregnant and lactating monkeys, resulting in a more metabolically efficient gut microbiome (Mallott et al., 2020; Mallott & Amato, 2018). Similarly, seasonal shifts in gut microbiome composition and SCFA production in black howler monkeys (Alouatta pigra) appear to compensate for decreases in energy intake, buffering against nutrient and energy shortfalls (Amato et al., 2015).

The extent to which changes in microbiome composition result in shifts in potential host‐relevant functions and microbial metabolism (beyond changes in SCFA production) remains unclear. The microbial functions actively being expressed vary more between individuals than either gene family or metabolic pathway presence in the microbiome, and variability in microbial gene expression appears to play a larger role in influencing host phenotypic plasticity than do changes in the taxonomic composition or functional potential of the gut microbiome (Barroso‐Batista et al., 2020; Heintz‐Buschart & Wilmes, 2018; Tanca et al., 2017). Thus, studies of gut microbiome composition and function using marker gene or metagenomic sequencing may miss important variation in actively expressed microbial functions. Metabolomics, identifying the small molecules produced during metabolism, allows us to examine how microbial metabolism responds to changes in nutrient and energy intake in the host (Bäckhed & Crawford, 2010; Ursell et al., 2014). Several studies have shown that microbially associated metabolites are strongly linked to host health and metabolism and respond to dietary change in much the same way as microbiome composition (Filippis et al., 2016; Maier et al., 2017; Mchardy et al., 2013; Sharon et al., 2014). Metabolomics offers greater insight into the fine‐scale plasticity in microbial metabolism that acts as a buffering agent against nutrient and energy shortfalls. Investigating how the metabolome responds to changing diet and/or energetic needs in wild systems helps us understand the role of the gut microbiome in buffering energetic or nutrient shortfalls in animals; however, few studies have used metabolomics to examine gut microbial functional shifts in wild animals (Garber et al., 2019; Gomez et al., 2015, 2016).

To determine how host–microbe cometabolic processes dynamically respond to seasonal changes in diet and nutrient intake, we examined the metabolome of a black howler monkey population experiencing seasonal changes in feeding behaviour, nutrient intake, gut microbiome composition and SCFA production, as reported in our previous research (Amato et al., 2014, 2015, 2017). This population experiences three distinct seasons—Wet, Fruit‐Dominated (WFD), Dry, Leaf‐Dominated (DLD), and Dry, Fruit‐Dominated (DFD)—each characterized by distinct energy and nutrient intake profiles (Table 1). The intake of protein, energy, lipids, neutral detergent fibre and nonstructural carbohydrates was highest in the WFD season and lower in the DLD and DFD seasons. Therefore, we found this population optimal for addressing three research questions: (i) Does the metabolome of black howler monkeys vary seasonally? (ii) If so, are seasonal changes in the metabolome associated with changes in gut microbial community composition and linked to changes in nutrient and energy intake? (iii) Do seasonal changes in the metabolome suggest that the gut microbiome is buffering shortfalls in energy or nutrient intake? If buffering is occurring, then when the consumption of a specific macronutrient declines, we expect to see an increase in metabolite concentrations related to those specific macronutrient metabolic pathways.

TABLE 1

Energy and nutrient intake profiles for each of the three seasons experienced by the study population of black howler monkeys. Energy consumption is expressed as energy per metabolic body weight (MBW) and nutrient intake values are expressed as grams per metabolic body weight

Season	Energy (kcal/MBW)	Protein (g/MBW)	Lipid (g/MBW)	Nonstructural carbohydrate (g/MBW)	Neutral detergent fibre (g/MBW)
Wet, Fruit‐Dominated (WFD)	High (177.4–182.9)	Intermediate to high (9.4–11.9)	High (3.2–4.2)	Intermediate to high (29.0–29.7)	High (42.3–50.3)
Dry, Leaf‐Dominated (DLD)	Intermediate (105.5–172.4)	Intermediate to high (7.7–10.0)	Low to intermediate (1.2–3.2)	Low (12.4–16.2)	Intermediate (27.2–39.9)
Dry, Fruit‐Dominated (DFD)	Low (106.6–114.5)	Low (4.7–5.3)	Low to intermediate (1.8–2.1)	Intermediate to high (16.3–17.0)	Low (21.9–24.7)

MATERIALS AND METHODS

Field data collection

Behavioural data and faecal samples for microbiome and metabolome analysis were collected from adult and juvenile black howler monkeys (N = 16 individuals) in Palenque National Park, Mexico. All juveniles included in this data set were over 1 year of age and foraging independently (Amato et al., 2014). Data collection occurred during three distinct 10‐week sampling blocks: September–November 2010 (WFD), January–March 2011 (DLD) and April–June 2011 (DFD). Activity and dietary data were collected during 20‐min instantaneous focal individual samples recorded 5 days per week. Nutrient consumption was estimated from observational feeding data and published plant nutritional content values. Plant material for metabolite analysis was sampled from the 10 most consumed food items in each season and preserved in 70% methanol. Faecal samples were collected every 2 weeks from all 16 individuals and stored in 96% ethanol. Once per month, an aliquot of one faecal sample from each individual was stored in 80% methanol for faecal metabolite analysis. In five instances, two samples for metabolite analysis were collected for an individual during a single month. At least one gut microbiome and gut metabolome sample was available for each season for all individuals included in the analysis. For more detailed data collection methods, see related publications (Amato et al., 2014, 2015, 2017; Amato & Garber, 2014).

16S rRNA gene sequencing and analysis

DNA was extracted from 115 biweekly faecal samples from eight of the 16 individuals (Mallott & Amato, 2021a), and the V1–V3 region of the 16S rRNA gene was amplified and sequenced as described previously (Amato et al., 2014, 2015). Sequences were re‐analysed in qiime2 (2019.4) using the dada2 algorithm, adjusted for use with 454 sequencing data, and taxonomy was assigned using a Naïve Bayesian classifier and the Greengenes 13.8 database.

Metabolite extraction and analysis

Polar and nonpolar metabolites were extracted separately from 81 monthly faecal samples from 16 adult and juvenile individuals (Mallott & Amato, 2021c) and 30 plant samples (Mallott & Amato, 2021b) following previously published protocols (Amato et al., 2017; Garber et al., 2019; Gomez et al., 2015). Briefly, metabolites were extracted in 1 ml of 70% methanol and sonicated. The resulting lysed cell pellets were fractionated with 5 ml 70% methanol and chloroform, centrifuged, and extracts of both fractions were evaporated under vacuum at −60°C. Untargeted metabolomics were performed on a GC/MS system. Spectra were compared to electron impact mass spectrum libraries NIST08 (NIST) and W8N08 (Palisade Corporation), and a library of 520 unique metabolites custom‐built by the University of Illinois Urbana‐Champaign Metabolomics Center. Data were normalized to an internal standard and sample dry weight, and all metabolite concentrations were reported relative to hentriacontanoic acid per gram dry weight. While the metabolomic methods used here may be biased against detecting lipids, we expect that bias to be consistent across the sampling period and therefore not affect seasonal comparisons in metabolite abundance. Metabolite extraction, metabolomics and chromatogram processing were carried out at the University of Illinois Urbana‐Champaign Metabolomics Center in the Roy J. Carver Biotechnology Center.

We examined differences in the metabolome across seasons using permutational multivariate analysis of variance (PERMANOVA) of a Bray–Curtis distance matrix of the entire set of metabolites detected in the vegan and pairwiseAdonis packages (Martinez Arbizu, 2017; Oksanen et al., 2019) in r (r‐project.org). Additionally, metaboanalyst 4.0 (www.metaboanalyst.ca) was used to examine variation in the metabolome across seasons. Only metabolites that were found in >50% of samples were analysed. After removing samples with >90% missing values, we calculated descriptive multivariate statistics and conducted a pathway enrichment analysis. For the pathway enrichment analysis, samples were normalized by sum and data mean‐centred and divided by the standard deviation of each variable. The SMPDB (https://smpdb.ca/) metabolite library was used to group metabolites into pathways. Pathways were considered to be more abundant in one season compared with another if two or more metabolites were detected in the pathway and all detected metabolites were more abundant. The p‐values were Holm‐adjusted to account for multiple comparisons for the pathway enrichment analysis.

Microbiome and metabolite interactions

Monthly averages for each individual were calculated for microbiome composition and matched with the monthly faecal metabolome data (eight individuals, 44 matched data points; for individuals with two faecal metabolome samples in a given month, an average value was calculated). Seasonal averages for each individual were calculated for comparisons between plant metabolome and faecal metabolome data (13 individuals, 39 matched data points). As seasonal averages were used for the comparison of plant and faecal metabolomes, this analysis was performed across the entire study period, not for each season. A Mantel test using a Spearman correlation method was used to compare mean monthly individual microbial community composition at the ASV (Amplicon Sequence Variant), genus, family and phylum levels with metabolite composition in the vegan package (Oksanen et al., 2019) in r. We used the ccrepe package (Bielski & Weingart, 2021) in r to construct co‐occurrence networks both within and between the microbiome and faecal metabolome, as well as between the plant metabolome and faecal metabolome. Compositionally corrected Spearman correlation matrices were constructed for microbial ASVs, genera, families and phyla present in >25% of samples and metabolites present in >50% of samples. Networks were graphed from significant positive correlations (q < 0.05 and rho >0.5) and network attributes were calculated in cytoscape 3.8.0 (cytoscape.org). Within cystoscape, the clusterMaker2 app was used to identify highly connected clusters. Linear models were used to examine if edge betweenness differed seasonally in r. Seasonal differences in edge betweenness for specific microbe–metabolite interactions were examined using false discovery rate (FDR)‐corrected Kruskal–Wallis tests in r.

RESULTS

Metabolome composition

We positively identified 282 compounds in the metabolome of black howler monkeys (Mallott & Amato, 2021d). Of the compounds found in >50% of samples, 38.9% were lipids, 14.6% were amino acids, 13.2% were carbohydrates, 11.1% were organic acids, 9.7% were benzenoids, 4.2% were phenylpropanoids and 8.3% were other classes of compounds. Palmitic acid was the most abundant metabolite (14.9 ± 9.5%), followed by stearic acid (14.4 ± 6.0%), beta‐amyrin (9.3 ± 7.2%), vaccenic acid (6.3 ± 3.4%) and glycerol (4.5 ± 2.5%).

Seasonal differences in metabolites

Faecal metabolite profiles, including rarely detected metabolites, differed significantly between seasons (PERMANOVA, F = 14.329, R 2 = .269, p < .001). Pairwise comparisons showed that each season was distinct (all p adj < .001). Partial least squares discriminant analysis (PLS‐DA) of a subset of frequently occurring metabolites (>50% of samples) did not find a significant separation between seasons (p = .106) (Figure 1).

FIGURE 1

Plot of the first and second components from the PLS‐DA analysis. Red triangles denote data points from the DFD season, green pluses denote data points from the DLD season and blue crosses denote data points from the WFD season. Ellipses encircle 95% confidence intervals

Pathway enrichment analysis found 48 pathways that were differentially enriched between the WFD and DLD seasons (39 with two or more metabolites detected in the pathway), 12 differentially enriched pathways between the DLD and DFD seasons (nine with two or more metabolites detected in the pathway), and 24 pathways that were differentially enriched between the WFD and DFD seasons (19 with two or more metabolites detected in the pathway) (Table 2). Histidine metabolism (both catabolism and anabolism), the malate–aspartate shuttle, nicotinate and nicotinamide metabolism and the degradation of several essential amino acids (lysine, valine, leucine, and isoleucine) were enriched during the DLD season compared with both the WFD and DFD seasons. While several pathways—aspartate metabolism, bile acid biosynthesis, citric acid cycle, mitochondrial electron transport chain, phospholipid biosynthesis, and the metabolism of essential and nonessential amino acids (arginine, proline, citrulline, cysteine, glutamate, phenylalanine, tyrosine, valine, leucine and isoleucine)—were enriched in the DFD season compared with the WFD season, no pathways were uniquely enriched in the DFD season compared to both the WFD or DLD seasons. Similarly, no pathways were uniquely enriched in the WFD season, though fructose and mannose degradation were enriched in the WFD season compared with the DLD season.

TABLE 2

Pathways that were differentially enriched between each pair of seasons and the associated holm‐corrected p‐value, and the ratio of detected metabolites in the pathway to total pathway metabolites

Metabolic pathway	WFD vs. DLD	DLD vs. DFD	DFD vs. WFD	No. of metabolites detected
Alanine metabolism	.005	.195	1.000	5/17
Alpha linolenic acid and linoleic acid metabolism	1.000	< .001	.187	1/19
Amino sugar metabolism	.002	.024	1.000	4/33
Ammonia recycling	.005	.123	1.000	6/32
Arachidonic acid metabolism	< .001	.005	1.000	1/69
Arginine and proline metabolism	< .001	.038	.004	8/53
Aspartate metabolism	< .001	.151	.002	4/35
Beta oxidation of very long chain fatty acids	1.000	.480	.933	5/17
Beta‐alanine metabolism	.005	.273	1.000	5/34
Betaine metabolism	.059	1.000	.216	1/21
Bile acid biosynthesis	.005	1.000	.001	4/65
Biotin metabolism	.039	1.000	1.000	1/8
Butyrate metabolism	1.000	.412	1.000	1/19
Cardiolipin biosynthesis	.003	1.000	.007	1/11
Carnitine synthesis	.134	.893	1.000	3/22
Catecholamine biosynthesis	1.000	1.000	1.000	1/20
Citric acid cycle	< .001	1.000	< .001	4/32
Cysteine metabolism	< .001	.139	1.000	3/26
d‐Arginine and d‐ornithine metabolism	1.000	1.000	.814	1/11
De novo triacylglycerol biosynthesis	.003	1.000	.007	1/9
Fatty acid biosynthesis	.057	1.000	< .001	6/35
Fatty acid elongation in mitochondria	.006	1.000	< .001	1/35
Fatty acid metabolism	.006	1.000	< .001	1/43
Folate metabolism	< .001	.005	1.000	1/29
Fructose and mannose degradation	.013	.363	1.000	3/32
Galactose metabolism	1.000	.005	.005	8/38
Gluconeogenesis	.003	1.000	1.000	5/35
Glucose–alanine cycle	.002	.146	1.000	4/13
Glutamate metabolism	.002	.129	1.000	8/49
Glutathione metabolism	.002	.141	1.000	5/21
Glycerol phosphate shuttle	.001	1.000	.008	3/11
Glycerolipid metabolism	< .001	.058	< .001	5/25
Glycine and serine metabolism	.020	.146	.995	10/59
Glycolysis	.006	1.000	1.000	3/25
Histidine metabolism	.002	.023	1.000	2/43
Homocysteine degradation	1.000	.224	1.000	2/9
Inositol metabolism	.005	1.000	< .001	3/33
Inositol phosphate metabolism	.005	1.000	.001	2/26
Ketone body metabolism	1.000	.412	1.000	1/13
Lactose degradation	1.000	1.000	1.000	2/9
Lactose synthesis	.027	1.000	1.000	2/20
Lysine degradation	.001	.027	1.000	3/30
Malate–aspartate shuttle	.001	.027	1.000	3/10
Methionine metabolism	.127	.347	.821	6/43
Mitochondrial beta‐oxidation of long chain saturated fatty acids	1.000	1.000	1.000	1/28
Mitochondrial beta‐oxidation of medium chain saturated fatty acids	1.000	1.000	1.000	1/27
Mitochondrial beta‐oxidation of short chain saturated fatty acids	1.000	1.000	1.000	1/27
Mitochondrial electron transport chain	< .001	1.000	< .001	3/19
Nicotinate and nicotinamide metabolism	.001	.011	1.000	4/37
Nucleotide sugar metabolism	1.000	1.000	1.000	1/20
Oxidation of branched chain fatty acids	1.000	.412	1.000	1/26
Pantothenate and CoA biosynthesis	.267	1.000	1.000	1/21
Pentose phosphate pathway	< .001	1.000	.624	1/29
Phenylacetate metabolism	.752	1.000	1.000	1/9
Phenylalanine and tyrosine metabolism	< .001	.180	.018	4/28
Phosphatidylcholine biosynthesis	1.000	.115	1.000	1/14
Phosphatidylethanolamine biosynthesis	1.000	.008	1.000	2/12
Phosphatidylinositol phosphate metabolism	.664	1.000	.002	1/17
Phospholipid biosynthesis	.026	.444	.019	2/29
Phytanic acid peroxisomal oxidation	1.000	.412	1.000	1/26
Plasmalogen synthesis	< .001	1.000	.001	3/26
Porphyrin metabolism	1.000	1.000	1.000	1/40
Propanoate metabolism	.005	.123	1.000	7/42
Purine metabolism	< .001	.128	.034	10/74
Pyrimidine metabolism	.059	1.000	.355	5/59
Pyruvaldehyde degradation	.063	1.000	1.000	1/10
Pyruvate metabolism	.006	1.000	1.000	3/48
Riboflavin metabolism	1.000	1.000	.784	1/20
Selenoamino acid metabolism	.010	.153	.207	5/28
Spermidine and spermine biosynthesis	1.000	1.000	1.000	1/18
Sphingolipid metabolism	1.000	.151	1.000	3/40
Starch and sucrose metabolism	1.000	.123	.187	4/31
Steroid biosynthesis	< .001	1.000	< .001	4/48
Steroidogenesis	1.000	1.000	1.000	1/43
Threonine and 2‐oxobutanoate degradation	1.000	1.000	1.000	2/20
Thyroid hormone synthesis	1.000	1.000	1.000	1/13
Transfer of acetyl groups into mitochondria	.315	1.000	1.000	3/22
Trehalose degradation	1.000	1.000	1.000	1/11
Tryptophan metabolism	.043	.141	1.000	3/60
Tyrosine metabolism	< .001	.061	.003	6/72
Ubiquinone biosynthesis	1.000	1.000	1.000	1/20
Urea cycle	< .001	.163	< .001	7/29
Valine, leucine and isoleucine degradation	.026	.021	1.000	5/60
Vitamin B6 metabolism	.007	1.000	.784	1/20
Warburg effect	< .001	.227	< .001	8/58

Notes: Cells highlighted in green are pathways with metabolites that are consistently higher in abundance in the first season listed, red cells are pathways with metabolites that are consistently lower in abundance in the first season listed and cells highlighted in yellow are pathways with no consistent direction of effect.

Metabolite and microbiome co‐occurrence networks

Few microbial genera were found to co‐occur. A small cluster of genera were significantly positively correlated when examining combined data from all three seasons—Desulfovibrio, Moglibacterium, Sutterella, the p.75a5 genus within Erysipelotrichaceae, an unclassified genus within Moglibacteriaceae, an unclassified genus within Synergistaceae and an unclassified genus within Veillonellaceae (Figure 2). Microbial association networks constructed using samples from each season yielded three (WFD), two (DLD) or no (DFD) significantly positively correlated taxa.

FIGURE 2

Faecal microbial interaction network calculated from all data. Nodes that were significantly positively correlated after FDR correction (q < .05 and rho > .5) are shown. Node colour denotes mean shortest path length (lowest = yellow, highest = purple) and node border width increases with higher values of betweenness centrality. Edge colour indicates edge betweenness values (lowest = purple, highest = blue) and edge width increases with higher values of rho

Faecal metabolite networks varied seasonally in both their connectivity and stability. The metabolite network during the WFD season had the lowest network connectivity—low mean shortest path length, edge betweenness, betweenness centrality and network density—which suggests that there are few interconnected metabolic processes during this season (Table 3). The WFD season metabolite network had a low mean clustering coefficient and a low network density, indicating intermediate network stability. The DLD season had intermediate to high network connectivity and high network stability, indicating that metabolite networks during this season might be robust to perturbations despite being somewhat interconnected. The DFD season metabolite networks had the highest network connectivity and lowest network stability, suggesting that many metabolic processes are intertwined during this season but that the metabolic processes present are ephemeral or easily perturbed.

TABLE 3

Network attributes for metabolites present in >50% of samples within each season

Season	Mean shortest path length	Mean edge betweenness	Mean betweenness centrality	Mean clustering coefficient	Network density	Number of clusters	Mean cluster size
WFD	1.382 ± 0.493	9.103 ± 7.731	0.077 ± 0.235	0.168 ± 0.290	0.058	4	6.25 ± 7.361
DLD	1.950 ± 0.789	18.826 ± 19.185	0.116 ± 0.223	0.167 ± 0.260	0.118	6	4.667 ± 5.086
DFD	2.839 ± 1.145	62.015 ± 65.232	0.064 ± 0.128	0.284 ± 0.352	0.059	8	8.375 ± 15.702

The metabolites that were most important in the network also varied seasonally. During the WFD season, metabolites related to ketone body metabolism and lipid metabolism had the highest betweenness centrality, suggesting that these metabolites are important to network functioning and communication (Figure 3). High edge betweenness values (important connections with high rates of information transfer) were found between several essential and nonessential amino acids, as well as between nucleic acids during the WFD season. During the DLD season, while some metabolites with high betweenness centrality were related to important metabolic pathways, such as propanoate metabolism and lipid metabolism, several highly connected metabolites with both high betweenness centrality and edge betweenness were tannins and other antifeedant compounds potentially harmful to the howlers (e.g., epicatechin and protocatechuic acid). Metabolites that played a central role in the DFD season were typically involved in lipid metabolism and the metabolism of sugars. High edge betweenness values were found between metabolites involved in glycolysis, the citric acid cycle and lipid metabolism.

FIGURE 3

Faecal metabolite interaction networks calculated from all data (a), samples from the WFD season (b), DLD season (c) and DFD season (d). Nodes that were significantly positively correlated after FDR correction (q < .05 and rho > .5) are shown. Node colour denotes mean shortest path length (lowest = yellow, highest = purple) and node border width increases with higher values of betweenness centrality. Edge colour indicates edge betweenness values (lowest = red, highest = blue) and edge width increases with higher values of rho

Associations between metabolites and microbial genera

Microbial community composition and faecal metabolite composition were not significantly associated (Mantel test; ASV: r = .082, p = .154; genus: r = .072, p = .196; family: r = .050, p = .265; phyla: r = −.045, p = .702). We did not find any significant correlations between individual microbial ASVs and metabolites, microbial genera and metabolites, microbial families and metabolites, or microbial phyla and metabolites after adjusting p‐values for multiple comparisons (all q > .05) (Tables S1–S4). Without FDR correction, we found one highly connected cluster containing Streptococcus, Desulfovibrio, an unclassified genus within Synergistaceae, an unclassified genus within Burkholderiales, an unclassified genus within Erysipelotrichaceae, an unclassified genus within Alphaproteobacteria, and several metabolites involved in lipid metabolism, fatty acid metabolism and nonessential amino acid biosynthesis (Figure 4).

FIGURE 4

Faecal metabolite–microbe interaction network calculated from all data points. Nodes that were significantly positively correlated (p < .05 and rho > .5) are shown. Node colour denotes mean shortest path length (lowest = yellow, highest = purple) and node border width increases with higher values of betweenness centrality. Edge colour indicates edge betweenness values (lowest = red, highest = blue) and edge width increases with higher values of rho

The strength of associations between microbial genera and metabolites varied seasonally. Edge betweenness differed significantly between seasons (Linear model; F = 1147.2, p < .001). The DLD season had significantly higher edge betweenness than either the DFD season (Tukey; t = 33.03, p < .001) or the WFD season (Tukey; t = 46.19, p < .001). The DFD season had significantly higher edge betweenness than the WFD season (Tukey; t = 17.35, p < .001).

Interactions between ingested plant metabolites and faecal metabolites

When examining associations between average seasonal consumption of plant metabolites and seasonal averages of faecal metabolites for each individual, we found significant positive associations (all q < .05 and rho > .5) between 64 pairs of ingested plant metabolites and faecal metabolites. Network analyses indicate that there were two highly connected clusters of plant and faecal metabolites, and several smaller clusters (Figure 5). The two highly connected clusters primarily contained sugars and fatty acids from both plant sources and faecal targets. This is probably linked to fatty acid biosynthesis and related metabolic functions. Here, we found plant secondary metabolites as faecal targets, including a plant diterpenoid (dehydroabietic acid) that significantly correlated with several sugars and fatty acids, and an organosulphur intermediate that also connected closely with sugars on the lower cluster. We also found 1,2,3‐trihydroxybutylbenzene as a faecal target linked to the plant source 1,3‐di‐tert‐butylbenzene, also known as pyrogallol. The other highly connected cluster contained plant sugars, long‐chain fatty acids, and faecal metabolites related to the metabolism of essential nutrients. Here, we found epicatechin as faecal targets, closely connected with fatty acids and a glycoside as plant sources.

FIGURE 5

Plant metabolite‐faecal metabolite interaction network calculated from all data points across all seasons. Nodes that were significantly positively correlated (p < .05 and rho > .5) are shown. Node colour denotes mean shortest path length (lowest = yellow, highest = purple) and node border width increases with higher values of betweenness centrality. Edge colour indicates edge betweenness values (lowest = red, highest = blue) and edge width increases with higher values of rho. Direction of arrow is from source (plant metabolite) to target (faecal metabolite)

DISCUSSION

In this study, we examined how seasonal variation in macronutrient intake influences the gut metabolome in wild black howler monkeys. We found that the gut metabolome varied seasonally in our study population, similar to previous studies (Amato et al., 2015). We found partial support for the hypothesis that changes in the gut microbiome and corresponding changes in the metabolome buffer seasonal energy and/or nutrient shortfalls. Additionally, we found strong associations between faecal metabolites and ingested plant metabolites. However, we did not find significant associations between the metabolome and microbiome composition.

The faecal metabolome of black howler monkeys is dominated by lipids, amino acids, carbohydrates and organic acids, similar to other species of nonhuman primates (Garber et al., 2019; Gomez et al., 2015, 2016; Ni et al., 2021). Some of the long‐chain fatty acids prominent in the black howler monkey metabolome—palmitic acid, stearic acid and vaccenic acid—are associated with high‐fat diets in mice (Daniel et al., 2014) and fruit‐dominated diets in lowland gorillas (Gomez et al., 2016). Palmitic acid and stearic acid are the most common saturated fatty acids in nature while vaccenic acid is one of the most common unsaturated fatty acids (Sommerfeld, 1983). For primates that principally consume plant foods, such as black howler monkeys, the majority of lipid intake comes from fruits and their seeds (Norconk et al., 2009). Given that this population of black howlers has not been observed to consume large amounts of seeds and commonly voids undigested seeds in their faeces, the fatty acid profiles of their faecal metabolomes are likely to be strongly influenced by a yearly diet composed of 57.3% fruits (per cent dry weight) (Amato & Garber, 2014).

We found significant seasonal differences in faecal metabolome profiles, both in composition and in pathway enrichment. Several metabolic pathways were differentially enriched between seasons, with lower concentrations of metabolites present during the WFD season when compared to both the DLD and DFD seasons. Additionally, metabolite network structure varied across seasons. The metabolite network was highly connected but less stable in the DFD season, when black howler monkey diets are calorically and nutritionally more constrained. In contrast, during the WFD season when black howler monkeys are consuming more energy‐ and nutrient‐rich diets, the metabolite network was diffuse. The most stable metabolite network occurred during the DLD, when most nutrients were being consumed at intermediate rates. In addition, edge betweenness in the microbial–metabolite interaction networks was highest in the DLD season, indicating that consistent, strong connections between individual bacteria and metabolic products may be contributing to metabolite network stability in this season. These patterns suggest that metabolic cross‐feeding may be more necessary for the gut microbial community when nutrients from the host diet are less readily available in the gut. Because cross‐feeding can increase the metabolic efficiency and/or ecological stability of the microbial community (Coyte et al., 2015; Coyte & Rakoff‐Nahoum, 2019; D'Souza et al., 2018; Evans et al., 2020; Goldford et al., 2018; Gudelj et al., 2016; Liu & Sumpter, 2017; Smith et al., 2019), our data provide preliminary evidence of improved nutritional buffering by the gut microbiome during the DLD. During the DFD season, when metabolite network connectivity was high but stability was low, metabolic cross‐feeding may be taking place, but the cross‐feeding relationships are not as consistant or stable over time.

Microbial responses to variation in the intake of specific nutrients over time also provide evidence of nutritional buffering by the gut microbiome. For example, given that leaves tend to be high in fibre and low in host‐metabolizable energy (Norconk et al., 2009), we expected to see enrichment of pathways related to structural carbohydrate and lipid metabolism during the DLD, as microbes degrade fibre to produce SCFAs. While the pathway enrichment analysis did not show an increase in carbohydrate and lipid metabolism, we found that metabolites related to lipid metabolism and SCFA metabolism became more important in metabolite networks during the DLD season. In addition, the metabolism of essential vitamins (B3) and amino acids were enriched in the DLD season, potentially compensating for a diet poor in specific nutrients and aiding in the digestion of protein‐rich leaves. Similarly, although fruits have more host‐metabolizable energy than leaves, they are lower in protein compared to leaves, and overall caloric intake during the DFD season was reduced by 37.4%–40.0% compared with the WFD season and 0%–33.6% compared with the DFD season (Table 1). Therefore, we expected microbial buffering to result in an enrichment of amino acid synthesis pathways and lipid metabolism pathways. We did see an enrichment of amino acids pathways during the DFD and lipid metabolites having a more central role in the metabolite network (Figure 3). In contrast, we observed a de‐enrichment of many pathways during the WFD season when the black howler monkey diet was least nutritionally constrained. These results confirm earlier work in this population (Amato et al., 2014, 2015), as well as other studies in mammals (Gomez et al., 2015; Koren et al., 2012; Mallott & Amato, 2018; Springer et al., 2017; Sun et al., 2016; Wu et al., 2017), that indicate the gut microbiome acts as a potential buffer limiting energy and nutrient shortfalls due to seasonal changes in diet or changes in host nutrient requirements.

The gut microbiome also provides nutritional benefits to hosts by processing plant secondary metabolites that otherwise act as toxins or digestive inhibitors. This relationship has been documented in desert woodrats (Dearing & Kohl, 2017; Kohl et al., 2014, 2016) and may be important to highly folivorous primates such as black howler monkeys, whose yearly diet contains 33.1% young and mature leaves (% dry weight) (Amato & Garber, 2014). Interestingly, our plant‐faecal metabolite interaction network from all seasons combined showed faecal target metabolites associated with essential nutrients, fatty acid metabolism, and the metabolism of chemical defensive compounds clustered with two groups of plant source metabolites: simple carbohydrates and fatty acids. Specifically, the plant secondary metabolites included dehydroabietic acid in the upper largest cluster, a diterpenoid for chemical defence commonly found in tree resin (Helfenstein et al., 2017), and epicatechin, a flavan‐3‐ol and major component of condensed tannin (Ferreira et al., 1999; Khanbabaee & Ree, 2001), in the right smaller cluster. This suggests that individuals ingested toxic plant secondary metabolites while consuming sugar and fatty acid metabolites used for fatty acid/glucose‐related biosyntheses. This relationship could indicate a potential trade‐off in that foods they consume with the highest protein and caloric value also contain high amounts of tannins and other difficult to digest plant secondary metabolites. While some mammals might avoid food items high in particular plant secondary metabolites, evidence suggests that other mammals readily consume foods high in secondary metabolites if those foods are also high in energy, protein or water (Felton et al., 2009; Lambert & Rothman, 2015; Remis et al., 2001; Simpson & Raubenheimer, 2001; Villalba & Provenza, 2005). The presence of tannin‐ or toxin‐degrading bacteria in the gut microbiome could facilitate this behaviour by allowing animals to tolerate higher concentrations of plant secondary metabolites in their diet. We also identified an aromatic thiol (organosulphur compound), positioned in the larger cluster of sugar/fatty acid metabolites along with dehydroabietic acid, that is either derived from plant leaves (Gonulalan et al., 2019), from soil (Shen et al., 2020), or is an intermediate byproduct of another aromatic hydrocarbon (such as butylbenzene) that was degraded by sulphur‐reducing gut bacteria (such as Desulfovibrio, which degrades hydrocarbons with a sulphate redox reaction) (Lyles et al., 2014; Widdel et al., 2006, 2010) after uptake in the diet. Either way, the presence of this aromatic thiol might be another indicator of how the black howler monkey gut microbiome facilitates a higher tolerance for increased amounts of plant secondary metabolites during dietary shifts. Mechanistic experiments testing the capacity of the howler monkey gut microbiome to degrade these potential plant toxins will be necessary to verify this relationship.

Although our data suggest that the gut microbiome buffers wild howler monkey hosts against nutritional challenges, they also indicate that these microbial services are likely to have limits. Our metabolite networks indicated increased cross‐feeding and community stability in response to nutritional constraints during the DLD compared to the WFD. However, during the DFD season, as black howler monkey diets became even more calorically and nutritionally constrained, the stability of the metabolite networks began to decline. This trend suggests there may be a threshold past which specific dietary changes alter the underlying structure of the microbiome in a way that compromises the nutritional services the microbiome provides to the host. For example, reductions in the intake of specific macronutrients consumed in large amounts such as fibre, or reductions that persist for extended time periods could lead to the loss of key microbial taxa (Sonnenburg et al., 2016). Studies of black howler monkeys in anthropogenically altered environments report reduced microbial diversity as well as reduced relative abundances of SCFA‐producing microbial taxa (Amato et al., 2013). These losses in microbial taxa are correlated with losses of plant species, and presumably particular macronutrients, from the howler diet, and are likely to inhibit nutritional buffering by the microbiome. Identifying the dietary thresholds past which the gut microbiome loses its ability to buffer hosts from nutritional shortfalls in black howler monkeys as well as other wild animals will provide important insight into host ecology and conservation as well as microbial community dynamics.

Data describing howler monkey physiology are critical for verifying the magnitude of impact of potentially beneficial microbiome functions on hosts in variable environments. While some changes in microbial metabolism could benefit the host by increasing the availability and subsequent absorption of nutrients or key vitamins lacking in the diet, they might also be harmful to the host if scarce nutrients are diverted into networks of gut microbial metabolism. Alternatively, the impacts of changes in microbial metabolism on the host might be negligible. These potential outcomes probably occur along a gradient and are context‐dependent. Evaluating host nutrient and energy balances using noninvasive markers as well as performing controlled laboratory experiments to assess microbial metabolism and interactions with host dietary substrates in real time will provide crucial insight.

Regardless of the extent to which the observed shifts in microbial metabolism buffer host nutrition, our data suggest they are strongly driven by host diet. Although previous research has shown that nutrient limitations in the large intestine resulting from dietary changes influence microbial community composition (Reese et al., 2018), wild animals systems such as this one make it difficult to test the extent to which host‐driven changes in the intestinal environment result in preferential recruitment of specific microbial taxa or genes. However, it is widely accepted that changes in host diet alter the nutritional environment in the gut such that microbes will differentially regulate their metabolic pathways (Fontaine & Kohl, 2020). This alters both competitive and mutualistic interactions between microbes, as suggested by the shifting networks of microbial taxa and metabolites that we observed across seasons. Furthermore, one of the strongest relationships we detected was between the faecal metabolome and the metabolite content of the plants consumed by black howler monkeys. For example, during the DLD season, we saw an increase in the centrality and importance of tannins and other plant secondary metabolites in the faecal metabolite networks that is probably related to the higher concentrations of these compounds in leaves compared with fruit. During the DLD, the consumption of mature leaves increased to 64% (per cent dry weight) compared with 37% during the WFD season and 43% during the DFD season (Amato et al., 2015; Amato & Garber, 2014). While mature leaves probably contribute more tannins to howler monkey diets compared with young leaves, rates of mature leaf consumption were found to be relatively low and constant in this population (1%–8% dry weight) (Amato & Garber, 2014). We also found a high number of significant direct associations between plant and faecal metabolites. Previous studies demonstrated similarly strong relationships between primate dietary intake and gut microbiome composition and function using DNA‐based approaches (Amato et al., 2015; Mallott et al., 2018; Orkin et al., 2019). Fewer have identified these relationships with faecal metabolites (Garber et al., 2019; Gomez et al., 2015, 2016).

In conclusion, we found strong relationships both between seasonal changes in diet and gut microbiome function and between the consumption of specific plant metabolites and faecal metabolite profiles. These patterns suggest that the gut microbiome might be buffering howler monkeys against seasonal variations in nutrient intake. However, we also identified evidence of a potential threshold in dietary intake past which the ability of the gut microbiome to buffer howler monkeys could be diminished. Moving forward, combining detailed studies of nutrient consumption and data on gut microbial community composition and function with biomarkers of host energy status and physiology will help to clarify the magnitude of these potential benefits and limitations and more precisely identify relevant mechanisms of interaction between both hosts and microbes as well as among different microbes.

AUTHOR CONTRIBUTIONS

E.K.M. analysed data and wrote the paper. L.J.S. analysed data and edited the paper. P.A.G. provided funding for metabolomics analysis and edited the paper. K.R.A. designed the study, provided funding for field data collection, performed data collection and edited the paper.

CONFLICT OF INTEREST

The authors have not conflict of interest to declare.

OPEN RESEARCH BADGES

This article has earned an Open Data Badge for making publicly available the digitally‐shareable data necessary to reproduce the reported results. The data is available at https://doi.org/10.6084/m9.figshare.14995407.v1

Figure S1 Scatterplots of pairwise Bray‐Curtis dissimilarities of microbial communities compared with those of metabolomes at each taxonomic level (ASVs, genera, families, and phyla).

Click here for additional data file.

Table S1 Compositionally corrected Spearman correlations between faecal metabolites and bacterial ASVs.

Table S2. Compositionally corrected Spearman correlations between faecal metabolites and bacterial genera.

Table S3. Compositionally corrected Spearman correlations between faecal metabolites and bacterial families.

Table S4. Compositionally corrected Spearman correlations between faecal metabolites and bacterial phyla.

Click here for additional data file.