Microbiota assembly, structure, and dynamics among Tsimane horticulturalists of the Bolivian Amazon

Samples and subjects

To explore whether ICABs are associated with patterns of microbial colonization in young children, stool samples (reflecting the distal gut) and swab samples from the dorsum of the tongue were collected from 47 Tsimane families living in six villages located along the Maniqui River in the Bolivian lowlands of the Amazon basin. Samples were collected from infants (0–2 years of age) and mothers (14 or more years of age) using a mixed longitudinal design (Supplementary Table 1, Fig. 1a).

Infant and maternal microbiota dynamics following birth

Overall, infant stool samples had high relative abundances of Bifidobacteriaceae, Enterobacteriaceae, and Veillonellaceae, while adult stool samples were dominated by Ruminococcaceae, Prevotellaceae, and Lachnospiraceae (Supplementary Fig. 1). On the tongue, both infants and adults had high relative abundances of Streptococcaceae, Veillonellaceae, and Micrococcaceae, although adult tongue samples also exhibited high relative abundances of Pasteurellaceae and Neisseriaceae (Supplementary Fig. 2). In addition, bacterial diversity increased with age in both stool and on the tongue of infants over the first 18 months of life (Fig. 1b, c). We used a linear mixed-effects model that accounted for the longitudinal structure of these data by treating the subject as a random effect. It demonstrated that diversity was positively correlated with age in both stool (P = 1.08 x 10⁻¹¹, conditional R² = 0.481) and tongue communities (P = .013, conditional R² = 0.397), although diversity increased at a faster rate in the stool samples (Fig. 1b, c; ANOVA, P = .004). The diversity of the maternal stool and dorsal tongue communities was stable during this time period (Fig. 1b, c, P = .065 and .099, respectively), as has been observed previously²⁸.

Variation in bacterial community composition among samples from mothers and infants was primarily explained by body site (13.7%, PERMANOVA, 1000 permutations, P < .001) and host age group (infant or adult, 10.5%, PERMANOVA, 1000 permutations, P < .001). Village type (proximal or distal to the regional market) was also a significant factor (PERMANOVA, 1000 permutations, P < .01), although it only explained 0.48% of the variation. Microbial communities colonizing the infant gut were distinct at all ages from those of the tongue and became increasingly differentiated with time (Fig. 1d), but were never as distinct as those of adults, owing largely to incomplete maturation of these infant stool microbiotas by 18 months of age. Despite their apparent immaturity, infant stool microbiotas shared significantly more bacterial taxa with those of the mother than with unrelated adults (Wilcoxon rank-sum test, P < .05, Supplementary Fig. 3). While the tongue swabs of 16–18 month olds harbored microbial communities similar in composition to adult tongue swabs (Fig. 1d), those communities were not more similar to those of their mothers than they were to those of unrelated adults (Supplementary Fig. 3), in contrast to the stool communities. This supports the findings of a previous cross-sectional study showing that the oral microbiotas of Tsimane infants (9–24 months old) were not more similar to those of their mothers than they were to unrelated adults, despite frequent maternal premastication²⁹. The number of amplicon sequence variants (ASVs) shared between body sites within infants was highest immediately after birth (33.3%), and then declined with time (Supplementary Fig. 4). However, the ASV overlap between body sites in 16–18 month olds was significantly higher than in adults (Supplementary Fig. 4, Wilcoxon rank-sum test, P = .03).

Growth faltering in the Tsimane relative to international standards might be due to pathogen burden, rather than lack of nutrition, as indicated by a much higher prevalence of height-for-age z-scores (HAZ) < −2 SD (53%), than for weight-for-age (WAZ, 16%) or weight-for-height z-scores (WHZ 9%) < −2 SD in children aged 2–5 years³⁰. Within the current study cohort of 47 Tsimane children, only one subject consistently exhibited HAZ and WAZ scores below −2 SD. Five other subjects had HAZ or WAZ scores that crossed from greater than to less than −2 SD after one year of age. All infants in the longitudinal cohort were breastfeeding at all study time points. However, among non-exclusively breastfed infants over 6 months of age, 63.8% of samples were collected at times when the infant had a minimum dietary diversity (MDD) score of < 4 (Supplementary Fig. 5), indicating potential micronutrient inadequacy in their complementary foods³¹ (see “Methods” for more details). A single linear mixed-effects model using MDD and age as predictors showed that infants’ MDD scores were positively correlated with infant stool microbiota diversity (Supplementary Fig. 5, P < .001, conditional R² = 0.374), and accounted for a larger effect than age (β_MDD = 1.2 vs. β_{Age (months)} = 0.6, P < .05 for both) on infant stool diversity. Chicha consumption, as well as vegetable and dairy consumption, was significantly associated with stool microbiota composition, even after controlling for age and village (Fig. 1e). Interestingly, consuming chicha made from manioc or corn had an opposite association with the infant stool microbiota than did consuming manioc or corn itself. Fecal neopterin, a marker of interferon-activated monocytes and macrophages, and of intestinal inflammation³², was associated with significantly lower infant stool diversity, even when controlling for age (Supplementary Fig. 5, P = .002, conditional R² = 0.333).

Potential sources of infant-colonizing microbes

To identify potential sources of infant-colonizing microbes, we first quantified the overlap between infant and adult microbiotas. 84.5% of the ASVs found in all infant stool samples at all time points were also found in at least one adult stool sample, while 78.3% of infant oral ASVs were also found in at least one adult oral sample. On average, 46.2% of the microbes colonizing an infant’s stool or tongue were also found in their mother’s samples from the same body site. This percentage increased in stool samples as the infant became older (Supplementary Fig. 6, blue line plus green line in panel a), increasing from 38.4% in 0–6 month olds to 50.5% in 12–18 month olds (Wilcoxon rank-sum test, P = .0013); in tongue swabs, the increase from 59.6% in 0–6 month olds to 63.8% in 12–18 month olds was not statistically significant (Wilcoxon rank-sum test, P = .52). However, the vast majority of these shared microbes were also found in other adults living in the same village as the infant (Supplementary Fig. 6, green line in panel a), reflecting a high degree of microbial sharing among nearby adults, similar to what has been observed in other indigenous populations³³. Mothers from the same village shared significantly more stool microbes, but not more tongue microbes, than did mothers from different villages (Supplementary Fig. 6, Wilcoxon rank-sum test, P < 0.001), suggesting that village ecology might play a role in structuring the stool microbiota. About half of the microbes found in an infant’s stool were found in at least one of their mother’s stool samples (Supplementary Fig. 6, purple line plus green line in panel c), and about half were not found in any of their mother’s samples (orange line), with very few infant stool microbes ever observed on their mother’s tongue (yellow line plus green line). The same was true for infant tongue microbes, of which the majority of infant tongue ASVs were observed on their mother’s tongue (Supplementary Fig. 6, yellow line plus green line in panel c) and very few were observed in their mother’s stool (purple line plus green line).

Even though they accounted for only about half of the taxa in the infant stool and on the tongue, the bacteria shared within a mother–infant dyad were found at higher relative abundances in both the infants and their mothers than non-shared taxa. This difference in average relative abundance was statistically significant in infants until they reached 12 months of age, after which the relative abundances of shared and non-shared taxa did not significantly differ (Fig. 2a, c). In mothers, shared taxa remained at higher abundance during the first 18 months of the infant’s life (Fig. 2a, c). At both body sites, the frequency at which a given microbe was shared within a dyad was largely a function of its relative abundance in the mother. Microbes with an average maternal abundance >0.1% were shared within a dyad at the same rate as its prevalence, while microbes at 0.01% or lower abundance in mothers were almost never shared within a dyad, even if the ASVs were highly prevalent in infants (Fig. 2b, d). Overall, microbes that were highly abundant in Tsimane mothers were more likely to be shared with their infants.

Microbiota assembly rules

Since many Tsimane ICABs may increase the opportunities for neutral dispersal of microbes from mothers to their children, we sought to quantify the relative contributions of deterministic and stochastic (i.e., neutral) processes during microbiota assembly. For this purpose, we used the neutral community model (NCM)³⁴, which has previously been applied to the assembly of many host-associated microbial communities^6,7. The NCM predicts the prevalence of each microbe given its average relative abundance in a regional species pool (RSP). Microbes that fit the prediction are inferred to have assembled neutrally from the RSP, while microbes at higher or lower prevalence are inferred to have been under local positive or negative ecological selection, respectively. However, species whose distributions fit the neutral model may still have experienced selection, if the selection were of similar direction and magnitude to the selection present in the RSP. In addition, the NCM assumes that all microbes in a RSP have an equal per capita ability to disperse to a local area, and once established, to have equal per capita fitness, growth, and death rates³⁴.

We defined the RSP as all of the microbes able to disperse to and colonize a local area, which in this case was the infant stool or tongue. In order to estimate the composition of this pool, we summed all of the local communities observed at a given body site in all of the infants in the cohort, and then assessed the degree to which the observations fit the model predictions by calculating the model’s root mean squared error (RMSE). RMSE is the square root of the mean of the squared differences between observed and predicted values, with 0 indicating that the data perfectly fit the model and higher values indicating a greater divergence from the model predictions. To account for the longitudinal sampling scheme, we randomly selected one sample per infant with 1000 permutations to calculate a bootstrapped estimate of model fit. Since the NCM does not account for deterministic factors and makes many simplifying assumptions regarding the functional equivalence of species, it should not be interpreted as a complete description of community assembly. Nevertheless, its ability to accurately predict a species’ prevalence suggests that neutral dispersal is an important force shaping microbial community structure. The NCM also serves as a useful mechanistic null model to help identify factors leading to a divergence from neutral dynamics or specific taxa that assemble in a nonneutral manner, such as strong competitors or active dispersers.

Of the 780 ASVs observed in infant stool, 659 were also found in at least one adult stool sample (Fig. 3a). The 484 ASVs observed in adult but not infant stool were inferred as being unable to colonize the infant gastrointestinal tract, either due to dispersal limitation or to negative selection. Microbes in infant stool fit the NCM with an average RMSE of 0.09 ± 0.006 (Figs. 3b, 1000 permutations, ±standard deviation). In total, 79.2% of infant stool ASVs were neutrally distributed in at least 90% of permutations, while only 4.0% were consistently under positive selection and 0.6% under negative selection. These low-prevalence or high-abundance taxa included obligate anaerobes such as Prevotella, Megasphaera, and Bifidobacterium, suggesting that oxygen exposure could be an impediment to transmission (Fig. 3c). The remaining 16.2% of bacterial taxa were variable in their predicted fit to the NCM across permutations. As a point of comparison, the NCM was also applied to adult stool samples. The RMSE of infant stool microbiotas was significantly lower than adult stool microbiotas (RMSE = 0.14 ± 0.003, P < 0.001), indicating that neutral processes play a more important role in shaping microbial community assembly in early life. Only 3 of the 31 consistently positively selected microbes were observed exclusively in infant stool (Fig. 3c), and of the remaining 28 ASVs, 21 were found at higher relative abundance in infant than in adult stools.

Similar to the communities found in stool, 166 of the 212 ASVs observed in infant tongue swabs were also found in at least one adult tongue swab (Supplementary Fig. 7). The 141 ASVs observed in adult but not infant tongue swabs were likely unable to colonize the infant tongue due to negative selection, especially since Tsimane ICABs may be facilitating dispersal. Microbes on the infant tongue fit the NCM with an average RMSE of 0.1 ± 0.003 (Supplementary Fig. 7, 1000 permutations, ±standard deviation). In total, 69.3% of infant tongue swab ASVs were neutrally distributed in at least 90% of permutations, while 10.8% were consistently under positive selection and 3.3% under negative selection. The remaining 16.6% of bacterial taxa were variable in their predicted fit to the NCM across permutations. Again, the NCM was also applied to adult tongue swab samples as a reference, and the RMSE of infant tongue swab microbiotas was significantly lower than adult tongue swab microbiotas (RMSE = 0.145 ± 0.004, P < 0.001), indicating that neutral processes play an important role in shaping microbial community assembly in early life in the oral cavity as well. Of the 23 consistently positively selected microbes, only 5 were observed exclusively in infant tongue swabs (Supplementary Fig. 7).

Infant microbiota assembly in a global context

We next sought to determine whether the assembly patterns that we observed in the Tsimane were conserved in infants from diverse cultures. We integrated our data with two publicly available datasets generated with comparable molecular methods for profiling the stool microbiota. The first of these datasets was created in a study of undernourished infants in Bangladesh²⁶ (we used the well-nourished control group only), and the second in a large observational study of type 1 diabetes in eastern Europe²⁷ (we used the Finnish infants only). The age ranges for each infant cohort were 0.2–19 months for Bolivia (Tsimane), 0–24 months for Bangladesh, and 1.2–37 months for Finland. Despite drastically different social and physical environments, access to health care, and diets, the stool microbiotas of the infants in these three studies were similar in composition at the earliest ages based on a range of distance metrics, indicating similar initial microbial colonists in societies across the globe (Fig. 4a, Supplementary Fig. 8). Based on the Jaccard similarity index and the Bray–Curtis dissimilarity, the microbiota trajectories diverged as the children aged, becoming more distinct as they approached a population-specific adult-like state (Fig. 4a, Supplementary Fig. 8). However, both weighted and unweighted Unifrac, which are phylogenetically informed and tend to emphasize similarities at higher taxonomic levels, revealed that microbiota developmental trajectories are much more consistent between populations (Supplementary Fig. 8).

In order to determine if neutral processes operate in a similar manner among infants in different cultural settings and geographies, we divided each dataset into three-month age intervals up to 1 year of age, and then reapplied the NCM model using only the samples from that interval as the estimate of the RSP, building a total of 12 additional models (Supplementary Fig. 9). These datasets exhibited striking consistency in the role of neutral and nonneutral processes across both age and geography, with the majority of microbes exhibiting patterns of neutral assembly across individuals (Fig. 4b). Furthermore, by partitioning ASVs according to country, age group, and how well they fit the neutral model, we found similar sets of neutral and nonneutral microbes across geography and ontogeny (Fig. 4c, Jaccard distance/PCoA). A heatmap displaying only positively selected microbes revealed Bifidobacterium (ASV 4271) as highly prevalent and abundant in all three cohorts (Fig. 4e). Country-specific patterns were also apparent, including four Bacteroides ASVs (5510, 5697, 5614, and 5651), Collinsella aerofaciens (4327), one Haemophilus ASV (6778), and two Veillonella ASVs (4061 and 4067) that were prevalent in Finnish and Bolivian (Tsimane), but not Bangladeshi infants, as well as two Prevotella ASVs (5359 and 5371), Lactobacillus ruminis (3618), and one Megasphaera ASV (4145) that increased in prevalence with age in both the Bolivian (Tsimane) and Bangladeshi infants, while remaining nearly absent in Finnish infants (Fig. 4e). Interestingly, Bacteroides and Prevotella have been identified as biomarkers for Western and non-Western lifestyles³⁵, respectively, yet both were prevalent in Tsimane infants.

The NCM also predicts the probability that when an individual microbe is removed from the local community, it is replaced by a microbe from the RSP rather than by reproduction from within the local community^7,34. This migration term, m, increased with age in all three populations, although it increased significantly more in the Tsimane infants (Fig. 4d). This significantly higher rate of microbial migration into the Tsimane infant gastrointestinal tract is consistent with their high exposure to maternal and environmental microbes.

Local ecology, market access, and Tsimane microbiota structure

In recent decades, the Tsimane have experienced profound changes to their culture and life history, including increased access to processed foods in markets and greater access to modern medical care. Shorter distance and ease of transport to local market towns is associated with lower infant mortality¹⁸, but also higher intake of processed food²⁴, increased obesity and metabolic syndrome³⁶, and higher levels of urinary phthalates³⁷, suggesting that increased market access can have important health impacts. In order to understand how these lifestyle shifts may be associated with the assembly of the gastrointestinal microbiota, we compared the stool microbiota structure of adults living in nine Tsimane villages. These villages belonged to two village ecotypes, ‘river’ and ‘forest’, which affect some seasonal subsistence activities (e.g., fishing) and frequency of access to market (Supplementary Table 1, Fig. 5a). During the sample collection period in river villages (September 2012–March 2013), proximal river villages (≤20 km to market) could regularly access the market via 1–3 h boat rides combined with approximately 1 h taxi trips during the dry season (approximately June–August). Distal river villages (>20 km to market) had access to market via 6–12 h trips on foot or by boat. Inhabitants of forest villages (>40 km from market) can only access market towns via 2–3 day trips on foot, or hitchhiking on infrequently passing vehicles (approximately 6 h). Samples from forest villages were collected during the dry season (July 2009), the only time roads are accessible (Supplementary Table 1).

A principal coordinates analysis on the weighted Unifrac distances among these microbiota data showed significant clustering of stool samples by ecotype of the village in which the subject lived (Supplementary Fig. 10, PERMANOVA, P = .01). Additional PCoAs based on unweighted Unifrac distance, Jaccard similarity index, and Bray–Curtis dissimilarity all showed similar results (Supplementary Fig. 10, PERMANOVA, P < .001 for all). PC1, which explained 25.5% of the variation in the data, was significantly correlated with the diversity of the microbiotas (Fig. 5b, P = .005, R² = 0.277). To evaluate whether these relationships might reflect an effect of season, local ecology, and/or more regular market access on assembly processes, we applied the NCM to stool samples from adults living in the two different village ecotypes with the regional species pool estimated from all of the samples from that age group and village ecotype, and found that the goodness of fit is lower in samples from river villages (t test, P < .003). This suggests that differences in local ecology and/or more regular market access, perhaps through increased consumption of processed foods (e.g., breads, pasta, and refined sugar) and/or medicine (e.g., antibiotics) may increase the role of ecological selection and diminish the role of stochastic dispersal during assembly of the gastrointestinal microbiota.

To identify how even subtle ecological and economic differences among villagers might be associated with stool microbiota structure, we used a phylogeny-based implementation of linear discriminant analysis. This approach uses both the leaves of the phylogenetic tree (i.e., the ASVs), as well as the nodes, to identify linear combinations of features that characterize two or more sample classes. The first model, built to discriminate between proximal and distal river villages, selected three predictors corresponding to ten ASVs. Nine of those ASVs, including two in the family Muribaculaceae (ASVs 5801 and 5805) and seven in the family Bacteroidaceae (ASVs 5614, 5639, 5697, 5510, 5620, 5651, and 5716) were discriminatory for the distal river samples, while only one ASV, in the family Prevotellaceae (ASV 5418), was discriminatory for the proximal river samples (Fig. 5c). The small number of proximal village-discriminatory features is consistent with studies suggesting that life in a high resource setting often leads to loss of bacterial diversity in the gastrointestinal tract^38,39. Despite the small number of optimum features following 10-fold cross validation, the model accurately predicted village type for 87.3% of samples (Supplementary Fig. 11).

The second model, built to discriminate between river and forest villages, identified 24 predictors corresponding to 56 ASVs, and correctly identified the village type for 91.8% of the samples (Supplementary Fig. 11). River villages were associated with ASVs in the bacterial families Prevotellaceae (18 ASVs), Lachnospiraceae (13 ASVs), Veillonellaceae (4 ASVs), Acidaminococcaceae (2 ASVs), and Ruminococcaceae (2 ASVs). Indicator taxa for Tsimane forest villages included ASVs in the bacterial families Ruminococcaceae (4 ASVs), Clostridiales Family XIII (2 ASVs), Lachnospiraceae (1 ASV), Muribaculaceae (1 ASV), Peptostreptococcaceae (1 ASV), and Christensenellaceae (7 ASVs), as well as 1 archaeal ASV in the family Methanobacteriaceae (1 ASV) (Fig. 5c).

Christensenellaceae and Methanobacteriaceae were previously identified in a study of adult twins as members of a co-occurrence network of heritable stool microbes that were also associated with body leanness^40,41. Since Tsimane have low but increasing rates of obesity²³, we sought to determine if a similar network were present in the Tsimane adults. The optimal network topology constructed from adult stool microbial communities contained a subnetwork involving the most highly connected node, ASV 12 (Methanobrevibacter), as well as five ASVs in family Christensenellaceae (Fig. 6a), indicating strong associations among these taxa that are driven by H₂ syntrophy⁴². However, we did not observe an association between the relative abundance of the taxa in this subnetwork and adult BMI. This subnetwork also included the second most highly connected node, Prevotella copri (ASV 5359). Although these two taxa frequently co-occurred, their abundances were inversely related (Fig. 6b, P = 1.04 × 10⁻¹⁰, marginal R² = 0.205). Methanobrevibacter and Prevotella have been shown to be positively correlated in the gastrointestinal tract of humans⁴³ and western lowland gorillas⁴⁴, especially in the setting of seasonal, carbohydrate-rich diets. These taxa also had distinct colonization dynamics in the gastrointestinal tract of infants from Bolivia (Tsimane), Bangladesh, and Finland. Prevotella copri, the same ASV that was abundant in Bolivian and Bangladeshi infants but essentially absent in Finnish infants, increased in abundance beginning at approximately 6 months of age, and remained at high levels throughout adulthood (Fig. 6c). Methanobrevibacter, on the other hand, was not detected in these stool microbiotas until later in adolescence (Fig. 6c). However, since Methanobrevibacter was more common in forest villages, and these villages did not include many samples from young children, this pattern may also have been a result of ecological, seasonal, or societal differences between villages.