A sparse covarying unit that describes healthy and impaired human gut microbiota development

Thirty-six members of a birth cohort with consistently healthy anthropometric scores living within the Mirpur district of Dhaka, Bangladesh, underwent monthly fecal sampling for the first 60 months of postnatal life [height-for-age Z score (HAZ), –0.92 ± 1.19 (mean ± SD); weight-for-height Z score (WHZ), –0.48 ± 1.33; n = 1961 fecal samples, 55 ± 4 samples collected per individual; table S1]. In Bangladesh, the median duration of breastfeeding is 4 months, whereas the weaning process is long, with a median of 25 months (14). Samples collected less frequently, or only after 36 months, from 19 other children from Mirpur were also included in our analysis (HAZ, –0.58 ± 1.12; WHZ, –0.25 ± 0.96; n = 25.7 ± 10.5 samples per child). Amplicons generated from variable region 4 (V4) of bacterial 16S rRNA genes present in these 2455 fecal samples were sequenced, and the resulting reads were assigned to operational taxonomic units with ≥97% nucleotide sequence identity (97%ID OTUs) (15, 16) (fig. S1). In total, 118 97%ID OTUs were represented at a relative fractional abundance of at least 0.001 (0.1%) in at least two of the samples collected over the 60-month period.

An initial broad description of microbiota development in this cohort was obtained by applying unweighted and weighted UniFrac to compute overall phylogenetic dissimilarity between gut communities from the 36 children sampled monthly from 1 to 60 months and 49 fecal samples collected in a previous study from 12 unrelated adults, aged 23 to 41 years, living in Mirpur (17). This metric indicated that the mean “infant/child-to-adult” distance decreases to “adult-to-adult” by 3 years of age (fig. S2, A and B). Alpha diversity also increased to adult-like levels during this time period (fig. S2, C and D). As an additional description of community development, we used the 16S rDNA dataset to construct a sparse Random Forests (RF)–derived model comprising age-discriminatory taxa (fig. S3, A to E). Microbiota “age” can be computed by noting the fractional abundances of these age-discriminatory taxa in a given sample obtained at a given time point (14). Applying the RF-generated model disclosed a high degree of correlation between microbiota age and chronologic age (R² = 0.8) (fig. S3C). Although these approaches provide measures of community development, they do not characterize interactions between community members during this process.

Principal components analysis (PCA) applied to taxa present in monthly fecal samples offers a way to mathematically characterize gut microbiota organization by defining principal components (eigenvectors). The result of PCA is a ranked list of principal components (principal component spectrum or “eigenspectrum”) where each principal component carries a percentage of data variance. Tracking the principal component spectrum through time offers a description of the evolving temporal organization of the gut microbiota. The approach we used, iterative PCA (iPCA), is described in fig. S4A. For each month, we created a matrix where the rows were fecal samples and the columns comprised the 118 taxa described above. In the example shown, time point 1 considers monthly fractional abundance data from month 1 and a reference time point. The dissimilarity between the two time points is reflected in the primary principal component (PC1). The system is considered to be “stable” at the time point where adding subsequent months’ data negligibly contributes to variance; mathematically, this is when the eigenvalue of PC1 reaches an asymptote.

We performed iPCA on sequentially joined monthly data with month 36 taken as a reference (fig. S4B). Month 36 was chosen on the basis of the results of phylogenetic dissimilarity and diversity measurements presented in fig. S2 [note that previous cross-sectional studies using these metrics had also indicated that an adult-like configuration was achieved by this time point; e.g., (18)]. iPCA revealed that month 20 and beyond signify a time period of minimal structural variation in the gut microbiota (fig. S4B). This conclusion was supported by using the very last time point in the 5-year longitudinal study as the reference (fig. S4C). Therefore, we were able to design a workflow to compute reproducible covariance (covariance conserved across time in a mature community assemblage, as opposed to transient covariance that may occur during community assembly) using months 20 to 60 without having to make any a priori assumptions about the importance of any taxon. For each month spanning postnatal months 20 to 60, we calculated the covariance between the 118 taxa over all individuals to generate monthly taxon-taxon covariance matrices (19) (see Fig. 1A, fig. S5, and table S2A). The matrices were averaged to a single taxon-taxon matrix (${\langle C_{\mathrm{bin}}^{i,j}\rangle }_t$) that represented a definition of consistent covariance where i and j are bacterial taxa and t designates the month (Fig. 1B and table S2B). PCA performed on this matrix revealed that PC1 encompassed 80% of the data variance (Fig. 1C; see supplementary text for a sensitivity analysis of the workflow). A group of 15 covarying taxa represented the top 20% of all taxa projections along PC1 (Fig. 1C; see table S3 for different threshold cutoffs). They include OTUs assigned to Bifidobacterium longum, another member of Bifidobacterium, Faecalibacterium prausnitzii, a member of Clostridiales, Prevotella copri, Streptococcus thermophilus, and Lactobacillus ruminis, all of which are age-discriminatory bacterial strains identified from RF-based analysis of bacterial 16S rDNA datasets generated from healthy members of this Bangladeshi cohort (fig. S3D).

The results of PCA performed on data generated from 478 samples collected from children sampled at postnatal months 50 to 60 provide an illustration of statistical covariation between these taxa: PC1 reveals that B. longum (OTU 559527) and L. ruminis (OTU 1107027) positively covary with one another across samples and negatively covary with two P. copri strains (OTUs 840914 and 588929); PC2 documents how two F. prausnitzii OTUs (514940 and 851865) positively covary with each other and negatively covary with S. gallolyticus (OTU 349024) and E. coli (OTU 1111294); PC3 discloses that the two P. copri OTUs negatively covary with the two F. prausnitzii OTUs (Fig. 1D).

Figure 1E provides a graphical depiction of this network of covarying taxa. Each green node represents one of the 15 OTUs that manifest a high degree of conserved covariance between months 20 and 60. Two nodes are connected by an edge if their temporally averaged covariance value (${\langle C_{\mathrm{bin}}^{i,j}\rangle }_t$ from Fig. 1B) is within the top 20% of all such values. Node size is proportional to the number of connections (edges) present. The green nodes collectively covary with one another. In contrast, gray nodes depict taxa that covary with green nodes but not with one another (Fig. 1E). The green nodes constitute an “insulated” ecostructure; its members exhibit significant intragroup covariation (fig. S6 and table S2C). We chose the term “ecogroup” to reflect the conserved collective statistical covariation of this sparse network of 15 organisms.

We asked whether components of the ecogroup provide a concise description of postnatal development of the microbiota in healthy members of the Bangladeshi cohort and, if so, whether changes in the representation of these taxa follow a pattern that is shared across other healthy birth cohorts representing distinct geographic locales and anthropologic features (20). Moreover, we postulated that if ecogroup taxa are informative biomarkers of normal community development, these taxa might be useful for characterizing impaired development and/or the extent to which community repair is achieved as a function of various therapeutic interventions (21).

Three different matrices were created where each row was a fecal sample collected from an individual at a particular month in the healthy Bangladeshi cohort and columns were either (i) all 118 taxa, (ii) the 15 ecogroup taxa, or (iii) the remaining 103 non-ecogroup taxa. PCA was performed on the rows of these matrices; fecal samples were plotted on the first three principal components. The left panel of Fig. 2A shows the results obtained when considering the fractional representation of all 118 taxa in fecal samples collected at postnatal months 4, 10, and 20. There is substantial interpersonal variation in gut community structure at postnatal month 1, as evidenced by the broad distribution along PC1, but this variation converges by month 4 (Fig. 2A, fig. S7, A and B, and table S2D). Thereafter, changes in the structure of the fecal microbiota are depicted by right-to-left movement along PC1, with minor variance observed along PC2 and PC3. Minimal movement along PC1 is observed after month 20 (fig. S7C), consistent with the results of iPCA in fig. S4, B and C. (Notably, children in this cohort had completed weaning by month 23; see fig. S8 for a description of the nature and timing of their dietary transitions.) Ecogroup taxa recapitulate the variance depicted by PC1, PC2, and PC3. Moreover, the ecogroup taxa capture (i) the significant interpersonal variation observed at postnatal month 1, (ii) the subsequent convergence to a B. longum–predominant microbiota at postnatal month 4, and (iii) temporal changes noted at postnatal months 10 and 20 (Fig. 2, A and B, fig. S7, A and B, middle panels, and table S2E). In contrast, the remaining 103 non-ecogroup taxa provide a less informative representation of developmental changes in the microbiota, as exemplified by the fact that PC1, PC2, and PC3 each capture ≤10% of the variance (Fig. 2A and fig. S7, A and B, right panels). The importance of taxa with low average fractional abundances and large standard deviations, such as P. copri (Fig. 2B, inset), is often overlooked when they are considered in isolation. However, analysis of taxon-taxon covariation can reveal relationships between member species, as illustrated by P. copri and B. longum (Fig. 1D, blue box).

To determine the extent to which the ecogroup is a generalizable descriptor of the microbiota in infants and children with healthy growth phenotypes, we turned to the MAL-ED network of study sites located in low- and middle-income countries (20, 21). Fecal samples had been collected monthly for the first two postnatal years, allowing sparse 30-taxon RF-generated models of normal community development to be generated from members of birth cohorts residing in Loreto, Peru (periurban area) and Vellore, India (urban area) (supplementary text, fig. S9, and table S4). Our ability to identify a network of covarying taxa in the Mirpur cohort depended on a high-resolution time-series study that extended well beyond the month at which the microbiota was determined to be “stable” (month 20). This duration of sampling did not occur at these other MAL-ED sites, obviating our ability to identify conserved covariance among taxa. However, to test how well the 15 ecogroup taxa identified in the Mirpur cohort could characterize the developing microbiota of children living in these countries, we created two matrices where each row was a fecal microbiota sample from the Indian or Peruvian cohorts and columns were either all taxa identified in the Peruvian and Indian samples or just the 15 ecogroup taxa identified from the Bangladeshi birth cohort. PCA was performed on the rows of these matrices and the same analysis performed as described for the healthy Bangladeshi birth cohort. The results show that the ecogroup taxa identified in members of the healthy Bangladeshi cohort also provide a concise description of community development in healthy members of these other two birth cohorts; that is, (i) they capture the variance depicted by PC1, PC2, and PC3 as compared to considering all taxa, and (ii) changes in their fractional abundances followed temporal patterns similar to those documented in the Bangladeshi cohort (fig. S10 and table S2F).

Bangladeshi children with acute malnutrition have perturbed microbiota development; their gut communities appear younger than those of chronologically age-matched individuals (14, 21). We examined whether ecogroup taxa provide a useful way to characterize the microbiota of children with moderate or severe acute malnutrition (MAM and SAM, respectively) prior to and after food-based therapeutic interventions. In the accompanying paper, Gehrig et al. describe 63 children from Mirpur diagnosed with MAM, aged 12 to 18 months, who were enrolled in a double-blind, randomized, controlled feeding trial of different microbiota-directed complementary foods (MDCFs) (21). Fecal samples were collected for 9 weeks at weekly intervals. The first 2 weeks comprised a pretreatment observation period. Over the next 4 weeks, children received either one of three MDCFs, or a ready-to-use supplementary food (RUSF) representing a form of conventional therapy that, unlike the MDCFs, was not designed to target specific members of the gut microbiota and repair community immaturity. The last 2 weeks represented the post-treatment observation period. In total, we identified 945 97%ID OTUs that had a fractional abundance of at least 0.001 (0.1%) in at least two fecal samples collected from one or more participants prior to, during, and after treatment (n = 531 samples). Gehrig et al. (21) also describe another trial involving 54 hospitalized Bangladeshi children with SAM, aged 6 to 36 months, where each participant was treated with one of three standard therapeutic foods and then followed over a 12-month period after discharge. In total, we identified 944 97%ID OTUs that had a fractional abundance of at least 0.001 in at least two fecal samples collected from one or more participants in this trial (n = 618 samples).

A matrix was created that included (i) all fecal samples from the SAM trial, (ii) pretreatment samples from children with MAM enrolled in all four arms of the MDCF trial, (iii) MAM samples obtained 2 weeks after treatment with one of the three MDCFs or the RUSF, and (iv) fecal samples from age-matched healthy Bangladeshi children (table S5). Each row of the matrix was a fecal sample, each column was an ecogroup taxon, and each element in the matrix was the fractional abundance of an ecogroup taxon within a particular fecal sample. PCA was performed on the rows of this matrix. Centroids for each cohort were computed and plotted on the PCA space (Fig. 3A). At the time of discharge, after receiving standard therapeutic foods, the microbiota of children with SAM remained in an incompletely repaired state: Although there was some improvement at 1 month after discharge, there was minimal additional improvement evident at 6 or 12 months, at which times their microbiota resembled that of untreated children with MAM (Fig. 3A). The microbiota of children with MAM that were treated with MDCF-1, MDCF-3, and RUSF clustered together, whereas the microbiota of those treated with MDCF-2 closely resembled that of healthy children. Notably, MDCF-2 was also distinct among the four treatment types in eliciting changes in the plasma proteome indicative of improved health status, including changes in biomarkers and mediators of metabolism, bone growth, central nervous system development, and immune function [see (21) for details].

PCA measures the effect of treatment on the gut microbiota by considering a constellation of changes in fractional abundance of ecogroup taxa, with the premise that the fractional abundances as well as the covariation of these taxa are important for characterizing community configuration. The left panel of Fig. 3B shows that the relationship between the fractional representations of B. longum (OTU 559527) and E. coli (OTU 1111294) determines microbiota position along PC1 in Fig. 3A. The center and right panels of Fig. 3B show that the relationship between the fractional representations of B. longum, E. coli, S. gallolyticus (OTU 349024), and P. copri (OTUs 588929 and 840914) determines position along PC2, whereas position along PC3 reflects the relationship between the abundances of S. gallolyticus and the two P. copri OTUs. B. longum, S. gallolyticus, and E. coli are the predominant ecogroup taxa represented in the microbiota of children with untreated SAM (Fig. 3C and table S2G). Treatment results in movement of their microbiota along PC1 and PC3 in Fig. 3A; this movement is associated with a decrease in B. longum, S. gallolyticus, and E. coli (Fig. 3C and table S2G). Differences between the microbiota of healthy children and those with SAM prior to and during the 12 months after treatment with standard therapeutic foods are manifest by differences in their respective positions along PC1 and PC3 (Fig. 3A). These differences signify incomplete repair to a “healthy” state and highlight the need to achieve further decreases in the fractional abundance of B. longum (associated with movement to the right of PC1) along with further decreases in the fractional abundance of S. gallolyticus and increases in P. copri (associated with positive movement along PC3). The representation of B. longum, P. copri, S. gallolyticus, and E. coli in the microbiota of 12- to 18-month-old children with untreated MAM accounts for their positive projection along PC1 and PC3 relative to the microbiota of children with untreated SAM (Fig. 3A). Among the tested therapeutic foods, MDCF-2 was uniquely associated with a positive movement along PC1 (Fig. 3A); this corresponds to decreased fractional abundance of B. longum (Fig. 3D and table S2G) and more complete community repair.

Two other methods, SparCC and SPIEC-EASI, have been used to describe microbiota organization (9, 10). As these methods were designed for cross-sectional studies, we adapted them (see supplementary text) so we could compare their ability to identify (i) temporally conserved aspects of community organization, and (ii) the degree to which SAM and MAM microbiota are repaired with different food-based interventions with the approach we had used to identify the ecogroup. SparCC identifies a subset of ecogroup taxa that describe healthy gut microbiota development in members of the 5-year healthy Bangladeshi cohort study (fig. S11, A and B). SparCC clearly separates the microbiota of children with untreated SAM from healthy controls and shows that treatment with standard therapeutic foods fails to repair their microbiota to a healthy state, or even to a state seen in children with untreated MAM. Compared to the approach described in Fig. 1A, SparCC does not as clearly separate MAM from healthy or (by extension) the differential effects of MDCF treatment, although it does place MDCF-2–treated microbiota closest to that of healthy children (fig. S11C). One explanation is that P. copri does not contribute as prominently to the collective group of correlated taxa identified by SparCC (fig. S11 and table S6, A and B). SPIEC-EASI identifies P. copri and other Prevotella OTUs as key microbes (fig. S12, A and B, and table S6, C to E). However, SPIEC-EASI does not provide as informative a description of the temporal pattern of healthy gut microbial development as does the ecogroup taxa [note the relative lack of movement over time of community configuration from right to left along PC1 in fig. S12C compared to Fig. 2A (ecogroup taxa) and fig. S11B (SparCC)]. The 15 interacting taxa identified by SPIEC-EASI separate untreated and treated SAM and MAM microbiota from one another and from healthy (fig. S12D). As with the two other approaches, although less clearly than with the ecogroup taxa, SPIEC-EASI shows that MDCF-2 is most effective in changing the configuration of the MAM-associated microbiota toward a healthy state relative to MDCF-1, MDCF-3, and RUSF. Together, these findings provide support for considering temporally conserved taxon-taxon covariance when characterizing the microbiota of children with undernutrition prior to and after various therapeutic interventions.

Our observations raise questions about the nature of the interactions among B. longum, P. copri, and other ecogroup taxa during postnatal development as a function of the dietary transitions that occur when children progress from exclusive milk feeding to complementary feeding to a fully weaned state. To address this issue, we colonized germ-free piglets with ecogroup taxa and tracked the dynamics of consortium members over time. We turned to gnotobiotic piglets rather than mice because the former have physiologic and metabolic qualities more similar to that of humans (22). Piglets were derived as germ-free at birth and were fed an irradiated sow’s-milk replacement (Soweena) for the first four postnatal days (fig. S13A). Piglets (n = 5) were then colonized, by oral gavage, with a consortium of seven cultured, sequenced B. longum strains recovered from the fecal microbiota of children living in Mirpur, Bangladesh as well as three other countries (Peru, Malawi, and the United States) (fig. S13A). On the basis of their genome sequences (table S7), six strains were classified as B. longum subspecies infantis and one as B. longum subspecies longum. The gavage mixture also contained two Bifidobacterium breve strains, which we used as comparators to delineate factors that contribute to the fitness of the B. longum strains, given the phylogenetic similarity of their genomes. Beginning on postnatal day 4, a diet representative of that consumed by 18-month-old children living in Mirpur [Mirpur-18 (21)] was added to food bowls containing Soweena. On postnatal day 7, piglets were gavaged with a second consortium consisting of 16 additional cultured sequenced ecogroup taxa (fig. S13A) representing 13 of the 15 species shown in Fig. 1C. During postnatal days 5 to 22, the amount of Mirpur-18 added to food bowls was progressively increased while the amount of Soweena was decreased; once a fully weaned state was achieved on day 22, animals were monotonously fed the Mirpur-18 diet until they were euthanized on postnatal day 29. Piglets increased their weight by 185 ± 31% (mean ± SD) between postnatal days 7 and 29.

To define features in ecogroup strains that relate to their fitness during the series of dietary transitions that mimic those experienced by children living in Mirpur, we performed short-read shotgun sequencing of community DNA prepared from rectal swabs obtained at 11 time points spanning experimental days 5 to 29 (fig. S13A) and along the length of the gut at the time of euthanasia. The results are presented in Fig. 4A and table S2H. After gavage of remaining ecogroup members, the representation of all B. longum strains diminished rapidly. From postnatal day 8 to day 22, as the animals were being weaned, S. gallolyticus, E. coli, E. avium, L. salivarius, and P. copri exhibited distinct patterns of temporal change in their representation. After the animals were fully weaned, there was a pronounced increase in P. copri, which became the dominant member of the cecal, colonic, and fecal microbiota (Fig. 4A and fig. S13B). The relationship between the abundances of P. copri and B. longum is comparable in these piglets to that observed in the healthy Bangladeshi children who were used to evaluate the microbiota configurations of untreated and treated children with MAM and SAM (Fig. 3, C and D).

The representations of 81 mcSEED metabolic modules (see methods) in strain genomes were used to make in silico predictions about their capacity to synthesize amino acids and B vitamins, utilize a variety of carbohydrates, and generate short-chain fatty acids. Predicted phenotypes were scored as either a “1” or a “0” signifying auxotrophy or prototrophy in the case of amino acid and B-vitamin biosynthesis, or the ability or inability to utilize various carbohydrates (table S8). PCA of a “binary phenotype matrix” of all strains present at a fractional representation of ≥0.001 in fecal samples collected from postnatal day 8 to day 18 identified 14 carbohydrate utilization pathways, plus the capacity to synthesize cysteine, folate, and pantothenate as genomic features that distinguish these strains from each other (table S9). Hierarchical clustering by these predicted metabolic phenotypes also grouped these strains by their fitness (Fig. 4, B and C).

We performed microbial RNA-seq using cecal contents to characterize the expression of genes encoding components of mcSEED metabolic modules present within the ecogroup strains. [The fractional representations of these strains in the cecum and feces at the time of euthanasia were highly correlated (r² = 0.98; table S10).] Figure S14A illustrates the workflow used to generate a mcSEED “enrichment matrix” (M_E) that signifies the extent to which the aggregate transcript levels of components of a given mcSEED metabolic module in a given bacterial strain quantitatively differ from that of a reference strain. Because P. copri had the highest fractional representation on postnatal day 29, it was used as the reference (fig. S14B and table S2I). PCA was performed on the mcSEED enrichment matrix (Fig. 5A and table S11A). The results revealed that the transcriptomes of Bifidobacterium strains cluster together and are distinct from those of P. copri, E. coli, B. luti, and E. avium. Moreover, the distribution of strains along PC1 based on their mcSEED enrichment profiles correlated with their fractional representation (fitness) in the cecal and fecal microbiota (Fig. 5A, inset).

To identify which expressed components of mcSEED metabolic modules contribute to the differences in the fractional representation, we required a way to relate the principal components of the rows (metabolic modules) and columns (strains) of the mcSEED enrichment matrix. To do so, we used singular value decomposition (SVD; fig. S14, C and D). Relative to P. copri, the most distinguishing features of the Bifidobacterium transcriptomes were markedly reduced or absent expression of pathways involved in (i) biosynthesis of cysteine, tyrosine, tryptophan, and asparagine; (ii) utilization of several carbohydrates (xylose and β-xylosides plus galacturonate/glucuronate/glucuronide); (iii) biosynthesis of queuosine; and (iv) uptake of cobalt related to cobalamin biosynthesis (Fig. 5B and tables S2J and S11B). Moreover, expression of four of these pathways (cysteine and asparagine biosynthesis; xylose/β-xyloside and galacturonate/glucuronate/glucuronide utilization) exclusively differentiate P. copri, B. luti, E. coli, and E. avium from all nine Bifidobacterium species and the other five strains whose transcripts were represented in the community metatranscriptome (Fig. 5B).

The biological significance of expression of these distinguishing mcSEED metabolic modules demanded a further contextualization based on whether these systems were complete or incompletely represented in the strain genomes. Figure 5C shows that all of the Bifidobacterium strains contain complete metabolic pathways for tyrosine, asparagine, and tryptophan biosynthesis but do not contain complete metabolic pathways for cysteine biosynthesis; utilization pathways for galactose, xylose, and glucuronides; and B-vitamin synthetic pathways for queuosine and cobalamin. In contrast, E. coli and B. luti have mcSEED binary phenotype profiles similar to that of P. copri and contain complete metabolic pathways for cysteine biosynthesis and xylose utilization (table S2J). These results indicate that genomic features of the Bifidobacterium strains examined limit their ability to thrive in the context of the Mirpur-18 diet and a community that contains the other ecogroup strains. In contrast, the fact that P. copri and other ecogroup strains contain and express these metabolic pathways provides support for their importance in maintaining their fitness under these conditions. As such, the feature-reduction approach used here provides a rationale for testing nutritional interventions that target these pathways in ecogroup members in children at risk for, or who already have, perturbed microbiota development.

We have developed a statistical approach to identify a group of 15 covarying bacterial taxa that we term an ecogroup. We found that the ecogroup is a conserved structural feature of the developing gut microbiota of healthy members of several birth cohorts residing in different countries. Moreover, the ecogroup can be used to distinguish the microbiota of children with different degrees of undernutrition (SAM, MAM) and to quantify the ability of their gut communities to be reconfigured toward a healthy state with a MDCF. Studies of gnotobiotic piglets subjected to a set of dietary transitions designed to model those experienced by members of the Bangladeshi healthy birth cohort demonstrate that temporal changes in the fitness of ecogroup taxa can occur in the absence of other gut community members. These observations suggest that the approach used to identify the ecogroup may be useful in characterizing microbial community organization in members of other longitudinally sampled (human) cohorts.

A critical feature of biological systems is that they function reliably, yet adapt when faced with environmental fluctuations (23, 24). An architecture of sparse but tight coupling enables rapid evolution to new functions in proteins (25, 26). Studies of macro-ecosystems such as ant colonies have argued that adaptive behaviors are dependent on proper network organization (27). The gut microbiota must satisfy the constraints of survival: namely, withstanding insult and maintaining functionality (robustness) while still having the capacity for plasticity. “Embedding” a sparse network of covarying taxa in a larger framework of independently varying organisms could represent an elegant architectural solution developed by nature to maintain robustness while enabling adaptation.

We are indebted to the families of study subjects for their active participation and assistance. We thank the staff and investigators at icddr,b for their contributions to the recruitment and enrollment of participants in the 5-year Bangladeshi birth cohort study plus the interventional studies of children with SAM and MAM, as well as the collection of biospecimens and data. We also thank the study team members and health care workers involved in the MAL-ED birth cohort studies; M. Gottlieb, D. Lang, K. Tountas, and M. McGrath, who provided invaluable assistance in coordinating the MAL-ED collaboration and providing access to key clinical datasets; M. Meier, S. Deng, and J. Hoisington-López for superb technical assistance; D. O’Donnell, J. Serugo, and M. Talcott for their indispensable help with gnotobiotic piglet husbandry; and R. Olson for technical support with the mcSEED-based genome analysis and subsystem curation. Funding: Supported by the Bill & Melinda Gates Foundation as part of the Breast Milk, Gut Microbiome, and Immunity (BMMI) Project. The 5-year birth cohort study of Bangladeshi children was funded by NIH grant AI043596 (W.A.P.). A.S.R. is a postdoctoral fellow supported by Washington University School of Medicine Physician Scientist Training Program and in part by NIH grant DK30292. D.A.R., A.A.A., and S.A.L. were supported by Russian Science Foundation grant 19-14-00305. J.I.G. is the recipient of a Thought Leader award from Agilent Technologies. Author contributions: R.H. and W.A.P. designed and oversaw the 5-year birth cohort study; they, together with T.A., were responsible for coordinating various aspects of biospecimen and metadata collection. S.H., M.M., R.H., W.A.P., and T.A. (Bangladesh), M.N.K. (Peru), G.K. (India), P.O.B. (South Africa), and A.A.M.L. (Brazil) oversaw the MAL-ED studies. S.H., I.M., M.I., M.M., and T.A. were responsible for studies involving the SAM and MAM cohorts. J.L.G. and S.S. generated 16S rDNA datasets from human fecal samples. M.J.B. managed the repository of biospecimens and associated clinical metadata used for the studies described above. H.-W.C. performed the experiments with gnotobiotic piglets with the assistance of A.S.R. S.V., and M.C.H. D.A.R., A.A.A., S.A.L., and A.L.O. performed in silico metabolic reconstructions based on the genome sequences of bacterial strains introduced into gnotobiotic piglets. A.S.R. conceived the mathematical approach and wrote all of the computational workflow for identifying ecogroup taxa, performed the sensitivity analysis of the workflow, compared the SparCC and SPIEC-EASI algorithms with the workflow, and undertook the analyses of gut microbial communities from subjects enrolled in the SAM, MDCF, Peruvian, and Indian cohort studies as well as the gnotobiotic piglet experiment, with J.L.G., S.V., M.J.B., and J.I.G. contributing in various supportive ways. A.S.R. and J.I.G. wrote the paper. Competing interests: J.I.G. is a co-founder of Matatu Inc., a company characterizing the role of diet-by-microbiota interactions in animal health. W.A.P. serves as a consultant to TechLab Inc., a company that makes diagnostic tests for enteric infections and has served as a consultant for Perrigo Nutritionals LLC, which produces infant formula. Data and materials availability: Bacterial V4-16S rDNA sequences in raw format (prior to postprocessing and data analysis), shotgun datasets generated from cultured bacterial strains, and COPRO-seq and microbial RNA-seq datasets obtained from gnotobiotic piglets have been deposited at the European Nucleotide Archive under study accession number PRJEB27068. Code has been archived at Zenodo (35). Fecal specimens from the MAL-ED birth cohorts in Bangladesh (icddr,b, Dhaka), Brazil (Federal University of Ceará, Fortaleza), India (Christian Medical College, Vellore), Peru (JHSPH/AB PRISMA), South Africa (University of Venda), and from the NIH birth cohort and SAM/MDCF studies at icddr,b, were provided to Washington University under material transfer agreements. This work is licensed under a Creative Commons Attribution 4.0 International (CC BY 4.0) license, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited. To view a copy of this license, visit http://creativecommons.org/licenses/by/4.0/. This license does not apply to figures/photos/artwork or other content included in the article that is credited to a third party; obtain authorization from the rights holder before using such material.