A synbiotics, long chain polyunsaturated fatty acids, and milk fat globule membranes supplemented formula modulates microbiota maturation and neurodevelopment

Original article A synbiotics, long chain polyunsaturated fatty acids, and milk fat globule membranes supplemented formula modulates microbiota maturation and neurodevelopment Tom as Cerd o a , b , c , 1 , Alicia Ruíz d , e , 2 , Inmaculada Acu~ na d , e , Ana Nieto-Ruiz a , b , f , Estefanía Di eguez a , b , Natalia Sepúlveda-Valbuena b , 3 , Mireia Escudero-Marín a , b , c , Jose Antonio García-Santos a , b , f , María García-Ricobaraza a , b , f , Florian Herrmann a , b , Jose Antonio Moreno-Mu~ noz g , Roser De Castellar g , Jesús Jim enez g , Antonio Su arez d , e , * , Cristina Campoy a , b , c , f , h a Department of Paediatrics, School of Medicine, University of Granada, Granada, Spain b EURISTIKOS Excellence Centre for Paediatric Research, Biomedical Research Centre, University of Granada, Granada, Spain c Neurociences Institute “Doctor Ol oriz”, Health Sciences Technological Park, Granada, Spain d Department of Biochemistry and Molecular Biology 2, University of Granada, Spain e Nutrition and Food Technology Institute “Jos e Mataix”, Biomedical Research Centre, University of Granada, Spain f Instituto de Investigaci on Biosanitaria de Granada (IBS-GRANADA), Health Sciences Technological Park, Granada, Spain g Ordesa Laboratories, S.L. Barcelona, Spain h Spanish Network of Biomedical Research in Epidemiology and Public Health (CIBERESP), Granada's node, Institute of Health Carlos III, Madrid, Spain article info Article history: Received 13 October 2021 Accepted 18 May 2022 Keywords: Infant formula Gut microbiota Enterotypes Neurodevelopment summary Background &aims: The critical window of concurrent developmental paths of the nervous system and gut microbiota in infancy provides an opportunity for nutritional interventions with potential health benefits later in life. Methods: We compared the dynamics of gut microbiota maturation and explored its association with neurodevelopment at 12 months and 4 years of age in 170 full-term healthy infants fed a standard formula (SF) or a new formula (EF) based on standard formula supplemented with synbiotics, long chain polyunsaturated fatty acids (LC-PUFA) and bovine milk fat globule membranes (MFGM), including a breastfed reference group (BF). Results: Using Dirichlet Multinomial Modelling, we characterized three microbial enterotypes (Mixed, anaerobic and aerobic profile; Bact, Bacteroides-dominant; Firm, Firmicutes-enriched) and identified a new enterotype dominated by an unidentified genus within Lachnospiraceae (U_Lach). Enterotypes were associated with age (Mixed with baseline, U_Lach with month 6, Bact and Firm with months 12 and 18). Trajectories or timely enterotype shifts in each infant were not random but strongly associated with type of feeding. Trajectories in SF shifted from initial Mixed to U_Lach, Bact or Firm at month. Microbiota maturation in EF split into a fast trajectory as in SF, and a slow trajectory with Mixed to U_Lach, Bact or Firm transitions at months 12 or 18, as in BF. EF infants with slow trajectories were more often inehome reared and born by vaginal delivery to mothers with pre-pregnancy lean BMI. At 12 months of age, language and expressive language scores were significantly higher in EF infants with fast trajectories than in BF. Neurodevelopmental outcomes were similar between EF infants with slow trajectories and BF at 12 months and 4 years of age. *Corresponding author. Lab 119, Nutrition and Food Technology Institute “Jos e Mataix”, Biomedical Research Centre University of Granada, Avd.Conocimiento sn, 18016, Armilla. Spain. E-mail address: [email protected] (A. Su arez). 1 Present address. Maim onides Institute for Research in Biomedicine Cordoba (IMIBIC), C ordoba, Spain. 2 Present address. Centre for Inflammation Research, Queen's Medical Institute, University of Edinburgh, Edinburgh, United Kingdom. 3 Present address. Nutrition and Biochemistry Department. School of Sciences, Pontificia Universidad Javeriana. Bogot a, Colombia. Contents lists available at ScienceDirect Clinical Nutrition journal homepage: http://www.elsevier.com/locate/clnu https://doi.org/10.1016/j.clnu.2022.05.013 0261-5614/©2022 The Author(s). Published by Elsevier Ltd. This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/licenses/by-nc-nd/4.0/). Clinical Nutrition 41 (2022) 1697e1711 Conclusions: Feeding a synbiotics, LC-PUFA and MFGM supplemented formula in a specific infant environment promoted probiotic growth and retarded gut microbiota maturation with similar neurodevelopment outcomes to breastfed infants. Clinical trial registry number: NTC02094547. ©2022 The Author(s). Published by Elsevier Ltd. This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/licenses/by-nc-nd/4.0/). 1. Introduction Microbial colonization of an infant's gastrointestinal tract after birth is a dynamic and non-stochastic process of crucial importance for nutrient processing, energy harvest, and maturation of host immune system [1]. At the same time, the formation and refinement of neural networks responsible for a vast repertoire of behaviours and learning processes is being stablished [2]. Physiological connections between both processes have been shown in animal models, suggesting a bidirectional microbiotaegutebrain axis that we are currently beginning to understand [3,4]. In healthy infants, shifts in gut microbial composition have been associated to fine motor skills at 18 months of life [5], to childhood temperament at 18e27 months of age [6], to cognition at 2 years of age [7], and to communication, motor, personal, and social skills at 3 years of age [8]. According to the WHO, UNICEF, and the ESPGHAN Committee on Nutrition recommendations, human milk is the gold standard for infant nutrition that has been safely replaced by formula in case of lactation failure or insufficient breast milk supply to meet infant's nutritional needs [9]. In an effort to match the composition of infant formulas to that of human milk, formulas are currently being supplemented with a variety of bioactive compounds and biotic components that influence gut microbial composition and neurodevelopment [10]. Preand post-natal long chain polyunsaturated fatty acids (LC-PUFA) supply influence gut microbial composition and improve visual acuity scores, psychomotor development scores and resting-state brain network functioning [11e14]. Prebiotic oligosaccharides, probiotics or their combination (synbiotics) have been shown to improve infant growth and influence cognitive and adaptive behaviour as well as attention deficit hyperactivity disorder and neurobehavioral outcomes [15e19]. Recently, commercial formulas are also being supplemented with milk fat globule membranes (MFGM), a tri-layered membrane rich in sphingolipids and proteins with beneficial effects on the microbial ecosystem and infant cognition [20e22]. The above studies lead to the hypothesis that an infant formula combining the beneficial effects of synbiotics, LC-PUFA, and MFGM may selectively promote gut microbial maturation and neurodevelopmentsimilarlytothoseofbreastfedinfants.Hereinweaimed 1) to compare the succession dynamics of microbial consortia and individual community members in the developing gut microbiota of infants fed a standard formula supplemented with synbiotics, LCPUFA and bovine MFGM with infants fed the un-supplemented standard formula, 2) to determine maternal, perinatal, anthropometric, environmental,lifestyle and dietary factorsthat influencethe maturation of the gut microbial community, and 3) to explore the association between gut microbiota maturation and neurodevelopment tested with the Bayley Scales of Infant Development III 3rd edition (BSID-III) [23] at 12 months of age, and language developmentwithOralLanguageTaskofNavarra-Revised(PLON-R)[24]at 4 years of age. Our group previously reported the beneficial effects of the combined supplementationwith synbiotics, LC-PUFA and bovine MFGM on infant plasma fatty acid levels, growth, visual potentials, childbehaviourat 2.5 years andlanguage skills at 4 years [11,25e28]. 2. Materials &methods 2.1. Subjects and experimental design TheCOGNIS study (ANeurocognitive and Immunological Studyof a New Formula for Healthy Infants) was designed as a prospective, double-blind randomized clinical trial witha nutritional intervention based on bioactive nutrients-enriched infant formula and long-term follow-up of recruited infants, registered at www.ClinicalTrials.gov identifier NTC02094547. Recruitment of participants was performed between 2010 and 2014 at the EURISTIKOS Excellence Centre for Pediatric Research, the School of Medicine, and the Mind, Brain and Behaviour Research Center (CIMCYC) at the University of Granada (Spain). Infants were also recruited from outpatient centres in the Granada province. The inclusion criteria of participants were: healthy terminfants(37weeksand41weeksofgestation); adequatebirth weight for gestational age (between 3 and 97 percentile); normal APGAR score at 1 0 and 5': 7e10; umbilical pH 7.10; age of inclusion: 0e2 months (60 days) in the formula fed groups, 0e6 months (180 days) in the breastfeeding group; maximum first 30 days of exclusive breastfeeding in the formula fed groups; after this 30 days, exclusive or majority infant formula intake (>70% or >4 doses/day of infant formula); exclusive breastfeeding minimum 2 months in the breastfeeding group; availability to continue throughout the study period; signature of informed consent by parents/guardians. The exclusion criteria of participants were: infants who were participating in other study; exclusively breastfed infants who received formula for more than 25% of milk intake after the initial month 2 till month 6; infants who suffered nervous system abnormalities (hydrocephalus, perinatal hypoxia, intraventricular haemorrhage, neonatal meningitis, Abbreviations AIC Akaike information criterion Bact Bacteroides enterotype BF exclusively breast-fed reference group for at least 2mo BSID-III Bayley Scales of Infant Development III 3 rd edition DMM Dirichlet Multinomial Mixture EF standard infant formula supplemented with synbiotics, LC-PUFA and MFGM FDR False Discovery Rate Firm strict Firmicutes anaerobes enterotype LC-PUFA long chain polyunsaturated fatty acids LOESS locally weighted regression spline smoothing Mixed mixed aerobic-anaerobic enterotype MFGM milk fat globule membranes OTUs operational taxonomic units PLON-R Oral Language Task of Navarra-Revised SF standard control infant formula U_Lach Unclass_Lachnospiraceae enterotype T. Cerd o, A. Ruíz, I. Acu~ na et al. Clinical Nutrition 41 (2022) 1697e1711 1698 septic shock, West’syndrome …); infants who suffered gastrointestinaldisorders(cow'smilkproteinallergyand/orlactoseintolerance); pathological background of the mother and/or history of mental illness during pregnancy (neurological diseases, metabolopathies, type 1 diabetes mellitus, chronic disease (hypothyroidism), maternal malnutrition, TORCHcomplex infections);motherstaking anxiolytics or antidepressants, and other treatments with drugs potentially affecting neurodevelopment; parents inability to follow the study. Eligible infants aging 0e2 months (initial visit from now on referred as to timepoint 1 month) were assigned to receive a standard infant formula (SF) or an experimental formula (EF) consisting of SF supplemented with bovine MFGM components (10% of total protein content (wt:wt)), LC-PUFA (arachidonic and docosahexaenoic acids) and synbiotics (mix of fructooligosaccharides:inulin (ratio 1:1) and sialic acid, B.infantis IM1 (Bifidobacterium longum subsp. infantis strain CECT 7210) and Lactobacillus rhamnosus LCS-742) (SITable1).A mathematical statistical method was applied to randomize infants intoSForEFgroups(ratio1:1).Formulaswerestoredatstudysitesina secure and limited access storage area protected from extremes of light, temperature, and humidity. Both formulas were delivered in identical containers with different colour labels and codes to parents in a box with 12 cans of 400 g of the corresponding infant formula, which covered infant feeding for approximately 1 month. Infants received initiation formula up to 6 months of age, and follow-on formula was given between 6 and 18 months of age. All infant formulas were provided by Laboratorios Ordesa, S.L. (Barcelona, Spain). Both infant formulas followed guidelines of the Committee on Nutrition of the European Society for Paediatric Gastroenterology, Hepatology and Nutrition [29] and the international and national recommendations for the composition of infant formulas. Fifty infants who were exclusively breastfed (BF) for at least 2 months were enrolled and included asa referencegroup. Data onanthropometrics, medication use, clinical symptoms and feeding pattern (type of feeding prior to the study, frequency offeeding inbreastfed infants or total milk volume intake per day) as well as data on both formula digestive tolerance (mainly stools, vomiting, regurgitation, and colic.) were collected by the pediatrician during follow-up visits, and subsequently included in the Pediatrician's Data Collection Questionnaire. In addition, a three-day dietary record was used to collect quantitative data about foodand drinks consumed duringthree days, including aweekend day and twoworking days.Double-blind design for parents, clinicians and researchers was continued throughout the study. 170 infants were enrolled in formula feeding groups whose baseline characteristics and study flowchart are shown in SI Table 1 and Fig. 1,respectively. 2.2. Ethics, consent, and permissions This study was carried out following the updated Declaration of Helsinki Principles [30], the Good Clinical Practice recommendations of the EEC (document 111/3976/88 July 1990), as well as the current Spanish legislation governing clinical research in humans (Royal Decree 561/1993 on clinical trials). The study protocols were also approved by the Research Bioethical Committee from the University of Granada, and the Bioethical Committees for Clinical Research of the Clinical University Hospital San Cecilio and the Mother-Infant University Hospital of Granada (Granada, Spain). All families were informed about procedures, and a signed written informed consent was obtained from each parent or legal guardian for their offspring. 2.3. Genomic DNA extraction Fresh stool samples were collected in sterile bottles at home by the parents, following provided instructions, and stored at 20  C for a maximum of 24 h until delivery to the laboratory, where they were stored to 80  C. Genomic DNA was extracted from faecal bacteria as previously described [31]. 2.4. 16S rRNA gene sequencing and data processing Genomic DNA from faecal bacteria collected at 1, 6, 12 and 18 months was used as template for 16S rRNA gene amplification using 27F and 338R universal primers for V1eV2 region as previously described (n ¼412) [31] on an Illumina MiSeq platform (University of Granada, Spain). Nucleic acid manipulations and reactions were performed by the same researcher and sequences were obtained with two independent arrays. Sequences were grouped in a single dataset that was quality-filtered. Operational taxonomic units (OTUs) were taxonomically classified using Ribosomal Database Project [32]. OTUs were considered unassigned when confidence value score was lower than 0.8, and were annotated using upper taxonomic ranks. 2.5. Assessment of infant neurodevelopment Infant neurodevelopment at 12 months of age using BSID-III [23] and at 4 years of age using PLON-R [24] were performed by trained psychologists in the presence of the mother of the child. Recorded scores were within the range for normal healthy infants of their corresponding age. 2.6. Statistical analysis To detect a minimum difference of 0.5 standard deviations in cortical visual evoked potentials (primary outcome) with a statistical power of 90% and a ¼0.05, the sample size was established in 71 infants/formulagroup.Thissamplesizealloweda statisticalpowerof 85% to detect a minimum difference of 0.6 SD in microbiota (secondary outcome). Due to long-term follow-up of CONGNIS study a potentialdropoutrate of20%wasconsideredand,asaresult,thefinal sample size was estimated in 85 infants/formula group. The effect size and significance of study variables on gut microbiota composition was determined using the envfitfunction in vegan package [33]. a -diversity was measured with a phylogenetic diversity measurement, Faith's phylogenetic diversity, the Shannon's Diversity Index, a non-phylogenetic measurement of bacterial abundance (richness) and Rao's quadratic entropy at the OTU level using phyloseq and picante packages in R [34,35]. a -diversity differences due to age and typeof feedingwere analysed using repeated-measures ANOVAwith Bonferroniadjustmentformultipletesting.Enterotypingof thestudy cohort was performed following the Dirichlet Multinomial Mixture (DMM) method on taxa classified at genus level present in 80% of samples [36]. Ordination of enterotypes was performed using NMDS basedonBrayeCurtisdissimilaritymetrics.Differencesin a -diversity between enterotypes were assessed using either ANOVA and post hoc Tukey test or KruskaleWallis and post hoc Wilcoxon test, with BenjaminieHochberg adjustment for multiple testing. Grouping of study variables is shown in SI Table 1. Significant differential phylotypeabundanceatseveraldifferenttaxonomylevelswasconstructed from non-normalized raw count tables with DESeq2 package adjusted for covariates using a two-sided Wald's test with multiple comparisons correction by the Benjamini-Hochberg method [37]. Probability distributions were analysed with ShapiroeWilk normality test, KolmogoroveSmirnov Goodness-of-Fit Test and Hartigans' Dip Test for Unimodality using stats,dgof and diptest packages. Logistic regression analyses were used to examine the association of study variables with enterotypes and trajectories after collinearity correction using nnet,car,survey and caret packages. For neurodevelopment, a multivariate ANCOVA model with covariate T. Cerd o, A. Ruíz, I. Acu~ na et al. Clinical Nutrition 41 (2022) 1697e1711 1699 correctionandpost hoc Benjamini-Hochbergadjustmentformultiple comparisons was performed. For all determinations, the significance cut-off was set at p0.05 or False Discovery Rate (FDR) 0.05 when multiple test correction was applied. 3. Results 3.1. Baseline characteristics of participants All children who participated in COGNIS study were born healthy at term, mostly by spontaneous vaginal delivery (SI Table 1). The median maternal age at delivery was 31.5 years old (interquartile range 27e35.5) and they were usually non-smokers during pregnancy (more than 86%). Median pre-pregnancy maternal BMI was 24.0 (interquartile range 21.64e26.91), and no differences in weight gain during pregnancy were observed between COGNIS study groups. Maternal IQ and parent's educational level were significantly higher in BF than in formula groups (p<0.01). Formulas were well tolerated by all infants and overall satisfaction of pediatricians and parents was excellent. No differences between study groups were found in infant anthropometric data at birth, including weight, length and head circumference. When considering drop-out subjects, EF infants showed higher weight and length at birth than SF ones (p¼0.012; p¼0.023, respectively). Due to the COGNIS study design, days of breastfed significantly differed between study groups (p<0.001), but not between infant formula groups [median duration of 2.5 days (interquartile range 1e25)]. No differences in terms of diarrhea episodes (p¼0.283), antibiotic treatment (p¼0.579) and fever episodes (p¼0.972) were observed between groups. In agreement with recommendations of Spanish Association of Paediatrics, solid food was introduced at a median age of 17 weeks (interquartile range 16e19) in infant formula groups. In BF, 68% of infants continued breastfeeding at 18 months of life, and solid food introduction was delayed till month 6 (p<0.001). 3.2. Taxonomic profiling and factors shaping the gut microbiota After quality filtering and high-confidence phylogenetic annotation, 18, 175, 653 16S rRNA sequence reads ranging from 8500 to 183,967 per sample (Mean ¼43,221; SD ¼26,837) rendered a gut microbial profile consisting of 686 OTUs assigned to 92 distinct genera belonging to 38 families. The phylogenetic composition and categorical breakdown of OTUs are shown in SI Table 1. The gut microbiota composition of each infant at genus level binned by age is shown in Fig. 2A. The most abundant taxa at genus level were Streptococcus, Bacteroides,Laschnospiracea_incertae_sedis and two unassigned genera within the families Enterobacteriaceae (unclass_Enterobacteriaceae) and Lachnospiraceae (unclass_Lachnospiraceae), accounting for 51.4% of total reads. We next sought to determine the significant factors associated with infant's gut microbial communities in our dataset as determined by EnvFit (Figs. 2B and SI Table 1). Twenty-nine microbiota covariates spanning maternal, nutritional, lifestyle and paediatric factors were tested after removing colinear variables. When all samples were analysed together, sampling age explained the greatest amount of variance (56.5%, p<0.001) in gut microbiota composition, followed by type of feeding (3.2%, p<0.004). Hierarchical clustering based on BrayeCurtis dissimilarity metrics illustrated the correlation of age and type of feeding with overall microbial community variation (SI Fig. 2A). Additionally, breastfeeding, antibiotic treatment, diarrhoea episodes and starting age of other dairy products and fish introduction together resulted in a total significant additive effect size of 10.2%. Breastfeeding was a better factor than mixed breast and formula-feeding, underscoring that duration of breastmilk feeding impacted gut microbial assemblages in our cohort rather than consumption of any breast Fig. 1. Participant flow chart from inclusion until 18 months of age per dietary group. D ¼drop outs and E ¼exclusions (1 infant to perinatal hypoxia, 1 infant to digestive surgical intervention, 1 infant to epileptic seizure, 1 infant to hydrocephalus, 3 infants to deficiency of growth, 25 infants to not ingestion of formula, 2 to colic of the infant, 5 infants to lactose intolerance and 1 infant of BF group to not is breastfeeding). BF: breastfed infants; SF: standard infant formula; EF: experimental infant formula; D: dropouts; E: exclusions. *SF and EF infants were randomized between 0 and 2 months of age; BF infants were randomized between 0 and 6 months of age. T. Cerd o, A. Ruíz, I. Acu~ na et al. Clinical Nutrition 41 (2022) 1697e1711 1700 Fig. 2. Microbiome profiles, variance explained by variables and ecological diversity along the chronosequence of the COGNIS study (A) Phylogeny of infants' gut microbiota at the genus level (top 30 genera) based on 16S rRNA gene sequencing data. (B) Horizontal bars show the amount of variance (r 2 ) explained by each variable in the model calculated by EnvFit in the dataset and at the specified time points. Significant variables (false discovery rate (FDR) p<0.05) are coloured based on overall metadata group. Colour bars of nonsignificant variables are translucent. (C) Gut microbial diversity was measured with a non-phylogenetic measurement of bacterial abundance (species richness), the Shannon's Diversity Index sensitive to species richness and evenness, and a phylogenetic diversity measurement, Faith's phylogenetic diversity. Curves show LOESS fit for the data per type of feeding, and shaded areas show permutation-based 95% confidence intervals for the fit. Samples are coloured by type of feeding: SF, Infants fed a non-supplemented infant formula; EF, Infants fed SF supplemented with synbiotics, LC-PUFA and milk fat globule membranes; BF, Infants fed with human milk. Repeated-measures ANOVA with Bonferroni adjustment was used. When statistical differences were observed, different letters indicate significant differences. T. Cerd o, A. Ruíz, I. Acu~ na et al. Clinical Nutrition 41 (2022) 1697e1711 1701 Fig. 3. Probabilistic modelling with Dirichlet Multinomial Mixtures based on lowest Laplace approximation of infant faecal samples revealed four enterotypes. (A) Heatmap showing the relative abundance of the 30 most important signature genera per enterotype. (B) Principal component analysis on genus-level data with samples coloured according to DMM cluster. Vertical and horizontal bar charts depict contribution of core genera to axis loadings. (C) Species richness, Shannon's diversity index and Faith's phylogenetic diversity of enterotypes. Different letters indicate significant differences. by the KruskaleWallis rank-sum test followed by Dunn's test with FDR correction p<0.05 (D) Importance to model prediction of enterotypes by microbial signature genera. (E) Contribution of microbial signature genera (mean relative abundances) to enterotype profiles. (F) Three-way T. Cerd o, A. Ruíz, I. Acu~ na et al. Clinical Nutrition 41 (2022) 1697e1711 1702 milk. Mode of delivery, birth weight, smoking during gestation, gestational age, maternal age, pre-pregnancy BMI and weight gain, gender, siblings, household pets or day care exposure showed no association with the overall microbiota phylogenetic makeup. For the temporal pattern of microbial maturity, covariates were also analysed by binning samples at each time point separately (Fig. 2B). At month 1, only father smoking habit influenced gut microbial assemblies. Breastfeeding and type of feeding had a significant effect on residual variance at month 6 and, specially, at months 12 and 18, underscoring the impact of the nutritional intervention on gut microbial patterns. Additionally, at month 6, day care exposure and age of meat and vegetable introduction explained 22.7% of the residual variation. The highest influence of complementary food was observed at month 12 when starting age of solid food consumption and of gluten-free cereals, fruits, other dairy product and fish intake accounted for a total additive effect size to 25.1%. Finally, antibiotic treatment had a significant effect (8.9%) on microbiota composition at 18 months of age. 3.3. Postnatal kinetics of gut microbiota diversity in experimental groups The above pattern of many rare compared to few abundant populations suggested strong inter-sample variation in infant's gut microbiota. Rao's quadratic entropy showed increased a -diversity but reduced b -diversity as a function of time, suggesting that the gut microbial ecosystem accumulated diversity into less heterogeneous configurations (SI Fig. 2B). Increasing average values of microbial richness, Shannon's diversity index and Faith's phylogenetic diversity characterized the temporal evolution of gut microbial communities in infants. When samples were stratified by time point (Fig. 2C), microbial richness was significantly higher in SF and EF compared to BF infants at months 6 and 18. Interestingly, microbial richness was significantly higher in SF than in EF and BF at month 12. At 18 months of age, BF had lower Faith's phylogenetic diversity than EF and SF. No differences were observed in Shannon's diversity index between experimental groups. 3.4. Enterotype stratification of infants’gut microbiota To investigate the potential association of age and feeding mode with prevalence of gut microbial community profiles, we identified community types or enterotypes in infant's samples by using Dirichlet Multinomial Mixtures (DMM) modelling that assigned them into clusters (lowest Laplace approximation) based on the relative abundance of the microbial groups at genus level of classification. Model fitting was set to k ¼8 and rendered an optimum number of 4 DMM community types, here referred as enterotypes (Figs. 3A and SI Fig. 3A). The four-gut microbial enterotypes had weights P ¼0.33,0.32,0.20and0.16(componentofweightforwhich smaller values correspond to less frequent communities). The enterotypes also differed in how variable their communities were with Q ¼6.22, 12.76, 36.43 and 26.08 (component of variability for which smaller values correspond to highly variable communities). Thus,gutmicrobialprofileswereorganizedintotwohighlyabundant and variable enterotypes, and two less abundant homogeneous enterotypes. Non-metric multidimensional scaling of BrayeCurtis distances illustrated these features of the clusters (Fig. 3B). We observed that enterotypes were classified by a core of nineteen genera that accounted for 86.8% of total reads with mean total reads from 0.45% (Escherichia/Shigella) to 16.1% (unclass_Lachnospiraceae) (SI Fig. 3B). Highly abundant genera (>1% of total reads) like Prevotella,Clostridium XVIII, Lactococcus and Erysipelotrichaceae_incertae_sedis did not drive enterotype classification. The diversity of the samples assigned to each enterotype indicated that microbial richness (Dunn's test FDR <0.001) and Faith's phylogenetic diversity (Dunn's test FDR <0.001) increased from the first to the fourth enterotype while Shannon's diversity index was not significantly different between the second and third enterotypes (Dunn's test FDR ¼1) (Fig. 3C). The firstenterotype contained samples whose top microbial signature genera were a mixed population of highly abundant facultative anaerobic microorganisms belonging to Actinobacteria,Firmicutes and Proteobacteria such as Streptococcus, Enterococcus,Klebsiella,Lactobacillus and unclass_Enterobacteriaceae, and the highest abundance of Bifidobacterium among all enterotypes (Figs. 3D and SI Fig. 3C). The second enterotype showed a high prevalence of unclass_Lachnospiraceae whereas the third enterotype was dominated by Bacteroides, these microbes being efficient degradersofdietaryfibers.Topmicrobialpredictivegeneraofthefourth enterotype were a group of strict Firmicutes anaerobes within Lachnospiracea_incertae_sedis,unclass_Lachnospiraceae,Fusicatenibacter, Blautia,Roseburia and Faecalibacterium, active producers of short chain fatty acids. On the basis of their respective genus-level dominance profiles, we referred to enterotypes as mixed aerobicanaerobic dominant type (Mixed), unclass_Lachnospiraceae dominant type (U_Lach), Bacteroides dominant type (Bact) and strict Firmicutes anaerobes dominant type (Firm). Core signature genera of enterotype classification accounted for 89.5% (Mixed; SD ¼6.2), 86.8% (U_Lach; SD ¼6.8), 86.8% (Bact; SD ¼1.6) and 77.9% (Firm; SD ¼4.3) of mean total abundances in their corresponding samples (Fig. 3E). These results showed that signature generawere dominant in their enterotypes where Bacteroides and unclass_Lachnospiraceae stoodoutascontributorstooverallenterotypeclassification(Fig.3D). At phylum level, the Mixed enterotype was characterized by the highest Firmicutes/Bacteroidetes ratio and relative abundances of Actinobacteria and Proteobacteria (SI Fig. 3D and E). The U_Lach and Firmenterotypessharedthe highestrelativeabundanceof Firmicutes while Bact enterotype was characterized by the lowest Firmicutes/ Bacteroidetes ratioandhighestrelativeabundanceofBacteroidetes.To test for associations between covariates and enterotypes, we fita multinomiallogisticregressionmodeladjustedforstudyvariables(SI Table 1). Multinomial logistic regression analyses indicated that belongingtoanyenterotypewasstronglyassociatedwithinfant'sage and feeding group. These significant associations were visualized with a three-way association plot and mosaic plot by time point, feeding group and enterotype (Figs. 3F and SI Fig. 3F). All but one samples at baseline belonged to Mixed enterotype and, therefore, represented the initial stage of gut microbial configurations in our dataset. As infants aged, infants' gut microbiota diversified and became progressively dominated by U_Lach, Bact and Firm communities. At 6 months of age, most samples in SF and EF groups belonged to U_Lach enterotype. Bact and Firm enterotypes were the mostdominantinSFandEFgroupsat12and18 months,respectively. In contrast, enterotypes in BF infants mostly belonged to Mixed and U_Lach assemblies at all time points, suggesting that driver genera of these enterotypes were determinant in BF microbial community configurations. In addition, regression models indicated that the Mixed enterotype was mainly associated to maternal pre-pregnancy BMI. The chance that an infant microbiota belonged to U_Lach enterotype was increased in mothers with pre-pregnancy obesity association and common angle plot on the prevalence of enterotypes and their association to age (months) and type of feeding. Colours represent the level of the residual for that cell/combination of levels. Blue means more observations and red fewer observations in that association than would be expected under the null model (independence). SF, Infants fed a non-supplemented infant formula; EF, Infants fed SF supplemented with synbiotics, LC-PUFA and milk fat globule membranes; BF, Infants fed with human milk. T. Cerd o, A. Ruíz, I. Acu~ na et al. Clinical Nutrition 41 (2022) 1697e1711 1703 andgestationalage>38weeks.TheBactenterotypewasassociatedto siblingswhilebelongingtoFirmenterotypewasstronglydetermined by gestational age and siblings (SI Table 1). No association was observed with gender, breastfeeding, mode of delivery, household pets or day care exposure. 3.5. Evolutive trajectories of enterotypes There were differences in the frequency of subject-independent enterotype transitions (changes in enterotype) between experimental groups. Transitions were quantified to create a general multi-state Markov Chain model where nodes were enterotypes and edges reflect transition rates by their weight, that is, probabilities of changing to another enterotype at any time point (Fig. 4A). Despite the high variability of transitions, the Markov chain model revealed preferences for transitions by type of feeding. Self-enterotype transitions were significantly more frequent in BF infants (62% of transitions) compared to SF (30%; Fisher's twotailed exact test, p¼0.01, OR ¼0.27) and EF (38%; Fisher's twotailed exact test, p¼0.03, OR ¼0.35) infants, suggesting that breastfeeding was associated with higher microbial community stability. The markov chain model showed that the most prevalent transition in SF infants was from Mixed to U_Lach followed by U_Lach to Bact and self-Bact transition. In EF infants, selftransitions were more prevalent compared to SF infants, except for Bact self-transition. Transitions in EF infants also involved Mixed to U_Lach, U_Lach to Bact and Bact to Firm enterotypes. In BF infants, the most prevalent transitions were self-Mixed and selfU_Lach transitions, followed by Mixed to U_Lach enterotype transition, showing scarce progress to either Bact or Firm enterotypes during their 18 months of life. We next modelled enterotype evolutive trajectories in each infant along the chronosequence. From months 1e18, the gut microbiota was highly dynamic and twenty-five distinct types of enterotype trajectories were identified (Fig. 4B) where the two most prevalent trajectories evolved from Mixed to U_Lach at month18 and step-by-step from Mixed, U_Lach, Bact to Firm enterotype. In SF and EF infants, trajectories from initial Mixed enterotype towards U_Lach, Bact and Firm enterotypes were temporarily unidirectional. This evolutive pattern was gradual with 77% and 80% of SF and EF infants, reaching enterotypes Bact or Firm at 18 months. In contrast, trajectories in BF infants only switched from Mixed to U_Lach enterotype at 12 but not at 6 months of age, indicating that introduction of complementary food (median ¼23.5 weeks, range 17e30) did not impact enterotype transitions in BF infants. Most BF infants stayed in U_Lach enterotype at month 18. Surprisingly, a few retrogressions in trajectories, that is, returning to a previous enterotype with time, were observed in 20.3% of formula-fed infants but not in BF infants. Retrogressions from Bact to U_Lach, Firm to U_Lach and Firm to Bact enterotype were observed between 12 and 18 months of age. Since retrogressions can be considered as a major disturbance in microbial ecosystem evolution, we wondered whether extrinsic factors like diarrhoea episodes or antibiotic treatment were associated to gut microbiota retrogressions. While the risk of retrogression associated to antibiotic treatment was not significant (Fisher's two-tailed exact test, p¼0.37, OR ¼0.50), 69% of infants with retrogressions reported diarrhoea episodes compared to 50% of infants without retrogression events (Fisher's two-tailed exact test, p¼0.02, OR ¼0.23). 3.6. Maternal and perinatal factors determine enterotype trajectories in EF infants When we ranked progressive trajectories by enterotype and time point in each individual, we observed that distributions of enterotype trajectories in intervention groups were not normal for SF and BF groups (ShapiroeWilk's normality test: SF, p¼0.01; BF, p¼0.001) and were significantly different between BF and formula-fed infants (KolmogoroveSmirnov test: BF vs SF D ¼0.66, p¼0.012; BF vs EF D ¼0.58, p¼0.04). Indeed, weighted density plots of trajectories showed two peaks of varying abundance in enterotype trajectories in SF and BF infants but skewed in opposite directions while the trajectory distribution in EF infants had two distinct equally-weighted abundance peaks (Fig. 5A). The bimodality of the distribution of enterotype trajectories in EF infants was confirmed by the Hartigans' Dip Test for Unimodality (D ¼0.11, MonteeCarlo test p¼0.007) but not for SF and BF distributions (D ¼0.087, MonteeCarlo test p¼0.156 for SF; D ¼0.0625, MonteeCarlo test p¼1 for BF). These results indicated that there were two different evolutive trajectories in EF infants with Mixed to U_Lach, Bact or Firm transition soon at 6 months (“fast”trajectory) or late at 12 or even 18 months (“slow”trajectory) (Fisher's twotailed exact test, p¼0.017, OR ¼0.32). The ecological richness and diversity during early life of “slow”and BF infants was lower compared to “fast”infants, and was significantly different at month6 (richness: “fast”vs “slow”p<0.0001, “fast”vs BF p<0.001; Shannon's diversity index: “fast”vs “slow”p<0.002, “fast”vs BF p<0.01) and month18 (richness: “fast”vs BF p<0.0001). We next questioned whether these distinct transition trajectories in EF infants were associated to gestational, anthropometric, nutritional and clinical factors. Multicollinearity checks between predictor variables revealed a significant relationship between siblings and day care exposure, gestational age and maternal prepregnancy BMI, and between gender and residence. Because of the large number of variables under consideration, the multivariable binomial logistic regression model was fit using variables with a starting p<0.20 value in univariate logistic regression analyses and a backwards elimination procedure to keep those with significance (p<0.05) was undertaken. The model discarded antibiotic treatment (p¼0.413), maternal age (p¼0.540) and IQ (p¼0.368), breastfeeding (p¼0.442), diarrhoea episodes (p¼0.371), vitamin D intake (p¼0.879), residence (p¼0.519), and age of solid food introduction (p¼0.572). Thus, the model was fit with mode of delivery, maternal pre-pregnancy BMI, smoking during gestation, day care exposure and household pets. The model fit revealed that “slow”and “fast”trajectories of EF infants were associated to mode of delivery (p¼0.049), day care exposure (p¼0.048) and strongly to maternal pre-pregnancy BMI (p¼0.0317) (Fig. 5B, C and E). The chance that an infant fed EF formula belonged to the “fast”trajectory was strongly increased among those born by C-type delivery [OR (95% CI) ¼12.4 (1.93e127.8)] to mothers with higher BMI (OR (95% CI) ¼0.74 (0.55e0.92)] and spending more time at day care centres [OR (95% CI) ¼0.18 (0.03e0.57)]. In return, “slow infants”in EF cohort were more frequently inehome cared infants born by vaginal delivery to mothers with pre-pregnancy lean BMI (AUC ¼0.84) (Fig. 5D). The logistic model was tested against the null hypothesis and was confirmed to be significant (Wald's test p<0.003). 3.7. Neurodevelopment outcomes and microbiota maturation at 12 months and 4 years of age We next explored whether distinct gut microbial maturation paces in EF infants associated with neurodevelopmental outcomes at 12 months of age. Bayley-III assessment showed that cognitive, motor and language scores were lower in BF and infants with “slow”gut microbial maturation compared with infants with “fast” gut microbial maturation though none reached statistical significance (Table 1). When controlling for covariates that affect infant development (gender, siblings, maternal and paternal ages, T. Cerd o, A. Ruíz, I. Acu~ na et al. Clinical Nutrition 41 (2022) 1697e1711 1704 maternal and paternal educational levels, maternal IQ and smoking during pregnancy), the language and expressive scores were significantly different in infants with “fast”gut microbial maturation compared with BF. At 4 years of age, PLON-R assessment showed no significant differences in language performance between maturation groups and BF, even after adjustment for confounders (SI Table 1). 3.8. Bifidogenic and Lactogenic effect of experimental formula To account for the effect of formula feeding on gut microbial compositional changes, we used the statistical software DESeq2 to identify taxa at the genus level with differential abundances between SF and EF infants at each time point. We identified sixteen genera overabundant in SF infants and eleven genera in EF infants along the chronosequence (Fig. 6A). These overabundant genera in SF and EF infants were phylogenetically different. Discriminating genera belonged to Firmicutes,Proteobacteria and unclassified_Bacteria in EF infants, and to Bacteroidetes,Firmicutes and Fusobacteria in SF infants. Focusing on the effect of formula feeding on signature genera of enterotype configurations, the abundance of Blautia, driver of U_Lach, Bact and Firm enterotypes, and Roseburia and Faecalibacterium, drivers of Bact and Firm enterotypes, was influenced by SF intake. Conversely, EF formula intake was associated with significantly higher abundances in Lactobacillus, driver of Mixed enterotype. At 6 months of age, the gut microbiota of SF infants showed a high prevalence of Hungatella,Roseburia,Anaerotruncus,Clostridium XI, Dorea,Faecalibacterium and Blautia within Clostridia. A higher phylogenetic diversity in differential genera was observed in EF infants. EF infants were enriched in genera within Erysipelotrichia (Catenibacterium), Clostridia (Coprococcus,Ruminococcus), Bacilli (Lactobacillus and an unassigned genus within Bacilli) and Betaproteobacteria (Parasuterella). At 12 months of age, an enrichment in Acidaminococcus,Holdemania,Megamonas,Roseburia,Paraprevotella, and an unassigned genus within Clostridiales_Incertae Sedis XIII was observed in SF infants while EF infants had higher abundance of Coprococcus,Parasutterella,an unassigned genus within Firmicutes and in three genera within Fig. 4. Dynamics of enterotype transitions and trajectories in infants fed a standard formula (SF), standard formula supplemented with synbiotics, LC-PUFA and milk fat globule membranes (EF) and human milk (BF). (A) Markov chain model depicting probabilities of subject-independent enterotype transitions in infant feeding groups. Nodes represent enterotypes whose size is proportional to their prevalence. Edges represent transition directions whose rates are shown numerically and by edge weight (thickness). Edges of most prevalent enterotype transitions within each intervention group are shadowed. Red arrows depict trajectory regressions. (B) Complete set of subject-dependent trajectories organized from Mixed to U_Lach, Bact and Firm enterotypes (bottom-top) and time point (left-right). Each row on the y axis represent a distinct trajectory. Trajectories within black boxes illustrate regressions. Prevalence of each trajectory is coloured by type of feeding. T. Cerd o, A. Ruíz, I. Acu~ na et al. Clinical Nutrition 41 (2022) 1697e1711 1705