DOI: 10.5061/dryad.brv15dvcj
title: The generality of cryptic dietary niche differences in diverse large-herbivore assemblages

description: Ecological niche differences are necessary for stable species coexistence
 but are often difficult to discern. Models of dietary niche
 differentiation in large mammalian herbivores invoke the quality,
 quantity, and spatiotemporal distribution of plant tissues and
 growth-forms but are agnostic towards food-plant species identity.
 Empirical support for these models is variable, suggesting that additional
 mechanisms of resource partitioning may be important in sustaining
 large-herbivore diversity in African savannas. We used DNA metabarcoding
 to conduct a taxonomically explicit analysis of large-herbivore diets
 across southeastern Africa, analyzing ~4,000 fecal samples of 30 species
 from 10 sites in 7 countries over 6 years. We detected 893 food-plant taxa
 from 124 families, but just two families—grasses and legumes—accounted for
 the majority of herbivore diets. Nonetheless, herbivore species almost
 invariably partitioned food-plant taxa; diet composition differed
 significantly in 97% of pairwise comparisons between sympatric species,
 and dissimilarity was pronounced even between the strictest grazers (grass
 eaters), strictest browsers (non-grass eaters), and closest relatives at
 each site. Niche differentiation was weakest in an ecosystem recovering
 from catastrophic defaunation, indicating that food-plant partitioning is
 driven by species interactions, and stronger at low rainfall, as expected
 if interspecific competition is a predominant driver. Diets differed more
 between browsers than grazers, which predictably shaped community
 organization: grazer-dominated trophic networks had higher nestedness and
 lower modularity. That dietary differentiation is structured along
 taxonomic lines complements prior work on how herbivores partition plant
 parts and patches and suggests that common mechanisms govern herbivore
 coexistence and community assembly in savannas.
description: This archive presents data on the diets of 30 species of large
 mammalian herbivores in eastern and southern Africa. Data were generated
 via DNA metabarcoding of fecal samples collected in 10 protected areas (7
 countries) across multiple seasons and years. The methods summary below is
 abstracted from Pansu et al. (2022); please see that paper for additional
 details, references, and context.
 <strong><em>Sample collection and
 analysis</em></strong> Fresh fecal samples
 were collected during vehicle-based road surveys. Methods were similar for
 all samples collected in Laikipia, Serengeti (Tanzania), Niassa
 (Mozambique), Gorongosa (Mozambique), Nyika (Malawi), Kafue (Zambia),
 Hwange (Zimbabwe), Hluhluwe-iMfolozi (South Africa) and Kruger (South
 Africa) are described here. Although there were subtle differences in the
 pipeline used for some sets of these samples (described below), all
 samples were collected and processed in a similar way. Samples from Addo
 Elephant National Park (South Africa) were collected as part of an
 independent project; as a result, there were more substantive
 methodological differences in the way those samples were processed
 (described under <strong><em>Addo
 samples</em></strong>, below). Fecal
 samples were collected and processed as described by Kartzinel et al.
 (2015) and Pansu et al. (2019). Fresh samples without any adhering plant
 tissue were collected in unused plastic bags and kept cool until returned
 to camp the same day, where they were pre-processed as follows. We
 homogenized samples by kneading the bag and then transferred ~200 mm3 of
 sample (avoiding plant macroremains) into tubes containing silica beads
 and a stabilization/lysis buffer (Zymo Xpedition Stabilization/Lysis
 Solution, Zymo Research); tubes were vortexed for 30 s to lyse cells and
 then frozen for transport to Princeton University. Before import into the
 United States, samples were subjected to one of two antiviral treatments,
 as mandated by the US Department of Agriculture’s Animal and Plant Health
 Inspection Service (permits 122489, 123156, 130123 to R.M.P.). Samples
 collected in Laikipia from 2013–2016 were treated with proteinase K,
 heated to 95°C for 15 min, and treated with RNase A. Following issuance of
 revised regulations, samples collected from 2016–2018 were subjected to
 heat- only treatment of 72°C for 30 min. On arrival at Princeton
 University, samples were frozen and stored at −80°C and later extracted in
 a facility dedicated to fecal DNA analysis, using Zymo Xpedition
 Soil/Fecal DNA MiniPrep kit according to manufacturer’s instructions. We
 performed one extraction control (sample-free extract) per extraction
 batch (~20 to 30 samples). We amplified a short and
 variable region of the chloroplast genome, the P6 loop of the
 <em>trn</em>L intron, using the universal primers
 Forward (g):
 5’-GGGCAATCCTGAGCCAA-3’ Reverse (h):
 5’-CCATTGAGTCTCTGCACCTATC-3’ Tags composed of 8 base
 pairs (bp), each differing by ≥4 nucleotides, were added to the 5’ end of
 each primer to enable the multiplexing of multiple PCR products per
 library before high-throughput sequencing. PCRs were carried out in a 20
 μL reaction volume including 2 μL template fecal DNA extract; 0.2 μM each
 primer; 0.2 mM each dNTP; 1X GenAmp PCR buffer II; 2.5 mM
 MgCl<sub>2</sub>; 0.5U AmpliTaq Gold DNA Polymerase (Applied
 Biosystems); 4% dimethyl sulfoxide (Sigma- Aldrich); and 0.1 mg ml−1 of
 Bovine Serum Albumin (New England Biolabs). Thermocycling followed a
 program of initial denaturing at 95°C for 10 min, followed by 35 cycles of
 95°C for 30 s, 55°C for 30 s, and 72°C for 30 s, with a 2-min final
 extension at 72°C. Extraction controls, PCR controls (using DNA-free water
 instead of DNA), and positive controls (made of DNA extract of known
 plants) were also amplified and later sequenced.
 Amplification, purification, and sequencing strategies differed
 slightly between samples from Laikipia (2013–2016) and those collected in
 other sites (2016–2018). First, for all non-Laikipia samples, we performed
 multiple PCR replicates (2 or 3 per extract in 2016, 3 per extract in 2017
 and 2018) to monitor reproducibility of results and stochasticity of the
 PCR and sequencing processes. Second, PCR products from Laikipia were
 purified using a SequalPrep Normalization Plate Kit (Applied Biosystems);
 for other sites, samples were pooled per plate and purified with a
 MinElute PCR Purification Kit (Qiagen). Finally, libraries for Laikipia
 samples were prepared using a PCR-based approach and sequenced in
 single-end (170 bp), whereas we used PCR-free library preparation and
 2×150 bp paired-end sequencing for other sites. Samples from Laikipia in
 2013 and 2016, Gorongosa in 2016, and Serengeti in 2017 and 2018 were all
 sequenced on separate sequencing runs. Samples from Laikipia in 2014 and
 2015 were sequenced together, as were the 2017 samples from all sites
 except Serengeti; in these cases, samples from different sites/years were
 placed in different libraries. All libraries were sequenced on Illumina
 HiSeq 2500 platform at Princeton’s Lewis-Sigler Institute for Integrative
 Genomics. Sequence data were curated using OBITools
 v.1.2. For libraries sequenced in paired-end (i.e., all but those from
 Laikipia), paired-end reads were first aligned and assembled using the
 Illuminapairedend command; sequences with a low alignment-quality score
 (&lt;40, the value corresponding to perfect alignment between the last
 10 bases of each read) were discarded. Consensus sequences were then
 assigned to their original sample from the tag information attached to the
 primers, using the <em>ngsfilter</em> command (with default
 parameters allowing zero errors on tags and a maximum of two errors on
 primers). Identical sequences were merged with the
 <em>obiuniq</em> command, which retains information about
 their occurrence in each sample. Low-quality sequences were filtered out
 with <em>obigrep</em>; these included sequences with ambiguous
 nucleotides, those with a size outside the expected length of the barcode
 (&lt;8 and &gt;180 bp) and those represented by only one read in
 the entire dataset. Taxonomic assignment was performed using the
 <em>ecotag</em> command against multiple reference databases:
 a comprehensive local database for Mpala Conservancy in Laikipia, a
 partial local database for Gorongosa, a grass-specific local database for
 Serengeti, and a global reference database generated by in silico PCR from
 the EMBL database (release 134) using ecoPCR. We used the
 <em>obiclean</em> command (with parameters d = 1 and r = 0.25)
 to detect sequences potentially resulting from PCR and/or sequencing
 errors. For each PCR product, this program determines if a sequence is
 more likely to be a true sequence (‘head’), a sequence derived from
 another one (‘internal’), or a sequence from which no other sequence is
 derived and is itself not derived from another (‘singleton’). This
 information was used later in the filtering process to remove probably
 erroneous sequences. The fasta file was then converted into a
 sequence-by-sample matrix using the <em>obitab</em> command.
 Additional filtering steps were conducted in R v.3.5.3.
 For sites with local reference databases, sequences were
 preferentially assigned to the local reference database. However, if the
 assignment score obtained with the local database was &lt;98% and
 lower than that obtained with the global database, then the sequence was
 reassigned to the global database. For sites without local databases,
 sequences were assigned to the global reference database unless the
 assignment score was higher with one of the local reference
 databases. To further curate the dataset, we first
 discarded PCR products with low numbers of reads. For this, we compared
 the density distribution of the log-transformed number of reads in
 controls and in true samples within each library, using the intersection
 of the two distributions as a threshold. We then removed sequences that
 were likely to have resulted from PCR or sequencing errors. For this, we
 first used the outputs of the <em>obiclean</em> analysis as
 follows: for each site and sampling bout, we discarded all sequences that
 were more frequently considered to be errors (internal) than true
 sequences (head or singleton) for that site in that year. For sites with
 local reference databases (Laikipia, Gorongosa, Serengeti), these
 sequences also needed not to match perfectly with any sequence from the
 local reference database, or else they were retained. We also filtered out
 putative contaminants by discarding any sequence that had its maximal
 average relative read abundance (RRA) in negative controls. Similarly,
 sequences that displayed low similarity (&lt;80% identity) with their
 closest match were considered likely to be chimaeras and/or highly
 degraded sequences and therefore filtered out. Next, we
 removed outlier PCR replicates (those with non-reproducible results). For
 each library, we iteratively determined the density distributions of
 within- and between-sample distances and discarded replicates that fell
 within the distribution of between-sample distances, the threshold being
 defined as the intersection of the two density distributions. This process
 was iterated until no further replicate was removed. Last, to reduce the
 impact of low-abundance false positives that can arise from tag- jumps
 during Illumina sequencing, we removed sequences representing &lt;1%
 of reads in each sample. Remaining sequences were considered molecular
 operational taxonomic units (mOTUs) in subsequent analyses.
 <strong><em>Addo
 Samples</em></strong>  The main
 methodological difference between samples collected in Addo and those from
 all other sites was the DNA extraction method: whereas we extracted total
 DNA at other sites, an extracellular DNA-extraction protocol was used for
 samples from Addo (Taberlet et al. 2012). Fresh samples from all
 ruminants, except buffalo, were collected in tubes containing silica gel
 and stored dry until DNA extraction in the field. Samples from hindgut
 fermenters and buffalo were collected in unused plastic bags and kept cool
 until DNA extraction in the field on the same day. Samples were extracted
 in the field from a much larger volume of initial material than used at
 other sites. For this, fecal material was mixed with an equivalent amount
 of saturated phosphate buffer (Na2HPO4; 0.12 M; pH = 8) for 15 minutes to
 desorb extracellular DNA. DNA extraction followed methods described in
 Kerley et al. (2018). Although the general
 metabarcoding approach based on the P6 loop of the trnL intron and
 Illumina sequencing was similar for all samples, laboratory steps and
 data-filtering protocols performed at Université Grenoble Alpes differed
 slightly from those used at Princeton University. Specifically, the primer
 pair used was identical to other sites, but the composition of the PCR mix
 and PCR conditions differed. PCRs were performed in a 20 μL reaction
 volume containing 10 μL of AmpliTaq Gold 360 master mix (Applied
 Biosystems), 0.5 μM of each primer, 0.16 μL (20 mg/mL) of bovine serum
 albumin (BSA, Roche Diagnostic), and 2 μL template fecal DNA extract
 (diluted 10 times). Polymerase activation was performed at 95°C for 10
 min, followed by 40 cycles of 95°C for 30 s (denaturation), 50oC for 30 s
 (primer annealing), and 72°C for 60 s (extension), with a final elongation
 for seven minutes at 72oC. Three technical PCR replicates were performed.
 All experiments included extraction controls, blanks, and negative and
 positive PCR controls. All PCR products (samples and controls alike) were
 mixed together and purified using the MinElute PCR Purification Kit.
 Libraries were prepared using the MetaFast protocol
 (https://www.fasteris.com/dna/?q=content/metafast-
 protocol-amplicon-metagenomic-analysis) and sequenced in paired-end on a
 HiSeq 2500 platform (2×150 bp) by Fasteris (Geneva,
 Switzerland). Sequence data were processed with
 OBITools software. Paired-end read alignment, assignment of sequences to
 their original samples, sequence dereplication, removal of sequences with
 ambiguous nucleotides, and selection of sequences based on their length
 (range 10–220 bp here) were conducted as described above for the non-Addo
 sites. In addition, sequences obtained in fewer than two different PCR
 replicates (with a minimum of 10 reads in at least one of them) were
 discarded, as were those represented by &lt;100 reads over the entire
 dataset. The ecotag program was used for taxonomic assignment using a
 comprehensive local reference database for Addo Elephant National Park,
 comprising 473 plant species. Methods employed to build the reference
 database followed methods described in Taberlet et al. (2018, p. 25),
 using PCR with the same P6-loop primers. Only sequences with an assignment
 score ≥0.99 identity were retained. Sequences with &gt;10% of their
 reads observed in controls were discarded, as were samples with a
 sequencing depth &lt;2,000 reads after filtering.
 <strong><em>Confirmation of sample
 identity</em></strong> Whenever possible,
 fecal samples were collected directly after observing animals defecate. In
 some cases, we were forced to collect samples without a direct observation
 (e.g., for animals that are most active at night and/or dangerous, such as
 elephant, rhinoceroses, buffalo, and hippopotamus). We are highly
 confident of the species assignments used in this study for multiple
 reasons. First, the collection team always included at least one member
 with considerable experience in identifying herbivores and their feces.
 Second, fresh dung is readily distinguishable from even hours-old dung.
 Third, relatively few sympatric species have dung similar enough to be
 confused. Fourth, we performed confirmatory analyses whenever we
 considered it possible that identifications might be mistaken, as detailed
 below. We used DNA analyses to confirm the species
 identity of 200 uncertain samples from 2017-2018. These DNA samples were
 amplified with tagged primers targeting a 16S metabarcode allowing
 identification of mammals to species (MamP007F:
 5’-CGAGAAGACCCTATGGAGCT-3’; MamP007R: 5’-CCGAGGTCRCCCCAACC-3’). To limit
 amplification of human DNA, a blocking primer (MamP007_B_Hum1,
 5’-GGAGCTTTAATTTATTAATGCAAACAGTACCC3-3’) was added to the mix. PCR
 amplifications were conducted in a final volume of 20 μL containing 2 μL
 of template DNA, 0.2 μM of each primer, 2 μM of human blocking primer, 0.2
 mM of each dNTP, 0.2 mg.mL-1 of BSA, 1X GenAmp Gold II buffer; 2 mM MgCl2,
 and 1U of AmpliTaq Gold DNA Polymerase. Thermocycling conditions included
 an initial denaturation at 95°C for 10 min, followed by 35 cycles of
 denaturing at 95°C for 30s, annealing at 50°C for 30s, extension at 72°C
 for 30s, and a final elongation step (7 min at 72°C). PCR products were
 pooled and purified using the MinElute PCR Purification Kit. The library
 was constructed using a PCR-free protocol and sequenced in paired-end
 (2×250 bp) on an Illumina MiSeq platform. Data processing was performed
 with OBITools. Paired-end alignment, assignment of sequences to their
 original sample, dereplication of identical sequences, and removal of
 low-quality sequences were conducted as described above. We generated a
 mammal DNA reference database by using in silico PCR from the EMBL
 database (release 134) to assign sequences to herbivore taxa. PCRs with a
 low number of reads (&lt;350) were discarded and non-mammalian reads
 were filtered out. The following criteria were used to define an
 acceptable assignment: (i) &gt;98% similarity (i.e., maximum of 1 bp
 difference) with a reference sequence, (ii) no multiple assignments, and
 (iii) the top sequence was &gt;50% of reads and was at least twice as
 abundant as the second. In addition, when the identity score was
 &lt;100%, sequences were manually inspected to check if the barcodes
 of the putative species were dissimilar enough to avoid
 misidentification. A comparable approach was applied to
 all samples from Addo (2013-2014) using the same primer pair as above. In
 Laikipia (2013–2016), a barcoding approach based on primer pairs targeting
 different regions of mitochondrial COI (cytochrome c oxidase subunit I)
 was used to confirm the identity of 217 of the samples used in the present
 study. In total, we performed confirmatory analyses for 1634 samples of 23
 species, encompassing 7 of the 10 sampling sites.
 <em><strong>References</strong></em> T. R. Kartzinel, et al., DNA metabarcoding illuminates dietary niche partitioning by African large herbivores. <em>Proc. Natl. Acad. Sci. U. S. A.</em> 112, 8019–8024 (2015). G. I. H. Kerley, et al., Diet shifts by adult flightless dung beetles <em>Circellium bacchus</em>, revealed using DNA metabarcoding, reflect complex life histories. <em>Oecologia</em> 188, 107–115 (2018). J. Pansu, et al., Trophic ecology of large herbivores in a reassembling African ecosystem. <em>J. Ecol.</em> 107, 1355–1376 (2019). J. Pansu, et al., Generality of cryptic dietary niche differentiation in diverse large-herbivore assemblages. <em>Proc. Natl. Acad. Sci. U. S. A.</em> doi: 10.1073/pnas.2204400119 (2022). P. Taberlet, et al., Soil sampling and isolation of extracellular DNA from large amount of starting material suitable for metabarcoding studies. <em>Mol. Ecol</em>. 21, 1816–1820 (2012). P. Taberlet, A. Bonin, L. Zinger, E. Coissac, <em>Environmental DNA for Biodiversity Research and Monitoring</em> (Oxford University Press, 2018).