Identifying Environmental Refuges (“Coldspots”) from Infection by Batrachochytrium Dendrobatidis of Amphibians in Eastern Europe ^†

Abstract

Amphibians are the most threatened group of vertebrates. While habitat loss poses the greatest threat to amphibians, a spreading fungal disease caused by Batrachochytrium dendrobatidis (Bd) is seriously affecting an increasing number of species. Although Bd is widely prevalent, there are identifiable heterogeneities in the pathogen’s distribution that are linked to environmental parameters. Our objective was to identify conditions that affect the geographic distribution of this pathogen using species distribution models (SDMs), with a special focus on Eastern Europe. SDMs can help identify hotspots for future outbreaks of Bd, but perhaps more importantly, they can identify locations that may be environmental refuges (“coldspots”) from infection. In general, climate is considered a major factor in driving amphibian disease dynamics, but temperature in particular has received increased attention. Here, 42 environmental raster layers containing data on climate, soil and human impacts were used. Mean annual temperature range (or ‘continentality’) was found to have the strongest constrain on the geographic distribution of this pathogen. Using the partial dependence visualization module in the R package ‘embarcadero’, a number of corresponding coldspots were identified.

1. Introduction

Amphibians are the most threatened group of vertebrates, with a third of currently known species in danger of extinction. Although habitat loss clearly poses the greatest threat to amphibians, a newly recognized fungal disease is seriously affecting an increasing number of species [1]. This disease, caused by the chytrid fungus Batrachochytrium dendrobatidis (Bd), has been linked to the global decline of amphibian species and represents the greatest documented loss of biodiversity attributable to a pathogen [2]. Although Bd is widely prevalent, there are identifiable heterogeneities in the pathogen’s distribution that are linked to environmental parameters [3]. In this respect, species distribution models (SDMs) have proven to be useful tools for predicting Bd distribution and elucidating the importance of a wide range of environmental covariates considered to affect Bd occurrence. The first Bd SDMs to be developed were global in scope [4,5]. Using SDMs, our objective was to identify conditions that constrained the geographic distribution of this pathogen in Eastern Europe in an aim to identify hotspots for future outbreaks of Bd but, perhaps more importantly, identify locations that may be environmental refuges (“coldspots”) from infection [6]. Undoubtedly, both aspects are essential for proactive conservation planning [7].

2. Materials and Methods

Localities for Bd were gathered from GBIF (https://www.gbif.org, accessed on 17 January 2021) (https://doi.org/10.15468/dl.bhqawb, accessed on 17 November 2021) and the literature [8,9], etc. Because many uncertainties are associated with SDM projections, particularly when it comes to building an SDM for a species expanding its home range into a new area, we used only records of European localities for the analysis. In total, there were 648 such records. These were filtered out by enforcing a distance of 50 km between records; we used this filtering process because ecological niche models are sensitive to sample bias [10]. In the end, the total number of records was reduced to 116. To build models, environmental values at localities of known occurrence were determined and then used to identify geographic regions that have similar combinations of environmental values. Several types of environmental variables at a geodetic resolution of 5 arc minutes were used as a proxy for the fundamental niche [11]: (1) the Bioclim dataset (https://www.worldclim.org/bioclim, accessed on 27 December 2020); (2) several eco-attributes, such as human fragmentation, accessibility and appropriation (https://databasin.org, accessed on 27 December 2020); (3) the ENVIREM dataset (https://envirem.github.io/, accessed on 17 November 2021); and (4) the Global Soil Dataset (http://globalchange.bnu.edu.cn/research/soilw, accessed on 27 December 2020) (see also: [12]). SDMs were generated by employing Bayesian additive regression trees (BART), a powerful machine learning approach. Running SDMs with BARTs has recently been greatly facilitated by the development of an R package, ‘embarcadero’ [13], including an automated variable selection procedure that is highly effective for identifying informative subsets of predictors. Additionally, the package includes methods for generating and plotting partial dependence curves and visualization called spatial partial dependence plots, which can reclassify predictor rasters based on their partial dependence plots and show the relative suitability of different regions for an individual covariate. Habitat suitability values ranged from 0 to 1. Model performance was assessed using measures of accuracy: the area under the receiver–operator curve (AUC, [14]) and the true skills statistic (TSS, [15]).

3. Results

Both measures of accuracy showed that the SDM performed very well (AUC = 0.92 and TSS = 0.73). The automated variable selection procedure identified informative subsets of predictors, for which continentality (°C × 10), the minimum temperature of the coldest month (°C × 10), Thornthwaite’s aridity index [16,17], pH (measured in a calcium chloride solution), and human appropriation were all significant [18,19]. Temperature is considered one of the most important environmental factors driving chytridiomycosis [20], with a lower thermal limit below 4 °C [21]. Amongst the used covariates, continentality is perhaps the most distinguishable dimension of the climatic niche of Bd in the study area, featuring the seasonal amplitude at an ambient temperature. As shown in Figure 1, it is clear that large differences between low temperatures in the cold season and high temperatures in the hot season are limiting factors for the pathogen, with its habitat suitability (the ‘response’) dropping from above 0.6 in the west of the continent to below 0.3 in the east (Figure 2). An indication that low temperatures limit the spread of Bd is the increasing trend of the partial dependence curve plotted for the minimum temperature of the coldest month, showing a steep increase in habitat suitability (to over 0.5) upon approaching the 0 °C mark. Earlier bioclimatic variables associated with precipitation were found to make high contributions to SDMs considering Bd [6]. In our case, Thornthwaite’s aridity index, commonly used for measuring the aridity of an area and based on both precipitation and temperature, showed a better performance than most other bioclimatic variables, perhaps because of its compound character. The partial dependence curve built for this index clearly highlights wet, humid and marginally moist sub-humid (values between 1 and 31, or slightly over 32) climates as suitable for pathogens. Of all soil features, pH turned out to make the highest contribution to the SDM. Outbreaks of chytridiomycosis may be affected by pH, but the pH optimum (pH 6–7) for B. dendrobatidis is not different from common pH values of freshwater systems [21]. Our model indicated an optimum of around 6.5, where the suitability was the highest. Finally, human appropriation, the only human-related covariate selected for its contribution to the SDM, provided a useful measure of human intervention in the biosphere through the appropriation of net primary production. Our model explicitly points towards areas of greater human intervention as areas more likely to be suitable for Bd.

4. Discussion

Using the SDM and an arbitrarily selected threshold of 0.2, we identified locations that may be environmental refuges (“coldspots”) for amphibians from infection by Bd in Eastern Europe (Figure 3). These particular areas are close to the dividing threshold in Poland and Moldova (also shared with Romania and Ukraine), where adequate management and conservation plans for protecting amphibians have to be designed to begin with. The results of the modelling assume that large portions of Latvia, Lithuania, Ukraine, and the Kaliningrad Province of the Russian Federation will not be favorable for the pathogen; niche case of Belarus, this seems to hold for the entire country

5. Conclusions

SDMs, on the one hand, predict the geographic extent of a species, and on the other hand, they can identify the contribution of habitat covariates to explaining this distribution. When applied to Bd, we find SDMs to be useful for identifying conditions that constrain the geographic distribution of this pathogen and identify locations that may be environmental refuges from infection.