Three sheep (one Najdi, one Noaimi, and one Harrei) were obtained at Jeddah's slaughterhouse market. The samples used in this study were taken immediately after animal sacrifice, following best practice veterinary care. The animals used in this study were 12–20 months old and were fed a dried mix of Ammophila arenaria, Medicago sativa, Hordeum vulgare, and Sorghum bicolor for 10–25 days at the slaughterhouse market before sacrifice. A schematic drawing of the sheep digestive system and of the sampling sites is presented in Figure 1. About 30–35 g sample from mid-small and mid-large intestine and from the rectum region were collected into sterile 50 ml Falcon tubes. The intestine samples were mostly semi liquid while the rectal samples were often solid. Immediately after collection the samples were kept refrigerated for 24–28 h then immersed in AquaStool kit solution (MultiTarget Pharmaceuticals LLC, Colorado, USA) according to the manufacturer's instructions.

Bacterial DNA was extracted and purified using the commercial QIAamp DNA Stool Mini Kit (Qiagen: https://www.qiagen.com/) according to the manufacturer's protocol. A quantity of 200 mg of each animal sample was used for DNA extraction. Quantification and assessment of extracted DNA was carried out using a NanoDrop ND-1000 spectrophotometer (NanoDrop, Inc.). The extracted materials were stored at −80°C until use. The concentration of the extracted DNA is reported in Supplementary Table 1.

The samples were submitted to LGC Genomics GmbH (Berlin, Germany) for all sequencing tasks. The sequencing technology chosen was Illumina MiSeq V3 delivering 300 bp paired-end reads. All sequence reads were deposited in the Sequence Read Archive (Leinonen et al., 2011).

The V3–V4 regions of the bacterial 16S rRNA gene were amplified using two universal primers 341F (CCTACGGGNGGCWGCAG) and 785R (GACTACHVGGGTATCTAAKCC) complemented with sample-specific barcodes for multiplexing during sequencing. The primary sequencing output data were demultiplexed using the Illumina bcl2fastq 1.8.4 tool to generate the FASTQ files corresponding to each library. Before clipping the barcode and the adapter sequences, reads with the following criteria were filtered out: (i) reads with more than one mismatch per barcode, (ii) reads with missing barcodes, one-sided barcodes, or conflicting barcode pairs and (iii) reads with a final length <100 bases.

Additional filtering criteria were considered before primer detection and clipping. Thus, only reads with pairs of primers (Fw-Rev or Rev-Fw) and reads with less than three mismatches per primer were kept for further analysis. Furthermore, all reads were turned into forward-reverse primer orientation after removing the primer sequences. The forward and reverse reads were combined using BBMerge v34.48 (http://bbmap.sourceforge.net/). Read quality was controlled using FastQC (http://www.bioinformatics.babraham.ac.uk/projects/fastqc/) for each sample. Reads with an average Phred quality score <33 or containing ambiguous bases were removed.

Operational Taxonomic Units (OTU) clustering has been performed using Mothur v1.35.1 using a 97% identity threshold (cluster.split method). The representative sequences of each OTU were queried against the Ribosomal Database Project release 11.4 (https://rdp.cme.msu.edu/) using BLAST with the following parameters: E-value of 1 or less and percent identity of 90% or higher. Next, the OTU diversity was analyzed using QIIME v1.9.0 (http://qiime.org/) following the standard workflow to determine the membership of all OTUs to different taxonomic levels and groups. The within-sample (alpha diversity and Shannon's diversity index H) as well as between-sample diversity (beta) were calculated using the R package “vegan” (Oksanen et al., 2017). Downstream analysis and result visualization were done with custom R scripts based on the “ggplot2” package (Wickham, 2009). Hierarchical clustering between samples was performed using “hclust” function of R (R Core Team, 2016) with “complete” method applied to a distance matrix derived from the beta diversity using the “betadiver” function from the “vegan” package.

The pre-processing of the shotgun sequencing output data was performed following the same workflow as described above. Assembly of the reads was performed using Trinity (Grabherr et al., 2011) with default parameters. Supplementary Table 2 reports the sequencing and assembly statistics for each sample.

The detection of the assembled DNA sequences that encode carbohydrate-active enzymes was done using FASTY (Pearson et al., 1997) against a custom sequence library made of the catalytic domains of ~150,000 sequences of glycoside hydrolase (GH), polysaccharide lyase (PL), glycosyl transferases (GT), carbohydrate esterases (CE), auxiliary activity (AA), and of their associated carbohydrate-binding modules (CBMs) derived from the carbohydrate-active enzymes database (www.cazy.org; Lombard et al., 2014). We retained as positive hits those that had 30% identity or more and an e-value of 10⁻⁶ or better with a sequence in the library.

A total of 625,937 16S rRNA sequences were obtained from the three intestinal sites sampled (small intestine, large intestine, and rectum) of three sheep. The obtained 16S rRNA sequences represented 7,147 unique OTUs which could be assigned to 24 phyla, 69 order, and 186 genera. Overall, the bacterial communities were mostly dominated by Firmicutes (42.5% of the total sequences and 72% of different OTUs) followed by Bacteroidetes (21.5% sequences and 13.9% OTUs), Proteobacteria (11.5% of total sequences and 2.6% of OTUs), Candidate Division TM7 (7.7% of sequences and 3.4% of OTUs), Actinobacteria (2.5% of total sequences and 2.8% of OTUs), and Planctomycetes (6% of sequences and 1.5% of OTUs). Other phyla were represented in minor proportions ranging from 2% to <0.1% of sequences and OTUs (Table 1).

The taxonomic profiles were globally similar between the three animals along the intestinal tract especially at the phylum level. The core intestinal tract microbiota of the three sheep was composed of the major phyla including Firmicutes, Bacteroidetes, Actinobacteria, Candidate division TM7 and Proteobacteria.

The large intestine and rectal microbiota profiles of the Najdi sheep were largely similar despite minor taxonomic differences (Supplementary Figure 2). In fact, the rectal microbiota of this animal was characterized by higher number of OTUs within Firmicutes especially Clostridiales order with 1,152 different OTUs vs. 845 OTUs in large intestine. Bacteroidales (main order of Bacteroidetes) were similarly represented in both large intestine and rectum while this phylum was more abundant in the large intestine (33.7% of sequences) than in the rectum (17.6% of sequences). The diversity within the three GIT subsites in the three animals is given in Supplementary Figure 4. The small intestine of Najdi sheep had a unique taxonomic profile compared to the large intestine and rectum, with a very low diversity (α-diversity index 39.5) and a particularly high number of sequences that belong to Proteobacteria with 66% of total sequences corresponding to Escherichia and Shigella genera. Lactobacillus genus (Firmicutes) was also highly represented in small intestine of Najdi sheep with 29% of total sequences (Supplementary Table 3).

Samples derived from the large intestine and rectum of the Noaimi sheep had very similar taxonomic profiles at both phylum and order level. They were characterized by the predominance of Firmicutes (mainly Clostridiales), Bacteroidetes (mainly Bacteroidales) and Planctomycetes (mainly Planctomycetales). The number of sequences and OTUs at the Order level are given in Supplementary Table 4.

Here again the small intestine had a distinct taxonomic profile compared to other intestinal subsites. It was characterized by a low diversity of different bacterial groups (only 555 OTUs in total and α-diversity index 92.1) and a particularly low relative abundance of Bacteroidetes members which represent <1.5% of sequences in total. In addition, the small intestine was enriched in other phyla such as Actinobacteria, Tenericutes, Proteobacteria, and Candidate division TM7 (new candidate phylum) compared to the other subsites (Figure 2 and Table 2).

The microbial communities within the intestinal tract of the Harrei sheep followed the same trend as with the two other sheep with the highest microbial diversity residing in the large intestine and rectum, and a very low abundance of Bacteroidetes in the small intestine. We observed the highest abundance of Candidate division TM7 and Cyanobacteria in the small intestine compared to the two other subsites (Figure 2 and Table 2).

The distribution of microbial communities at a given intestinal subsite thus appears to be stable across the three animals studied here. The microbial diversity and the relative abundance were correlated to the body subsite. Overall, the small intestine was characterized by low microbial diversity and abundance compared to the large intestine and rectum (Supplementary Figure 4). At the phylum level, the most striking observation is the very low abundance of Bacteroidetes in the small intestine of the three sheep compared to the other intestinal subsites. The large intestine and rectal bacterial communities were characterized by a high predominance of members of Firmicutes (Clostridiales) and Bacteroidetes (Bacteroidales).

Hierarchical clustering based on beta diversity between different sheep and intestinal subsites confirmed the results described above. Supplementary Figure 3 shows that sheep microbiota diversity in large intestine and rectum were similar to each other (especially within the same animal) while all samples from small intestine showed a particular profile and clustered together.

Finally we note that a large portion of identified OTUs from different samples did not belong to any known or candidate phylum and was flagged as “unclassified.” The proportion of unclassified OTUs varied between ~8% in the small intestine to ~18% in the large intestine of Harrei sheep (Figure 2 and Table 2).

The reads of each sample were assembled into contigs. While the final number of reads obtained for each sample was of the same order of magnitude (between 1.5 and 2.9 million reads), the number of assembled contigs was lower and with a larger average length for the three small intestine samples (Supplementary Table 2), suggesting that the small intestine of sheep may be characterized by a lower organismal concentration than the other two subsites. By contrast there was no major difference in the number of assembled contigs and in the contig length for the large intestine or rectal samples.

Because carbohydrate-active enzymes often have a modular structure where the catalytic domain is appended to a variable number of other domains that may be catalytic or carbohydrate-binding, all contigs were analyzed by pairwise alignment to a library of ~150,000 individual functional domains (modules) derived from the five main categories of modules described in the carbohydrate-active enzymes database, namely glycoside hydrolases (GH), polysaccharide lyases (PL), carbohydrate esterases (CE), glycosyltransferases (GT), auxiliary activity (AA), and CBMs (Lombard et al., 2014). From a total number of ~1.5 million contigs, 59,337 genes encoding carbohydrate-active enzymes were detected (average ~6,600 per body site), representing ~4% of the total of contigs that were assembled. The count of CAZy families for each body site and each animal is given in Supplementary Table 5. No similarity to enzymes of the AA category could be detected, a finding not unexpected given that AA enzymes are oxygen-utilizing oxidoreductases and the three samples subsites are strictly anaerobic. The relative abundance of the families of CAZymes involved in glycosidic bond cleavage [the glycoside hydrolases (GHs) and polysaccharide lyases (PLs)] was computed and the resulting profiles are shown in Figure 3 for the 32 most abundant families. Regardless of the animal, the profile of the GHs and PLs of the small intestine samples is much lower than that in the other sites likely due to the lower concentration of bacteria, a situation similar to that encountered in the human duodenum, where only a very limited number of CAZymes can be found due to the lower concentration of resident bacteria (Angelakis et al., 2015). On the other hand, when the 32 most abundant families were ranked by their average abundance across subsites, the same families tended to be found at the same rank, regardless of the animal (Figure 3) for the large intestine and rectal subsites. By contrast the family profile of the three small intestine samples was markedly different, with a much lower average value and different rank for the most abundant families. This observation mirrors what was observed upon 16S analysis, with the large intestine and rectal subsites being much richer and diverse than the small intestine.

We have then examined the putative function encoded by the most abundant carbohydrate-active enzymes families in the large intestine and rectal samples, viz. the sites with the most abundant profiles (Figure 3). In these sites, regardless of the animal and for both sites, GH13 is consistently the most abundant family, reflecting the universal role of starch/glycogen as a central nutrient source and reserve macromolecule. The other most abundant families include families with miscellaneous substrates (GH2 and GH3), and more specific families like GH43 (α-L-arabinofuranosidases), GH20 (N-acetyl β-glucosaminidases), GH78 (α-L-rhamnosidases), GH97 (starch), GH92 and GH38 (α-mannosidases), GH77 (amylomaltase), GH29 and GH95 (α-L-fucosidases), GH36 (α-galactosidases), CE4 and CE1 (deacetylases), GH33 (sialidases and other ulosonic acid hydrolases), GH25 (lysozymes). Interestingly, the relative abundance of enzymes that target cellulose and xylan is lower than what was found in a large scale metagenomics investigation of the cow rumen (Hess et al., 2011; Table 3).

The handling Editor declared a shared affiliation, though no other collaboration, with several of the authors ED, VL, BH and states that the process nevertheless met the standards of a fair and objective review. The other authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.