The gut microbiome is suggested to play a role in the pathogenesis of autoimmune disorders such as type 1 diabetes. Evidence of anti-islet cell autoimmunity in type 1 diabetes appears in the first years of life; however, little is known regarding the establishment of the gut microbiome in early infancy. Here, we sought to determine whether differences were present in early composition of the gut microbiome in children in whom anti-islet cell autoimmunity developed. We investigated the microbiome of 298 stool samples prospectively taken up to age 3 years from 22 case children in whom anti-islet cell autoantibodies developed, and 22 matched control children who remained islet cell autoantibody–negative in follow-up. The microbiome changed markedly during the first year of life, and was further affected by breast-feeding, food introduction, and birth delivery mode. No differences between anti-islet cell autoantibody–positive and –negative children were found in bacterial diversity, microbial composition, or single-genus abundances. However, substantial alterations in microbial interaction networks were observed at age 0.5 and 2 years in the children in whom anti-islet cell autoantibodies developed. The findings underscore a role of the microbiome in the pathogenesis of anti-islet cell autoimmunity and type 1 diabetes.
Type 1 diabetes is the result of a complex interplay of genetic susceptibility (1) and environmental determinants leading to anti-islet cell autoimmunity against pancreatic islet β-cells and autoimmune β-cell destruction (2). Anti-islet cell autoimmunity precedes the clinical onset of type 1 diabetes and often develops within the first years of life (3). This suggests that early shaping of the immune system in children is critical for the initiation of autoimmunity (3). There is increasing evidence that the immune response is shaped by factors that include how the host establishes a stable ecosystem with a large cohort of accompanying bacteria (4–7). With this, the role of microbiota in type 1 diabetes pathogenesis has become an important subject of investigation (8–12). The largest community of bacteria is established in the gastrointestinal tract (13,14) where beneficial host-bacteria interactions have been demonstrated for food degradation or pathogen defense (14–16).
Relatively few studies of the human gut microbiome have been performed in children <5 years old. These studies suggest that the phylogenetic composition of the bacterial communities evolves toward an adult-like configuration within the 3-year period after birth (14,17–19). Hence, it is conceivable that the evolution of the microbiome in infancy could influence the risk of anti-islet cell autoimmunity in susceptible children. Indeed, studies from Finland have provided evidence for this hypothesis (10,20). The aim of our study was to investigate gut bacterial community structures during the early period from birth to the age of 3 years from the perspective of complex interaction networks. We estimated interaction on the basis of covariation of bacterial abundances to compare children in whom anti-islet cell autoantibodies developed with children in whom such autoantibodies did not develop. We took advantage of the prospective BABYDIET study (21), where infants at increased risk of type 1 diabetes were monitored from birth for the development of anti-islet cell autoantibodies and type 1 diabetes. The gut microbiome composition was estimated based on measurements of 16S rRNA gene sequences from fecal samples that were obtained at 3-month intervals up to the age of 3 years. Analyses were focused on bacterial diversity, community composition, individual bacterial species, and microbial interaction networks. Results show that complex bacterial interaction networks, rather than single genera, appear to be relevant to early preclinical type 1 diabetes.
Research Design and Methods
BABYDIET Study Material
Analysis of microbiota was performed on 298 stool samples from 44 children participating in the BABYDIET study (21). The BABYDIET study randomized 150 infants with a first-degree relative with type 1 diabetes and with the type 1 diabetes risk HLA genotypes DR3/4-DQ8 or DR4/4-DQ8, or DR3/3 to gluten exposure at 6 or 12 months of age. The intervention had no effect on anti-islet cell autoimmunity outcome. Blood and stool samples were collected at 3-month intervals from age 3 to 36 months and subsequently at 6-month intervals. Anti-islet cell autoantibodies (i.e., autoantibodies to insulin, glutamic acid decarboxylase, insulinoma-associated antigen-2, and zinc transporter 8) were measured at each study visit. Written informed consent was obtained from the parents. The study was approved by the ethics committee of the Ludwig-Maximilians-University, Munich, Germany (Ethikkommission der Medizinischen Fakultät der Ludwig-Maximilians-Universität no. 329/00).
Stool samples chosen for the study included 147 samples from the 22 BABYDIET cohort children in whom persistent anti-islet cell autoantibodies developed at a median age of 1.54 years (interquartile range 0.90 years; maximum age 2.45 years) and 151 samples from 22 children who remained anti-islet cell autoantibody–negative, and were matched for date of birth. Of the 22 children with persistent islet cell autoantibodies, persistent multiple islet autoantibodies had developed in 15 children, and diabetes had developed in 10 children after a median follow-up time of 5.3 years. For the 44 children, stool samples were taken from age 0.24 to 3.2 years with an average of 6.8 probes per child (Supplementary Table 1).
Sample Processing and Deep Sequencing
Stool samples were collected at home and shipped by express courier overnight to the clinical study center, where they were processed and immediately frozen at −80°C. DNA was extracted from the stool samples as described previously (10). Bacterial 16S rRNA genes present within fecal samples were amplified using the primers 515F and 806R (22) modified with a sample-specific barcode sequence and Illumina adapter sequences.
PCR was performed at an initial denaturation temperature of 94°C for 3 min, followed by 20 cycles of 94°C for 45 s, 50°C for 30 s, and 65°C for 90 s. A final elongation step at 65°C was run for 10 min. PCR products were purified using the Qiagen PCR purification kit following the protocol of the manufacturer. Illumina high-throughput sequencing of 16S rRNA genes was conducted as described (23). Illumina sequencing was performed with 101 cycles each. Sequences were trimmed based on quality scores using a modified version of Trim2 (24), and the first 11 bases of each paired read were removed to eliminate degenerate bases derived from primer sequences. The prokaryotic database (25) used for 16S rRNA gene analysis was formatted using TaxCollector (26). Sequences were compared with the TaxCollector-modified RDP database using CLC Assembly Cell version 3.11 using the paired reads and global alignment parameters. Two parameters were used in this step: a 98% length fraction and similarity values dependent on the desired taxonomic level (i.e., 80% at domain/phylum, 90% to class/order/family, 95% to genus levels) (27). Pairs that matched different references at the species level were classified at the lowest common taxonomic level. Unresolved pairs were discarded. Henceforth, successfully paired reads are referred to as reads.
Data on breast-feeding (yes, no), the duration of breastfeeding (weeks), and the introduction of solid food (gluten-free and gluten-containing cereals, vegetables, fruits) were analyzed from daily food records as previously described (21). Data on caesarean section (yes, no) were obtained from obstetric records.
Shannon evenness and Chao richness indices were estimated at genus level, as previously described (28,29). To correct bacterial diversity for the influence of confounding factors, stepwise multiple regression was performed with diversity as a dependent variable. Age, breast-feeding at sampling time, introduction of solid food, first gluten exposure, and delivery by caesarean section were used as confounding variables. Akaike information criterion (30) was used in the stepwise regression procedure to select confounding factors associated with diversity. To avoid bias due to violations of normality, rank regression (31) was used to estimate P values of the regression coefficients corresponding to confounding factors associated with diversity. The R Package fields (32) was used for cubic spline regression of age versus Shannon evenness. Chao richness was corrected for the influence of caesarean section by using the residuals of a regression model with richness as a dependent variable and caesarean section as an independent variable. Diversity analyses were performed on the entire age range, and after grouping reads into three age classes of 0.5 ± 0.25, 1.0 ± 0.25, and 2.0 ± 0.5 years. At most, one single probe closest to 0.5, 1, and 2 years was used, respectively, for each child.
For further analyses, phyla and genera with <0.01% abundance within the total number of reads were neglected. This reduced the number of genera from 452 to 75, and the number of phyla from 21 to 8. For the analysis of bacterial community compositions, Bray-Curtis distances (33) were estimated on Hellinger transformed data (34). Differences in community compositions were tested with the nonparametric MANOVA (npMANOVA) (35) method. To visualize the results, Principal Coordinate Analysis (PCoA) was performed. Relative abundances of individual phyla and genera were compared by Wilcoxon-Mann-Whitney tests. To account for heterogeneity in variances, Brunner-Munzel tests (36) were used for bacteria where a Bartlett test (37) showed evidence for unequal variances (P < 0.05). Second, we adjusted for confounding factors by using the residuals of stepwise Akaike information criterion models with bacterial abundance as the dependent variable and the confounding variables as independent variables. For all independent variables with a P value <0.1, the model was again estimated and the resulting residuals were used as adjusted abundance values. P values were corrected for multiple testing with the Benjamini and Hochberg method (38).
Correlation-based networks reflecting covariations of bacterial abundances were used as a surrogate for bacterial interaction networks. To construct interaction networks, we computed the Spearman rank correlation (ρ) for all possible pairs of genera. We used 1,000 random permutations and set an edge, if P < 0.05 and ρ > 0.3, considering positively correlated genera. Networks were plotted with the Fruchterman and Reingold (39) algorithm. Eigenvector centrality (EC) was estimated as described (40), and Kolmogorov-Smirnov tests were used to test for differences in distributions (41,42). Differences in the number of isolated nodes were analyzed by comparing the number of nodes of degree ≤1. To test whether the observed differences were due to lower sample size in children who became anti-islet cell autoantibody–positive, networks were estimated with all combinations of N autoantibody-negative samples, where N denotes the number of autoantibody-positive samples. All statistical analyses were performed with the statistical software R, version 2.15.2.
Anti-Islet Cell Autoantibodies and Diversity of the Gut Microbiome
Bacterial diversity was analyzed for 452 bacterial genera. Diversity can be described by its evenness and richness. Evenness measures the similarity of proportions of taxa within a community, while richness represents the number of taxa in the community. We first analyzed covariates for Shannon evenness and Chao richness via stepwise regression models for all 298 stool probes. We observed an association of evenness with age (P = 0.025), caesarean section (P = 0.0026), gluten exposure (P = 0.0095), and an association of richness with caesarean section (P = 0.0025). Cubic spline regression of evenness with age showed that evenness increased until 2 years of age and saturated for older children (Fig. 1A). In contrast, richness remained constant over time (Fig. 1B). We found no association of anti-islet cell autoantibody–positive/negative status with evenness (P = 0.27; Fig. 1C) and richness (P = 0.56; Fig. 1D). Richness and evenness were also not different after adjustment for associated covariates (Shannon evenness, P = 0.62; Chao richness, P = 0.40).
Supporting the association of bacterial evenness with increasing age, the analysis of the distribution of all 21 phyla revealed a considerable shift between ages 6 months and 1 year (Fig. 1E). This shift was primarily due to an increase in the relative abundance of Firmicutes in both the autoantibody-positive (P = 8.7 × 10−6) and autoantibody-negative (P = 0.016) groups. The distribution of phyla remained relatively constant between ages 1 and 2 years.
To account for the effect of age- and sample-dependent collinearities, the data were grouped into the following three age intervals: 0.5 ± 0.25 years (anti-islet cell autoantibody–positive N = 19; anti-islet cell autoantibody–negative N = 21); 1 ± 0.25 years (anti-islet cell autoantibody–positive N = 16; anti-islet cell autoantibody–negative N = 19), and 2 ± 0.5 years (anti-islet cell autoantibody–positive N = 18; anti-islet cell autoantibody–negative N = 20). At each age interval, no differences between anti-islet cell autoantibody–positive and anti-islet cell autoantibody–negative children were observed for evenness (P0.5 = 0.22, P1 = 0.83, P2 = 0.29; Supplementary Fig. 1A–C) and richness (P0.5 = 0.12, P1 = 0.1, P2 = 0.63; Supplementary Fig. 1D–F).
Anti-Islet Cell Autoantibodies and Bacterial Community Composition
Differences in bacterial community composition were tested by comparing the intragroup distances of bacterial abundances between case and control children based on Bray-Curtis distances (33) estimated on the 75 genera that remained after filtering bacteria with low abundances (<0.01%). Single-variable and multivariable npMANOVA (35) models, including the covariates age at sampling time, caesarean section, duration since solid food introduction, duration since introduction of gluten, and breast-feeding at sampling time were applied for each of the three age intervals.
The covariates with the strongest effects on bacterial community composition at age 0.5 years were breast-feeding (P = 0.002, Supplementary Table 2) and the duration since first solid food introduction (P = 0.001, Supplementary Table 2). At age 1 year, only duration since first solid food introduction (P = 0.049, Supplementary Table 2) was associated with bacterial community composition and at age 2 years the effects of nutrition vanished. Children who became anti-islet cell autoantibody–positive did not show significant differences in community composition in univariable (P0.5 = 0.52, P1 = 0.36, and P2 = 0.35) and multivariable npMANOVA (P0.5 = 0.20, P1 = 0.38, and P2 = 0.18; Supplementary Table 2) analyses for all of the three analyzed age intervals. PCoA plots did not show a noticeable clustering of anti-islet cell autoantibody–positive and –negative children (Fig. 2).
Anti-Islet Cell Autoantibodies and Bacterial Abundances
Abundances at the phylum and genus levels were assessed at ages 0.5, 1, and 2 years. None of the eight analyzed phyla differed in bacterial abundances between anti-islet cell autoantibody–positive and –negative children. Of the 75 analyzed genera, Dorea and Barnesiella abundances at age 0.5 years (P = 0.003 and P = 0.035), Candidatus Nardonella abundances at age 1 year (P = 0.031), and Erwinia and Enterobacter abundances at age 2 years (P = 0.024 and P = 0.045) differed between anti-islet cell autoantibody–positive and –negative children (Supplementary Tables 3–5). None of these were significant after correction for multiple testing. None of the genera showed a persistent difference between anti-islet cell autoantibody–positive and –negative children across all three age groups.
Since nutrition affected bacterial composition in our cohort, we also compared bacterial abundances between children in whom anti-islet cell autoantibodies developed and children who remained autoantibody-negative after adjustment for confounding factors. None of the phyla abundances were significantly different after the adjustment. Some additional genera showed differences after adjusting for confounding factors (Supplementary Tables 3–5). These include Veillonella abundances, which were lower in children in whom anti-islet cell autoantibodies developed (average 3%) than in children who remained autoantibody-negative (average 10%; P = 0.0098) at age 0.5 years, and Enterococcus abundances, which were higher in children in whom anti-islet cell autoantibodies developed (average 3%) than in children who remained autoantibody-negative (average 0.8%, P = 0.00011) at age 0.5 years. Again, these differences were not significant after correction for multiple comparisons.
Anti-Islet Cell Autoantibodies and Bacterial Interaction Networks
Since the gut microbiome constitutes an ecosystem where bacteria depend on each other, and in particular compete for nutrition, we hypothesized that a functional interplay of bacteria is crucial for the development of the gut microbiome and that differences in the interaction between bacteria might be associated with the development of anti-islet cell autoimmunity. We therefore used microbial correlation networks at the genus level (N = 75) as an approximation for bacterial interactions using two different scores: EC and the number of isolated nodes. Correlation-based bacterial interaction networks were estimated at ages 0.5, 1, and 2 years for the anti-islet cell autoantibody–positive and –negative groups.
EC is a measure of the relative importance and connectivity of each node in a network (40). EC of a node accounts for the centrality of its neighbors, assuming that a node is more central if the surrounding neighbors also have high centrality (43). Differences in EC indicate that the information flow varies throughout the network. The networks of children who became anti-islet cell autoantibody–positive showed significantly different centrality distributions at age 0.5 years (P = 0.0024; Fig. 3A, D, and G) and again at age 2 years (P = 0.013; Fig. 3C, F, and I). Most of the genera that had high centrality at age 0.5 years had also high centrality at age 2 years for children who became autoantibody-positive (88%) and children who remained autoantibody-negative (77%). No differences were observed between the two groups at age 1 year (Fig. 3B, E, and H).
In both groups, a cluster of nodes had high ECs (>0.5). In contrast to the autoantibody-negative group, more nodes with intermediate levels of EC (>0.05 and <0.5) were found in the anti-islet cell autoantibody–positive group (Fig. 3). At age 0.5 years, the bacterial genera Enterococcus, Sarcina, Prevotella, and Corynebacterium showed high EC in networks of children who became anti-islet cell autoantibody–positive and low EC (<0.05) in networks of children who remained autoantibody-negative (Supplementary Fig. 2A). A detailed overview of EC values at age 1 year can be found in Supplementary Fig. 2B. At age 2 years, Barnesiella and Candidatus Nardonella showed high EC in the autoantibody-positive and low EC in the autoantibody-negative group (Supplementary Fig. 2C). In contrast, Staphylococcus and Nocardioides had high EC in the autoantibody-negative group and low EC in the autoantibody-positive group (Supplementary Fig. 2C).
While EC measures the capacity of overall information flow, node degree measures the connectivity in a topological sense. In the following, we refer to genera with node degree ≤1 as isolated nodes. More isolated bacterial genera (Fig. 3A–F) were found in children who developed anti-islet cell autoantibodies at ages 0.5 years (P = 0.00012) and 2 years (P = 0.0044; Fig. 4). A detailed overview of node degrees of all genera in the bacterial networks for anti-islet cell autoantibody–positive and –negative children is shown in Supplementary Fig. 3. Sample sizes differed slightly between anti-islet cell autoantibody–positive children (N0.5 = 19, N1 = 16, N2 = 18) and anti-islet cell autoantibody–negative children (N0.5 = 21, N1 = 19, N2 = 20) in the three age groups. We therefore performed interaction network estimates for all possible equal number subsets of children who remained autoantibody-negative. At 0.5 and 2 years of age, there was no single combination of children who remained autoantibody-negative that resulted in a similar high number of isolated bacterial genera as that observed for children who became anti-islet cell autoantibody-positive. No differences in the number of isolated bacterial genera between the two groups were observed at 1 year of age.
In murine models, associations between gut microbiome composition and type 1 diabetes or anti-islet cell autoimmunity have been found (11,44–46). Little is known regarding the early establishment of the gut microbiome in children in whom anti-islet cell autoimmunity develops. In this study, stool samples from children participating in the prospective BABYDIET cohort were analyzed within the first 3 years of life to compare bacterial diversity, composition, individual phyla, and genera abundances, and interaction networks for children who became anti-islet cell autoantibody–positive to those of children who remained anti-islet cell autoantibody–negative. No differences in bacterial diversity or community composition were observed between autoantibody-positive and -negative children. After correction for multiple testing, there were no individual bacterial genera that showed significantly different abundances between anti-islet cell autoantibody–positive and anti-islet cell autoantibody–negative children. However, children who became anti-islet cell autoantibody–positive showed significantly different distributions of EC in correlation-based interaction networks of bacterial communities, and their networks consisted of more isolated nodes than those of children who remained autoantibody-negative.
The strengths of our study lie in the relatively large number of stool samples analyzed, the early and sequential sample collection, and the homogenous cohort of first-degree relatives of individuals with type 1 diabetes with similar type 1 diabetes high-risk HLA DR-DQ genotypes. Furthermore, case and control children were matched by date of birth so that stool samples were collected within the same year and season, and under similar study conditions between groups. Although the control islet autoantibody–negative children in our study are well-matched to the case children, they are enriched for type 1 diabetes susceptibility genes and may therefore not be the most suitable control subjects. We have not examined the microbiome in children without an enriched genetic susceptibility, and it is possible that our findings may not represent the microbiome status of children from the general population. Development of multiple islet autoantibodies, a status that confers extreme risk for diabetes (47), occurred in the majority of the islet cell autoantibody–positive children in our study. However, findings may differ if only case children in whom diabetes subsequently developed were analyzed. Further limitations of the study include the collection of samples at home with overnight transport at room temperature, and that the data on other potential confounding factors, such as infections or antibiotic therapy, were not available for this analysis.
Two other studies (10,20) from Finland have investigated the role of the gut microbiome in children in whom anti-islet cell autoantibodies developed. One study (10) investigating four children who seroconverted and an equivalent number who remained anti-islet cell autoantibody–negative found a lower bacterial diversity and differences in the abundances of the phyla Bacteroidetes and Firmicutes in children with anti-islet cell autoantibodies before, at, and after islet autoantibody seroconversion. At the genus level, the same study found differences in Eubacterium and Faecalibacterium abundances, and reported that community composition was more similar in children who remained anti-islet cell autoantibody–negative compared with children who became anti-islet cell autoantibody–positive (10). A second study (20) from Finland compared the gut microbiome of 18 children in whom anti-islet cell autoimmunity developed with 18 children in whom it did not develop. In contrast to our study, the anti-islet cell autoantibody–positive children in the study by de Goffau et al. (20) were already autoantibody-positive at the time of sampling, and the probands were older. The authors reported significant differences in the abundance of the phylum Bacteroidetes, the genus Bacteroides, and several bacteria on the species level (20). In addition, a trend of increased bacterial diversity in anti-islet cell autoantibody–negative children was reported (20). We did not find these reported differences between anti-islet cell antibody–positive and anti-islet cell antibody–negative children in our cohort. There was no single phylum that showed differences between anti-islet cell antibody–positive and anti-islet cell antibody–negative children. The reported differences in the abundances of the genera Faecalibacterium, Eubacterium (10), and Bacteroides (20) were also not observed in the data presented here. The deviating results may be explained by differences in sample sizes (10), in the different study design used by de Goffau et al. (20), and/or geographical differences between the German and Finnish children. Finally, a recent study found differences in the abundances of several bacteria, including Prevotella, Clostridium, Veillonella, Bifidobacterium, Lactobacillus, and Bacteroides in the microbiome of children with established diabetes compared with healthy control children (48).
Nutrition has an important effect on the early microbial community. For example, exclusively breast-fed infants have different distributions of bacteria than formula-fed infants (49). In line with these data, our analysis of microbial community composition revealed that early microbiome composition is associated with breast-feeding duration and the age at which solid food was introduced. Caesarean section was also found to be associated with bacterial abundances at the genus level. These and other confounders should be considered when analyzing bacterial abundances in young children. We analyzed abundances with and without adjustment for such confounders. Although some differences in the genera were observed between children in whom anti-islet cell autoantibodies did and did not develop, most of the significant genera had low abundances, and none of the genera was significant after multiple testing corrections.
We hypothesized that, instead of individual bacterial abundances, the interplay between bacteria might be compromised in children who became anti-islet cell autoantibody–positive, and that differences in the interaction between bacteria might be associated with the development of anti-islet cell autoimmunity. The estimation of microbial co-occurrence networks was recently successfully applied to a large microbiome data set from different body sites (50). This encouraged us to use covariation of microbial abundances as a surrogate for bacterial interaction and to compare the networks of children who became anti-islet cell autoantibody–positive with the networks of children who did not. An increased number of isolated nodes can cause a reduced number of possible communication paths, and therefore impair the flexibility of the network and the adaptability of the bacterial community. Interestingly, we found significantly higher numbers of isolated nodes in children who became anti-islet cell autoantibody–positive. In addition, we observed significant differences in the distributions of EC, suggesting differences in the information flow between bacteria. Differences were observed at ages 6 months and 2 years but not at age 1 year. Since many of the children changed from breast-feeding to solid food and we found the most pronounced shift of large-scale bacterial distributions between 6 months and 1 year, we suspect that strong nutritional effects may mask anti-islet cell autoimmunity associations with bacterial networks at ∼1 year of age.
To our knowledge, this is the largest study of the early gut microbiome in children in whom autoimmunity is developing. Potentially relevant findings in relation to anti-islet cell autoantibodies did not appear to be focused on individual microbiota, but on their connectivity. Moreover, the gut microbiome at an early age was strongly influenced by factors such as delivery mode, fundamental changes in nutrition, and the shift to an adult like microbiome. We suggest that a systemic view is necessary to understand the complex relationship between the development of type 1 diabetes, the environment, and the gut microbiome.
Acknowledgments. The authors thank Sandra Hummel (Helmholtz Zentrum München), Christiane Winkler (Helmholtz Zentrum München), Annette Knopff (Helmholtz Zentrum München), and Melanie Bunk (Helmholtz Zentrum München) for data management and clinical assistance.
Funding. This work was supported by grants from the Juvenile Diabetes Research Foundation (17-2012-16 and 17-2011-266), Deutsche Forschungsgemeinschaft (DFG; grants ZI-310/14-1 to -4), and the Children With Type 1 Diabetes Foundation (Stiftung Das Zuckerkranke Kind). E.B. is supported by the DFG Research Center and Cluster of Excellence, Center for Regenerative Therapies Dresden (grant FZ 111).
Duality of Interest. No potential conflicts of interest relevant to this article were reported.
Author Contributions. D.E., W.z.C., and M.H. performed the data analysis, drafted the manuscript, and critically reviewed and approved the manuscript. A.A., A.G.D.-R., K.A.G., J.R.F., J.C.D., C.T.B., B.K., and E.W.T. performed the 16S sequencing and data analysis, and critically reviewed and approved the manuscript. P.A. and M.P. contributed to data collection, islet autoantibody measurement, and interpretation and analysis of the data, and critically reviewed and approved the manuscript. M.A. and D.S. contributed to the design of this study, data analysis, data interpretation, and manuscript writing, and critically reviewed and approved the manuscript. E.B. contributed to the design of this study and data analysis, and critically reviewed and approved the manuscript. E.W.T. contributed to the design of this study, performed the 16S sequencing and data analysis, and critically reviewed and approved the manuscript. A.-G.Z. is principal investigator of the BABYDIET study; contributed to the design of this study, data collection, data analysis, and manuscript writing; and critically reviewed and approved the manuscript. A.-G.Z. is the guarantor of this work and, as such, had full access to all the data in the study and takes responsibility for the integrity of the data and the accuracy of the data analysis.