dataesr/scientific-paragraphs-categorization huggingface dataset:
3b ). Of the 11 P. copri MAGs in our dataset, these two were the only MAGs assigned to this species of which the abundances were significantly positively correlated with WLZ. Both MAGs are members of a P. copri clade (clade A) that is broadly distributed geographically  19, 20  (Supplementary Fig. 4a ; see Supplementary Fig. 4b  for the predicted carbohydrate-utilization pathways represented in all 51 MAGs assigned to the genus Prevotella that were identified in our 1,000-MAG dataset). Although P. copri MAGs were the greatest source of leading-edge transcripts related to arabinose and aAOS utilization, other MAGs in the microbiome display expression responses consistent with their participation in metabolizing MDCF-2 glycans (or their breakdown products); these include MAGs that are negatively correlated with WLZ.
The representation of these 10 PULs varied among the other nine P. copri MAGs, which span three of the four principal clades of this organism (Fig. 4a  and Supplementary Table  15 ). Notably, the degree of genomic conservation of these PULs is significantly associated with the strength of WLZ association for each of the 11 P. copri MAGs in our MAG dataset across both treatment groups (Pearson r between the Euclidean distance from Bg0019 PUL profile and β 1 (MAG) = -0. 79 (P = 0. 0035); Fig. 4b  and Supplementary Table  3  (WLZ associations)). Five of the seven highly conserved PULs are related to utilization of mannan and galactan-glycans that are significantly more abundant in MDCF-2 than RUSF. Expression of three of these seven PULs, as well as two of the conserved but structurally distinct PULs, is also related to the enrichment of transcripts in carbohydrate-utilization pathways that distinguish upper-from lower-quartile WLZ responders ('WLZ-response quartile' or 'WLZ quartile × study week' terms in Fig.
Collection and handling of biospecimens obtained from participants in the randomized controlled clinical study of the efficacy of MDCF-2 The human study entitled 'Community-based clinical trial with microbiota-directed complementary foods (MDCFs) made of locally available food ingredients for the management of children with primary moderate acute malnutrition (MAM)' was approved by the Ethical Review Committee at the icddr,b (protocol PR-18073; ClinicalTrials. gov: NCT04015999)  4. Informed consent was obtained for all of the participants. The objective of the study was to determine whether twice daily, controlled administration of a locally produced, microbiota-directed complementary food (MDCF-2  3, 4  ) for 3 months to children with MAM provided superior improvements in weight gain, microbiota repair and improvements in the levels of key plasma biomarkers/mediators of healthy growth, compared with a commonly used rice-and lentil-based RUSF formulation.
Preprocessed, short-read shotgun data were aggregated from each participant's faecal sample set (n = 7-8 samples per participant; 118 participants) before MAG assembly. This strategy was adopted to enable the contig abundance calculations required by the MAG assembly algorithms used below, while at the same time mitigating the risk of chimeric assemblies inherent to co-assembly across individuals. Assemblies were generated for all 118 datasets using MegaHit  33  (v. 1. 1. 4), and the resulting contigs were quantified in each assembly by mapping preprocessed reads to the assembled contigs using kallisto  34. Contigs were assembled into MAGs using MaxBin2  35  (v. 2. 2. 5) and MetaBAT2 36 (v. 2. 12. 1). The parallel results of both binning strategies were merged and dereplicated using DAS Tool 37 (v.
9', '-sa 0. 99' in dRep  42  (v. 2. 6. 2)). The final dataset contained 681 ± 99 (mean ± s. d. ) MAGs per participant. MAG assembly summary statistics were collected from CheckM  40  and quast 43 analyses (v. 4. 5) and aggregated (Supplementary Table  1  and Extended Data Fig. 10 ). Initial MAG annotations were performed using prokka  44 (v. 1. 14. 6). Further details about benchmarking the methods we used for assembly are presented in the Supplementary Discussion. To quantify the abundance of each MAG in each sample, MAGs were processed to create a single kallisto quantification index  45. Reads from each faecal DNA sample were then mapped to this index. MAG taxonomy. Taxonomic assignments were initially made using the Genome Taxonomy Database Toolkit  13  (GTDB-Tk) and the corresponding database (release 95).
1. 1. 4). Certain Bifidobacterium species consist of multiple closely related subspecies (such as B. longum). We therefore calculated a pan-genome for 34 Bifidobacterium MAGs in our dataset (Supplementary Table  1 ), plus 14 reference isolate genomes (Supplementary Fig. 5 ), using Roary 50 (v. 3. 12. 0) and a 60% minimum sequence identity threshold for BLASTp  51. The reference isolate genomes included ten Bifidobacterium species and three subspecies of B. longum (ssp. longum, infantis and suis). Concatenated nucleotide sequences of 142 identified core genes were aligned using MAFFT  52  (v. 7. 313). The resulting alignment was trimmed (microseq  53  R package (v. 2. 1. 4)) and then used to construct a maximum-likelihood phylogenetic tree (IQ-TREE  54  (v. 1. 6. 12)). The B. gallicum DSM 20093 genome was selected as an outgroup.
We then returned to the unfiltered count output from kallisto, performed a variance-stabilizing transformation (VST, DESeq2  55  (v. 1. 34. 0)) to control for heteroskedasticity, and filtered the dataset to the same 837 MAGs. We subsequently fit linear mixed-effects models to the transformed abundances of each MAG across all 707 faecal samples (lme4 56 , v. 1. 1-27. 1; lmerTest  57  , v. 3. 1-3). We used analysis of variance to determine the statistical significance of the fixed effects in our model-specifically, the relationship between MAG abundance and WLZ. WLZ-associated MAGs were defined as those with false-discovery-rate-adjusted P (q) < 0. 05.
As an initial pre-processing step, raw reads were aggregated by sample across the two NovaSeq runs, resulting in a total of 5. 0 × 10 7 ± 4. 7 × 10 6 paired-end 150-nucleotide reads per sample (mean ± s. d. ). Adapter sequences and low-quality bases were removed from raw reads (Trim Galore  31  , v. 0. 6. 4), and pairs of trimmed reads were filtered out if either one of the paired reads was less than 100 nucleotides long. Pre-and post-trimmed sequence quality and adapter contamination were assessed using FastQC 59 (v. 0. 11. 7). Filtered reads were pseudoaligned to the set of 1,000 annotated, dereplicated high-quality MAGs to quantify transcripts with kallisto  34. Reads that pseudoaligned to rRNA genes were excluded, leaving an average of 7. 1 × 10 6 ± 3. 9 × 10 6 bacterial mRNA reads (mean ± s.
d. ) per sample. Count tables were further filtered to retain only transcripts that pseudoaligned to the 837 MAGs that passed the abundance and prevalence thresholds described above. To minimize inconsistently quantified counts related to low-abundance MAGs, we assigned a transcript count of zero, on a per-sample basis, to any MAG with a DNA abundance < 0. 5 TPM in that sample. Differential expression analysis (edgeR  60  , v. 3. 32. 1) was conducted using the following steps: (1) transcript filtering for presence/absence and prevalence; (2) library-size normalization using trimmed mean of M values (TMM); (3) estimating per-gene count dispersions; and (4) testing for differentially expressed genes. Transcripts were first filtered using the edgeR default parameters, followed by a parameter sweep of transcript abundance and prevalence threshold combinations.
On the basis of this analysis, transcripts with ≥5 counts per million mapped reads in ≥35% of samples were retained for differential expression analysis. The transcripts that passed this filtering were normalized using a TMM-based scaling factor. We next estimated negative binomial dispersions and fit trended per-gene dispersions (using the power method) to negative binomial generalized linear models. These models were used to characterize (1) the effect of treatment group and study week among all participants and (2) the effect of WLZ quartile and study week among MDCF-2 participants in the upper and lower quartiles of WLZ response using the following model formula: From these models, we identified genes that exhibited significant differential expression using the quasi-likelihood F-test (edgeR, function glmQLFTest), which accounts for the uncertainty in estimating the dispersion for each gene.
PCA PCA was performed on VST-transformed DNA or transcript counts for the 837 MAGs passing the filter described in the 'Defining the relationship between MAG abundances and WLZ' section above. The PCA performed on transcript abundances encompassed 27,518 genes expressed by these MAGs at thresholds for levels and prevalence that are described in the 'Microbial RNA-seq analysis of MAG gene expression' section above. PCA was performed in R using the 'prcomp' function, with each data type centred but not scaled as the dataset was already VST-normalized. The functions 'get_eigenvalues,' 'get_pca_ind' and 'get_pca_var' from the factoextra  62  (v. 1. 0. 7) package were used to extract (1) the variance explained by each PC; (2) the coordinates for each sample along PCs; and (3) the contributions of each variable to PC1-3.
Purified DNA was prepared for long-read sequencing using the SMRTbell Express Template Prep Kit 2. 0 from Pacific Biosciences (Pacific Biosciences; PacBio); we followed the manufacturer's instructions for creating HiFi Libraries from low DNA input, with adjustments made to accommodate a 96-well plate configuration  69. The DNA concentration and fragment-size distribution of resulting libraries were evaluated (genomic DNA ScreenTape; TapeStation instrument) and the libraries were pooled at equimolar concentrations after normalizing for expected genome size. Pooled libraries were sequenced on the Sequel long-read DNA sequencer (PacBio) using the Sequel Binding Kit 3. 0 and Sequencing Primer v4, with 24 h of data collection. The samples were demultiplexed, and Q20 circular consensus sequencing reads were generated using the Cromwell workflow from PacBio.
Constructing a marker-gene-based phylogeny for P. copri MAGs and cultured P. copri isolates We used CheckM  40 (v. 1. 1. 3)  to extract and align the amino acid sequences of 43 single-copy marker genes in each of the 11 P. copri MAGs, each of the 6 cultured P. copri strains, plus the type strain of Bacteroides thetaiotaomicron VPI-5482 (accession number: 226186. 12). Concatenated, aligned marker gene sequences were analysed (fasttree  71  ; v. 2. 1. 10) using the Jones-Taylor-Thornton model and 'CAT' evolution rate to create a phylogenetic tree, which was subsequently rescaled using the 'Gamma20' optimization. The resulting tree was rooted to the B. thetaiotaomicron genome ('ape', v. 5. 6-2  72  ) before extracting the phylogenetic distances between MAGs and isolates. The tree was visualized using ggtree  73  (v.
3. 2. 1). Subsystem-based annotation and prediction of functional capabilities (inferred metabolic phenotypes) of MAGs and cultured P. copri strains MAG and isolate genes were assigned functions, and metabolic pathways were reconstructed using a combination of (1) public-domain tools for sequence alignment and clustering; (2) custom scripts to process the results of sequence alignments (for example, for domain annotation in multifunctional proteins); and (3) a reference collection of 2,856 human gut bacterial genomes for which we have reconstructed and manually curated metabolic pathways related to 98 distinct metabolites and 106 metabolic phenotypes  16. These annotations are captured in the mcSEED database, a microbial community-centred adaptation of the SEED genomic platform  14, 15  , featuring subsystem-based annotation and pathway reconstruction applied to representative human gut bacterial genomes that were initially automatically annotated by RAST or downloaded from the PATRIC (recently renamed Bacterial and Viral Bioinformatic Resource Center, BV-BRC) database  74.
Each mcSEED subsystem includes a set of functional roles (for example, enzymes, transporters, transcriptional regulators) that contribute to the prediction of functional metabolic pathways and pathway variants  75  that are involved in the utilization and catabolism of carbohydrates and amino acids, biosynthesis of vitamins/cofactors and amino acids, and generation of fermentation end-products such as short-chain fatty acids. A complete list of MAG genes comprising these pathways, their abbreviations and functions is provided in Supplementary Table  5  (annotations are based on the January 2021 version of the mcSEED database). Our annotation workflow is shown in Supplementary Fig. 1. In brief, we constructed a reference database containing 995,591 functionally annotated proteins comprising the entire set of curated metabolic subsystems from the 2,856 reference genomes plus an additional 2,988,751 proteins (outgroup is not included in these metabolic subsystems), clustered at 90% amino acid identity ('cluster' command, MMSeqs  76  , v.
1-c7a89). We aligned the predicted protein sequences from the set of 1,000 high-quality MAGs against this reference protein database (DIAMOND  77  , v. 2. 0. 0). To account for any influence of MAG fragmentation on metabolic reconstruction, we also identified gene fragments using prodigal  78  (v. 2. 6. 3) and annotated them in parallel. We implemented the following method to account for instances of multidomain structure that require multiple annotations. For each MAG query protein, we used the top 50 hits based on the bitscore, clustered the start and end position coordinates of the corresponding alignments (DBSCAN function, Scikit-learn  79, 80  , v. 0. 22. 1), used the centre of each clustered start and end position as potential domain boundary coordinates, and split query proteins into domains with database hits attributed to the corresponding domains.
Next, for each domain of ≥35 amino acids, we used Gaussian kernel density modelling (KernelDensity function, neighbours module, Scikit-learn  80  ) of the sequence identity distribution of each set of hits to that domain. A highest local minimum (argrelextrema function, signal module, Scikit-learn  80  ) was used as a threshold to remove low-confidence hits. Finally, functional annotations were applied from the reference database to each query protein or domain by majority rule within each set of high-scoring, domain-specific reference hits. High-identity hits to proteins from the outgroup of the reference database were used as criteria to vote against applying annotation to each query. This procedure yielded a set of 199,334 annotated MAG proteins, representing 1,308 unique protein products across a set of 80 mcSEED subsystems (Supplementary Table  5 ).
Phenotype prediction strategies. We integrated the results of gene-level functional annotation into in silico predictions of the presence or absence (denoted as binary: '1' for presence or '0' for absence) of 106 functional metabolic pathways (Supplementary Table  7 ) using a semi-automated process based on a combination of the following three approaches: Pathway-rules-based phenotype predictions. This approach uses explicit logic-based 'pathway rules' to assign binary phenotypes. These rules combine (1) expert curators' knowledge regarding the gene composition of various metabolic pathway variants contained in the mcSEED database with (2) a decision tree method to identify patterns of gene representation in reference genomes corresponding to an intact functional pathway variant (and a respective binary phenotype value denoted as '1').
A total of 106 functional pathway-specific decision trees was generated (Rpart  81 , v. 4. 1. 15) , where the presence or absence of a particular phenotype was the response variable, and the presence or absence of functional roles (encoded by genes) in each reference pathway were predictor variables. The resulting pathway rules were formally encoded into a custom R script that enabled us to process MAG gene data and assign values (1 or 0) for each of the 106 functional metabolic pathways. Machine-learning-based phenotype predictions. We compared >30 ml methods (Caret 82 , v. 6. 0. 86), using a leave-one-out cross-validation approach in which we removed a single reference genome from the set of 2,856 reference genomes, trained machine-learning models on the remaining genomes, then applied the models to the 'test' genome to predict phenotypes.
This procedure was then repeated for each genome and each metabolic phenotype. The results of this analysis identified random forest as the best-performing method (that is, it produced the greatest number of correctly predicted phenotypes in our reference training dataset). We then built random-forest models for each phenotype based on the reference dataset, optimized model parameters using a grid search and used these models to predict binary (1/0) values for the same set of 106 phenotypes for all MAGs. Neighbour-group-based phenotype predictions. This approach identifies reference bacteria that are closely related to the MAGs in this study and uses these high-quality reference genomes for phenotype predictions that are robust to variation in MAG quality. Examination of groups of closely related reference organisms suggested that close phylogenetic neighbour genomes tend to either possess or lack an entire pathway variant, whereas more distant neighbours (such as other neighbour groups) often carry more divergent pathway variants that specify the same phenotype.
Consensus phenotype predictions. We established a procedure to reconcile inconsistent phenotype predictions between the three strategies described above, based on observing discordant gene patterns and/or discordant predicted phenotypes within a given group of neighbour genomes. In the rare case of an irreconcilable disagreement between the prediction methods, assignment of a consensus phenotype defaulted to that produced by the machine-learning method. We assigned consensus confidence scores to each prediction on the basis of the degree of concordance between the three techniques and our confidence in the accuracy of each (Supplementary Discussion and Supplementary Table  7a ). The complete phenotype prediction process was validated using the 2,856 reference genomes in the mcSEED database, their functionally annotated genes and the accompanying patterns of presence/absence of functional metabolic pathways (curator-inferred binary phenotypes).
Annotation of CAZymes and PULs in P. copri MAGs and cultured P. copri strains CAZymes were annotated according to the CAZy classification scheme  84. Amino acid sequences from MAGs and isolate genomes were analysed using a bioinformatics workflow that performs homology searches against the CAZy database 85 and specifically accommodates the modular structure of CAZymes (which often carry a variable number of ancillary modules in addition to their catalytic domain). Details of this workflow and its application are provided in a previous publication  86. PULs were identified/predicted by combining information from marker genes (SusC/SusD pairs), operon structure and CAZyme annotation information, as previously described  24. The experimentally validated substrate specificities of CAZyme homologues contained in each PUL were used to infer the carbohydrate substrate(s) of each PUL.
Acknowledgements We thank  M. Meier ,  J. Serugo ,  S. Deng ,  K. Ahsan ,  J. Hoisington-Lopez  and  M. Crosby  for technical assistance; members of the  Genome Technology Access Center at Washington University School of Medicine  for Illumina NovaSeq-based generation of faecal DNA shotgun sequencing and microbial RNA-seq datasets; and  E. Martin  and  B. Koebbe  for high-performance computing support. This work was supported by the  Bill & Melinda Gates Foundation , the  National Institutes of Health  ( DK30292 ) and the  Washington University-Centene Personalized Medicine Initiative. J. I. G. is the recipient of an  Agilent Thought Leader Award. Article Extended Data Fig. 3 | LC-MS analysis of glycans present in MDCF-2, RUSF and their component ingredients. (a, b) Analysis of monosaccharides (panel a) and glycosidic linkages (panel b) liberated by hydrolysis of glycans present in MDCF-2 and RUSF, and in the food ingredients used to formulate them.
Data availability Shotgun DNA sequencing and microbial RNA-seq datasets generated from faecal samples, plus annotated P. copri isolate genome sequences are available at the European Nucleotide Archive (PRJEB45356). These data were analysed to generate Figs. 1, 3 and 4, Extended Data Figs. 1, 2, 4, 5, 7, 8 and 10 and Supplementary Figs. 4  and 5. Anthropometry data from a previous publication  4  were used in combination with data from the current study to generate Figs. 1b  and 4b. The mcSEED database (https://doi. org/10. 5281/zenodo. 10041395) that was used to predict the presence or absence of metabolic pathways shown in Fig. 3 , Extended Data Fig. 2  and Supplementary Figs. 1, 4 and 5 was described in a previous publication  16. The CAZy  84, 85  and PULdb 25 databases were used to identify and analyse genes encoding carbohydrate-active enzymes and PULs in Figs.
4  and 5  and Extended Data Figs. 5  and 7-9. LC-MS datasets of monosaccharide, glycoside linkage and polysaccharide data are deposited in GlycoPOST (GPST000244) and were used to generate Figs. 2  and 5  and Extended Data Figs. 3, 6  and 7. All other relevant and/or supporting data are available in the Supplementary Information.
Code availability Code detailing the steps in the MAG assembly workflow and analyses of microbial RNA-seq and glycan datasets are available at GitLab (https://gitlab. com/hibberdm/hibberd_webber_et_al_mdcf_poc_mags) and have been accessioned at Zenodo (https://doi. org/10. 5281/ zenodo. 8000098). Code for annotation of bacterial genes and prediction of metabolic pathway presence/absence is available at GitHub (https://github. com/rodionovdima/PhenotypePredictor) and has been accessioned at Zenodo (https://doi. org/10. 5281/zenodo. 10049439).
Immunofluorescence staining RFe and MMy cells were plated on m-slides (Ibidi, 80286). After 4 days, cells were washed once with DPBS and fixed with Cytofix/ Cytoperm solution (Becton Dickinson, BDB554714) for 20 minutes at 4 C. Cells were rinsed with Perm/Wash buffer (Becton Dickinson, BDB554714) and then incubated overnight at 4 C in Perm/Wash buffer containing primary anti-Afp (R&D Systems, AF1369), anti-Pax6 (BioLegend, 901301), J2 anti-dsRNA (Scicons, 10010200), anti-(gag/pol)HERVK (Austrial Biological, HERM18315), FIPV3-70 anti-Pan Corona (Life Technologies, MA1-82189) or directly conjugated anti-Oct3/4-AF488 (Santa Cruz, sc-5279-AF488), anti-Brachyury (R&D Systems, IC2085G), anti-Otx2 (R&D Systems, AF1979), anti-Zic2 (Abcam, ab150404), anti-Tfe3 (Sigma Aldrich, HPA023881), anti-Tfcp2l1 (R&D Systems, AF5726) or anti-Cdx2 (Life Technologies, PA5-20891) in a 1:50 (anti-Oct3/4) or 1:100 dilution (all others).
11. 9, 99 and visualized using MultiQC v1. 9. 120 The mean phred score was around Q35 across each base position in the RFe and BEF and BiPS samples and no filtering or processing was necessary. In the differentiated samples (EB, endoderm, mesoderm, ectoderm differentiation), the quality of the first nucleotide was less than Q20 in many cases and the reads were processed with Trimmomatic v0. 39. 100  To carry out the differential expression analysis, the genome of Rhinolophus ferrumequinum was used as reference genome (RefSeq accession GCF_004115265. 1 assembled and annotated by the Vertebrate Genomes Project (https:// vertebrategenomesproject. org) or GenBank accession GCA_014108255. 1 assembled and annotated by the Bat1K project (https://bat1k. com) as indicated. The reads were mapped with HISAT2 v2.