This study was conducted in the TNNR, the area where South China sika deer is distributed. The TNNR is located on the south bank of the middle and lower reaches of the Yangtze River, Pengze, Jiangxi Province, China. The total area of TNNR is 12,500 hm², the core area is 2,670 hm², the experimental area is 1,830 hm², and the buffer zone is 8,000 hm². The TNNR mainly consists of low mountains and hills (Wang et al. 2021). The TNNR lies in a climatic zone transitional from tropical to middle subtropical and has transitional climate characteristics that are warm, with a humid monsoon, and four distinct seasons. The frost-free period is up to 247 days with little snow cover (Zhou 2019). The vegetation type is mainly composed of mixed evergreen–deciduous broad-leaved forest, coniferous forest, mixed coniferous–broad-leaved forest, broad-leaved forest, and bamboo (Zhou 2019). From December 2020 to February 2021, 90 fecal samples were collected from sika deer, Reeve’s muntjac, and Chinese hare in the TNNR and stored at −80 °C. Sampling sites mainly focused on Nursery bases (MP), Fir forests (SS), NieJiashan (NJS), XianLingAn (XLA), WuGuiShi (WGS), and Bamboo Garden (ZY) (Fig. 1).

Figure 1. 

Fecal samples sites of three sympatric species at Taohongling Sika Deer Nature Reserve (MP: Nursery bases; SS: Fir forests; NJS: NieJiashan; XLA: XianLingAn; WGS: WuGuiShi; ZY: Bamboo Garden).

In our study, to minimize possible bias caused by variation in individual digestibility, five fecal pellets were randomly taken from each fecal sample and mixed to form a single composite sample. Total DNA was extracted using the DNA extraction kit (TIANGEN, Beijing) following the liquid nitrogen grinding method. The final DNA concentration and purification were determined by NanoDrop 2000 UV-vis spectrophotometer (Thermo Scientific, Wilmington, USA), and DNA quality was checked by 1% agarose gel electrophoresis. The P6 loop region of the trnL(UAA)intron region was amplified with universal primers g (5՚-GGGCAATC CTGAGCCAA-3՚) and h (5՚-CCATTGAGTCTCTGCACCTATC-3՚) by thermocycler PCR system (Gene Amp 9700, ABI, USA). PCR amplifications were carried out in a total volume of 25 μl containing 12.5 μl PCR mix (Tiangen, Beijing, China), 1 μl DNA, 1 μl of each primer, and 9.5 μl H₂O. The reaction conditions were as follows: denaturation at 95 °C for 3 min followed by 35 cycles at 95 °C for 30 sec, 56 °C for 30 sec, and 72 °C for 45 sec, with a final 10 min at 72 °C and storage at 4 °C for 10 h. The PCR products were detected by Agarose gel electrophoresis and sequenced by Shanghai Personal Biotechnology Co., Ltd.

Purified amplicons were pooled in equimolar and paired end sequenced (2 × 300) on an Illumina MiSeq platform (Illumina, San Diego, USA) according to the standard protocols. The analysis was conducted by following the tutorial of QIIME2 docs along with customized program scripts (https://docs.qiime2.org/2019.1/). Briefly, raw FASTQ files were demultiplexed using the QIIME2 v. 2019.4 demux plugin based on their unique barcodes (Caporaso et al. 2010). Demultiplexed sequences from each sample were quality filtered and trimmed, denoised, and merged, and then the chimeric sequences were identified and removed using the QIIME2 dada2 plugin to obtain the feature table of operational taxonomic units (OTUs). The QIIME2 feature-classifier plugin was then used to align OTUs sequences to the National Center for Biotechnology Information (NCBI) database to generate the taxonomy table. Diversity metrics were calculated using the core-diversity plugin within QIIME2. Feature level alpha diversity indices, such as observed OTUs, Chao1 richness estimator, and Goods coverage index were calculated to estimate the diet diversity within an individual sample. A difference significance test was performed by R v. 4.1.3. Beta diversity distance measurements using Bray–Curtis were performed to investigate the structural variation of fecal plant communities across samples and then visualized via principal coordinate analysis (PCoA) and nonmetric multidimensional scaling (NMDS) (Vázquez-Baeza et al. 2013). Venn diagram analysis was performed to explore the common and special OTUs types among the three herbivores.

Formulas of forage plants diversity and niche analysis were conducted as follows:

The Shannon–Wiener diversity index (H′) was calculated to explore diet diversity (Shannon 1948), according to the following formula:

$H^{\text{'}}=-\sum _{i=1}^n{P}_i\ln P_i$ (1)

Where P_i is the proportion of food item i out of all foods, and n is the total number of food items.

Pielou evenness index (Pielou 1969), according to the following formula:

J′ = H′ / H_max (2)

H _max = ln n (3)

Where n is the number of plant species in the fecal sample, and the number of plant species is represented by the number of plant OTUs types.

The Levin index (Smith 1982) was applied to standardize the trophic niche measure, and the formula was as follows:

$\mathrm{B}=1/\sum _{i=1}^S{\mathrm{p}}_i^2$ (4)

The niche overlap index was obtained using the Pianka index (Pianka 1973; Hou et al. 2021), and the formula was as follows:

$Q_{jk}=\frac{{\sum }_{i=1}^S{P}_{ij}{P}_{ik}}{\sqrt{{\sum }_{i=1}^S{P}_{ij}^2{\sum }_{i=1}^S{P}_{ik}^2}}$ (5)

Where Q_jk is Pianka’s niche overlap index between species j and species k; P_ij is the proportion of resource i out of all resources used by species j, and P_ik is the proportion of resource i out of all resources used by species k. The values range from 0 (no food item in common) to 1 (complete overlap in resource use).

After processing the raw reads, a total of 11,411,958 counts were obtained from 90 fecal samples. The mean OTUs length was 67.96 bp with a range from 32 bp to 189 bp. Venn diagram showed that the OTUs in the overlap were commonly shared, and those in the nonoverlapping parts were special OTUs. In total, 764 OTUs, 833 OTUs, and 843 OTUs were obtained from sika deer, Reeve’s muntjac, and Chinese hare samples, respectively. The number of OTUs among the three herbivore groups was 373, and the specific OTUs in sika deer, Reeve’s muntjac, and Chinese hare were 391, 449, and 470, respectively (Fig. 2).

Figure 2. 

Venn analysis of OTUs in three herbivores of Taohongling nature reserve.

Based on OTUs sequences alignment in the NCBI database, the diets of sika deer, Reeve’s muntjac, and Chinese hare includes 203 genera in 90 families, 203 genera in 95 families, and 163 genera in 75 families, respectively (Editorial Committee of FRPS 2004; Jiang 2009; see Suppl. material 1). The three species consumed common and specific forage plants but varied greatly in their use of the available forages. On the whole, the 10 most abundant families in the diet of sika deer include Rosaceae (36.94%), Hamamelidaceae (25.75%), Pentaphylacaceae (13.41%), Theaceae (3.89%), Celastraceae (3.01%), Poaceae (2.98%), Ericaceae (2.43%), Moraceae (2.09%), Cupressaceae (1.68%), and Cannabaceae (1.16%). For the Reeve’s muntjac, the 10 most abundant families include Rosaceae (52.41%), Pentagliaceae (11.01%), Celastraceae (6.12%), Poaceae (4.52%), Cannabaceae (3.27%), Moraceae (2.64%), Sabiaceae (2.14%), Asteraceae (2.10%), Oleaceae (1.50%), and Smilacaceae (1.26%), and the 10 most abundant families in the Chinese hare consist of Rosaceae (16.35%), Poaceae (16.33%), Smilacaceae (15.58%), Anacardiaceae (10.71%), Fabaceae (9.74%), Asteraceae (6.77%), Rubiaceae (5.56%), and Cupressaceae (4.49%).

The dominant genera foraged by sika deer were Rubus (36.49%) and Loropetalum (25.52%), followed by Eurya (13.41%), Camellia (3.89%), Euonymus (2.94%), Phyllostachys (2.46%), Maclura (1.88%), Cunninghamia (1.67%), Rhododendron (1.45%), and Celtis (0.99%), and others (9.30%). For Reeve’s muntjac the diet was strongly dominated by Rubus (51.77%), other genera high abundance were Eurya (11.01%), Euonymus (6.11%), Celtis (3.21%), Arrhenatherum (3.02%), Sabia (2.14%), Maclura (1.67%), Ligustrum (1.49%), Phyllostachys (1.31%), and Smilax (1.26%). Chinese hare consumed almost equal proportions of Rubus (15.81%) and Smilax (15.58%), followed by Rhus (10.64%), Campylotropis (9.52%), Bidens (5.37%), Hedyotis (5.02%), Cunninghamia (4.47%), Eleusine (3.57%), Digitaria (3.09%) and Miscanthus (2.55%). To sum up, these three herbivores all mostly feed on Rubus in winter. The species composition heatmap was drawn from the species and sample levels. The 20 genera with the highest abundance were selected based on species annotation information for 90 samples of the three species. The clustering results show the differences in the relative abundances of sika deer, Chinese hare, and Reeve’s muntjac (Fig. 3).

Figure 3. 

Species composition heat map at the genus level. (Ch: Chinese hare; Sd: Sika deer; Rm: Reeve’s muntjac).

High-throughput sequencing can be used to detect the diet on species levels in most samples combined with a background survey of TNNR, except for OTUs that were undetectable. Sika deer food items included Rubus chingii, Loropetalum chinense, Eurya japonica, Camellia japonica, Euonymus grandiflorus, etc. The forage plants of Reeve’s muntjac consisted of Rubus chingii, Eurya japonica, Euonymus grandiflorus, Arrhenatherum elatius, Celtis sinensis, etc. The diet of Chinese hares mainly focused on Rubus chingii, Smilax china, Rhus chinensis, Campylotropis sp., Hedyotis diffusa, etc. The detailed top 30 forage species in three species are shown in Table 1.

Alpha diversity reflects the abundance and diversity of species communities. The Chao1 and Observed species indices showed the highest community richness was Reeve’s muntjac (Chao1 index; Reeve’s muntjac = 242.46, Sika deer = 236.52, Chinese hare = 192.03, on average). The Shannon and Simpson indices showed the highest community diversity was Chinese hare (Shannon index; Chinese hare = 2.36, Reeve’s muntjac = 2.21, Sika deer = 1.82, on average), with no significant differences (P > 0.05). The goods coverage of 0.998 indicated that an average of 99% of the species were annotated (Fig. 4a; Table 2). Rarefaction curves describe the increase in species diversity as the sample size increases. It is crucial to point out that the characterization of species diversity was considered very reliable since the depth of rarefaction applied (35,000) was found to be sufficiently satisfactory (e.g., rarefaction curves had already reached a plateau at ~35,000 sequences in all samples) (Fig. 4b). The rank abundance curve reflects the richness and evenness of species in the sample through the flatness. The evenness of community composition of sika deer and Reeve’s muntjac was higher, while the lowest of Chinese hare (Fig. 4c).

Figure 4. 

a box-plot of the alpha diversity index. In each panel, the abscissa is the group, and the ordinate is the value of the corresponding alpha diversity index b sample rare faction curves c rank abundance curve. The abscissa is the sequence number of OTUs arranged according to the Abundance size. The ordinate is the abundance value of each OTU in this grouping by Log2 log transformation (Ch: Chinese hare; Sd: Sika deer; Rm: Reeve’s muntjac).

We assessed the beta diversity using the Bray–Curtis distance. When the distance between the samples was smaller, the species-composition structure was more similar, and the PCoA diagram and NMDS analysis revealed the similarity of the composition of the diet between sika deer and Reeve’s muntjac (Fig. 5a, b), which was consistent with our expected results. In addition, a hierarchical clustering heat map was convenient for the intuitive identification of the species present in corresponding samples; Rubus had the highest abundance among the three species. The tree plot indicated that there was much similar distance in the samples between sika deer and Reeve’s muntjac but less similar with Chinese hare (Fig. 5c).

Figure 5. 

a PCoA analysis chart, in which each point represents a sample b NMDS analysis chart. Diagram analysis with 95% confidence ellipse c hierarchical clustering diagram. Analysis of the hierarchical clustering tree diagram and the stacked bar diagram of the top 10 genera in abundance (Ch: Chinese hare; Sd: Sika deer; Rm: Reeve’s muntjac).

The intergroup difference analysis shows the difference between the intragroup and intergroup sample distances. Compared with Reeve’s muntjac, the intragroup distance of Reeve’s muntjac was smaller than the intergroup distance of sika deer and Chinese hare (Fig. 6a). While compared with sika deer, the intragroup distance of sika deer was slightly higher than the intergroup distance of Reeve’s muntjac (Fig. 6b). With different species composition, the difference between the intragroup should be smaller than the intergroup. In comparison with the NMDS, this phenomenon was speculated to be the fact that the diet composition of Reeve’s muntjac was much similar to those of sika deer, and the Reeve’s muntjac has more forage plant diversity.

Figure 6. 

Intergroup difference analysis a shows the boxplots of the distances between samples in the sika deer group and the distances between samples in this group and samples in other groups b shows the boxplots of the distances between samples in Reeve’s muntjac group and the distances between samples in this group and samples in other groups.

By analyzing the niches of the three sympatric herbivorous animals, we found that the highest niche breadth was Chinese hare (~7.78), followed by sika deer (~4.53) and Reeve’s muntjac (~3.44). The niche overlap index between sika deer and Reeve’s muntjac was 0.83, sika deer and Chinese hare was 0.62, and Reeve’s muntjac and Chinese hare was 0.69. The overlap index ranges from 0 to 1, where 0 indicates that the food ranges do not overlap at all, and 1 indicates that the food ranges overlap entirely. Our results suggested that Reeve’s muntjac and sika deer have the highest diet overlap (Table 3).

Quantitative analysis is of great significance for the families, genera, and species of the herbivores’ diet. Liu (2012) investigated the diet of wild northeast sika deer in winter and found the sik deer feed on 35 plant species belonging to 25 genera in 17 families. Liu (2020) showed that 131 plants, including nine herbs, 31 shrubs, and eight trees, were foraged by South China sika deer in the Tianmu Mountains. One hundred and thirty-nine food items were identified in the feces of sika deer in our study. Comparative analysis shows that the metabarcoding method based on high-throughput sequencing provides more detailed forage plant information.

Yao et al. (2010) found that the sika deer in TNNR mainly eat Homonoia riparia, Fallopia multiflora, Lespedeza bicolor, Pueraria montana, and Vicia faba, and other species, and then especially fed on young leaves and shoots (Li et al. 2014). Compared with our study, the dominant forages in the feces of sika deer were Rubus chingii, Loropetalum chinense, Eurya japonica, Camellia japonica, Euonymus grandiflorus, and Phyllostachys edulis in winter; these are all Chinese medicinal herbs, which ensure the prevention of various diseases throughout the life cycle of sika deer (Yao et al. 2010; Ye 2015). The differences may be attributed to geographical and seasonal differences. The vegetation condition varied according to the different regions, which means that diverse forage plants are available for the sika deer.

In seasons when plant resources are scarce, sika deer will choose to eat non-favorable plants or the available food resources at the moment. Studies on Japanese sika deer showed that they mainly choose their favorite deciduous species from summer to autumn, such as Cornus controversa, and Quercus sp., but from early winter to spring, non-favored herbaceous and tree species, such as Juncus decipiens and Cryptomeria japonica, will be foraged (Nakahama et al. 2020). Therefore, sika deer have a diverse diet and feed on various plants in different seasons (Yao et al. 2010; Wang et al. 2019). In the TNNR, oaks included Quercus acutissima, Q. aliena, Q. chenii, Q. fabri, and Q. serrata (Jiang 2009). Significantly, oak leaves, which are rich in tannins and toxic to most mammals, including cattle (Doce et al. 2009), are conversely found to increase the reproductive rate and fawn survival rate of sika deer in captive breeding by some farmers (Xing et al. 2022). Oak leaves are essential for maintaining healthy sika deer in wild and farmed populations. However, we found a lower abundance in Quercus for three herbivorous mammals’ diets. We speculated this might be related to deciduous Quercus in winter or the gestation period of female sika deer, even though some studies claimed that tannins were not toxic to sika deer because of their rumen microbes and fermentation processes (Li et al. 2013). There are differences in the forage plants of South China sika deer in different seasons and regions.

In our study, we identified most of the forage species; however, someof the forage species that we identified were not previously known from Jiangxi Province. It may be difficult to identify all forage plants using a single gene fragment, and continuous succession of plant communities caused by invasive species adds a level complication to this. Therefore, it is necessary to add auxiliary barcodes as well as strengthen the overall investigation of potential food resources in the reserve. The construction of a database of plant species barcodes for the Taohongling Sika Deer Reserve would provide a reference and source of sequence alignments. Such as database would allow for a more accurate determination of the diets of herbivores and allow for better comparisons of the diets of sympatric herbivores.

Competition theory indicates that the greater the overlap of resources between species, the greater the competition coefficient because of the widespread use of niche overlap to estimate competition for resources (Colwell and Futuyma 1971). Tibetan red deer have a similar diet to sympatric ungulates, which inevitably leads to interspecific conflicts in food use (Lu et al. 2020). However, spatial-temporal variations in dietary consumption of the two dominant rodent species on Mount Kilimanjaro, Tanzania, have been found serve as a mechanism of resource portioning that enable these species to coexist with a niche overlap (Mulungu et al. 2011; Thomas et al. 2022); this contradicts the key assumption of competition theory. Thus, we cannot limit the research on interspecific competition in sympatric species to trophic ecological niches, and spatial-temporal dimensions should also be considered. For sika deer in the TNNR, we do not know the exact reasons for spatial variations in dietary overlaps with its two sympatric herbivores. We suggest that a cautious approach is required to interpret the high dietary overlaps and their implications for competitive interactions among the three studied herbivores in the TNNR. As pointed out by others (Belovsky 1986; Jenkins and Wright 1988; Gordon and Illius 1989), high dietary overlaps may not necessarily imply competition (Andersen et al. 2017) and may simply indicate that the food item is sufficient, permitting sympatric species to share resources. In the TNNR there are increasing populations of sympatric wild boar and Reeve’s muntjac. Sika deer face competition stress in space and food resources, especially when food is scarce in winter (Li et al. 2014). Thus, even when the determinants of competition mechanisms are uncertain, competition does exist and its important role among species cannot be denied.

Competition among sympatric species is mostly expressed as a compensatory mechanism in ecological niches when species are similar in one dimension, they differ on another. Food resources, habitat, and temporal partitioning are the most common dimension partitioned (Bagchi et al. 2003). For example, high dietary overlap among the species may result in niche differentiation (Reitz and Trumble 2002; Yin et al. 2007; Cao et al. 2009); large herbivores are forced to expand their food range to avoid competition during periods of food scarcity (Noor et al. 2013). In northeastern China, red deer tend to increase their browsing intensity to maintain their high food intake, but sika deer meet their relatively constant food intake and potential nutritional requirements by increasing their bite diameter in winter. This reflects the short-term foraging strategies by sharing similar foods with the sympatric ungulates (Zhong et al. 2020). Currently, we do not know the competition mechanism of sika deer and further studies are needed to determine the coexistence mechanism with its sympatric species.

In our study, we found the niche breadth of the sika deer was higher than the Reeve’s muntjac. Optimal forage theory suggests that preference and palatability will be selected for the animals in abundant food periods. While in a period of scarce food resources, feeding generalization will occur by selecting different forage plants (Belovsky 1978). It remains to be studied whether the higher dietary niche breadth of sika deer results from the physiological characteristics of digesting a wide range of foods or to avoid competition. Moreover, the niche breadth of Chinese hare was larger than sika deer and Reeve’s muntjac. We speculate that Chinese hare, as opportunistic feeders, have a broader range of forage plants but consume less due to their smaller body size. The diet composition of Chinese hares includes trees, shrubs, and herbs, but this may be due to a passive and random proximity foraging strategy or even the indirect ingestion from the process of grinding teeth.