The sacred deer conflict of management after a 1000‐year history: Hunting in the name of conservation or loss of their genetic identity

Human–wildlife conflict involves diverse stakeholders with conflicting values. Resolving such conflicts necessitates the development of management plans rooted in scientific knowledge and establishing of social consensus. In Nara City, situated within the Japanese Archipelago, wild sika deer (Cervus nippon) hold sacred significance due to religious reasons and have been protected for over a millennium, resulting in a distinct genetic identity. However, the escalating deer population has caused significant agricultural damage in the areas surrounding the sanctuary. Consequently, a debate has arisen regarding the advantages and disadvantages of implementing lethal measures to address individuals that might be considered sacred deer but are regarded as pests in the vicinity of the sanctuary. Here, we used mitochondrial DNA (mtDNA) control region and nuclear simple sequence repeat (SSR) markers to detect the origin of deer in the areas neighboring the sanctuary (management areas). As a result, two genetic populations of deer were detected in Nara City. In the sanctuary, we detected only one specific mtDNA haplotype (S4). On the other hand, seven haplotypes, including S4, were detected in the management area. SSR analysis also suggested that the sika deer in the management area may be an admixed population of multiple origins from the sanctuary and out of Nara City. Interbreeding populations may expand into the sanctuary, and unique genetic populations for more than 1000 years may disappear. This study suggests that ordinary deer could soon replace the deer revered and protected by the people of Nara. Additionally, the proximity of sanctuary deer to tourists worldwide and the interaction with wild deer in other areas pose a potential risk of spreading zoonotic diseases. Urgent decisions are required to determine whether to advocate for extermination in managed areas for ‘conservation’ purposes or risk losing the genetic identity of the sacred deer.


| INTRODUCTION
Since ancient times, humans have competed with other species for habitats and resources through innovation and adaptation (Waters et al., 2016).Human activities have contributed to the extinction of many species (Cowie et al., 2022;Diamond, 1989) and the modification of ecosystems (Estes et al., 2011).While the negative impacts of humans on ecosystems have been enormous, taboos based on traditional human religious sentiments and belief systems can positively affect the conservation of certain species and ecosystems (Colding & Folke, 2001;Saj et al., 2006).Furthermore, sacred natural sites established by sociocultural systems are equipped with biological diversity, including plant and animal species, habitats and ecosystems, and ecosystem dynamics and functions, which support life within and beyond the site (Verschuuren et al., 2010).Research on animals and plants in such sacred sites may affirm their significance in natural or social science.However, such studies might also reveal the losing process of their essence or facing significant conflicts between the value of human society and nature.
Human attitudes and traditional religious feelings towards the natural environment and wildlife are changing in the context of modernization.Sacred natural environments can contribute to biodiversity conservation strategies, but their effective utilization also requires the involvement of all relevant stakeholders (Dudley et al., 2009).For example, in some parts of Asia, macaque monkeys (genus Macaca) have been protected by humans due to religious sentiment.Monkeys have adapted to urban environments such as temples (Aggimarangsee, 1992;Fooden, 1971;Watanabe et al., 2007) and these religiously protected areas provide a refuge for wildlife and bring economic value as ecotourism destinations (Priston & McLennan, 2013).As a result, monkeys are increasingly valued as economic animals rather than for their religious significance.At the same time, protected wildlife populations that have increased in size may be repelled by farmers, as they cause crop damage, or they may be sold illegally by animal traders (Loudon et al., 2006;Schilaci et al., 2010).Human-wildlife conflict is a complex interplay of ecological, economic, and socio-political factors.Wherefore, human-wildlife conflicts are some of the most pervasive and intractable problems facing conservation biologists (Dickman & Hazzah, 2016).
The conflict between Japanese sika deer (Cervus nippon) and humans is complicated by factors such as religious sentiment, the economic value of tourism, the history of sacred deer indicated by recent scientific analysis, and damage to farming in Nara, Japan.Around the 8th century the capital of Japan (Heijo-kyo) was located in Nara City, on the northern part of the Kii Peninsula, Honshu, Japan.In the city there is a sacred site where religious buildings such as the Kasuga Taisha Shrine and Todaiji Temple stand.Within this sanctuary, wild Japanese sika deer are considered sacred and have been protected as an isolated population for over 1000 years (Torii & Tatsuzawa, 2009).In particular, Tokugawa Ieyasu, the first shogun of the Edo period (1603-1867), established severe penalties for deer killing to flaunt the ruler's power (Nara Park History Editorial Committee., 1982;Torii & Tatsuzawa, 2009).This sanctuary was named Nara Park in 1880 and was managed as an urban park (Nara Park History Editorial Committee., 1982).Since 1957, the sika deer in Nara City have been protected as a 'natural monument', and today more than 1000 deer inhabit Nara Park alone (Nara Park History Editorial Committee., 1982;Watanabe, 2012).A landscape where tourists can interact with wild deer within an urban area is a rarity in the world, with approximately 17 million tourists visiting Nara Park recorded in FY2019 (Nara City, 2021).Recent genetic analysis found a unique mitochondrial DNA (mtDNA) haplotype and nuclear simple sequence repeat (SSR) alleles in the sacred deer population (Takagi et al., 2023).Estimated divergence times between the sacred deer population and the other populations were 1400 years or more.These results suggest the extinction or fragmentation of deer populations in the large area of the Kii Peninsula due to historical human activity.
However, the sacred sika deer population, with a history of approximately 1400 years, has become troublesome for the surrounding present-day farmers.While encountering numerous wild deer up close is initially a delightful experience for tourists visiting from distant places, it can quickly shift to the distressing sight of pest being left unchecked for local farmers.Since the end of the 20th century, the increase in the number of deer individuals, due to strict protection within the sanctuary, has caused agricultural damage in the surrounding area.The farmers around Nara Park filed two lawsuits in 1979 and 1981 against the national government (Agency for Cultural Affairs), the Kasuga Taisha Shrine, and others, demanding compensation for damage to their crops caused by sika deer (Watanabe, 2001).In 1985, a settlement was reached between the farmers and the national government, and it was decided to divide Nara City into two areas for the management of sika deer: the deer protected area and the deer management area.In the management area, farmland was fenced to prevent deer intrusion, and population management, including the culling of deer, was permitted (Figure 1, Watanabe, 2001Watanabe, , 2012)).However, even after the settlement, the administration did not cull the sika deer because they were worried about the negative impact on tourism.
A critical focus in wildlife management is to increase the scientific basis of management decisions and plans rather than relying solely on opinion and experience (Apollonio et al., 2017).Unfortunately, the problem of deer management in Nara City has become more complex over the last 40 years because there has been a lack of scientific basis and accurate solutions, including the absence of a master plan supported by legislation and stakeholder consensus.Over the past 40 years, the distribution of sika deer in the Japanese Archipelago has increased 2.7 times, with about 600,000 deer captured in 2019 (Oka et al., 2022).Nara City is no exception to this trend as over the past 20 years the number of deer of unknown origin has increased in the eastern part of the city (Figure 2, in the management area), which was suspected to have been a vacuum in the distribution of sika deer (Watanabe, 2017).Local government, Nara Prefecture and Nara City used tax revenues from citizens to compensate for agricultural damage caused by these deer of unknown origin (Watanabe, 2017).In response to this changing situation, the Nara Prefectural Government formulated a new policy for deer management in Nara City, "The Deer of Nara Management Plan" in 2016.Based on this new policy, a culling program for sika deer began in July 2017 (Nara Prefecture, 2019a).However, the origin of the sika deer in the management area is yet to be determined and institutional contradictions regarding management remain.The minimum distance from the protected area's center to the management area's deer cull point is only about 4 km.In other words, if a sacred sika deer in Nara Park were to leave the sanctuary, it would be exterminated.Conversely, sika deer of the 'common' lineage from outside the sanctuary are not exterminated when they disperse into the sanctuary.The sacred sika deer, protected as a natural monument in Nara City, have created a conflict in deer management, as they can be both of economic value for tourism and a pest for agriculture.Wildlife management decision-making is a process of linking scientific knowledge with solutions to social problems such as human-wildlife conflict.The problem of conservation, especially of cultural animals or plants deeply connected to human society, ultimately might become a value system issue.Nevertheless, disciplinary, and interdisciplinary sciences provide stakeholders with knowledge-based options for solving a social problem (Yahara et al., 2021).For example, wildlife management units of a variety of taxa have been studied based on population genetic analysis (Budd et al., 2018;Edelhoff et al., 2020;Hohenlohe et al., 2021;Robinson et al., 2012).In this study, the local population genetic structure of Japanese sika deer in Nara City was analyzed using two genetic markers, mtDNA control region and nuclear SSR.These results allowed us to determine the origin of the sika deer that appeared in the management area and to re-examine the validity of the zoning in the current wildlife management plan for sika deer.This study aimed to publish these deer population genetics findings and stimulate discussion on the prospects for deer management in Nara City, considering various competing values such as human-wildlife conflict.

| Study area
Nara City (34 41 0 N, 135 50 0 E) is in the northern part of Nara Prefecture on central Honshu Island, Japanese archipelago (Figure 2).The total area of the city is approximately 276.9 km 2 , and the total human population is about 353,000 (Nara City., 2021).Nara Park (6.6 km 2 ), located in the center of Nara City, contains historical buildings, such as Todaiji Temple and Kasuga Taisha Shrine, which are World Heritage Sites of the Ancient Capital of Nara, and the virgin forest of Mount Kasuga (UNESCO, World Heritage Conservation, 2023).As mentioned above, Nara City is zoned according to the Japanese sika deer management plan.This plan has defined three areas: the protected area (districts A and B), the management area (district D), and the buffer area (district C) (Figure 2, Table S1) (Nara Prefecture, 2019bPrefecture, , 2022)).

| Sample collection and genetic analysis
Muscle or fecal samples were collected from 167 deer between 2017 and 2019 at nine sites in districts A, B, and D in Nara City (Figure 2, Table S1).In district A, fecal samples were collected simultaneously at two sites (Todaiji, NP1 and Tobihino, NP2) in Nara Park on September 6, 2018.In district B, fecal samples were collected at Nara University of Education (NUE) on December 29, 2018.One sample of one individual fecal pellet group was taken and genotyping further confirmed that there were no duplicates of individuals.After collection, fecal samples were immediately chilled on dry ice and stored at À30 C after delivery to the laboratory.In district D, muscle samples were collected from six sites (Ohyagyuu, OY; Tahara, TH; Touri, TU; Seika, SI; Semakawa, SM; and Yagyuu, YG) between 2017 and 2019.The samples from district D were collected from individuals captured in accordance with the sika deer management plan for Nara City (Nara Prefecture, 2019a, 2022).Deer was treated following the guidelines of the management plan and laws on hunting in Japan.The culling points of sika deer in district D may contain slight deviations due to the hunters' reporting system.In this study, we consider the nine sampling sites as sub-populations for our analysis.Muscle samples were stored in 70%-90% ethanol at room temperature.DNA was extracted from feces using DNeasy Stool Mini Kit (Qiagen, Hilden, Germany), and from muscles using DNeasy Blood & Tissue Kit (Qiagen, Hilden, Germany).

| Genetic variation in Nara City
Sequencing data for each primer and each individual were aligned using GeneStudio Pro 2.2 (GeneStudio, Inc.) with default settings.Tandem repeats (TRs) were excluded from the sequence data, as per a previous study (Nagata et al., 1999) and used for subsequent phylogenetic analysis because the mtDNA-CR of the Japanese sika deer has multiple TRs of 37-40 bp, with different lengths for each haplotype.All sequences were aligned using Clustal W (Thompson et al., 1994) in MEGA 11 (Tamura et al., 2021).The haplotypes of individuals were determined using DnaSP ver.6.12.03 (Rozas et al., 2017).The naming convention for mtDNA-CR haplotypes follows that of (Takagi et al., 2023).The haplotype diversity (h) and nucleotide diversity (π) for each subpopulation or area were calculated using Arlequin ver.3.5.2.2 (Excoffier & Lischer, 2010).Genetic differentiation between the protected and management areas was determined in Arlequin ver.3.5.2.2 (Excoffier & Lischer, 2010) by calculating F ST values based on haplotype frequencies.In addition, pairwise F ST values for the nine sub-populations were calculated using Arlequin ver.3.5.2.2 (Excoffier & Lischer, 2010) and p-values for significance tests were obtained by permutation procedures (1000 times).
We calculated the genetic variation of nuclear SSR genotype data in terms of the average number of alleles (N a ), observed heterozygosity (H O ), and expected heterozygosity (H E ) using GenAlEx 6.503 (Peakall & Smouse, 2012), and calculated allelic richness (Rs) and inbreeding coefficients (F IS ) in FSTAT version 2.9.4 (Goudet, 1995(Goudet, , 2003)).We tested the significance of deviations from Hardy-Weinberg equilibrium and heterozygous deficiencies at a significance level of 5% after Bonferroni correction using FSTAT version 2.9.4 (Goudet, 1995(Goudet, , 2003)).The pairwise F ST values of sub-populations based on SSR genotypes were calculated using GenAlEx 6.503 (Peakall & Smouse, 2012) and their significance tested using permutation tests.

| Individual-based genetic analyses
Population genetic structure was inferred by using STRUCTURE ver.2.3.4 (Falush et al., 2003;Pritchard et al., 2000) implementing the admixture model and correlated allele frequencies without prior sampling location information.Twenty independent runs were performed for each K from 1 to 8 with a burn-in of 500,000 and 1,000,000 Markov chain Monte Carlo iterations (MCMC).The delta-K method (Evanno et al., 2005) implemented in STRUCTURE Harvester (Earl & VonHoldt, 2012) was used to evaluate the probability of the data (Ln P(D)) for each K value and choose the best-fit K value.Outputs of 20 independent runs around optimal K values were integrated and visualized using CLUMPAK software (Kopelman et al., 2015).To confirm the geographical location of haplotypes and genetic clusters, the results of the mtDNA and STRUCTURE analyses were mapped using QGIS 3.26.3(QGIS.org,2022.).
The STRUCTURE analysis showed that ΔK was highest at K = 2, with one genetic cluster identified in the protected area (cluster 1) and one identified in the managed area (cluster 2) (Figure S1).To illustrate the genetic characteristics of individuals across geographic space, we plotted individuals (coded by mtDNA haplotype) by using their membership value (Q 1 ) for cluster 1 (protected area) at K = 2 for the y-axis, and the distance from the center of Nara Park to their sampling point for the x-axis.The proportion of each individual belonging to cluster 1 (Q 1 ) was calculated using the STRUCTURE results described above, and the distance from the center of Nara Park (34 41 0 00.4 00 N, 135 50 0 23.6 00 E) to the sampling point of each individual was calculated using GenAlEx 6.503 (Peakall & Smouse, 2012).See Table S3 for the genotypes of each sample and the membership values (Q-values) for each individual to each cluster.In addition, we used BAYESASS 3.0.4(Wilson & Rannala, 2003) to confirm the gene flow between the management unit of deer in Nara City.BAYESASS analysis units were the following three, based on the STRUCTURE analysis and mtDNA haplotype results: the protected area within 4 km from the center of Nara Park, the management area near the protected area between 4 and 10 km from the center of Nara Park, and the management area is not near the protected area more than 10 km from the center of Nara Park.We conducted a run for 1 Â 10 7 iterations for MCMC, with sampling chains for every 1000 iterations, discarding the first 1 Â 10 6 iterations as burn-in.We increased the mixing parameters of the allele frequency (ΔA = 0.3) and inbreeding coefficient (ΔF = 0.3), leaving the other parameters at default.Mixing and convergence of the MCMC were visually inspected by using Tracer ver.1.7.2 program (https://beast.community/tracer).

| Genetic variation of mtDNA in Nara City
The mtDNA haplotype analysis revealed single haplotype fixation in the protected area, and multiple haplotypes in the management area, including the haplotype from the protected area.The mtDNA-CR sequences excluding TRs (683 bp) were determined from 135 individuals, and a total of seven haplotypes (M1, M2, M4, S4, S5, S7, and S11) were identified across all areas (Table 1).These haplotypes were also detected in the comprehensive analysis of Kii Peninsula populations by (Takagi et al., 2023).In the protected area, all individuals sampled had only one haplotype, S4.No differences in genetic variation were detected among the three sites in the protected area (both h and π were 0).In the management area, all seven haplotypes, M1, M2, M4, S4, S5, S7, and S11, were found (Table 1).In the management area h was 0.60 and π was 0.0048.Nucleotide sequence data of mtDNA reported herein are available in the DDBJ/EMBL/Gen-Bank databases (The International Nucleotide Sequence Database Collaboration, INSDC) under the accession numbers LC730667-LC730801.
Despite the small spatial scale of the study area (less than 15 km), the distribution of haplotypes showed a clear trend: S4 haplotypes in and around the protected area and other haplotypes in the external boundaries of the management area.In the management area, the individuals with the S4 haplotype tended to be caught relatively close to the protected area (Figure 3a).The other six haplotypes in the management area were detected far from the protected area (Figure 3a).The average F ST value between sub-populations in the protected and management areas was 0.152 ( p < 0.01).The mtDNA had pairwise F ST values ranging from À0.020 to 0.596 (Table S2).Significant differentiation was detected between the protected area and the management area: NP1, NP2, and NUE were significantly different from SM and YG.Significant genetic differentiation was also detected within the management area between OY and TH, and SM and YG (Table S2).

| Genetic diversity of nuclear SSR in Nara City
The genetic diversity of nuclear SSR also differed between the protected and management areas.The nuclear SSR genotypes of 163 individuals could be T A B L E 1 The genetic diversity of 9 deer sub-populations in Nara City based on the analysis of mtDNA control region and 14 nuclear simple sequence repeat (SSR).
determined.The genotype data were confirmed by two or more independent PCRs in all samples, with an overall missing value of 1%.In the protected area, N a was 3.7, R S was 3.64, H O was 0.45, H E was 0.52, and F IS was 0.16 (Table 1).Two fecal samples in NP2 had 100% matching genotypes.In the management area, all measures of genetic variation were higher than in the protected area: N a = 6.0,R S = 4.91, H O = 0.62, and H E = 0.64 (Table 1).The F IS value was 0.04 (Table 1).The average F ST value between sub-populations in the protected and management areas was F ST = 0.057 (p = 0.01).The pairwise F ST values based on the SSR genotypes ranged from À0.010 to 0.159 and were lower than those of mtDNA (Table S2).

| Genetic structure and admixture among sika deer in two areas of Nara City
In the STRUCTURE analysis, although the highest mean posterior probability value was detected for K = 2, a higher variance of the values among runs was apparent when K ≥ 4 (Figure S1).However, the results for K = 2 and K = 3 were stable, and the highest ΔK was detected when K = 2 (Figure S1).Therefore, the results for K = 2 were used and are discussed in this study (Figure 4, Figure S1).
Cluster 1 was the dominant gene pool in the protected area and the F ST value, indicating the magnitude of genetic drift of the cluster, was 0.17 (Figure 4).Individuals assigned to cluster 2 were observed at high frequency in the management area sub-populations, SI, SM, and YG (Figure 4).These sub-populations are located relatively far from the protected area.The F ST value of cluster 2 was 0.01, suggesting that genetic drift from the putative common ancestral gene pool has had little influence.The individuals in subpopulations TH, TU, and OY, located relatively close to the protected area, were a mixture of cluster 1 and 2 (Figures 3b, 4, and S1).
In the area between 0 and 2 km from the center of Nara Park, which is a protected area, most individuals had a high membership rate to cluster 1 (Q 1 ), and all individuals had haplotype S4 (Figure 5).However, in the management area (4-14 km from the center of Nara Park), the individuals had very diverse Q 1 values: the Q 1 values of individuals with the S4 haplotype ranged from 0.028 to 0.967, whereas the values for individuals with other haplotypes ranged from 0.032 to 0.585 (Figure 5).Individuals with the S4 haplotype in the management area were only found at distances of 4-10 km from the center of Nara Park, whereas those with other haplotypes were found at distances of 4-14 km, which is relatively far from the center of Nara Park.BAYESASS analysis, assessing gene flow between areas, confirmed a migration between the protected area and the management area.The migration rate per generation from the protected area to the management area within 4-10 km from the center of Nara Park was determined to be 0.143 (Figure S2).The migration rate from F I G U R E 4 Results from genetic clustering analysis using STRUCTURE (K = 2) based on 14 loci for sika deer in Nara City.Full names of abbreviated sub-populations are provided in Figure 2 and Table S1.
the management area close to the protected area within 4-10 km from the center of Nara Park to the more distant management area more than 10 km from the center of Nara Park was 0.293 (Figure S2).The migration rate from the management to the protected areas was also measured at 0.012 and 0.019 (Figure S2).

| Bottleneck and loss of genetic variation in the sacred deer population
The present analysis revealed the existence of genetically heterogeneous populations in the protected area of Nara Park and the management area surrounding it.In the protected area, only the S4 mtDNA haplotype was found, a haplotype identified as unique to the sacred deer population in previous study (Table 1) (Takagi et al., 2023).The nuclear SSR results indicated lower genetic variation in the protected area compared to management areas (Table 1).The low level of genetic diversity in the protected area compared to other populations on the Kii Peninsula is likely due to the spatial isolation of the population and repeated bottlenecks, discussed below (Takagi et al., 2023).In protected areas, all individuals had a high proportion of assignment to cluster 1 (Figures 4  and 5).Furthermore, the F-value of cluster 1 (F 1 = 0.17) was higher than that of cluster 2 (F 2 = 0.01), indicating that the protected population was more strongly affected by genetic drift (Figure 4).Behavioral observations have indicated that female deer in Nara Park have a narrow home range (Kawamura, 1957;Miura, 1976;Nara Prefecture, 2019b).
It is also assumed that it is difficult for deer from the outside of cities to enter habitats and populations that exist in isolation within cities.Moreover, in Nara Park, the population density of deer is very high in the 1.2 km 2 flat area, with approximately 1000 individuals per km 2 (Tatsuzawa et al., 2002;Torii & Tatsuzawa, 2009), which may also make it difficult for deer from other areas to invade.
The low genetic diversity of sika deer in the protected area of Nara City probably results from multiple past population declines.Declines in the genetic diversity of Japanese sika deer in several areas of the Japanese archipelago have been attributed to population declines around the 19th century (Goodman et al., 2001;Iijima et al., 2023;Konishi et al., 2017;Nabata et al., 2004).The historical protected area in Nara City is no exception.In the 19th century, with the modernization of Japan, the sika deer in Nara City temporarily lost their special religious status, and the number of deer decreased to several dozen (Torii & Tatsuzawa, 2009).The sika deer population in the protected area recovered to about 1000 individuals at one point.However, the political turmoil and food shortages caused by World War II are said to have reduced the deer population to 79 individuals again (Nara Park History Editorial Committee., 1982).Conservation efforts have since restored the deer population in Nara Park, which exceeded 1000 in the 1960s and which has maintained a population of around 1000-1300 since then (Nara City, 2021;Torii & Tatsuzawa, 2009).These extreme population fluctuations over a period of less than 200 years have undoubtedly contributed to the loss of genetic diversity in the population of the current protected area.Furthermore, the Nara populations may have experienced similar additional bottlenecks in their more than long history (Takagi et al., 2023).The population in the management area was characterized by higher genetic variation and included individuals that had the S4 mtDNA haplotype, which characterizes the protected area.Only in the past 20 years has deer been steadily observed in the management area.Sika deer was unconfirmed in the management area in the 1998-1999 survey, and the density of sika deer rapidly increased to 22.8 deer/ km 2 at the time of the 2018 survey (Nara Prefecture, 2022).The rapid increase in genetic diversity and density in the management area population can be explained by the expansion of individuals from the protected area and the invasion of individuals from surrounding areas.Previously, these populations had little gene flow.We detected seven mtDNA haplotypes of sika deer in the management area (Table 1), including the S4 haplotype that is observed in the protected area.The other six haplotypes (M1, M2, M4, S5, S7, and S11) have been observed in various parts of the Kii Peninsula (M1); in other municipalities in Nara Prefecture adjacent to Nara City and Mie Prefecture (M2, M4, S5, and S11); and in Kyoto Prefecture (S7) (Takagi et al., 2023).The high mtDNA diversity in the management area suggests that deer have migrated from multiple regions or are descendants of these migrants.The establishment of new populations of deer that have multiple lineages of mtDNA has been observed throughout Japan and genetic diversity tends to be higher in these new populations (Sato et al., 2013;Toma et al., 2021;Yamazaki, 2018).
The spatial distribution of mtDNA haplotypes and nuclear SSR genotypes indicates that the sacred deer, which have the unique S4 haplotype, are now becoming admixed with deer targeted for extermination in other areas.In previous study (Takagi et al., 2023), samples were collected from Nara Park in 2003-2005, whereas the samples used in this study were collected in 2018.A similar composition of haplotypes was observed in the protected area in both studies, suggesting that the composition of haplotypes has stayed the same over the past 10-15 years, and that female deer from outside Nara Park have not immigrated to the protected area.Individuals in the protected area had the S4 mtDNA haplotype and were assigned with high confidence to the cluster 1 gene pool by STRUCTURE analysis of the nuclear SSR genotypes.These genetic characteristics differ from those of individuals sampled in the management area only a few kilometers away from the center of the protected area (Figures 3 and 5).In the management area, some deer with the S4 haplotype also had high Q 1 values, but notably there were a considerable number of individuals that had low to moderate Q 1 values, indicating admixture.In contrast, the Q 1 values of all individuals with mtDNA haplotypes from outside Nara City (M1, M2, M4, S5, S7, and S11) were relatively low (Figure 5).These results show how two genetically distinct lineages are becoming mixed in the management area; female deer from the protected area are mating with male deer from other areas.Nearly 10 years have passed since deer have been present in the management area (Nara Prefecture, 2022).Thus, breeding between individuals from the protected area and those outside Nara City has existed for several generations.This may explain the existence of individuals with the S4 haplotype and low Q 1 values: they are likely due to further mating, for example, F2 or backcrosses, between females descended from the S4 haplotype and males from other areas.In the management area, individuals with mtDNA haplotypes originating from outside Nara City (not S4) show low Q 1 values, suggesting less admixture with males from the protected area (Figure 5); thus, breeding between males from the protected area and females from outside Nara City might be rare.
Management practices in the protected area may influence such sexual differences in admixture patterns in the management areas.In Nara Park, sika deer are strictly protected, injured individuals are treated, and pregnant female deer are kept temporarily to provide a safe place to give birth.In addition, for the safety of tourists, males in the protected area are captured with tranquilizing guns and their antlers are cut off before being released back into the wild (Watanabe, 2014).Removing antlers reduces the fighting instinct of individuals (Takaragawa & Kawamichi, 1977) and may reduce the reproductive success rate of males that move from the protected area to the management area.Antlerless males immigrating outside Nara Park have been observed (Watanabe, 2012), suggesting that deer born in the protected area are migrating to management areas.The results of the BAYESASS analysis also indicate deer migration between the protected area and the management area (Figure S2).Indeed Nara City's management area, where there is a high-density of females with relatively less competitive males, would be an ideal environment for 'outside' male deer to invade; thus promoting the admixture of 'the sacred deer' lineage.The level of population admixture is likely to increase in the protected area in the future.

| For creative dialogue for the future 'conservation' of sacred deer
This study paints a picture of the human-wildlife and human-human conflicts surrounding the sacred deer that have been protected for over 1000 years.The results of this analysis lead to two contradictory scientifically optimal management policies.If we allow the disappearance of sacred deer lineage, we should continue the current conservation.Conversely, to maintain the lineage of sacred deer, we must restore habitat gaps through intensive eradication around the sanctuary.Over the past several 100 years, hunting pressure has resulted in the decline and fragmentation of the distribution range of sika deer in the Japanese archipelago (Iijima et al., 2023;Tsujino et al., 2010).This has also affected the population dynamics of other species, for example, on the Kii Peninsula, the speciation of dung beetles (Phelotrupes auratus), which feed heavily on deer dung (Araki & Sota, 2023), coincided with the fragmentation of deer populations (Takagi et al., 2023).The recent re-expansion of deer populations in the area surrounding Nara and their contact with the Nara Park population may indicate that deer are experiencing a release from the negative influence of human hunting of the past 2000 years.Gene flow into the Nara City population may represent an opportunity to enhance the genetic variation of the protected area population, whose diversity has been lost due to isolation and bottlenecks.
However, the reinstatement of deer will bring new problems and anguish to various stakeholders in human society.The damage to agriculture in the management area, which is already a problem, will become more serious.The stakeholders promoting deer conservation in protected areas will also have difficult decisions.Ongoing admixture of the sacred deer with other deer could lead to a loss of the sacred deer's genetic specificity that has long history.It might be unwelcome for those who value the unique genetic identity of the deer in Nara Park.In order to preserve this genetic identity, deer population densities around the protected area need to be reduced.For this purpose, capturing deer (both sacred and others) inhabiting management areas is required in the name of 'conservation'.The economic or cultural standpoints and values of stakeholders will inevitably influence decisions on these issues.Therefore, it will be necessary to reconcile and agree on the values of various stakeholders based on the results of scientific analysis.
While Nara Park is a precious environment where wild deer can be seen up close, the proximity of humans and wildlife is causing various public health problems.Nara Park is an internationally recognized tourist attraction.If Nara Park were to become an epicenter of infectious disease due to zoonosis, it could cause widespread problems.Several tick species have been identified in Nara Park (Fujimoto & Yamaguchi, 1990;Sakai & Torii, 2014), indicating a concern for the spread of tickborne diseases (TBD), such as Severe Fever with Thrombocytopenia Syndrome (SFTS) and Japanese spotted fever.Exposure to the SFTS virus has been recorded in some sika deer populations in the Japanese archipelago (Lundu et al., 2018); thus, it is possible that an influx of individuals from other populations may carry SFTS virus into the Nara Park population.Also, Ikushima et al. (2021Ikushima et al. ( , 2023) ) showed that Quinolone-resistant Escherichia coli (QREC) is prevalent among sika deer in Nara Park and highlighted the risk of QREC transmission to humans.Deer species have also been reported to be capable of contracting Mycobacterium (Böhm et al., 2007;Ghielmetti et al., 2021;Itoh et al., 1992), and care should also be taken in Nara Park.Finally, it is predicted that cervids have a high risk of SARS-CoV-2 infection (Damas et al., 2020) and widespread infection of white-tailed deer (Odocoileus virginianus) and human-to-deer spillover events have been reported (Caserta et al., 2023).Thus, sika deer infection with SARS-CoV-2 and deer-to-human transmission are important issues for the management of Nara Park, an international tourist destination with the closest proximity between people and deer in the world.
One of the goals of this study was to provide scientific knowledge that could help resolve human-wildlife and/or human-human conflicts that have been created due to an increase in the number of wildlife populations that are protected for religious reasons.This research has allowed us to identify the genetic population structure of one such wildlife population, but this is only one piece of the puzzle towards resolving the conflict.We still need a clearer understanding of the factors contributing to conflict on the human side.For example, rural communities are often particularly frustrated by the urban elite's perceived protection or imposition of wildlife damage (Zscheischler & Friedrich, 2022).We must not focus only on living organisms, but also consider the more comprehensive socio-economic, ecological, and cultural conditions in which conflicts occur (Dickman, 2010).In the future, we must publish such broad findings and provide decision-making criteria for stakeholders from various positions.The resolution of human-wildlife conflicts often takes a long time.However, in the case of sika deer in Japan, time is of the essence: various problems caused by deer population and distribution expansion are becoming increasingly serious and little time may be left to resolve conflicts and appropriately manage the sacred deer population with a 1000-year history.cooperation was provided by Ms. Ayako Takano, Nara University of Education, Mrs. Yumena Nakamura-Kojo, Mr. Nobuaki Kojo, and Miss Rio Kojo, the Wildlife Partnership Office Yamagata Dormouse Research group.We would like to thank Dr. L. Faulks, Scientific English Editor (https://scientificenglisheditor.com/), for English editing the manuscript.

F
I G U R E 1 Differences in the way people deal with Japanese sika deer in Nara City.(a) Deer are a tourist resource in the protected areas and tourists feed them.(b) In the management area, deer are vermin and are prevented from feeding on crops by fences and other means.Photos by Harumi Torii.U R E 2 (a) The location of Nara City in Nara Prefecture, Japanese Archipelago.(b) The sampling areas and sub-populations examined in this study.
Distribution map of sika deer samples in Nara City (with known location information) and their (a) mtDNA haplotype, and (b) the estimated proportion of the membership coefficient (Q 1 ) obtained by STRUCTURE analysis where a Q-value of zero indicates cluster 2 and a Q-value of 1 indicates cluster 1.For STRUCTURE results, see also Figure4.See Figure1for the color coding of districts.
The estimated proportion of the membership coefficient (Q 1 ) to cluster 1 plotted against distance from the center of Nara Park.Q values indicate attribution rate to cluster 1.The mtDNA haplotypes are indicated by pink circles for S4 and black triangles for the other haplotypes (M1, M2, M4, S5, S7, and S11).