RES E A R C H Open Access
Evolution of gut microbiota composition
from birth to 24 weeks in the INFANTMET
Cohort
Cian J. Hill
1,2
, Denise B. Lynch
1,2
, Kiera Murphy
1,2,3
, Marynka Ulaszewska
5
, Ian B. Jeffery
1
, Carol Anne O’Shea
4
,
Claire Watkins
3
, Eugene Dempsey
4
, Fulvio Mattivi
5
, Kieran Tuohy
5
, R. Paul Ross
1,2
, C. Anthony Ryan
2,4
,
Paul W. O’ Toole
1,2
and Catherine Stanton
2,3*
Abstract
Background: The gut is the most extensively studied niche of the human microbiome. The aim of this study was to
characterise the initial gut microbiota development of a cohort of breastfed infants (n = 192) from 1 to 24 weeks of age.
Methods: V4-V5 region 16S rRNA amplicon Illumina sequencing and, in parallel, bacteriological culture. The
metabolomic profile of infant urine at 4 weeks of age was also examined by LC-MS.
Results: Full-term (FT), spontaneous vaginally delivered (SVD) infants’ microbiota remained stable at both phylum and
genus levels during the 24-week period examined. FT Caesarean section (CS) infants displayed an increased faecal
abundance of Firmicutes (p < 0.01) and lower abundance of Actinobacteria (p < 0.001) after the first week of life
compared to FT-SVD infants. FT-CS infants gradually progressed to harbouring a microbiota closely resembling FT-SVD
(which remained stable) by week 8 of life, which was maintained at week 24. The gut microbiota of preterm (PT)
infants displayed a significantly greater abundance of Proteobacteria compared to FT infants (p <0.001)atweek1.
Metabolomic analysis of urine at week 4 indicated PT-CS infants have a functionally different metabolite profile than FT
(both CS and SVD) infants. Co-inertia analysis showed co-variation between the urine metabolome and the faecal
microbiota of the infants. Tryptophan and tyrosine metabolic pathways, as well as fatty acid and bile acid metabolism,
were found to be affected by delivery mode and gestational age.
Conclusions: These findings confirm that mode of delivery and gestational age both have significant effects on early
neonatal microbiota composition. There is also a significant difference between the metabolite profile of FT and PT infants.
Prolonged breastfeeding was shown to have a significant effect on the microbiota composition of FT-CS infants at 24 weeks
of age, but interestingly not on that of FT-SVD infants. Twins had more similar microbiota to one another than between two
random infants, reflecting the influence of similarities in both host genetics and the environment on the microbiota.
Background
The gut microbiota is increasingly regarded as an ‘invis-
ible organ’ of the human body and considered an im-
portant factor for host health. This dynamic microbial
population develops rapidly from birth until 2 to 3 years
of age, when adult-like composition and stability is
established [1, 2]. If the establishment of the stable adult
microbiota is programmed in infancy, it may lead to a
lifelong signature with significant effects on health.
Bacterial colonisation begin s at birth, although recent
papers have suggested microbiota acquisition occurs in
utero [3], challenging the traditional dogma of uterine
sterility. The developing gut microbiota of neonates dif-
fers widely between individuals [2] and both internal
host properties and external factors inf luence the estab-
lishment of the microbiota [4]. At birth, the infant mi-
crobial population resembles the maternal vagina or skin
microbiota depending on mode of delivery, i.e. by spon-
taneous vaginal delivery (SVD) or Caesarean section
(CS), respectively [5]. Birth mode has a significant effect
on the nascent neonatal gut microbiota after these initial
* Correspondence: Catherine.stanton@teagasc.ie
2
APC Microbiome Institute, University College Cork, Cork, Ireland
3
Teagasc Moorepark Food Research Centre, Fermoy, Co. Cork, Ireland
Full list of author information is available at the end of the article
© The Author(s). 2017 Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0
International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and
reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to
the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver
(http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.
Hill et al. Microbiome (2017) 5:4
DOI 10.1186/s40168-016-0213-y
founder populations have been replaced [6–9]. At 1 week
of age, the microbiota of the SVD infant gut is charac-
terised by high levels of Bifidobacterium and Bacteroides,
while Clostridium is more abundant in CS neonates [10].
Numerous other factors have been sh own to exert an
influence on this de velopment , including antibiotic ex-
posure [11] and breastfeeding [12, 13]. Development of
the microbiota occurs as bacteria are replaced in a dy-
namic, n on-random pattern [14, 1 5]. The use of infant
milk formula (IMF ) impac ts on metabolism [16] and
development of the neonatal immune system [17]. This
introduction of IMF or solid food perturbs bacterial
colonisation [18, 19] and may reduce the benefit s of ex-
clusive breastfeeding [17].
Preterm (PT) neonates experience a number of unique
challenges to the establishment of their microbiota. CS
delivery, maternal and neonatal exposure to antibiotics
and the sterile environment of the neonatal intensive
care unit (NICU) may all alter the natural pattern of ac-
quisition of microbiota. A few published studies with
high subject numbers examining the P T gut microbiota
mainly focus on the initial hospitalised period [15, 20]. A
knowledge gap surrounding PT gut microbiota develop-
ment was recently highlighted [21], and to our know-
ledge, the current study is the largest, well-phenotyped
analysis of the longitudinal microbiota development of
PT infants after leaving the hospital environment. It has
previously been suggested that post-conceptional age, ra-
ther than post-birth age, is the main determinant of the
bacterial community profile in preterm infants [15]; the
aforementioned factors were found to influence the
pace, but not the sequence, of microbial acquis ition.
Metabolites have been shown to influence regulatory T
cells in the gut [22], with changes posited to contribute to
autoimmune diseases including inflammatory bowel disease,
asthma, allergies, arthritis and multiple sclerosis [23, 24].
These conditions ha ve also been linked to CS and PT birth.
In this prospective study, we compared the gut micro-
biota of initially breastfed infants from a single geo-
graphical area (Cork, Ireland) who were born under
different birth modes (SVD or CS) and different gesta-
tional ages (FT or PT), in the same maternity hospital.
We investigated the effect of both of these factors on
the establishment of the nascent gut microbiota of
breast fed infants. We also examined the link between
the microbiota and the metabolome in early life through
comparison of urine metabolomic data with 16S gut
microbiota data.
Methods
Participants and sample collection
The infants included in this study are part of the
INFANTMET study cohort. Mothers were approached
for consent between February 2012 and May 2014 at the
Cork Univ ersity Maternity Hospital, with ethical ap-
proval provided by the Cork University Hospital Re-
search Ethics Committee (ethical approval reference:
ECM (w) 07/02/2012). The study design was to recruit
groups of infants according to birth mode and gestation:
FT-SVD, FT-CS, PT-SVD and PT-CS infants (PT; less
than 35 weeks gestation). Information about the in-
fants was collected at delivery using medical records.
Further data were collected using detailed question-
naires given to the mothers when the infants were
1 year old (Additional file 1: Table S1). Faecal samples
were collected from the infant s at 1, 4, 8 and
24 week s of age (Table 1). P T infants were sampled
at 1 week of age and the same time point s (i.e. weeks
4, 8 and 24) after the due delivery date . Samples were
collected and placed at 4 °C by the mother, before
collection in a temperature -controlled transport col-
lection case by the re search nurse for transport to the
lab for DNA extraction. An additional sample was ac-
quired at the due date of delivery for PT infants.
Urine samples were also collected at 4 weeks of age
for metabolomic analysis using Sterisets Uricol Urine
Collection Pack (Medguard, Ireland). A pad was placed
in the diaper and used to collect an unsoiled urine sam-
ple from the infant. The pad was then plac ed in a bio-
hazard bag and frozen immediately by the mothers. This
frozen sample was collected in conjunction with the
week 4 faecal sample and placed in a −80 °C freezer
upon arrival at the lab prior to processing.
The P T infant s in the study had a mean gestational
age of 31 weeks and 6 days (SD ± 2 weeks 5 days) and
mean birth weight of 1715 g (SD ± 5 64 g). Twenty six
of the P T infant s were born between 32 and 35 week s,
while the remaining infa nt s were less than 32 weeks
gestation (range 24–32 weeks). There were 10 multiple
births (9 twin and 1 triplet set) and 20 singleton births;
two thirds were male and one third was female. All but
four PT infants were born by CS (emergency 73% and
elective 12%). The average length of stay in the neo-
natal u nit was 39 days (SD ± 3 9.14, range 4–190 days).
All infants under 32 weeks gestation received o ne
course of antibiotics , with a third re ceiving at lea st one
additional course. In comparison, only one third of
Table 1 Breakdown of total number of faecal samples collected
in the study
FT-CS PT-CS FT-SVD PT-SVD Total
Week 1 70 35 83 4 192
Week 4 56 30 63 3 152
Week 8 62 27 74 4 167
Week 24 62 30 74 4 170
Due date N/A 30 N/A 4 34
Total 250 152 294 19 715
Hill et al. Microbiome (2017) 5:4 Page 2 of 18
infants born between 32 a nd 35 week s gestation re-
ceived a course of antibiotics and only 4% re ceived a
second course. See Additional file 1: Table S2 for fur-
ther details o n PT infant s.
Sample extraction and processing
Faecal samples were processed within 24 h of collection
after storage at 4 °C, without freezing. Microbial DNA
was extracted from 0.2-g stool samples using the repeat
bead beating (RBB) method described by Yu and Morri-
son [25], with some modifications. A 0.2-g stool sample
was incubated with 1 ml RBB lysis buffer (500 mM
NaCl, 50 mM tris-HCL, pH 8.0, 50 mM EDTA and 4%
sodium dodecyl sulphate (SDS)) in a 2-ml screw cap
tube with 0.5 g sterile zirconia beads (A single 3.0 mm
bead, 0.1 g of 0.5 mm beads and 0.3 g of 0.1 mm beads).
It was homogenised for 90 s (Mini-Beadbeater™, BioSpec
Products, Bartlesville , OK, USA), with the tubes cooled
on ice for 60 s before another 90 s of homogenisation.
Samples were incubated at 70 °C for 15 min to further
lyse the cells. Samples were centrifuged (16,000g), the
supernatant was removed, and the RBB steps were re-
peated with 0.3 ml of RBB lysis buffer. The supernatants
were pooled and incubated with 350 ml of 7.5 M ammo-
nium acetate (SIGMA). The DNA was pre cipitated by
isopropanol, centrifuged at 16,000g into a nuclear pellet
which was washed with 70% (v/v) ethanol. The pellet
was allowed to dry, then re-suspended in TE buffer, and
treated with RNAse and Proteinase K. It was cleaned
with QIAGEN buffers AW1 and AW2 using a QIAGEN
column and eluted in 200 μl of AE buffer (QIAamp
DNA Stool Mini Kit, QIAGEN, UK). DNA wa s visua-
lised on a 0.8% agarose gel and quantified using the
Nanodrop 1000 (Thermo Scientific, Ireland). DNA was
then stored at −80 °C.
Primers used for PCR amplification were the V4 –V5
region primers 520F (AYTGGGYDTAAAGNG) and
926R (CCGTCAATTYYTTTRAGTTT) (Additional file 1:
Table S3). Initial primers for Illumina sequencing contain
the sequencing primer binding sites, forward or reverse
16S rRNA gene specific primer and a 10nt in-line multi-
plexing identifier (MID). Dual separate MIDs were at-
tached to both ends of the PCR product (Additional file 1:
Table S3).
The V4–V5 amplicons for Illumina seque ncing were
generated using a two-step amplification procedure. The
first step reaction mix contained 50 μl BIO-X-ACT™
Short Mix (BIOLINE), 10 μl of 2 nM forward and re-
verse primers, 50 ng genomic DNA and ddH
2
0 to give a
final volume of 100 μl. Cycling conditions were the fol-
lowing: an initial 95 °C, 5-min denaturation step; 30 cy-
cles of 95 °C for 15 s, 42 °C for 15 s and 72 °C for 30 s;
and a final 10-min extension at 72 °C. The products
were purified using SPRIselect beads (Beckman Coulter,
Indianapolis, IN) as per manufacturer’s instructions,
using a 0.9:1 volum e ratio of beads to product. The puri-
fied PCR products were eluted in 40 μl of ddH
2
O. DNA
quantity was assessed via Quant-iT™ PicoGreen® dsDNA
Assay Kit (Invitrogen™). The samples were pooled in
equimolar amounts (20 ng DNA per sample) and se-
quenced at the University of Exeter (UK) using Illumina
MiSeq 2 × 300 bp paired-end sequencing, on multiple
sequencing runs. Nextflex Rapid library preparation was
carried out by the sequencing laboratory to attach bridge
adaptors necessary for clustering.
LC-MS metabolomic analysis of urine
Urine samples were extracted as previously described
[26]. A 100-μl urine sample was placed on a 96-well
plate with PVDF filter 0.45 μm, together with 200 μlof
internal standard in methanol (see Additional file 2:
Supplementary materials for details). Samples were
then filtered using a positive pressure-96 manifold
(Waters, USA). The eluate was diluted with 200 μlof
MiliQ water containing cinnamic acid standard.
Untargeted LC-MS assays were performed with a hy-
brid linear ion trap Fourier Transform (LTQ F T )
Orbitrap mass spectrometer (Thermo Fisher, Bremen,
Germany), in positive and negative ionisation modes.
The XCMS Online portal (https: //xcmsonline.scr ipp-
s.edu/) was used for data processing (alignment, peak
picking, zero peak re-integrations , features grouping
and assessment of q uality control samples); p leas e see
Additional file 2: Supplementary materials for details.
Data obtained from this p rocessing consisted of a list
of m/z features and its relative intensities , which vary
between sample groups. Such matrix file, with infor-
mation about sample codes , m/z feature and its in-
tensity, wa s used for statistical analysis. In positive
ionisation mode, 2380 statistically significant features
were found. In neg ative ion isation mode, there were
3832 statistica lly significant features. To annotate
compounds, a selection strategy was used based on
the most abundant and the most statistically signifi-
cant fea tures. The procedure for annotation of com-
pounds was adapted from standard metabolomic
initiatives (see Additional file 2: Supplementary mate-
rials for details). Le vels of identification were a s fol-
lows: level I corresponds to compounds identified by
matching masses and retention times with authentic
standards in the laboratory, or by matching w ith
high-resolution LC-MS and LCMSn spectra of stan-
dards reported in the literature; and level II corre-
sponds to compounds identified by matching with
high- and low-resolution LC-MS and LC-MSn spectra
from databases and literat ure. Compounds identified
only by spectral similarities to a similar compound
Hill et al. Microbiome (2017) 5:4 Page 3 of 18
class and lit erature knowledge are reported a s level
III. Unknown compounds are reported a s level I V.
Bioinformatic analysis
The Illumina MiSeq 2 × 300 bp paired-end sequencing
reads were joined using the Fast Length Adjustment of
SHort reads to improve genome assemblies (FLASH)
programme [27]. MIDs were extracted and sequences
were assigned to their corresponding individual samples
by QIIME’s split_libraries_fastq.py, permitting two am-
biguous bases per MID (Ns), and using QIIME’s default
quality settings. The USEARCH sequence analysis tool
[28] was used for further quality filtering. Sequences
were filtered by length, retaining sequences with lengths
of 350–370 bp. This range was used to select the most
abundant sequences for the base of each operational
taxonomic unit (OTU) with reads of all lengths then
aligned to the OTU sequences. Single unique reads were
removed, and the remaining reads were clustered into
OTUs. Chimaeras were removed with UCHIME, using
the GOLD refere nce databa se. The original input se-
quences were mapped o nto the OT Us with 97% simi-
larity. All reads were taxonomically cla ssified by the
classify. seqs command within the mothur su ite of tools
(v1.31.2), using the RDP reference database (training
set 14) [29]. OTUs were cla ssified from these when
>50% of the re ads agreed on a classification at each
phylogenetic level. The returned read numbers varied
greatly from 129 to 815,400 reads (average = 69,410
reads per sa mple). To adjust for the influence of the
number of s equences in a sample on diversity and other
statistical te st s, any sample with less than 10,000 se-
quence reads was eliminated from the study. This re-
sulte d in the loss of e ight samples from th e data set.
Fifteen samples had been sequenced in duplicate, so the
samples with the lower read numbers of duplicated
pairs were removed, as we belie ved that these may not
be the best representations of those samples due to the
lower read count s. The OTU table containing the
remaining 715 samples wa s rarefied to 10,000 reads , to
remove any bias from variation in sample read num-
bers. The remaining samp les w ere from variable mode s
of delivery and time points (data not shown).
Culture-dependent analysis
One gramme of fresh faecal sample per infant was serially
diluted in maximum recovery diluent (Fluka, Sigma
Aldrich, Ireland). Enumeration of bifidobacteria was per-
formed by spread-plating serial dilutions onto de Man,
Rogosa, Sharpe agar (Difco, Becton-Dickinson Ltd.,
Ireland) supplemented with 0.05% L-cysteine hydrochlor-
ide (Sigma Aldrich), 100 μg/ml mupirocin (Oxoid, Fannin,
Ireland ) and 50 units nystatin suspension (Sigma
Aldrich). Agar plates were incubated anaerobically at 37 °C
for 72 h (Anaerocult A gas packs, Merck, Ocon Chemicals,
Ireland). Enumeration of lactobacilli was determined by
plating samples onto Lactobacillus selective agar (Difco)
with 50 units nystatin and incubated anaerobically at 37 °C
for 5 days. Bacterial counts were recorded as colony form-
ing units (CFU) per gram of faeces and were log10 trans-
formed prior to statistical analyses.
Statistical analysis
Statistical analysis was performed using the R statistical
framework, using a number of software packages or li-
braries including, made4, vegan, DESeq2, car, nlme and
lme4. Relative abundance bar charts were generated with
Microsoft Excel. Whe re possible, statistical analyses of
changes over time take the subject numbers into ac-
count, such as the alpha diversity linear modelling, and
DESeq2 tests for differential abundance.
To assess alpha diversity, we calculated the Shannon
Diversity Index with the diversity function from the R
vegan package. After fitting Shannon Diversity to mul-
tiple distributions and performing Shapiro-Wilk normal-
ity tests, we found that it best approximated a normal
distribution, as determined by Quantile-Quantile plots
(qqplots; not shown). Therefore, differences of alpha di-
versity between infants of different modes of delivery at
a given time were detected using mixed effect linear
models (R package nlme), which allow for the adjust-
ment of sequencing run (random effect), while testing
for differences due to mode of delivery (mixed effect). In
order to compare alpha diversity over time, mixed effe ct
linear models were applied (R package lme4, and Ana-
lysis of Deviance using the ANOVA command from the
Car package to test for significance), which allow for
controlling for the subjects and the age of the infants,
along with sequencing run.
Multiple beta diversity metrics were also calculated, in-
cluding weighted and unweighted UniFrac and Spearman
distance ((1 – Spearman Correlation)/2). Principal coordi-
nates analysis was performed on each beta diversity metric
to highlight the separation of infants based on mode of
delivery and sampling time point. Differences between
groups were tested for using permutational multivariate
analysis of variance (PerMANOVA) on beta diversity
matrices, adjusting for sequencing run. False discovery
rate was adjusted for with Benjamini-Hochberg [30].
To identify taxa (phyla and genera) that may be driv-
ing the significant differences detected between time
points and mode of delivery, differential abundance ana-
lysis was determined using DESeq2 on raw phylum- and
genus-level count data. We determined that DESeq2 was
an appropriate tool for differential abundance analysis as
the negative binomial model best fit all genera, deter-
mined by the “goodfit” command from the “vcd” R pack-
age. A heatplot was generated to highlight the major
Hill et al. Microbiome (2017) 5:4 Page 4 of 18
genera driving clustering of samples from different modes
of delivery at different time points and to identify bacterial
co-abundance. We used only genera that were found in at
least 10% of the samples, and utilised Spearman correl-
ation and Ward clustering on log10 of the rarefied genus
count matrix.
We determined significant differences of culture-
dependent count data betw een time points or mode of
delivery using the Wilcoxon rank sum test, and adjusted
for false discovery rate with Benjamini-Hochberg. Corre-
lations between culture-dependent (plate count) and
culture-independent (16S sequencing) data were deter-
mined using Pearson’s product-moment correlation.
Pearson’s product-moment correlation wa s also used to
determine if abundance of any genera correlated with
that of any other genera, and the false discovery rate was
adjusted with Benjamini-Hochberg. To determine if
twins were more closely related to ea ch other than ran-
dom infants, we performed t tests with Monte-Carlo
simulations on beta diversity between samples.
The urine metabolomics dataset was unit scaled before
significant features were identified using the ANOVA
statistical test with term and delivery mode as explana-
tory variables. This analysis gave consistent results when
compared to pareto scaled data and ANCOVA as the
statistical test with 84% of the identified metabolites be-
ing return ed (data not shown). The logged fold change
and the mean value for each variable were calculated
and the results were filtered using the false discovery
rate (FDR) cal culated from the raw p values. To aid the
identification of metabolites, an additional clustering
analysis was performed by WGCNA cluster analysis
using the Spearman correlation and a soft threshold of
nine [31].
Results
Drivers of infant gut microbiota
Gut microbiota is influenced by mode of delivery and
gestational age
The structure of the infant gut microbiota is clearly af-
fected by mod e of delivery (Fig. 1, Additional file 1:
Table S4). The results demon strate that there was a sig-
nificant difference in microbiota composition at genus
level across the four different groups from week 1 to
week 24, when analysed by Spearman distance matrix
and visualised by principal coordinates analysis (PcoA).
At 1 week of age, the microbiota composition of the
FT-CS group was significantly different from that of
both PT-CS and FT-SVD groups (p va lues <0.001).
PT-CS and F T-SVD were also distinct from one an-
other (p < 0.001). The low number of P T-SVD infant s
(n = 4) did not permit significant t esting at this time
point, but it is worth noting that this microbiota clus-
ter is situated between the P T-CS and FT-SVD
groups. At 4 weeks of age, FT-CS microbiota was sig-
nificantly differe nt to all other groups (p < 0.001). The
PT infant microbiota mainly separated across the x-
axis. At 8 weeks of age, the F T-CS group is distinct
from both PT-CS and FT-SVD (p < 0.001), separated
on both a xes. The FT-SVD and P T-CS a re also dis-
tinct (p < 0.001); all three groups have significantly
different microbiota composition a t 8 weeks of age.
By 24 weeks, there were no significant differences between
PT-CS and FT-CS microbiota, while FT-CS and FT-SVD
microbiota were still significantly different (p <0.01). At
this time point, mode of delivery remains influential while
differences due to gestational age have been elimi-
nated. At all time points , there was wide diversity of
individual population structures within each group,
showing the heterogeneous composition of the de vel-
oping infant gut microbiota.
Distinctive metabolomic profiles are associated with
microbiota profiles
Co-inertia analysis of the week 4 microbiota data at the
OTU level and the metabolomic dataset showed that
there was a significant (p < 0.05) amount of co-variation
in the two datasets (Fig. 2). There was little separation
observed between FT birth modes (FT-CS and FT-SVD);
however, the PT-CS samples separated distinctly from
the FT samples. The co-inertia analysis showed that
there were greater differences between the group micro-
biota profiles than between the group metabolomic pro-
files. These differences are evident where the FT
metabolomic baricentres are overlaid while there is a
separation between the microbiota baricentres. The PT-
SVD metabolomic and microbiota baricentres were rela-
tively distant from one an other, but this separation may
be due to the low number of samples in each of these
groups. The compounds associated with the PT-FT split
are from multiple different sources and were represented
by a diverse selection of metabolites (Additional file 3:
Figure S1, Additional file 1: Table S5). Annotated metab-
olites were grouped based on their origin and chemical
character: (i) amino acids and metabolites; (ii) carboxylic
acids and phenolic acids and their metabolites; (iii) vita-
mins and their me tabolites; (i v) drugs and their me -
tabolites; (v) carnitines; (vi) indole metabolites; and
(vii) fatty ac ids a nd their metabolites (see Additional
file 1: Tables S5 and S6). Urea and it s associated me-
tabolite derivatives were situate d in the centre of the
metabolite cluster, suggesting it is abundant in both
groups, providing confidence in our classifications.
We found a number of paracetamol metabolites t o be
significantly higher in P T infant s, as well as several dif-
ferent vitamins and their metab olites s uch as riboflavin,
CECH—a tocopherol metabolite or p yridoxic acid
(Additional f ile 1: Table S5). These metabolites may be
Hill et al. Microbiome (2017) 5:4 Page 5 of 18