Assessment of seasonal variation of diet composition in rodents using DNA barcoding and Real-Time PCR

The study of animal diet and feeding behaviour is a fundamental tool for the illustration of the ecological role of species in the ecosystem. However, size and quality of food intake samples make it hard for researchers to describe the diet composition of many small species. In our study, we exploited genomic tools for the analysis of the diet composition of the Savi’s pine vole (Microtus savii) using DNA barcoding and qPCR techniques for the identification of ingested plant species retrieved from stomach contents. In contrast with previous studies, we found that, despite being a fossorial species, the Savi’s pine vole is a selective feeder that undergoes intense superficial activity in search for food. In addition, our study shows that with a a priori knowledge of the candidate plant species included in animal diet, qPCR is a powerful tool to assess presence/absence, frequency of occurrence and electivity of ingested species. We conclude that this approach offers new opportunities to implement the analysis of food selection in small animals, thereby revealing a detailed picture of plant-animal interactions.


Materials and Methods
Study area. This study was conducted in a 1 ha peach orchard located in an agricultural area in Emilia Romagna, northern Italy (44°21′N, 11°42′E) from November 2014 to September 2015. Average annual rainfall was 750 mm and temperatures varied between +2.6 °C and +23.7 °C. The orchard had trees between 5 and 15 years old planted in rows 4.5 m apart at a distance between 1.5 m and 3 m from each other. The area was cultivated following traditional practices and periodically treated with insecticides, fungicides and herbicides. No rodenticides were used.
Vegetation sampling. Sampling  Food availability was evaluated by sampling vegetation using the quadrat method. We established a sampling grid consisting of 2,500 2 × 2 m quadrats. Each quadrat was then partitioned into 100 20 × 20 cm sub-quadrats. We randomly selected 40 quadrats by simple random sampling with replacement and we sampled 10 out of each of the 100 sub-quadrats by random sampling without replacement. We assessed species composition and richness for each quadrat, which was then rated using percent vegetation cover. Herbaceous plants were collected and placed between two sheets of blotted paper, gently patted to absorb moisture and subsequently wrapped in folded paper. Identification of plant material was conducted by means of Identification of DNA barcode sequence. We searched the BOLD and Genbank sequence database for mat K and rcbL gene sequences for each of the 45 species of plants characterized during our survey that were potentially part of the diet of Savi's pine vole. We identified a 417 bp region of the rcbL gene available on Genbank for 40 candidate species (Supplementary Table 1). The partial rcbL gene sequence was located between position 455 and 872 of the Arabidopsis thaliana reference rcbL complete gene sequence (GenBank accession: U91966.1). Species-specific segregating sites were identified by comparing barcode sequences for 30 plants species, while a genus-specific segregating site was characterized for five groups of species, with each group consisting of two species belonging to the same genus (Supplementary Table 2). The segregating sites were used as target positions to design 35 unique Taqman assays using the Custom Taqman Assay Design Tool for gene expression (Thermo Fischer Scientific).
DnA extraction. DNA was extracted from 30 plant species for which species-specific segregating sites were identified, and from one species for each of the five groups with genus-specific segregating sites. Extractions were conducted using a protocol modified from Doyle & Dickinson 63 by incubating 200 mg of homogenized plant sample in 1 ml lysis buffer containing 200 mM Tris-HCl, 1.4 M NaCl, 20 mM EDTA, 20 mg cetyl trimethylammonium bromide (CTAB), 0.2% 2-mercaptoethanol and 10 mg silica powder for 2 h at 65 °C. Samples were centrifuged for 15 min at 13,000 rpm. One volume chloroform:iso-amyl alcohol (24:1) was added to the supernatant. The mixture was centrifuged for 20 min at 13,000 rpm and DNA precipitated by first incubating the supernatant with 1 volume isopropanol for 30 min at −80 °C. After the second round of centrifugation, the pellet was washed with 500 μl 70% ethanol, centrifuged for 15 min at 13,000 rpm and resuspended in DNase-free water. DNA was accurately quantified with a Qubit dsDNA BR assay kit in a Qubit 4.0 fluorometer and used as a reference for quantification of DNA from stomach contents.
DNA was then extracted from plants ingested by Savi's pine vole via incubation of the entire stomach content (average weight: 443.5 ± 44.8SE mg) in a 2 ml microcentrifuge tube with 1 ml lysis buffer containing 0.1 mM Tris-HCl, 1.4 M NaCl, 20 mM EDTA and 20 mg CTAB for 3 h at 65 °C. DNA isolation was then conducted as described in the CTAB method by Mafra et al. 64 . Analysis of the whole stomach content ensures that all plant species contained in the stomach are sampled for DNA extraction.
Real time PCR assay, conditions and thermal profiles. Because of DNA degradation in stomach contents, the length of barcoding regions that can be successfully amplified by PCR is generally limited to 100-250 bp fragments (see 26 and references therein). In our study, we used real-time qPCR for dietary analysis by designing species-and genus-specific Taqman assays. Each assay included two PCR primers and a target-specific oligonucleotide probe labelled with flourescin (FAM) reporter and non-fluorescent quencher. Forward and reverse primers were designed over a 200 bp region stretching 100 bp in the 3′-5′ direction and 100 bp in the 5′-3′ direction from the target site, respectively. Taqman DNA probes had a conjugated minor groove binding (MGB) moiety attached to the 3′ end. The conjugated MGB folds into a minor groove formed in the DNA when the terminal 5-6 bp of the probe binds to the template. This provides the probe with a higher melting temperature (T m ), close to the T m of the primers, and an increased specificity for single base mismatches at elevated hybridization temperatures, thus strengthening probe binding. PCR products ranged from 59 bp to 102 bp (Table 1).
A total of 97 qPCRs were performed for each candidate plant species (or genus) to quantify and assess plant species presence/absence in Savi's pine vole stomach contents. Amplification reactions included a negative control, a standard dilution series made of three replicates of each of four 10-fold serial dilutions of candidate plant species (or genus) DNA of known concentration, and DNA samples with target-specific Taqman assay. Samples also included an exogenous internal positive control (Thermo Fisher Scientific). The internal positive control (IPC) is a single-stranded, short synthetic DNA template which is added to each amplification reaction along with a pair of specific primers and a Taqman probe labelled with a VIC fluorescent reporter. The IPC was used to distinguish between true negative results and negative results caused by PCR inhibitors, incorrect assay setup, or reagent or thermocycler failure.
Real www.nature.com/scientificreports www.nature.com/scientificreports/ fluorescent dye included in the Master Mix was used as an internal passive reference to normalized PCR fluorescent dye signals. We used a passive reference to correct for possible fluorescent fluctuations including well-to-well volume or light source intensity variations, minor changes in concentration, and non-PCR related fluctuations caused, for instance, by pipetting errors. Thermal-cycling profiles consisted of a denaturation step at 95 °C for 20 s, followed by 40 cycles of 1 s at 95 °C and annealing/extension of 20 s at 60 °C. Filter sets were x1-m1, x1-m2 and x4-m4 for FAM, VIC and ROX, respectively. The number of initial cycles of the PCR during which background fluorescent signal is produced (baseline) and the threshold value whereby enough amplified product has accumulated to yield a detectable fluorescent signal were set automatically by the QuantStudio Real-Time PCR software 1.0.
Standard curve. Amplification of four 10-fold serial dilutions of DNA template of known concentration was used to generate a standard curve by plotting the log-scaled starting quantity of DNA template (N 0 ) against the threshold cycle (C T ) value obtained during amplification of each dilution. The C T values of the samples of unknown concentration were compared to the standard curve to derive the quantity of starting DNA concentration. Three replicates of each dilution point in the standard curve were performed to ensure statistical significance.
Performance Statistical analysis. The mean percentage cover of each plant species and bare ground, along with associated variance, were calculated by averaging the Horvitz-Thompson estimates of percent coverage obtained from the 40 vegetation survey quadrats 66 . As we were only interested in the proportional abundance of plant species, we excluded bare ground data and re-scaled plant cover data to between 0-100. The sampling variances of the scaled estimates were calculated using the delta method (e.g. 67 ). The amount of plant DNA recovered from each stomach content by qPCR was used as a proxy for the proportion of each plant species ingested by an individual vole. As the peach, Prunus persica, is an arboreal species and the Savi's pine vole is known to feed upon roots rather than aerial parts of the plant, we were not able to estimate its availability. Therefore, no analysis on food selection could be performed and the amount of P. persica DNA recovered in the stomach contents was not considered when quantifying the proportion of plant species ingested.
For each trapping session a sign test was performed to assess the selection of plant species by the Savi's pine vole 68 . The test statistic was based on the number of animals with a percent of plant DNA in the stomach higher than the plant percent availability in the study area. Plant availability was estimated using the previously described sampling strategy and could therefore be equal to zero for some species, because either a species was not available in the study area or it was not detected in the sampled sub-quadrats. Presence of DNA in at least one stomach content for a species of plant estimated as not available in the study area highlighted that the zero estimate was due to a sampling error. Then, the availability of that species was considered greater than its proportional use also for those individuals that did not feed on that plant and therefore had a percentage of use equal to zero. For each plant species, p-values of sign tests were derived by means of the binomial probability distribution and subsequently combined in a test statistic to assess the overall null hypothesis of no plant selection by Savi's pine voles. Statistical significance of the overall null hypothesis was determined by permuting sample observations 68 . The hypothesis of no plant selection was rejected for all six trapping sessions. The p-values of the tests performed for each plant species were used to partition the set of available plant species into preferred, avoided and proportionally used www.nature.com/scientificreports www.nature.com/scientificreports/ food items. Significance level of the tests for each plant species was set equal to 0.05. Analyses were performed in R 69 using the "phuassess" package 70 , available from the Comprehensive R Archive Network (CRAN). Table 2 displays the percentages of estimated plant cover in the study area collected during each sampling session. Approximately 50% were perennial species while the other half were annual plants. Couch grass (Elytrigia repens), ribwort plantain (Plantago lanceolata), broadleaf plantain (Plantago major) and common dandelion (Taraxacum officinale) were the dominant species, which, despite seasonal variation of vegetation cover, accounted for the majority of plants available to the Savi's pine vole.

Results
We trapped 20 voles in November, 15 in January, 10 in March, 15 in May, 13 in July and 11 in September. Each stomach sample contained an average of 17.5 ± 1.67SE species of plants (range: [8][9][10][11][12][13][14][15][16][17][18][19][20][21][22][23][24]. Results of the permutation test showed that the proportion of plant species found in the voles' stomachs did not mirror their availability in the study area (P < 0.001). The majority of plant species were found in the voles stomachs in a greater proportion with respect to their percent availability (Tables 3 and 4). In addition, we found seasonal variations in the Savi's pine vole diet, a periodic selection of 6-8 species of plants and avoidance of between 14 and 20 other species. Seasonal selection of plants included rare species such as Amaranthus retroflexus, Avena barbata, Lolium sp. between November and May, and Cardamine hirsuta and Geranium dissectum from July to September. Although plant species selection by the Savi's pine vole changed across the entire sampling period, we found that voles never fed on T. officinale, Bellis perennis, Geranium pusillum, P. lanceolate, P. major, Setaria verticillata and Trifolium pratense regardless of season and relative abundance. In March, no samples of Setaria verticillata were found in any of the sampling plots or SAvi's pine vole stomachs, and it was therefore not considered in the analysis. The average proportion of the peach, P. persica, in the stomach contents was low, from 0.12% in September to 5.52% in January (Table 5). Nevertheless, this species was contained in all stomach samples during every single sampling session.  Table 2. Estimated percent of plant availability for each trapping session.

Discussion
The aim of many dietary studies is not simply to assess food item diversity but to acquire quantitative data on the relative amounts of plant species or preys ingested by an organism 26,27 (and references therein). Our study shows that with a relatively comprehensive a priori knowledge of the candidate plant species an animal can possibly feed upon, the employment of qPCR can provide a good estimate of presence/absence, frequency of occurrence and electivity of each ingested species. Under this assumptions, qPCR can offer either an alternative or a complementary method to HTS, in which even a well designed dietary barcoding study is likely to provide semi-quantitative estimates of the diet of a species or frequencies of sequencing reads as a proxy of the relative abundance of dietary items 36,37,43,71 . Moreover, our results indicate that Taqman assays based on short fragments of the rcbL gene can perform relatively well as DNA barcodes even in significantly degraded samples such as those found in stomach contents 28,72 .
Although plant species availability in our study area was, to some extent, affected by anthropogenic disturbance (e.g. mowing and plowing) which may have altered the natural phenological cycle of plants, we found significant seasonal variability in the diet composition of the Savi's pine vole. Indeed, our results show high levels of selectivity for some species of herbaceous plants, including A. retroflexus, A. barbata, C. arvensis, Portulaca oleracea, Senecio vulgaris and Soncus sp., the latter being almost always selected throughout the year. On the other hand, Savi's pine vole appears to avoid other species such as T. officinale, B. perennis, G. pusillum, P. lanceolate, P. major, S. verticillata, and T. pratense. Although P. persica averaged only 5.5% of the overall food intake, with peaks of up to 20% in a few samples, this species was found in all stomach samples suggesting that the peach was likely consumed throughout the year. The seasonal presence of rare species of plants found in stomach contents, including A. retroflexus, A. barbata, Lolium sp. between November and May, and C. hirsuta and G. dissectum from July to September, suggests that the Savi's pine vole actively selects the plant species to include in its diet. This  Table 3. Percent of Savi's pine voles feeding on a plant species in a greater proportion than its estimated availability for each trapping session.
www.nature.com/scientificreports www.nature.com/scientificreports/ implies intense search activities and specific behavioral and ecological patterns that may have been so far widely overlooked 73 . Based on our results, Savi's pine vole can be indeed regarded as a pest to agroecosystems particularly at high population densities.
Interestingly, we found that despite the number of individuals varies between 2 and 32 per hectare 74 , and therefore regardless of the competition for food resources, the Savi's pine vole feeds upon S. vulgaris, a poisonous plant which contains high concentrations of secondary toxic compounds such as pyrrolizidine alkaloids that were showed to cause liver damages and even lead to the death of a number of other herbivore species 75,76 . We suppose that the Savi's pine vole have evolved specific physiological mechanisms that allow them to metabolise these toxins. However, because detoxification processes of chemical compounds are energy-consuming, further investigations on the factors affecting food selection would greatly contribute to the understanding of the species ecology. Particular attention should be draw on the use of pesticides in agroecosystems and the understanding of their role in mediating diet selection in the Savi's pine vole.    Table 4. Food preference of Savi's pine voles. For each trapping session, plant species are ranked according to the proportion of voles feeding on a palnt in a greater proportion than its estimated availability in the field. Samples size (n) and significance values of the null hypothesis of proportional vegetation use (p) are also reported. Symbols to the right of each plant species represent avoidance (−), proportional use (▯) or preference (+). Mean ± SE 1.90 ± 0.72 5.52 ± 1.14 1.96 ± 1.26 5.50 ± 1.68 1.64 ± 1.04 0.12 ± 0.07 Table 5. Mean percentage of Prunus persica found in the stomach contents of Savi's pine voles. n: sample size.