The economic loss caused by herbivore browsing in forest plantations is a concerning problem in many areas around the world. Information on the spatial distribution of browsing damage is important for forest owners when selecting locations for new plantations, because planting trees in areas of high browsing pressure increases economic losses. Although it is difficult to survey browsing damage across large areas, sporadic sampling data on browsing damage are often collected by foresters, governments, and researchers. Thus, in this study, we applied a generalized additive model (GAM) for analysis of sporadic data to reveal large-scale spatial variation in deer (Cervus nippon) browsing damage. A map of browsing pressure produced by a GAM that used years after planting (plantation age) and location as independent factors showed a few areas of high browsing pressure. In addition, browsing damage increased with increasing plantation age, and plantation stands aged 2+ years showed more browsing damage. Areas with high browsing damage estimated based on plantation stands aged 2+ years generally coincided with areas of high deer abundance, with some exceptions. Thus, this model reflects large-scale browsing damage relatively well and will help forest owners to avoid locating new plantations in areas of high browsing pressure.
Damage to forest plantations due to herbivore browsing causes economic losses (Conover 1997, Reimoser 2003) and is a problem in many areas (Gill 1992, Côté et al. 2004, Brousseau et al. 2017). In particular, tree seedlings in clearcuts face high browsing pressure by terrestrial animals (Tremblay et al. 2007), and, when browsing pressure is high, forest plantations can fail to reach maturity (Takatsuki & Gorai 1994, Bergström & Bergqvist 1997, Borkowski et al. 2017, Uno et al. 2019). This means that planting trees in areas with a potential high risk of browsing damage might increase economic losses. Thus, it is important to determine the spatial distribution of browsing pressure to identify the most suitable areas for tree planting taking into account that damages show a spatially heterogeneous distribution (Bergquist & Orlander 1998, Bulinski & McArthur 2000). Nonetheless, as it is often difficult to survey the damages across large areas, the large-scale spatial distribution (e.g., at the landscape level) of the damages remains poorly understood. However, it has often been surveyed in sporadic sampling sites by foresters, governments, and researchers. We hypothesized that the large-scale spatial distribution of the damage, as estimated by using these sporadic sampling data, could help in siting new plantations in areas where damages would be minimized.
Large-scale spatial distributions are often estimated from sporadic data in studies on marine animals. For example, Yano et al. (2017) estimated the spatial distribution of the North Pacific spiny dogfish (Squalus suckleyi) across the North Pacific region from catch records. Also, Suzuki et al. (2018) and Furuichi et al. (2020) estimated the location of spawning grounds of small pelagic fishes in the East China Sea and the Sea of Japan from fixed-point observation data. Yamamoto et al. (2015) clarified the density distribution pattern of streaked shearwaters (Calonectris leucomelas) in the Pacific Ocean based on observations collected by vessel surveys and individual tracking. In these studies, distributions were estimated in a single plane by smoothing sporadic data using a habitat-based additive model (Forney et al. 2012) and producing a distribution map (Grüss et al. 2014).
In Japan, the sika deer (Cervus nippon) population increased rapidly in recent decades (Kaji et al. 2000) and deer damage to forests is a serious problem (Takatsuki 2009). On Kyushu Island in southwestern Japan, deer browsing of forest plantations rapidly increased after 1990, and recent years’ economic losses in forest plantations reached the equivalent of approximately US$6 million per year (Yabe & Yayota 2012). To understand local economic losses, data on browsing damage have been collected sporadically from many areas of Kyushu Island. We produced a map of browsing damage distribution based on this dataset by using a generalized additive model (GAM), with the aim of identifying areas of low browsing damage on Kyushu Island.
In addition, deer density is considered to be an essential factor in determining browsing pressure, because removal of deer decreases browsing damage in small-scale (few square kilometers) experimental sites (Enoki et al. 2016). However, it remains unclear whether deer density affects browsing pressure at a larger scale. On Kyushu Island, deer density at the prefectural scale has been surveyed sporadically for wildlife management by prefectural governments where sika deer are abundant. We produced a map of deer distribution from available density data by using GAM and performed a preliminary analysis of the relationship between browsing pressure and deer density at the large scale.
Material and methodsStudy area
Kyushu Island (Fig. 1) is the third largest (approximately 36,750 km2) of the four main islands of Japan. Forest vegetation on Kyushu Island is divided into two zones: the evergreen broad-leaf forest zone at low altitude and the deciduous broad-leaf forest zone at high altitude. The border between the two vegetation zones is at approximately 800-1000 m a.s.l. Fifty-six percent of forests on Kyushu Island are plantations comprised of Japanese cedar (Cryptomeria japonica D.Don) and Japanese cypress (Chamaecyparis obtusa Endl.). In recent years, harvest of these plantations has increased, and Japanese cedar cuttings and/or Japanese cypress seedlings have been planted in clearcuts.
On Kyushu Island, deer browsing is a major factor that inhibits regeneration of plantations. On Hokkaido Island, with its cooler temperatures, not only deer but also Sasa senanensis (Fr. & Sav.) Rehd. (dwarf bamboo) prevent the establishment of young trees (Uno et al. 2019). In contrast, on the warmer Kyushu Island its distribution was limited (Tsuyama et al. 2011), though Sasamorpha borealis (Hack.) Nakai (evergreen bamboo) can also prevent tree seedling establishment. However, heavy browsing of the bamboo caused its decline (Saruki et al. 2004) and currently it is not a major factor inhibiting tree growth in plantations.
Sika deer inhabit six of the seven prefectures on Kyushu Island (Fukuoka, Kumamoto, Oita, Miyazaki, Kagoshima, and Nagasaki, but not Saga). Nagasaki Prefecture was excluded from the analysis in this study because its sika deer population is totally isolated from the other populations on Kyushu.
Browsing damage
In 2018, browsing damage was surveyed in 338 recently planted forest stands (Tab. 1) in the five focal prefectures on Kyushu Island by the Forest Management Center (Fig. 1). The plantation stands were composed of Japanese cedar and/or Japanese cypress, and the size of the plantation stands was 4.9 ± 4.1 ha (mean ± standard deviation). In each plantation, the area damaged (ha) by sika deer and damage levels were recorded. The damage levels, averaged among the trees within each selected site, were ranked into four main categories (<25%, ≥25% to <50%, ≥50% to <75%, and ≥75%). As years after planting (plantation age) varied from 1 to 5 years, the damage is cumulated from planting to the survey day.
To clarify the distribution of browsing damage, we constructed a GAM with a binomial distribution with a logit-link function. The proportion of damaged area to plantation size, produced by using the cbind function in R ver. 3.5.2 (R Core Team 2018), was used as the dependent variable. Spatial effect, the tensor product of latitude and longitude, was used as an independent variable. In addition, we used a smoothed term for plantation age and the interaction between location and plantation age as independent variables, because we considered that browsing damage accumulates over time. We used the following GAM formula (eqn. 1):
logit(y) = a + s(age) + te(lat, lon) + ti(lat, lon, by = age)
where y is the proportion of damaged area produced by cbind (damage area, planted area - damage area), s is the smoothed term, age is the plantation age, te is the tensor product, lat is the latitude, lon is the longitude, ti is the interaction including tensor production, and a is intercept. We evaluated the significance of the independent variables using a χ2 goodness-of-fit test. In addition, to produce the browsing distribution map, we extrapolated browsing pressure 100 km2 around each surveyed plantation stand based on the results of the GAM. The extrapolation area was limited to forested areas only. The accuracy of the model was evaluated using a four-fold cross-validation (Picard & Cook 1984) and a Spearman’s correlation test.
In addition, we clarified the relationship between the proportion of damaged area and the damage level, by using a generalized linear model (GLM) with a binomial distribution. We used the following GLM formula (eqn. 2):
y = a + DL
where y is the proportion of damaged area produced by cbind (damage area, planted area - damage area), a is the intercept, and DL are the damage levels. All calculations were run in R ver. 3.5.2 (R Core Team 2018).
Deer density
To estimate deer density, we used the fecal pellet count method (Iwamoto et al. 2000, Kondoh et al. 2007). Pellet surveys were performed in Fukuoka, Kumamoto, and Miyazaki prefectures by each prefecture’s forestry experimental station. Pellets were counted a total of 1880 times at 480 sites from 1995 to 2018 (Tab. 2). Based on the pellet counts and average temperature in the survey month, deer density per square kilometer was estimated by using the FUNRYU software (Iwamoto et al. 2000) in Microsoft Excel® (Microsoft Inc., Redmond, WA, USA). By using this method, however, deer densities in a few sites were estimated to be > 700 deer km-2. In Japan, even in areas with the highest deer density, such as Kinkazan (Ito 1986) and Nakanoshima Islands (Kaji et al. 1988), deer density ranged from approximately 30 to 60 deer km-2. Such outliers can increase dispersion. Thus, to avoid large dispersion due to a few outliers, we transformed the deer density (D) to a density index by using the following formula (eqn. 3):
y = log(D+ 1)
where y is the density index. The deer density data contained many zeros. In these cases, addition of a small value such as 10% of the mean or 1 is traditionally used (ICCAT 1997). Because mean deer density in this study was 20.5, we considered that adding 1 in the transformation was acceptable.
To estimate deer density distribution, we constructed a GAM; because approximately 17% of the survey area showed no deer, we resolved the zero-data problem with a Tweedie distribution (Candy 2004, Shono 2008). The density index was treated as a dependent variable. Location, the tensor product of latitude and longitude, was used as the independent variable. Also, because deer density fluctuates year by year, we included the survey year and the interaction between location and year in the model. We used the following GAM formula (eqn. 4):
y = a + s(year) + te(lat, lon) + ti(lat, lon, by = year)
where y is the density index, s is the smoothed term, year is the surveyed year, te is the tensor product, ti is the interaction including tensor production, and a is the intercept. We evaluated the significance of the independent variables using an F test. The accuracy of the model was evaluated using a four-fold cross validation (Picard & Cook 1984) and a Spearman’s correlation test. In addition, to produce the deer density distribution map, we extrapolated the density index 100 km2 around each surveyed plot based on the results of the GAM. The extrapolation area was limited to forested areas only.
ResultsBrowsing damage
Deer browsing was recorded at 76 of the 338 plantation stands, and deer browsing damage by plantation age is shown in Tab. 1. The results of the GAM for browsing damage (r2 = 0.555, deviance explained = 60.5%) demonstrated significant effects of plantation age and location as well as interactions between these variables on browsing damage (Tab. 3). Plantation age was shown to have a positive effect on browsing damage, with the damaged area increasing with age, and the relationship was relatively linear (Fig. 2).
We compared actual damage areas with damage areas estimated using by four-fold cross validation of the GAM. Spearman’s correlation test showed a significant positive correlation (ρ = 0.479, P < 0.001). Thus, the estimation accuracy of the model was relatively high.
We produced damage distribution maps (Fig. 3) based on the model. A few core areas showed high browsing damage, and those locations differed in plantation age, as the effect of the interaction showed. In the northern part of the study area (>33.0° N), one core area containing plantation stands aged 2+ years showed serious damage. In contrast, plantation stands with an age of 1 year showed slight damage. In the central (>32.5 to ≤33.0° N) and southern (≤32.5° N) parts of the survey area, even 1-year-old plantation stands faced relatively high damage, and plantation stands with trees aged 3 or 4 years suffered extensive damage.
The GLM showed a strong positive correlation between the proportion of damaged area and the damage level (n = 76, χ2 = 88.59, p < 0.001). We layered damage level plots on the map showing mean proportion of damage area across all plantation ages produced by the GAM (Fig. 4). Plantations with higher damage levels tended to be concentrated around the extensive damaged areas. In contrast, some plantations with a low damage level were located in areas with a high proportion of damaged area.
Relationship between browsing damage and deer density
Results of the GAM for sika deer density index (r2 = 0.363, deviance explained = 30.3%) showed that survey year and location as well as the interaction between the two variables affected deer density (Tab. 4). We compared actual deer index with deer index estimated by using four-fold cross validation of the GAM. Spearman’s correlation test showed a significant positive correlation (ρ = 0.548, P < 0.001). Thus, the estimation accuracy of this model was relatively high. Although density index tended to increase gradually on Kyushu Island (Fig. 5), the increase rates differed among locations (Fig. 6). In the northern area, deer density gradually increased year by year, and density index exceeded 3.5. The central to southern areas had a few core areas of high deer density. During the study period, the high-density area around 32.5° N expanded and the density increased. In the other core areas, deer density did not change markedly, but the size of the core area expanded slightly.
Locations of core areas with high browsing damage appeared to be a relatively good fit with areas of high deer abundance (Fig. 3 and Fig. 6). For example, the highest browsing pressure occurred in the north, where deer density indexes were consistently ≥3.5. In the central area, very high browsing pressure occurred around the core area with highest deer density. In the south, density indexes ≥3.5 were observed in two of the three areas with high browsing pressure. However, when we viewed the data at a finer scale, few discrepancies between damage and deer density became apparent. For example, in the north, although deer density was relatively high in 2017, browsing damage was low in plantation stands aged 1 year, planted in 2017. In addition, although deer densities in the north were lower than those in the central area or the south during the study period, damage was higher in the north than in the southern part of the central area. In addition, in the southwest deer density was low but browsing pressure was high.
Discussion
Because the estimation accuracy of this model in four-fold cross validation was relatively high, we consider that the GAM represents large-scale browsing damage relatively well. In addition, areas of high browsing damage generally coincided with areas of high deer abundance. This means that damage can actually occur.
We identified a few core areas with high browsing pressure. Also, areas with high damage level matched areas with a high proportion of damaged area relatively well. Thus, we suggest that planting should be avoided in these areas, because it will likely meet with low success. However, there is a caveat. The model showed that browsing damage increases as plantation age increases, which indicates that browsing damage accumulates year on year and/or the amount of forage increases with plantation age. Thus, low browsing pressure in the first year of planting may not be a recommendation to plant more trees in the area. As in the northern part of our study area, browsing pressure may markedly increase in plantations aged 2 years and older even if browsing pressure is low in the first year. Thus, plantation stands aged at least 2 years may be better indicators of where browsing pressure will be high.
A few studies have predicted the distribution of forest damage by using additive models (Schmidt et al. 2010). These studies predicted damage from physical factors such as storms and wind. In contrast, this study indicates that a GAM can be used to predict the damage caused by biological factors. In this study, we considered the association between browsing damage and sika deer density. However, deer species other than sika deer also damage forests (Partl et al. 2002, Balazy et al. 2016). We consider that this model could also be used to predict browsing damages by other deer species.
We studied the patterns of browsing damage caused by sika deer in tree plantations on Kyushu Island, Japan. Our comparison of the distribution of browsing damage with deer density at the landscape level is novel, because many previous studies examining the effect of deer density on browsing were performed at a small scale, in single plantations or site unit level (Daubenmire 1940, Horsley et al. 2003, Verheyden et al. 2006, Koda & Fujita 2011, Iijima & Nagaike 2015, Enoki et al. 2016). We used data from 1995 to 2018 in the model and estimate deer density from 2013 to 2017, because estimate accuracy usually increases as the size of the dataset increase. Although core areas of browsing damage and high deer density generally coincided, there were discrepancies between browsing pressure and deer density in some areas. Browsing pressure by deer is determined not only by density but also by the habitat, especially food abundance (Moore et al. 1999, Partl et al. 2002, Iijima & Nagaike 2017, Takarabe & Iijima 2019). We did not evaluate the effect of the local environment because we wanted to focus on the large-scale comparison between deer density and browsing pressure. To clarify the factors causing the discrepancies, the environmental factors affecting browsing should be studied at a smaller scale.
Conclusion
Reducing browsing damage is an important factor for sustainable forestry management. In this study, we modelled a large-scale spatial variation in deer browsing damage by using a GAM. Because high browsing pressure on young trees is an important factor that reduces the success of mature trees, using this model to predict browsing pressure could help to reduce economic losses caused by deer browsing. In addition, because our model is very simple, consisting of only two independent factors (plantation age and location) it can easily be reproduced in other areas. We expect that the model will help to select locations for new plantations in deer habitats.
Acknowledgements
Data from pellet count surveys were provided by the Kumamoto, Miyazaki, and Fukuoka prefecture governments. We thank Dr. N. Kotaka for helpful comments concerning this study. We acknowledge the invaluable comments of the two anonymous reviewers on earlier versions of this manuscript.
Funding
This study was funded by the Forestry and Forest Products Research Institute (grant no. 201703).
ReferencesBalazy R, Ciesielski M, Sterenczak K, Borowski ZThe role of topography in the distribution and intensity of damage caused by deer in Polish mountain forests. PLoS One 11: e0165 967.2016Bergquist J, Orlander GBrowsing damage by roe deer on Norway spruce seedlings planted on clearcuts of different ages: 1. Effect of slash removal, vegetation development, and roe deer density. Forest Ecology and Management 105: 283-293.1998Bergström R, Bergqvist GFrequencies and patterns of browsing by large herbivores on conifer seedlings. Scandinavian Journal of Forest Research 12: 288-294.1997Borkowski J, Dobrowolska D, Dabrowski W, Banul R, Zaluski DYoung conifer stands form a deer browsing refuge for an oak admixture: silvicultural implications for forest regeneration under herbivore pressure. European Journal of Forest Research 136: 787-800.2017Brousseau M, Thiffault N, Beguin J, Roy V, Tremblay J-PDeer browsing outweighs the effects of site preparation and mechanical release on balsam fir seedlings performance: implications to forest management. Forest Ecology and Management 405: 360-366.2017Bulinski J, McArthur CSpatial distribution of browsing damage and mammalian herbivores in Tasmanian eucalypt plantations. Australian Forestry 63: 27-33.2000Candy SGModelling catch and effort data using generalised linear models, the Tweedie distribution, random vessel effects and random stratum-by-year effects. CCAMLR Science 11: 59-80.2004Conover MRMonetary and intangible valuation of deer in the United States. Wildlife Society Bulletin 25: 298-305.1997Côté SD, Rooney TP, Tremblay J-P, Dussault C, Waller DMEcological impacts of deer overabundance. Annual Review of Ecology, Evolution, and Systematics 35: 113-147.2004Daubenmire RFExclosure technique in ecology. Ecology 21: 514-515.1940Enoki T, Yabe T, Koizumi TChanges in spatial patterns of sika deer distribution and herbivory of planted seedlings: a comparison before and after deer population control by culling. Journal of Forest Research 21: 84-91.2016Forney KA, Ferguson MC, Becker EA, Fiedler PC, Redfern JV, Barlow J, Vilchis IL, Ballance LTHabitat-based spatial models of cetacean density in the eastern Pacific Ocean. Endangered Species Research 16: 113-133.2012Furuichi S, Yasuda T, Kurota H, Yoda M, Suzuki K, Takahashi M, Fukuwaka MDisentangling the effects of climate and density-dependent factors on spatiotemporal dynamics of Japanese sardine spawning. Marine Ecology Progress Series 633: 157-168.2020Gill RMAA review of damage by mammals in north temperate forests: 1. Deer. Forestry 65: 145-169.1992Grüss A, Drexler M, Ainsworth CHUsing delta generalized additive models to produce distribution maps for spatially explicit ecosystem models. Fisheries Research 16: 113-133.2014Horsley SB, Stout SL, De Calesta DSWhite-tailed deer impact on the vegetation dynamics of a northern hardwood forest. Ecological Applications 13: 98-118.2003ICCATReport of the bluefin tuna methodology session. International Commission for the Conservation of Atlantic Tunas - ICCAT, Collective Volume of Scientific Papers 46: 187-268.1997Iijima H, Nagaike TThe factors that determine the intensities of deer browsing and debarking on broadleaf tree around artificial grasslands. Journal of Forest Research 22: 199-203.2015Iijima H, Nagaike TAppropriate vegetation indices for measuring the impacts of deer on forest ecosystems. Ecological Indicators 48: 457-463.2017Ito TPopulation trend of Sika deer on Kinkazan Island. Mammalian Science 26: 29-31. [in Japanese]1986Iwamoto T, Sakata T, Nakazono T, Utaoka H, Ikeda K, Nishishita Y, Tsuneda K, Doi TImprovement of the pellet count method for the estimation of Sika deer density. Mammalian Science 40: 1-17. [in Japanese with English abstract]2000Kaji K, Koizumi T, Ohtaishi NEffects of resource limitation on the physical and reproductive condition of Sika deer on Nakanoshima Island, Hokkaido. Acta Theriologica 33: 187-208.1988Kaji K, Miyaki M, Saitof T, Ono S, Kaneko MSpatial distribution of an expanding sika deer population on Hokkaido Island, Japan. Wildlife Society Bulletin 28: 699-707.2000Koda R, Fujita NIs deer herbivory directly proportional to deer population density? Comparison of deer feeding frequencies among six forests with different deer density. Forest Ecology and Management 262: 432-439.2011Kondoh H, Ikeda K, Koizumi TSpatial estimation of Sika deer population density distribution. Journal of Forest Planning 13: 1-14.2007Moore NP, Hart JD, Langton SDFactors influencing browsing by fallow deer Dama dama in young broad-leaved plantations. Biological Conservation 87: 255-260.1999Partl E, Szinovatz V, Reimoser F, Schweiger-Adler JForest restoration and browsing impact by roe deer. Forest Ecology and Management 159: 87-100.2002Picard RR, Cook RDCross-validation of regression models. Journal of the American Statistical Association 79: 575-583.1984R Core TeamR: a language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria.2018Reimoser FSteering the impacts of ungulates on temperate forests. Journal for Nature Conservation 10: 243-252.2003Saruki S, Inoue S, Shiiba Y, Nagasawa H, Ohosaki S, Kubota KDistribution and growth situation of Suzutake (Sasamorpha borealis) damaged by grazing of Shika deer (Cervus nippon nippon) in Miyazaki Forest of Kyushu University: case study in 2003. Bulletin of Kyushu University Forest 85: 47-54.2004Schmidt M, Hanewinkel M, Kändler G, Kublin E, Kohnle UAn inventory-based approach for modeling single-tree storm damage experiences with the winter storm of 1999 in southwestern Germany. Canadian Journal of Forest Research 40: 1636-1652.2010Shono HApplication of the Tweedie distribution to zero-catch data in CPUE analysis. Fisheries Research 93: 154-162.2008Suzuki KK, Yasuda T, Kurota H, Yoda M, Hayashi A, Muko S, Takahashi MSpatiotemporal variations in the distribution of round herring eggs in the East China and Japan Seas during 1997-2013. Journal of Sea Research 142: 1-10.2018Takarabe K, Iijima HContrasting effect of artificial grasslands on the intensity of deer browsing and debarking in forests. Mammal Study 44: 173-181.2019Takatsuki SEffects of sika deer on vegetation in Japan: a review. Biological Conservation 142: 1922-1929.2009Takatsuki S, Gorai TEffects of Sika deer on the regeneration of a Fagus crenata forest on Kinkazan Island, northern Japan. Ecological Research 9: 115-120.1994Tremblay J-P, Hout J, Potvin FDensity-related effects of deer browsing on the regeneration dynamics of boreal forests. Journal of Applied Ecology 44: 552-562.2007Tsuyama I, Nakao K, Matsui T, Higa M, Horikawa M, Kominai Y, Tanaka NClimatic controls of a keystone understory species, Sasamorpha borealis, and an impact assessment of climate change in Japan. Annals of Forest Science 68: 689-699.2011Uno H, Inatomi Y, Ueno M, Iijima HEffects of sika deer (Cervus nippon) and dwarf bamboo (Sasa senanensis) on tree seedlings in a cool-temperate mixed forest on Hokkaido Island, Japan. European Journal of Forest Research 138: 929-938.2019Verheyden H, Ballon P, Bernard V, Saint-Andrieus CVariations in bark-stripping by red deer Cervus elaphus across Europe. Mammal Review 36: 217-234.2006Yabe T, Yayota CManagement techniques for effective control of overabundant sika deer populations in Kyushu. Mammalian Science 52: 114-115. [in Japanese]2012Yamamoto T, Watanuki Y, Hazen EL, Nishizawa B, Sasaki H, Takahashi AStatistical integration of tracking and vessel survey data to incorporate life history differences in habitat models. Ecological Applications 25: 2394-2406.2015Yano T, Ohshimo S, Kanaiwa M, Hattori T, Fukuwaka M-A, Nagasawa T, Tanaka SSpatial distribution analysis of the North Pacific spiny dogfish, Squalus suckleyi, in the North Pacific using generalized additive models. Fisheries Oceanography 26: 668-679.2017
Effect of plantation age on browsing damage. The continuous line is the best-fit curve generated by the generalized additive model; dashed lines show 80% confidence intervals.
Distribution map of areas of browsing damage on Kyushu Island, Japan. Browsing damage measured at each plantation was extrapolated to the surrounding 100 km2. Grey and colored areas indicate forest; white areas are not forested. Age in each panel indicates the plantation age.
Comparison of proportion of damaged area (contour) estimated by generalize additive model and damage levels in each site (circle).
Temporal change of estimate deer density index. The continuous line is the best-fit curve generated by the generalized additive model, dashed lines show 80% confidence intervals.
Estimated deer density on Kyushu Island, Japan, from 2013 to 2017. Deer density was extrapolated to the 100 km2 surrounding density survey points. Grey and colored areas indicate forest; white areas are not forested.
Summary of surveyed plantations. (SD): standard deviation.
Parameter
Year planted
2017
2016
2015
2014
2013
Plantation age (years)
1
2
3
4
5
No. of plantations
47
78
78
64
71
No. of damaged plantations
3
13
20
15
25
Mean size of plantation (ha)
4.12
4.04
4.49
6.07
5.52
SD of plantation size
2.67
3.76
3.06
5.95
4.08
Mean size of damaged area (ha)
0.24
0.61
0.8
1.51
1.49
SD of damaged area
1.27
1.72
1.83
4.8
3.01
Summary of fecal pellet count surveys.
Prefecture
Survey years
Number of survey plots
Number ofsurveys
Fukuoka
1995-2018
105
878
Kumamoto
2002, 2003, 2006, 2009, 2010, 2014
247
627
Miyazaki
2015, 2017
128
375
Summary of the generalized additive model used to examine browsing damage to plantations. (age): plantation age; (lat): latitude; (lon): longitude.
Independent variables
edf
χ2
P value
s(age)
3.35
48.38
<0.001
te(lat, lon)
8.69
40.59
<0.001
ti(lat, lon, by = age)
33.97
191.38
<0.001
Summary of the generalized additive model of density of sika deer. (year): survey year; (lat): latitude; (lon): longitude.