Not annotated metagenome-assembled genomes recovered from rumen samples from goats

Protozoa comprise a major fraction of the microbial biomass in the rumen microbiome, of which the genera <em>Entodinium</em> has been consistently observed to be dominant across a diverse genetic and geographical range of ruminant hosts. Despite the apparent core role that species such as <em>Entodinium caudatum</em> exert, their greater biological and metabolic contributions to rumen function remain largely undescribed. Here, we have leveraged (meta)genome-centric metaproteome datasets from rumen fluid samples originating from both cows and goats fed contrasting diets, to detail the specific metabolic niches that <em>E. caudatum</em> occupies in the context of its bacterial and archaeal co-habitants. Initial proteome estimations via total protein counts and label free quantification highlight that <em>E. caudatum</em> populations comprise an extensive fraction of the total rumen metaproteome. Our analysis also suggested increased microbial predation and volatile fatty acid metabolism by <em>E. caudatum</em> to occur in high methane emitting animals, although with no apparent direct metabolic link to methanogenesis. Despite <em>E. caudatum</em> having a well-established reputation of digesting starch, it was unexpectedly less active in low methane emitting animals fed high starch diets, which were instead dominated by propionate/succinate-producing bacterial populations suspected of being resistant to predation. Collectively, our results illuminate the substantial metabolic influence under-explored eukaryotic populations have in the rumen, with greater implications towards both digestion and methane metabolism. The experimental procedures were approved by the Auvergne-Rhône-Alpes Ethics Committee for Experiments on Animals (France; DGRI agreement APAFIS#3277–2015121411432527 v5) and complied with the European Union Directive 2010/63/EU guidelines. Experiments were performed at the animal experimental facilities of HerbiPôle site de Theix at the Institut National de la Recherche pour l’Agriculture, l’Alimentation l’Environnement (INRAE, Saint-Genès-Champanelle, France) from February to July 2016. Experimental design, animals and diets were previously described by Fougère <em>et al.</em> and Martin <em>et al</em>. Briefly, 4 Holstein cows and 4 Alpine goats, all lactating, were enrolled in respectively two 4 x 4 Latin square design trials. The experimental diets were CTL (first-cut grassland hay and concentrate) and COS (CTL diet with an additional supply of corn oil and cracked-wheat starch -5.0 % of total DMI (Table 1). Corn oil (Olvea, Saint Léonard, France) was added to the concentrate, at 5% of total DMI and contained (g/kg of total FA): 16:0 (114), 18:0 (16.4), cis-9 18:1 (297), cis-11 18:1 (6.30), 18:2n-6 (535), 18:3n-3 (7.57), 20:0 (3.48), 22:0 (1.0), 24:0 (1.5), and total FA (1000 g/kg). Detailed diet composition is available in Martin <em>et al</em>. Each experimental period lasted for 28 days. Rumen fluid was collected through stomach-tubing before the morning feeding on day 27 of each experimental period. The stomach tube consisted of a flexible 23 mm external diameter PVC hose fitted to a 10 cm-strainer at the head of the probe for cows, and a flexible 15 mm PVC hose with a 12 cm-strainer for goats. Samples were filtered through a polyester monofilament fabric (280 μm pore size), dispatched in 2-ml screw-cap tubes, centrifuged at 15000g for 10 mins and the pellet snap-frozen in liquid nitrogen. Samples were stored at -80°C until DNA extraction using the Yu and Morrison bead-beating procedure. In total, 32 rumen fluid samples (8 animals fed four diets) were sent to the Norwegian University of Life Sciences (NMBU) for metagenomic and metaproteomic analysis. Respiration chambers were used to measurements methane emissions over a 5-day period, while VFA and NH3 concentrations were determined by gas chromatography using flame ionization detector. Metagenomic shotgun sequencing was performed at the Norwegian Sequencing Centre on 2 lanes of the Illumina HiSeq 3/4000 generating 150 bp paired end reads in both lanes. Samples were prepared using TruSeq PCR-free library prep prior to sequencing. All 32 samples were run on both lanes to prevent potential lane-to-lane sequencing bias. FASTQ files were quality filtered and Illumina adapters removed using Trimmomatic 28 (v. 0.36) with parameters -phred33 for base quality encoding, leading and trailing base threshold set to 20, sequences with average quality score below 15 in a 4-base sliding window were trimmed and the minimum length of reads was set to 36 bp. MEGAHIT 29 (v.1.2.9) was used to co-assemble reads originating from samples collected from cow and goats separately, with options --kmin-1pass, --k-list 27,37,47,57,67,77,87, --min-contig-len 1000 in accordance with 30. Bowtie231 (v. 2.3.4.1) was used to map reads back to the assemblies and SAMtools32 (v. 1.3.1) was used to convert SAM files to BAM format and index sorted BAM files. The two co-assemblies (one co-assembly of reads from samples originating from cow and one co-assembly of reads originating from goat) were binned using Maxbin2, MetaBAT2 and CONCOCT. MetaBAT2 (v. 2.12.1) was run using parameters --minContig 2000 and --numThreads 4, Maxbin2 (v. 2.2.7) ran with default parameters and -thread 4, min_contig_length 2000, and CONCOCT (v. 1.1.0) with default parameters and –length_threshold 2000. Further, resulting bins from the three mentioned binners were filtered, dereplicated and aggregated using DASTool(v. 1.1.2 ) with the parameters –write_bins 1, --threads 2 and BLAST37 as search engine. This resulted in a total of 244 dereplicated metagenome-assembled genomes (MAGs) across the two host species (104 originating from cow and 140 from goat). CheckM40(v. 1.1.3) lineage workflow was used to determine completeness and contamination for each MAG, with parameters --threads 8, --extension fa, and CoverM (v. 0.5.0) (https://github.com/wwood/CoverM) was used to calculate relative abundance of each MAG, while GTDB-tk (v. 1.3.0) was used for taxonomic annotation. 90% (219 of 244) of the recovered MAGs were considered high or medium quality MAGs according to MIMAGs threshold for completeness and contamination for genome reporting standards44. Gene calling and functional annotation of the final MAGs were performed using DRAM45 with the databases dbCAN, Pfam, Uniref, Merops, VOGdb and KOfam.