Next Article in Journal
Wildlife Tourism and Climate Change: Perspectives on Maasai Mara National Reserve
Previous Article in Journal
Impact of Climate Change on the Bioclimatological Conditions Evolution of Peninsular and Balearic Spain During the 1953–2022 Period
 
 
Search for Articles:
Title / Keyword
Author / Affiliation / Email
Journal
Article Type
 
 
Advanced Search
 
Section
Special Issue
Volume
Issue
Number
Page
 
Journals
Climate
Volume 12
Issue 11
10.3390/cli12110184
climate-logo
Submit to this Journal Review for this Journal Propose a Special Issue
► Article Menu

Article Menu

share Share announcement Help format_quote Cite question_answer Discuss in SciProfiles thumb_up
...
Endorse textsms
...
Comment
first_page
settings
Order Article Reprints
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Species Distribution Modeling of Ixodes ricinus (Linnaeus, 1758) Under Current and Future Climates, with a Special Focus on Latvia and Ukraine

by
Volodymyr Tytar
1,
Iryna Kozynenko
1,
Mihails Pupins
2,
Arturs Škute
2,
Andris Čeirāns
2,
Jean-Yves Georges
3 and
Oksana Nekrasova
1,2,3,*
1
I. I. Schmalhausen Institute of Zoology, National Academy of Sciences of Ukraine, 01030 Kyiv, Ukraine
2
Institute of Life Sciences and Technologies, Daugavpils University, LV5400 Daugavpils, Latvia
3
Université de Strasbourg, CNRS, IPHC UMR7178, F-67000 Strasbourg, France
*
Author to whom correspondence should be addressed.
Climate 2024, 12(11), 184; https://doi.org/10.3390/cli12110184
Submission received: 27 July 2024 / Revised: 4 November 2024 / Accepted: 7 November 2024 / Published: 11 November 2024
(This article belongs to the Special Issue Ecological Modeling for Adaptation to Climate Change)
Browse Figures
Versions Notes

Abstract

:
This study assesses the impact of climate change on the distribution of Ixodes ricinus, which transmits Lyme disease, a growing public health concern. Utilizing ensemble models from the R package ‘flexsdm’ and climate data from WorldClim, ENVIREM, and CliMond, we project habitat suitability changes for the focus species. The models, validated against Lyme disease incidence rates, predict a 1.5-fold increase in suitable habitats in Latvia, contrasted with a 4.5-fold decrease in suitable habitats within Ukraine over the coming decades. SHAP values are analyzed to determine the most influential climatic features affecting tick distribution, providing insights for future vector control and disease prevention strategies. The optimal bioclimatic environment for I. ricinus seems to be an intricate balance of moderate temperatures, high humidity, and sufficient rainfall (bio7, 14, 18, 29). Also, radiation during the wettest quarter (bio24) significantly influences tick distribution in northern countries. This implies an increased presence of ticks in Scandinavian countries, Baltic states, etc. These findings largely coincide with our projections regarding bioclimatic suitability for ticks in Latvia and Ukraine. These shifts reflect broader patterns of vector redistribution driven by global warming, highlighting the urgent need to adapt public health planning to the evolving landscape of vector-borne diseases under climate change.

1. Introduction

Climate change is widely recognized as a key factor influencing biodiversity [1,2]. There is a specific interest in evaluating the potential shifts in habitat suitability for medically or economically important species as a result of climate change. According to reports from the IPCC, climate change is expected to alter the distribution of certain vectors and vector-borne diseases [3,4]. Tick-borne illnesses have emerged as a significant health concern in recent years. Ticks belonging to the suborder Ixodida are found across the globe and are highly sensitive to climate conditions due to their reliance on a complex interplay of climate factors for their development and survival [5]. Ticks serve as carriers for various pathogens, such as protozoa, viruses, bacteria, and nematodes [6,7], which are known to cause several diseases, including Lyme borreliosis. In Europe, Lyme borreliosis is among the most prevalent tick-borne diseases in humans [8], with Ixodes ricinus (Linnaeus, 1758) being the primary vector species. The disease is particularly widespread in Central Europe, with the highest infection rates recorded in ixodid ticks [9]. A recent study reported that Central Europe has the highest rate of residents with Lyme disease compared to North America (21% and 9%, respectively) [10].
Ixodes ricinus is the most prevalent and widely distributed tick species in Europe [11]. It parasitizes a diverse range of wild vertebrates, as well as incidental hosts like humans, livestock, and pets [12]. In Europe, I. ricinus is found across a broad geographical area including Scandinavia, the British Isles, Central Europe, France, Spain, Italy, the Balkans, Eastern Europe, and Ukraine [13] (ECDC 2020). Global average temperatures have increased by 0.7 °C over the past century, with a projected additional increase of 1.1 °C by the end of the 21st century [14]. Climate change can have a dual impact on tick populations and the transmission of tick-borne diseases by influencing vector biology. On the one hand, studies have demonstrated that climate change contributes to the expansion of I. ricinus in northern regions and/or increases their local abundance [15,16,17]. Conversely, high temperatures can lead to increased mortality in ticks due to desiccation [6]. Additionally, ticks are less active in seeking hosts during periods of high temperatures, which can result in higher mortality rates by reducing their success in finding hosts [18]. Therefore, global warming may negatively affect the species’ habitat suitability in the southern parts of its range of distribution [19] and may potentially lead to geographical contractions of its range in areas where the combination of climate variables for development and survival turns out to become less or entirely not favorable. Under these circumstances, there is a need to construct the current climatic niche of I. ricinus and project it into future conditions to assess how climatic change will affect the geographic distribution of the vector. Ultimately, climate-induced changes in vector distribution may affect the epidemiology of vector-borne diseases [20] and are expected to be of profound character [21].
With the global incidence of Lyme disease and other tickborne diseases on the rise, there is a critical need for data-driven tools to assess the magnitude of this problem and provide scientifically grounded guidelines for public health decision makers. Moreover, there is an interest in identifying potential areas of distribution of vector species in order to assess the future infection risks with vector-borne diseases and improve surveillance efforts [22]. In this respect, species distribution modeling (SDM) can be considered as a data-driven numerical tool that combines field observations of species occurrence or abundance with environmental estimates (usually climatic) [23]. Modeling the ecological requirements of species to explore future disease transmission patterns has been found to be challenging [24]; however, it is essential to assess the potential impacts of climate change on vector distribution, which can provide a theoretical basis for the prevention and control of vector-borne diseases [25,26].
In this study, we used an SDM approach to investigate the distribution of the Lyme disease vector I. ricinus under the current and forecasted (by 2030 and 2050) climate and to rank habitat suitability areas for the tick at a country scale with a focus on Latvia and Ukraine. Of particular concern is Russia’s military action in eastern Ukraine since 2022, which has resulted in the substantial relocation of Ukrainian citizens to mainly western and central regions of the country [27,28]. Therefore, this raises the question of how safe these regions are in terms of disease transmission. This could substantially help health officials in managing and monitoring tick distributions, assessing their overlap and potential contact with humans because it is vital to decrease the risk of zoonotic disease transmission.

2. Materials and Methods

2.1. Tick Input Data

As input, SDMs require georeferenced biodiversity observations. Faced with limited data for Ukraine, we built species distribution models using European and adjacent occurrences. Presence data were retrieved from online public databases, which were accessed using the R package ‘rgbif’, an interface to the Global Biodiversity Information Facility [29,30] (GBIF.org, 2022), supplemented from the literature [31,32,33,34], including Ukrainian and Latvian sources [35,36,37].
The conducted search for occurrence data yielded a total of 12,095 non-duplicate georeferenced records of I. ricinus across Europe and adjacent areas. Using the pre-modeling functions from the ‘flexsdm’ R package, the number of presence records was reduced to 1827 [38].

2.2. Climate Data

SDMs are predominantly influenced by climate, with the variables used in their development typically reflecting climatic factors [39,40]. This is logical as climate plays a key role in determining environmental suitability [41]. Bioclimatic parameters were gathered as raster layers from three climatic databases and utilized separately to construct the projected SDMs and assess their performance.
Nineteen bioclimatic variables indicating general patterns in precipitation and temperature, including extremes and temperature seasonality, were obtained at a resolution of 2.5′ from the WorldClim website (https://www.worldclim.org/; accessed on 21 November 2022; [42]; Table S1; Miroc6 ssp126).
For the modeling in this study, a set of 16 climatic and 2 topographic variables (from the ENVIREM dataset, downloaded from http://envirem.github.io; accessed on 26 November 2022; Table S2) at a resolution of 2.5′ was utilized. These variables were selected based on a recent reassessment of their biological significance, with many likely to have direct relevance to ecological or physiological processes that influence species distributions [43]. These variables are important to consider in species distribution modeling applications, particularly as some (such as potential evapotranspiration) are directly linked to processes crucial for species ecology. The inclusion of topographic variables is also significant as they can alter the impacts of climate descriptors.
Datasets from CliMond v.1.2 were downloaded from https://www.climond.org/ at a resolution of 10′ (accessed on 28 November 2022; Table S1; A1B). These datasets comprised the core set of 19 bioclimatic variables (temperature and precipitation) and an additional set of 16 variables (solar radiation and soil moisture) [44,45].
In previous studies, we used various climate change scenarios to predict future conditions [3,40]. We believe that scenarios utilizing average or minimum climate change values are the most optimal for our purposes. After reviewing the literature, we concluded that scenarios such as A1B (CliMond) and Miroc6 (ssp126—WorldClim) are the most suitable [33,34,39,40,42,44,45]. These scenarios assume a balanced emphasis on all factors and represent a middle path between various possible climate changes, making them useful for our forecasts and analysis.
Climate variables often exhibit high collinearity, and most SDM approaches necessitate selecting one variable from highly correlated ones [46]. To address this, the “removeCollinearity” function in the “virtualspecies” R package was utilized [47]. This function evaluates the correlation among variables in the provided stack of environmental variables, identifying variables that are not collinear and grouping them based on their collinearity levels. Collinearity was assessed across Europe and neighboring regions to avoid shifts in collinearity when projecting training extents [48]. Given that climate variables frequently have skewed distributions or outliers, the Spearman correlation method was applied with a threshold of 0.8 [49].
Removing high collinearity and the selection of one amongst strongly correlated predictors gave the following results, presented in Table 1.

2.3. Modeling Procedure

To prepare the input data, pre-modeling functions from the ‘flexsdm’ R package were utilized [38]. The calibration area was delineated by creating 500 km buffers around presence points. Occurrence data filtering was implemented to mitigate sample bias by randomly eliminating points that were densely clustered (oversampling) in both environmental and geographical space. The modeling functions within the package facilitated the fitting and evaluation of various modeling approaches, encompassing individual algorithms, tuned models, and ensemble models. Seven machine learning (ML) SDM (species distribution models) methods were applied, including ‘Generalized Additive Models’, ‘Gaussian Process Models’, ‘Generalized Linear Models’, ‘Maximum Entropy’, ‘Artificial Neural Networks’, ‘Random Forest’, and ‘Support Vector Machine’. This led to the development of a comprehensive ensemble model, where predictions from individual algorithms were amalgamated to generate a consensus distribution, thereby reducing model uncertainty and enhancing model transferability [50]. To create a global model across Europe using SDMs, we applied the Maxent method based on Maximum Entropy (see Supplementary Materials, Figure S3, and Table S1, https://biodiversityinformatics.amnh.org/open_source/maxent/; accessed on 28 November 2022) [51,52].
Initially, the models were assessed using the area under the receiver operating characteristic curve (AUC) and the true skill statistic (TSS) [53,54]. AUC scores range from 0 to 1, with 0 indicating consistently incorrect model predictions and 1 signifying consistently accurate model predictions. AUC values between 0.7 and 0.8 are deemed acceptable, while values exceeding 0.8 are considered good to excellent. TSS values range from −1 to +1, with −1 representing systematic inaccuracies and +1 indicating systematic correctness. TSS values below 0.4 are considered poor, those between 0.4 and 0.8 are considered useful, and values above 0.8 are classified as good to excellent [25]. In addition to AUC, the continuous Boyce index from the ‘flexsdm’ package was employed to address the limitations of the AUC. This index, which only requires presence data, evaluates how closely model predictions align with the observed distribution of presences across prediction gradients [55]. It ranges from −1 to +1, with positive values indicating consistent predictions, values near zero suggesting similarity to a random model, and negative values indicating contradictory predictions [56].
Subsequently, it was assumed that maps depicting habitat suitability for the Lyme disease vector could serve as proxies for the potential distribution of the disease agent itself. The performance of each ensemble model was assessed by correlating the derived habitat suitability values with the reported incidence of Lyme disease (expressed as the number of reported cases per 100,000 people) sourced from the Johns Hopkins Lyme and Tick-Borne Disease Dashboard (https://www.hopkinslymetracker.org; accessed on 23 December 2022). To address spatial dependencies, a modified t-test was employed to calculate the statistical significance of the correlation coefficient (a corrected Pearson’s correlation) based on geographically effective degrees of freedom, as implemented in the ‘SpatialPack’ package [57,58].
Finally, the ‘best’ ensemble SDM was selected based on the aforementioned performance criteria and the correlation between I. ricinus habitat suitability values and reported Lyme disease incidence. The top-performing SDM was categorized into areas of low, medium, and high potential habitat suitability. These thresholds were defined using Jenks natural breaks, a method that optimizes the grouping of numerical variables by minimizing the deviation within each class while maximizing the deviation from the means of other groups. Jenks natural breaks provide a standardized approach to determining class intervals for continuous numerical variables [59].

2.4. Feature Importance

A critical challenge in the adoption of ML models is their inherent lack of interpretability, which is often referred to as the “black box” problem. Shapley Additive exPlanations (SHAP) [60,61,62] represents a significant advancement in addressing this concern by providing a framework for understanding how ML models arrive at their predictions. SHAP offers various advantages over other feature importance methods, including its model-agnostic nature. Leveraging game theory concepts, SHAP provides a robust framework for feature importance attribution in ML models, despite the number of used variables. The SHAP value quantifies the magnitude and direction (positive or negative) of the feature’s influence on the prediction [60,61]. In our case, the relative importance of predictors was measured with mean absolute Shapley values, which can be interpreted as the magnitude of the relative contribution of a predictor (feature) towards a model output [63]. The R package ‘shap-values’ (https://github.com/pablo14/; author Pablo Casas; accessed on 10 December 2023) was used in a modified version to calculate SHAP values for a selected model. A summary plot offers a comprehensive view of the most influential features in a model and ranks features based on their effect on the model’s predictions. In other words, SHAP values can be integrated into global explanations such as variable importance, but using a completely different method from the various feature selection methods widely applied so far [64]. Until only now has the use of SHAP for explaining which environmental factors (including climatic) influence the spatial distribution of species occurrences started to be explored more widely [65,66,67,68,69].
Maps of habitat suitability in the GeoTIFF format were processed and visualized in SAGA GIS v.9.3.0 [70]. Statistical data were analyzed using the PAST v.4.04 software package [71,72] and/or the R environment [73].

3. Results

As a result of the modeling, we obtained models with a sufficiently high evaluation (Table 2). SDMs built using various selected subsets of predictors from the WorldClim v.2, ENVIREM, and CliMond v.1.2 datasets performed satisfyingly, with the AUC values considered to mostly correspond to good and excellent, and sound TSS values. In this respect, more weight should be given to the continuous Boyce index, which is considered one of the most appropriate metrics for assessing model predictions applied to presence-only datasets and which is a more reliable metric than the AUC [74,75,76,77]. Also, the continuous Boyce index shows a high uniformity of performance of the obtained SDMs and it now reveals a challenge to picking the ‘best’ model. Importantly, values reported in Table 2 indicate models in which present predictions are highly consistent (>0.9) with the distribution of presences in the datasets.
The ‘flexsdm’ package supports model evaluation based on a number of the abovementioned performance metrics. Amongst the threshold criteria, maximum training sensitivity plus specificity was considered a reasonable choice [78,79].
Further, the ensemble models for current and projected climate scenarios were clipped to Latvia and Ukraine and classified using the corresponding module in SAGA GIS, in accordance with Jenks natural breaks, into three categories, i.e., “high”, “medium” and “low” habitat suitability (Figure 1, Figure 2 and Figure 3).
To gain a comprehensive understanding of tick distribution in Europe, we developed a global European model of tick distribution. Our analysis revealed that factors such as the annual temperature range (bio7) and the precipitation of the driest week (bio14) also influence tick distribution (Table 3, Figure 4, Figure 5, Figures S3 and S4).
For the spread of ticks in Ukraine, the most significant parameters are bio18—precipitation of the warmest quarter—and bio7—annual temperature range. In Latvia, the most important factors are bio29—the highest weekly moisture index—and bio24—the radiation of the wettest quarter (Figure 5, Figures S1, S2 and S4).

4. Discussion

The distribution of ticks, especially Ixodes ticks, depends on numerous environmental factors that influence their life cycle, survival (including host seeking), and distribution. Key factors include temperature and humidity, which play a critical role in the physiological processes of ticks. These parasites can exist outside of hosts for extended periods, and their sensitivity to climate changes makes them particularly vulnerable to fluctuations in temperature and humidity levels. To better understand pan-European trends, we conducted a study within Europe. Our analysis revealed that factors such as the annual temperature range (bio7) and the precipitation of the driest month (bio14) also influence tick distribution. Bio14, the precipitation of the driest month, is determined by identifying the driest month for each year in the period, i.e., the month with the least precipitation, and then averaging the precipitation across all of the driest months. It is important to note that the driest month may vary from year to year. The annual temperature range (bio7) measures the temperature range between the most extreme (warmest and coldest) months (Table 3, Figure 4). Sharp temperature fluctuations and insufficient humidity can also affect the viability and distribution of ticks [80,81,82,83].
In our study, we analyzed SHAP values to understand feature importance governing current tick distribution in Ukraine and Latvia by calculating the average absolute values for each feature. This permitted us to provide an idea of which features are generally most important in influencing predictions within the area of our countries of interest (Figure 5). Given their extended off-host existence as generalist parasites, ixodid ticks, including the castor bean tick, exhibit a marked sensitivity to both temperature and humidity [84]. This sensitivity arises from the critical influence these environmental factors exert on physiological processes and desiccation rates, both of which are essential for tick survival and life cycle completion [33]. In our study, SHAP analysis suggests that rainfalls during the warmest quarter (bio18 in Figure 5a) have a particularly strong influence on predictions for Ukraine. Checking the correlation between bio18 and the ensemble SDM for the current climate demonstrated a strong relationship too, which could be explained by the fact that ‘precipitation during the warmest quarter’ acts as a trade-off factor for other environmental variables, including temperature and humidity. Overall, the optimal bioclimatic environment for the tick species seems to be an intricate balance of moderate temperatures, high humidity, and sufficient rainfall to maintain that humidity, despite the fact that the correlation between humidity and precipitation is not straightforward [82,83].
In addition to climate, the distribution of I. ricinus may also be affected by factors such as land cover, particularly vegetation. Ticks require moisture and humidity to survive and dense vegetation provides these conditions by shading the ground and preventing moisture loss [83]. Additionally, ticks in some stages of the tick lifecycle feed on the blood of various animals that frequent forested and other vegetated areas [12]. Because in our modeling we focused on bioclimatic predictors, strictly speaking, vegetation features were not used; however, we tested the connection between bio18 and NDVI (Normalized Difference Vegetation Index). NDVI is a measure of vegetation health and density. Corresponding rasters for the warmest quarter (June–August) were downloaded from the EDIT geoplatform (http://edit.csic.es/GISdownloads.html; accessed on 4 April 2023). Previous studies have shown that areas with higher NDVI values tend to be more suitable habitats for I. ricinus ticks [81]. Indeed, we found fairly strong correlations for summer months (June, July, and August). In our opinion, this once again highlights the importance of the ‘precipitation of warmest quarter’ factor and/or, associated with it, factors in shaping the tick’s ecological niche and distribution, specifically in Ukraine.
In fact, for Latvia we have an analogous situation, where the most important features in influencing predictions within the country are suggested to be bio29 (‘highest weekly moisture index’) and bio24 (‘radiation of wettest quarter’) (Figure 5b). Both characterize the moisture regime, with a direct or indirect emphasis on the warm time of the year. The study area experiences its highest weekly moisture index in July. Secondly, the wettest quarter occurs from July to September, with peak rainfall in July.
As already implied, there may be another approach, tentatively tagged “biological”, to finding a model that would meet our objective aiming to aid health officials in managing and monitoring tick and tick-borne disease distributions. We assume this can be accomplished by correlating the obtained habitat suitability values of each SDM with the reported incidence of Lyme disease and hypothesize that a closer correlation would point towards a model that fits our needs the most. Existing studies in Europe mainly focus on acarological risk assessment, with very limited investigations exploring Lyme disease occurrence in humans [84]. As an example, incidences of Lyme disease in Ukraine are scarcely reported, so proper statistical analyses cannot be achieved.
Therefore, we used data from Romania since this neighboring country presents a fairly similar climate compared to Ukraine, especially if accounting for warm-season temperatures when tick activity is the highest [85]. Annual incidences of Lyme disease were downloaded for 41 counties and the Bucharest Municipality from the Johns Hopkins Lyme and Tick-Borne Disease Dashboard for a 10-year period (from 2012 to 2021) and averaged. Corresponding averaged habitat suitability values for each administrative entity, represented by a polygon shapefile, were obtained using the ‘grid statistics for polygons’ module in SAGA GIS. Firstly, using the R package ‘trafo’, suitable transformations depending on statistical requirements and the data being analyzed were assessed by checking assumptions of normality, homoscedasticity, and linearity. By log-transforming the annual incidences of Lyme disease using the formula ln(x + 1), all specified assumptions were shown to be acceptable (p > 0.05). Secondly, using log-transformed annual incidences of Lyme disease and the averaged habitat suitability values regarding the considered administrative units in Romania, corrected Pearson’s correlation coefficient (r) controlling for spatial non-independence were calculated. Results concerning each of the SDMs reveal the following: for data extracted from the SDM based on the WorldClim v.2 subset, r = 0.3815 and p-value = 0.0584, meaning no confirmation of statistical significance; in the case of the SDM based on the ENVIREM subset, r = 0.3989 and the p-value is 0.0422; and in the case of the SDM based on the CliMond v.1.2 subset, r = 0.59150 and the p-value is 0.0015. For the last two subsequent cases, conclusions regarding correlations are of sound statistical significance (p < 0.05); however, a stronger correlation coefficient was found for the SDM built on the CliMond v.1.2 subset. Therefore, we can conclude that this particular model is the best fit for our purpose.
In this study, besides the investigation of the present distribution of the tick I. ricinus, we aimed to assess how future climatic changes predicted for 2030 and 2050 could affect tick distribution in order to be prepared for proactive management actions involving both the vector and spread of Lyme disease. For this, we employed the CliMond v.1.2 dataset projection of future climate for 2030 and 2050 generated for the A1B scenario for emissions of greenhouse gases and sulfate aerosols (Figure 4, Figure 5 and Figure S4 [39]), which envisages a balanced emphasis on all energy sources. The produced ensemble models using the ‘flexsdm’ algorithm performed fairly well (Table 2), with acceptable AUC values, consistent TSS, and a high value of the continuous Boyce index. Somewhat poorer was the performance of the 2050 model, but it was within acceptable limits.
In terms of the current bioclimate, such as for Latvia, areas of high habitat suitability for I. ricinus occur on 22% of the territory of the country, mainly in the west and, to a lesser extent, in the north (Figure 1b). In Ukraine, areas of high habitat suitability for I. ricinus occupy around 20% of the country, with locations found predominantly in the west (Figure 1a). Other smaller patches of habitat suitability were found around Kyiv and in the south of Crimea, where ticks are the most abundant representative of the family [86]. Russia’s military action in Ukraine since 2022 has resulted in the internal displacement of over 3.5 million people (The UN Refugee Agency; https://www.unhcr.org/emergencies/ukraine-emergency; accessed on 4 April 2023), with potential health consequences [87,88,89,90]. Moreover, while refugee movements are associated with increased risks of infectious disease transmission and are likely to affect zoonotic disease risks, it remains unclear how forced migrations affect disease dynamics [91]. Human susceptibility to disease during forced migration might increase due to exhaustion, malnutrition, and stress arising from displacement, magnified by crowded and substandard living conditions [71,72].
In Ukraine, by 2030, our models predict a reduction in surface area of around two-fold of high habitat suitability for I. ricinus, decreasing to <10% of the country area (Figure 2a). These are expected to be located exclusively in the west of Ukraine, whereas patches of habitat suitability around Kyiv and in the south of Crimea are predicted to most likely disappear. By 2050, highly suitable habitats for I. ricinus are predicted to continue their contraction (to only 4.5% of the country area), with potential reduction to Lvivska and Ivano-Frankivska oblasts, despite a patch of high habitat suitability that is predicted to appear in the southern region of Odesa (Figure 3a).
On the contrary, in Latvia, areas of high habitat suitability for I. ricinus are predicted to expand and occupy up to 27% of the territory of the country by 2030, and up to 33% in 2050, engulfing the entire western region of Kurzeme (Figure 2b).
Upon analyzing the models, we can predict that in the north and south of its range, the potential distribution of I. ricinus may take different directions: for instance, in Ukraine, in the south of its range, predicted increased fragmentation may lead to a significant decrease in tick populations. On the contrary, in the northern part of Ukraine, the models predict an increase in the potential suitable territory for the spread of the species.
In this context, it is important to consider changes in the distribution of host species, which can significantly affect the spread of ticks such as I. ricinus. These ticks feed on the blood of various animals, including rodents, birds, foxes, hares, and larger animals such as deer, moose, humans, and domestic animals. Changes in the distribution of host species, especially in the context of the large-scale migration of people and domestic animals associated with military actions, can be a significant factor influencing the spread of ticks. This is an important topic for future research, especially for those dealing with the medical aspects of tick distribution and related diseases. Furthermore, understanding the ecological niche of I. ricinus is crucial for predicting its distribution under current and future climatic conditions. The ecological niche encompasses the range of environmental conditions that support the survival and reproduction of the species. Climate change scenarios, such as those used in our study, help model potential changes in these niches, allowing us to anticipate shifts in tick distribution. Species distribution modeling (SDM) primarily focuses on the theoretical or potential niche, also known as the fundamental niche [40]. This is because SDM typically uses environmental and species presence data to predict where a species could potentially exist in the absence of constraints such as competition or predation.
Our results regarding the prospects of bioclimatic suitability for the castor bean tick in Ukraine and Latvia are consistent with findings from previous studies. Indeed, previously published SDMs have predicted that Ixodes ticks will increase their distribution and abundance under global warming, with the most generally envisaged scenario being an expansion of the geographical range for ticks, including I. ricinus, with potential increases in abundance in some regions [92,93,94]. Bioclimatic models also predict a shift in the distribution of the tick species, particularly at higher altitudes and latitudes, like Scandinavia, the Baltic states, and Belarus. Conversely, tick populations are expected to decline in areas like the Alps, central and western Italy, and northwestern Poland. This means an increase in tick presence for Scandinavian countries (Sweden, Norway, and Finland), the Baltic states (Estonia, Latvia, and Lithuania), Denmark, and Belarus [95,96,97].

5. Conclusions

Global warming is one of the most important factors in the redistribution of tick populations and the outbreak of tick-borne diseases, and this is supported by recent studies [96,97,98]. SDMs built for estimating such redistribution using various selected subsets of predictors from the WorldClim v.2, ENVIREM, and CliMond v.1.2 datasets performed satisfyingly. In particular, the continuous Boyce index shows a high uniformity of performance of the obtained SDMs. In another approach correlating the obtained habitat suitability values of each SDM with the reported incidence of Lyme disease, we found a stronger correlation coefficient for the SDM built on the CliMond v.1.2 subset; therefore, we concluded that this particular model is the best fit for our purpose.
SDMs employing the CliMond v.1.2 dataset projection for future climates (2030 and 2050) generated for the A1B scenario for emissions of greenhouse gases and sulfate aerosols performed fairly well; somewhat poorer was the performance of the 2050 model, but this was within acceptable limits. In terms of the current bioclimate, areas of high habitat suitability for I. ricinus in Ukraine are predicted to occupy around one-fifth of the country, with locations found predominantly in the west, thus highlighting here the elevated risk of disease transmission. Regarding Latvia, areas of high habitat suitability for the tick are predicted to occur mainly in the west and, to a lesser extent, in the north of the country. Modeling suggests a significant reduction in highly suitable habitats for I. ricinus within Ukraine over the coming decades, whereas for Latvia the expansion of such suitable habitats is expected to continue, aggravating future health concerns accordingly.
We report that the relative importance of predictors in our case study was successfully measured with mean absolute Shapley values, underlining the pivotal significance of moisture conditions in the warm months of the year for creating highly suitable habitats for the tick species. Therefore, we believe that the use of multiple types of climatic factors can significantly improve the quality of predictions.

Supplementary Materials

The following supporting information can be downloaded at https://www.mdpi.com/article/10.3390/cli12110184/s1. Table S1: Bioclimatic variables (35) from the WorldClim v.2 (http://www.worldclim.org/; accessed on 21 November 2022; [42]) and the CliMond (https://www.climond.org/ at a resolution of 10′ (accessed on 28 November 2022; [44,45]), Table S2: Environmental variables (18) from ENVIREM (https://envirem.github.io/; accessed on 26 November 2022; [43])) datasets. Figure S1: Graphs of the most significant bioclimatic variables (CliMond) from different parts of the tick range: a—Latvia; b—Ukraine. Figure S2: Results of niche space exploration (CliMond): 3D ecological envelope: a—Latvia; b—Ukraine. Figure S3: Results of the analysis of SDM Maxent for I. ricinus in Europe: (A) current; (B) 2050 (2041–2060) (green-yellow color > 0.5, WorldClim dataset). The average training AUC for the replicate runs is 0.824, and the standard deviation is 0.001. Figure S4: Changes in key bioclimatic variables important for the spread of ticks in Ukraine and Latvia (CliMond, Figure 5 and Figure S1).

Author Contributions

Conceptualization, V.T., I.K., M.P. and O.N.; methodology, V.T. and O.N.; software, V.T. and O.N.; validation, V.T. and O.N.; formal analysis, V.T. and O.N.; investigation, V.T., O.N., M.P. and A.Č.; resources, V.T., I.K., O.N. and M.P.; data curation, V.T., I.K., O.N., M.P. and A.Š.; writing—original draft preparation, V.T. and O.N.; writing—review and editing, V.T., O.N., I.K., M.P., A.Č., J.-Y.G. and A.Š.; visualization, V.T.; supervision, V.T., O.N., I.K., M.P., A.Š. and J.-Y.G.; project administration, V.T., O.N., M.P., J.-Y.G. and A.Š. All authors have read and agreed to the published version of the manuscript.

Funding

This research was funded through the Emys-R project (https://emysr.cnrs.fr), funded through the 2020–2021 Biodiversa+ and Water JPI joint call for research projects, under the BiodivRestore ERA-NET Cofund (GA N°101003777), with the EU and the funding organizations Agence Nationale de la Recherche (ANR, France, grant ANR-21-BIRE-0005), Bundesministerium für Bildung und Forschung (BMBF, Germany, grant 16LW015), the State Education Development Agency (VIAA, Latvia, grant ES RTD/2022/2), and the National Science Center (NSC, Poland, grant 2021/03/Y/NZ8/00101), and by the project “Ecological and socioeconomic thresholds as a basis for defining adaptive management triggers in Latvian pond aquaculture” (lzp-2021/1-0247). O.N. was supported by Collège de France and Agence Nationale de la Recherche (ANR) through the PAUSE ANR Ukraine program (grant ANR-23-PAUK-0074).

Data Availability Statement

Data derived from public domain resources: 1. Species presence records Ixodes ricinus (Linnaeus, 1758): GBIF.org (22 May 2024) GBIF Occurrence Download, https://doi.org/10.15468/dl.bf6zjp (accessed on 22 May 2024). 2. Environmental data: The WorldClim v.2 website (https://www.worldclim.org/; accessed on 21 November 2022). The ENVIREM dataset was downloaded from http://envirem.github.io (accessed on 26 November 2022). The CliMond v.1.2 (current, 2030, 2050) datasets were downloaded from https://www.climond.org/ at 10′ resolution (accessed on 28 November 2022). 3. Software: Maxent, https://biodiversityinformatics.amnh.org/open_source/maxent/ (accessed on 28 November 2022); the R package ‘rgbif’, https://github.com/ropensci/rgbif (accessed on 28 November 2022); the R package ‘flexsdm’, https://sjevelazco.github.io/flexsdm/ (accessed on 28 November 2022); https://besjournals.onlinelibrary.wiley.com/doi/full/10.1111/2041-210X.13874 (accessed on 28 November 2022); the R package “virtualspecies”, https://borisleroy.com/virtualspecies/ (accessed on 28 November 2022); https://cran.r-project.org/web/packages/trafo/index.html (accessed on 28 November 2022); the R package ‘shap-values’, https://github.com/pablo14/ (accessed on 28 November 2022); the R-package ‘trafo’, https://github.com/cran/trafo (accessed on 28 November 2022); the R-package ‘SpatialPack’, http://spatialpack.mat.utfsm.cl/ (accessed on 28 November 2022); the PAST software (v.4.04; accessed on 28 November 2022) package and/or the R environment [72,73,74], https://www.nhm.uio.no/english/research/resources/past/ (accessed on 28 November 2022); the Saga GIS software, https://saga-gis.sourceforge.io/en/index.html (v.9.3.0; accessed on 28 November 2022).

Acknowledgments

We are grateful to project 16-00-F02201-000002 for the use of the mobile complex of scientific laboratories for research purposes.

Conflicts of Interest

The authors declare no conflicts of interest.

References

  1. Román-Palacios, C.; Wiens, J.J. Recent responses to climate change reveal the drivers of species extinction and survival. Proc. Natl. Acad. Sci. USA 2020, 117, 4211–4217. [Google Scholar] [CrossRef] [PubMed]
  2. Weiskopf, S.R.; Rubenstein, M.A.; Crozier, L.G.; Gaichas, S.; Griffis, R.; Halofsky, J.E.; Hyde, K.J.; Morelli, T.L.; Morisette, J.T.; Muñoz, R.C.; et al. Climate change effects on biodiversity, ecosystems, ecosystem services, and natural resource management in the United States. Sci. Total. Environ. 2020, 733, 137782. [Google Scholar] [CrossRef] [PubMed]
  3. Parry, M.L.; Canziani, O.F.; Palutikof, J.P.; van der Linden, P.J.; Hanson, C.E. Contribution of Working Group II to the Fourth Assessment Report of the Intergovernmental; Cambridge University Press: Cambridge, UK; New York, NY, USA, 2007. [Google Scholar]
  4. Tytar, V.; Nekrasova, O.; Pupins, M.; Skute, A.; Kirjušina, M.; Gravele, E.; Mezaraupe, L.; Marushchak, O.; Čeirāns, A.; Kozynenko, I.; et al. Modeling the Distribution of the Chytrid Fungus Batrachochytrium dendrobatidis with Special Reference to Ukraine. J. Fungi 2023, 9, 607. [Google Scholar] [CrossRef] [PubMed]
  5. Cumming, G. Comparing climate and vegetation as limiting factors for species ranges of African ticks. Ecology 2002, 83, 255–268. [Google Scholar] [CrossRef]
  6. Sonenshine, D.; Roe, R.M. Biology of Ticks, 2nd ed.; Oxford University Press: Oxford, UK, 2013; 496p. [Google Scholar]
  7. Čeirāns, A.; Pupins, M.; Kirjusina, M.; Gravele, E.; Mezaraupe, L.; Nekrasova, O.; Tytar, V.; Marushchak, O.; Garkajs, A.; Petrov, I.; et al. Top-down and bottom-up effects and relationships with local environmental factors in the water frog–helminth systems in Latvia. Sci. Rep. 2023, 13, 8621. [Google Scholar] [CrossRef]
  8. Wormser, G.P.; Dattwyler, R.J.; Shapiro, E.D.; Halperin, J.J.; Steere, A.C.; Klempner, M.S.; Krause, P.J.; Bakken, J.S.; Strle, F.; Stanek, G.; et al. The clinical assessment, treatment, and prevention of Lyme disease, human granulocytic anaplasmosis, and babesiosis: Clinical practice guidelines by the Infectious Diseases Society of America. Clin. Infect. Dis. 2006, 43, 1089–1134. [Google Scholar] [CrossRef]
  9. Derdáková, M.; Lencakova, D. Association of genetic variability within the Borrelia burgdorferi sensu lato with the ecology, epidemiology of Lyme borreliosis in Europe. Ann. Agric. Environ. Med. 2005, 12, 165. [Google Scholar]
  10. Dong, Y.; Zhou, G.; Cao, W.; Xu, X.; Zhang, Y.; Ji, Z.; Yang, J.; Chen, J.; Liu, M.; Fan, X.; et al. Global seroprevalence and sociodemographic characteristics of Borrelia burgdorferi sensu lato in human populations: A systematic review and meta-analysis. BMJ Glob. Health 2022, 7, e007744. [Google Scholar] [CrossRef]
  11. Gray, J.S.; Dautel, H.; Estrada-Peña, A.; Kahl, O.; Lindgren, E. Effects of climate change on ticks and tick-borne diseases in Europe. Interdiscip. Perspect. Infect. Dis. 2009, 593232. [Google Scholar] [CrossRef]
  12. Estrada-Peña, A.; de la Fuente, J. The ecology of ticks and epidemiology of tick-borne viral diseases. Antivir. Res. 2014, 108, 104–128. [Google Scholar] [CrossRef]
  13. ECDC. Ixodes Ricinus—Factsheet for Experts. 2020. Available online: https://www.ecdc.europa.eu/en/disease-vectors/facts/tick-factsheets/ixodesricinus (accessed on 23 November 2022).
  14. Intergovernmental Panel on Climate Change (IPCC). Summary for Policymakers. In Global Warming of 1.5 °C: IPCC Special Report on Impacts of Global Warming of 1.5°C above Pre-Industrial Levels in Context of Strengthening Response to Climate Change, Sustainable Development, and Efforts to Eradicate Poverty; Cambridge University Press: Cambridge, UK, 2022; pp. 1–24. [Google Scholar] [CrossRef]
  15. Jore, S.; Vanwambeke, S.O.; Viljugrein, H.; Isaksen, K.; Kristoffersen, A.B.; Woldehiwet, Z.; Johansen, B.; Brun, E.; Brun-Hansen, H.; Westermann, S.; et al. Climate and environmental change drives Ixodes ricinus geographical expansion at the northern range margin. Parasites Vectors 2014, 7, 11. [Google Scholar] [CrossRef] [PubMed]
  16. Lindgren, E.; Talleklint, L.; Polfeldt, T. Impact of climatic change on the northern latitude limit and population density of the disease-transmitting European tick Ixodes ricinus. Environ. Health Perspect. 2000, 108, 119. [Google Scholar] [CrossRef] [PubMed]
  17. Tomkins, J.L.; Aungier, J.; Hazel, W.; Gilbert, L. Towards an evolutionary understanding of questing behaviour in the tick Ixodes ricinus. PLoS ONE 2014, 9, e110028. [Google Scholar] [CrossRef] [PubMed]
  18. Randolph, S.E.; Storey, K. Impact of microclimate on immature tick-rodent host interactions (Acari: Ixodidae): Implications for parasite transmission. J. Med. Èntomol. 1999, 36, 741–748. [Google Scholar] [CrossRef]
  19. Porretta, D.; Mastrantonio, V.; Amendolia, S.; Gaiarsa, S.; Epis, S.; Genchi, C.; Bandi, C.; Otranto, D.; Urbanelli, S. Effects of global changes on the climatic niche of the tick Ixodes ricinus inferred by species distribution modelling. Parasites Vectors 2013, 6, 271. [Google Scholar] [CrossRef]
  20. Gage, K.L.; Burkot, T.R.; Eisen, R.J.; Hayes, E.B. Climate and vector borne disease. Am. J. Prev. Med. 2008, 35, 436–450. [Google Scholar] [CrossRef]
  21. Levy, S. Ticking time bomb? Climate change and Ixodes scapularis. Environ. Health Perspect. 2014, 122, A168. [Google Scholar] [CrossRef]
  22. Curriero, F.C.; Wychgram, C.; Rebman, A.W.; Corrigan, A.E.; Kvit, A.; Shields, T.; Aucott, J.N. The Lyme and Tickborne Disease Dashboard: A map-based resource to promote public health awareness and research collaboration. PLoS ONE 2021, 16, e0260122. [Google Scholar] [CrossRef]
  23. Franklin, J. Mapping Species Distributions: Spatial Inference and Prediction; Cambridge University Press: Cambridge, UK, 2009; pp. 1–340. [Google Scholar]
  24. Patz, J.A.; Graczyk, T.K.; Geller, N.; Vittor, A.Y. Effects of environmental change on emerging parasitic diseases. Int. J. Parasitol. 2000, 30, 1395–1405. [Google Scholar] [CrossRef]
  25. Zhang, L.; Liu, S.; Sun, P.; Wang, T.; Wang, G.; Zhang, X.; Wang, L. Consensus forecasting of species distributions: The effects of niche model performance and niche properties. PLoS ONE 2015, 10, e0120056. [Google Scholar] [CrossRef]
  26. Zhang, L.; Ma, D.; Li, C.; Zhou, R.; Wang, J.; Liu, Q. Projecting the potential distribution areas of Ixodes scapularis (Acari: Ixodidae) driven by climate change. Biology 2022, 11, 107. [Google Scholar] [CrossRef] [PubMed]
  27. Ukrainian Crisis: Situational Analysis, 5 March 2024. Data Friendly Space. 2024. 49p. Available online: https://reliefweb.int/report/ukraine/ukrainian-crisis-situational-analysis-05-march-2024 (accessed on 5 March 2024).
  28. Ukrainian Refugee Crisis (2022–Present). Available online: https://w.wiki/A2Hi (accessed on 9 May 2024).
  29. Chamberlain, S.; Oldoni, D.; Barve, V.; Desmet, P.; Geffert, L.; Mcglinn, D.; Ram, K.; rOpenSci; Waller, J. rgbif: Interface to the Global Biodiversity Information Facility API. 2022. R Package Version 3.7.4. Available online: https://CRAN.R-project.org/package=rgbif (accessed on 9 March 2023).
  30. GBIF.org. GBIF Occurrence Download. Available online: https://www.gbif.org/occurrence/download/0031420-240506114902167 (accessed on 22 May 2024).
  31. Estrada-Peña, A. Distribution, Abundance, and Habitat Preferences of Ixodes ricinus (Acari: Ixodidae) in Northern Spain. J. Med. Èntomol. 2001, 38, 361–370. [Google Scholar] [CrossRef] [PubMed]
  32. Estrada-Pena, A.; Venzal, J.M.; Acedo, C.S. The tick Ixodes ricinus: Distribution and climate preferences in the western Palaearctic. Med. Vet. Èntomol. 2006, 20, 189–197. [Google Scholar] [CrossRef] [PubMed]
  33. Alkishe, A.A.; Peterson, A.T.; Samy, A.M. Climate change influences on the potential geographic distribution of the disease vector tick Ixodes ricinus. PLoS ONE 2017, 12, e0189092. [Google Scholar] [CrossRef] [PubMed]
  34. Noll, M.; Wall, R.; Makepeace, B.L.; Newbury, H.; Adaszek, L.; Bødker, R.; Estrada-Peña, A.; Guillot, J.; da Fonseca, I.P.; Probst, J.; et al. Predicting the distribution of Ixodes ricinus and Dermacentor reticulatus in Europe: A comparison of climate niche modelling approaches. Parasites Vectors 2023, 16, 384. [Google Scholar] [CrossRef]
  35. Yemchuk, Y.M. Fauna of Ukraine Ixodid Ticks; Academy of Sciences of the Ukrainian SSR Press: Kyiv, Ukraine, 1960; Volume 25, 168p. (In Ukrainian) [Google Scholar]
  36. Akimov, I.A.; Nebogatkin, I.V. Distribution of Ixodes ricinus (Arachnida, Ixodidae) in Ukraine in the context of tick hazard, and factors favoring its persistence in conditions of fast-going environmental change. Zoodiversity 2022, 56, 429–434. [Google Scholar] [CrossRef]
  37. Capligina, V.; Seleznova, M.; Akopjana, S.; Freimane, L.; Lazovska, M.; Krumins, R.; Kivrane, A.; Namina, A.; Aleinikova, D.; Kimsis, J.; et al. Large-scale countrywide screening for tick-borne pathogens in field-collected ticks in Latvia during 2017–2019. Parasites Vectors 2020, 13, 351. [Google Scholar] [CrossRef]
  38. Velazco, S.J.E.; Rose, M.B.; de Andrade, A.F.A.; Minoli, I.; Franklin, J. flexsdm: An r package for supporting a comprehensive and flexible species distribution modelling workflow. Methods Ecol. Evol. 2022, 13, 1661–1669. [Google Scholar] [CrossRef]
  39. Kriticos, D.J. Regional climate-matching to estimate current and future sources of biosecurity threats. Biol. Invasions 2012, 14, 1533–1544. [Google Scholar] [CrossRef]
  40. Nekrasova, O.D.; Tytar, V.M.; Kuybida, V.V. GIS Modeling of Climate Change Vulnerability of Amphibians and Reptiles in Ukraine. In NAS of Ukraine; Schmalhausen Institute of Zoology NAS: Kyiv, Ukraine, 2019; 204p, ISBN 978-966-02-8956-7. (In Ukrainian) [Google Scholar]
  41. Schrodt, F.; Santos, M.J.; Bailey, J.J.; Field, R. Challenges and opportunities for biogeography—What can we still learn from von Humboldt? J. Biogeogr. 2019, 46, 1631–1642. [Google Scholar] [CrossRef]
  42. Fick, S.E.; Hijmans, R.J. WorldClim 2: New 1-km spatial resolution climate surfaces for global land areas. Int. J. Climatol. 2017, 37, 4302–4315. [Google Scholar] [CrossRef]
  43. Title, P.O.; Bemmels, J.B. ENVIREM: An expanded set of bioclimatic and topographic variables increases flexibility and improves performance of ecological niche modeling. Ecography 2018, 41, 291–307. [Google Scholar] [CrossRef]
  44. Kriticos, D.J.; Webber, B.L.; Leriche, A.; Ota, N.; Macadam, I.; Bathols, J.; Scott, J.K. CliMond: Global high-resolution historical and future scenario climate surfaces for bioclimatic modelling. Methods Ecol. Evol. 2012, 3, 53–64. [Google Scholar] [CrossRef]
  45. Kriticos, D.J.; Jarošik, V.; Ota, N. Extending the suite of bioclim variables: A proposed registry system and case study using principal components analysis. Methods Ecol. Evol. 2014, 5, 956–960. [Google Scholar] [CrossRef]
  46. Braunisch, V.; Coppes, J.; Arlettaz, R.; Suchant, R.; Schmid, H.; Bollmann, K. Selecting from correlated climate variables: A major source of uncertainty for predicting species distributions under climate change. Ecography 2013, 36, 971–983. [Google Scholar] [CrossRef]
  47. Leroy, B.; Meynard, C.N.; Bellard, C.; Courchamp, F. virtualspecies, an R package to generate virtual species distributions. Ecography 2016, 39, 599–607. [Google Scholar] [CrossRef]
  48. Feng, X.; Park, D.S.; Liang, Y.; Pandey, R.; Papeş, M. Collinearity in ecological niche modeling: Confusions and challenges. Ecol. Evol. 2019, 9, 10365–10376. [Google Scholar] [CrossRef]
  49. Shrestha, N. Detecting multicollinearity in regression analysis. Am. J. Appl. Math. Stat. 2020, 8, 39–42. [Google Scholar] [CrossRef]
  50. Araujo, M.; New, M. Ensemble forecasting of species distributions. Trends Ecol. Evol. 2007, 22, 42–47. [Google Scholar] [CrossRef]
  51. Peterson, A.T. Uses and requirements of ecological niche models and related distributional models. Biodivers. Informatics 2006, 3. [Google Scholar] [CrossRef]
  52. Peterson, A.T.; Soberón, J.; Pearson, R.G.; Anderson, R.P.; Martínez-Meyer, E.; Nakamura, M.; Araújo, M.B. Ecological Niches and Geographic Distributions; Princeton University Press: Oxford, UK, 2011; 314p. [Google Scholar]
  53. Allouche, O.; Tsoar, A.; Kadmon, R. Assessing the accuracy of species distribution models: Prevalence, kappa and the true skill statistic (TSS). J. Appl. Ecol. 2006, 43, 1223–1232. [Google Scholar] [CrossRef]
  54. Lobo, J.M.; Jimenez-Valverde, A.; Real, R. AUC: A misleading measure of the performance of predictive distribution models. Glob. Ecol. Biogeogr. 2008, 17, 145–151. [Google Scholar] [CrossRef]
  55. Boyce, M.S.; Vernier, P.R.; Nielsen, S.E.; Schmiegelow, F.K. Evaluating resource selection functions. Ecol. Model. 2002, 157, 281–300. [Google Scholar] [CrossRef]
  56. Hirzel, A.H.; Le Lay, G.; Helfer, V.; Randin, C.; Guisan, A. Evaluating the ability of habitat suitability models to predict species presences. Ecol. Model. 2006, 199, 142–152. [Google Scholar] [CrossRef]
  57. Medina, L.; Kreutzmann, A.-K.; Rojas-Perilla, N.; Castro, P. The R Package trafo for transforming linear regression models. R J. 2019, 11, 99–123. [Google Scholar] [CrossRef]
  58. Osorio, F.; Vallejos, R.; Cuevas, F. SpatialPack: Computing the association between two spatial processes. arXiv 2016, arXiv:1611.05289. [Google Scholar]
  59. Jenks, G.F.; Caspall, F.C. Error on choroplethic maps: Defnition, measurement, reduction. Ann. Assoc. Am. Geogr. 1971, 61, 217–244. [Google Scholar] [CrossRef]
  60. Lundberg, S.M.; Lee, S.I. A unified approach to interpreting model predictions. In Proceedings of the 31st International Conference on Machine Learning 2017, Long Beach, CA, USA, 4–9 December 2017; pp. 3547–3555. [Google Scholar]
  61. Lundberg, S.M.; Nair, B.; Vavilala, M.S.; Horibe, M.; Eisses, M.J.; Adams, T.; Liston, D.E.; Low, D.K.-W.; Newman, S.-F.; Kim, J.; et al. Explainable machine-learning predictions for the prevention of hypoxaemia during surgery. Nat. Biomed. Eng. 2018, 2, 749–760. [Google Scholar] [CrossRef]
  62. Lin, K.; Gao, Y. Model interpretability of financial fraud detection by group SHAP. Exp. Syst. Appl. 2022, 210, 118354. [Google Scholar] [CrossRef]
  63. Bhattacharyay, S.; Milosevic, I.; Wilson, L.; Menon, D.K.; Stevens, R.D.; Steyerberg, E.W.; Nelson, D.W.; Ercole, A.; the CENTER-TBI investigators participants. The leap to ordinal: Detailed functional prognosis after traumatic brain injury with a flexible modelling approach. PLoS ONE 2022, 17, e0270973. [Google Scholar] [CrossRef]
  64. Song, L.; Estes, L. itsdm: Isolation forest-based presence-only species distribution modelling and explanation in r. Methods Ecol. Evol. 2023, 14, 831–840. [Google Scholar] [CrossRef]
  65. Cha, Y.; Shin, J.; Go, B.; Lee, D.-S.; Kim, Y.; Kim, T.; Park, Y.-S. An interpretable machine learning method for supporting ecosystem management: Application to species distribution models of freshwater macroinvertebrates. J. Environ. Manag. 2021, 291, 112719. [Google Scholar] [CrossRef] [PubMed]
  66. Ghafarian, F.; Wieland, R.; Lüttschwager, D.; Nendel, C. Application of extreme gradient boosting and Shapley Additive explanations to predict temperature regimes inside forests from standard open-field meteorological data. Environ. Model. Softw. 2022, 156, 105466. [Google Scholar] [CrossRef]
  67. Bourhis, Y.; Bell, J.R.; Shortall, C.R.; Kunin, W.E.; Milne, A.E. Explainable neural networks for trait-based multispecies distribution modelling—A case study with butterflies and moths. Methods Ecol. Evol. 2023, 14, 1531–1542. [Google Scholar] [CrossRef]
  68. Lehnen, S.E.; Lombardi, J.V. Climate envelope modeling for ocelot conservation planning: Peering inside the black box. Ecosphere 2023, 14, e4477. [Google Scholar] [CrossRef]
  69. Scavuzzo, C.M.; Scavuzzo, J.M.; Campero, M.N.; Anegagrie, M.; Aramendia, A.A.; Benito, A.; Periago, V. Feature importance: Opening a soil-transmitted helminth machine learning model via SHAP. Infect. Dis. Model. 2022, 7, 262–276. [Google Scholar] [CrossRef]
  70. Conrad, O.; Bechtel, B.; Bock, M.; Dietrich, H.; Fischer, E.; Gerlitz, L.; Wehberg, J.; Wichmann, V.; Böhner, J. System for Automated Geoscientific Analyses (SAGA) v. 2.1.4. Geosci. Model Dev. 2015, 8, 1991–2007. [Google Scholar] [CrossRef]
  71. Hammer, C.C.; Brainard, J.; Hunter, P.R. Risk factors and risk factor cascades for communicable disease outbreaks in complex humanitarian emergencies: A qualitative systematic review. BMJ Glob. Health 2018, 3, e000647. [Google Scholar] [CrossRef]
  72. Hammer, Ø.; Harper, D.A.T.; Ryan, P.D. PAST: Paleontological Statistics software package for education and data analysis. Palaeontol. Electron. 2001, 4, 9. Available online: http://palaeo-electronica.org/2001_1/past/issue1_01.htm (accessed on 28 November 2022).
  73. R Core Team 2020. R: A Language and Environment for Statistical Computing; R Foundation for Statistical Computing: Vienna, Austria, 2020; Available online: https://www.R-project.org/ (accessed on 28 November 2022).
  74. Cianfrani, C.; Le Lay, G.; Hirzel, A.H.; Loy, A. Do habitat suitability models reliably predict the recovery areas of threatened species? J. Appl. Ecol. 2010, 47, 421–430. [Google Scholar] [CrossRef]
  75. Manzoor, S.A.; Griffiths, G.; Lukac, M.; Manzoor, S.A.; Griffiths, G.; Lukac, M. Species distribution model transferability and model grain size—Finer may not always be better. Sci. Rep. 2018, 8, 7168. [Google Scholar] [CrossRef] [PubMed]
  76. Liu, C.; Berry, P.M.; Dawson, T.P.; Pearson, R.G. Selecting thresholds of occurrence in the prediction of species distributions. Ecography 2005, 28, 385–393. [Google Scholar] [CrossRef]
  77. Shabani, F.; Kumar, L.; Ahmadi, M. Assessing accuracy methods of species distribution models: AUC, specificity, sensitivity and the true skill statistic. Glob. J. Hum. -Soc. Sci. B Geogr. Geo-Sci. Environ. Sci. Disaster Manag. 2018, 18, 1–13. Available online: https://globaljournals.org/GJHSS_Volume18/2-Assessing-Accuracy-Methods.pdf (accessed on 28 November 2022).
  78. Macleod, J. Ixodes ricinus in relation to its physical environment: II. The factors governing survival and activity. Parasitology 1935, 27, 123–144. [Google Scholar] [CrossRef]
  79. Alasmari, S.; Wall, R. Metabolic rate and resource depletion in the tick Ixodes ricinus in response to temperature. Exp. Appl. Acarol. 2021, 83, 81–93. [Google Scholar] [CrossRef]
  80. Gethmann, J.; Hoffmann, B.; Kasbohm, E.; Süss, J.; Habedank, B.; Conraths, F.J.; Beer, M.; Klaus, C. Research paper on abiotic factors and their influence on Ixodes ricinus activity—Observations over a two-year period at several tick collection sites in Germany. Parasitol. Res. 2020, 119, 1455–1466. [Google Scholar] [CrossRef]
  81. Zając, Z.; Kulisz, J.; Bartosik, K.; Woźniak, A.; Dzierżak, M.; Khan, A. Environmental determinants of the occurrence and activity of Ixodes ricinus ticks and the prevalence of tick-borne diseases in eastern Poland. Sci. Rep. 2021, 11, 15472. [Google Scholar] [CrossRef]
  82. Wongnak, P.; Bord, S.; Jacquot, M.; Agoulon, A.; Beugnet, F.; Bournez, L.; Cèbe, N.; Chevalier, A.; Cosson, J.-F.; Dambrine, N.; et al. Meteorological and climatic variables predict the phenology of Ixodes ricinus nymph activity in France, accounting for habitat heterogeneity. Sci. Rep. 2022, 12, 7833. [Google Scholar] [CrossRef]
  83. Medlock, J.M.; Hansford, K.M.; Bormane, A.; Derdakova, M.; Estrada-Peña, A.; George, J.-C.; Golovljova, I.; Jaenson, T.G.T.; Jensen, J.-K.; Jensen, P.M.; et al. Driving forces for changes in geographical distribution of Ixodes ricinus ticks in Europe. Parasites Vectors 2013, 6, 1. [Google Scholar] [CrossRef]
  84. Fu, W.; Bonnet, C.; Septfons, A.; Figoni, J.; Durand, J.; Frey-Klett, P.; Rustand, D.; Jaulhac, B.; Métras, R. Spatial and seasonal determinants of Lyme borreliosis incidence in France, 2016 to 2021. Eurosurveillance 2023, 28, 2200581. [Google Scholar] [CrossRef]
  85. Qviller, L.; Grøva, L.; Viljugrein, H.; Klingen, I.; Mysterud, A. Temporal pattern of questing tick Ixodes ricinus density at differing elevations in the coastal region of western Norway. Parasites Vectors 2014, 7, 179. [Google Scholar] [CrossRef] [PubMed]
  86. Evstafiev, I.; Simferopol, U.C.S.-E.S. Results of a 30-years-long investigation of small mammals in Crimea. Part 3. Parasites and epizootiology. Proc. Thériol. Sch. 2017, 2017, 111–135. [Google Scholar] [CrossRef]
  87. Vasyliuk, O.V.; Nekrasova, O.D.; Shyriaieva, D.V.; Kolomytsev, G.O. A review of major impact factors of hostilities influencing biodiversity in the eastern Ukraine (modeled on selected animal species). Vestn. Zool. 2015, 49, 145–158. [Google Scholar] [CrossRef]
  88. Nekrasova, O.; Marushchak, O.; Redinov, K.; Pupins, M.; Čeirāns, A.; Skute, A.; Tytar, V.; Moysiyenko, I.; Theissinger, K.; Georges, J.-Y. Assessing ecocide impacts for developing a conservation strategy in Ukraine. In Proceedings of the World Biodiversity Forum, Davos, Switzerland, 16–21 June 2024; pp. 95–96. [Google Scholar]
  89. Marushchak, O.; Nekrasova, O.; Zinenko, O.; Drohvalenko, M.; Mykytynets, H.; Suriadna, N.; Kotserzhynska, I.; Kotserzhynska, S.; Brusentsova, N.; Kuzmenko, Y.; et al. Herpetofauna at the frontline: So many ways to die. Responsible Herpetoculture J. 2024, 14, 114–128. [Google Scholar]
  90. Moysiyenko, I.I.; Khodosovtsev, O.Y.; Vasyliuk, O.V.; Parkhomenko, V.V.; Rusin, M.Y.; Viter, S.H.; Kuzemko, A.A.; Drapalyuk, A.M.; Biatov, A.P.; Sadogurska, S.S.; et al. Consequences of the Russian terrorist attack on the Kakhovka HPP for wildlife. In Traditions of Reserve Management, Modern Problems of Conservation and Post-War Restoration of Protected Areas: Collection of Scientific Papers Based on the Materials of the All-Ukrainian Round Table Dedicated to the 160th Anniversary of Friedrich Falz-Fein; Shapoval, V.V., Ed.; “Hlyboki Balyky”, Balyko-Shchuchynka village, Ukraine; Series: “Conservation Biology in Ukraine”; Issue 32; Druk Art: Chernivtsi, Ukraine, 2023; pp. 158–163. (In Ukrainian) [Google Scholar]
  91. Tarnas, M.C.; Desai, A.N.; Lassmann, B.; Abbara, A. Increase in vector-borne disease reporting affecting humans and animals in Syria and neighboring countries after the onset of conflict: A ProMED analysis 2003–2018. Int. J. Infect. Dis. 2021, 102, 103–109. [Google Scholar] [CrossRef]
  92. Boeckmann, M.; Joyner, T.A. Old health risks in new places? An ecological niche model for I. ricinus tick distribution in Europe under a changing climate. Health Place 2014, 30, 70–77. [Google Scholar] [CrossRef]
  93. Li, S.; Gilbert, L.; Harrison, P.A.; Rounsevell, M.D.A. Modelling the seasonality of Lyme disease risk and the potential impacts of a warming climate within the heterogeneous landscapes of Scotland. J. R. Soc. Interface 2016, 13, 20160140. [Google Scholar] [CrossRef]
  94. Ogden, N.H.; Ben Beard, C.; Ginsberg, H.S.; Tsao, J.I. Possible effects of climate change on Ixodid ticks and the pathogens they transmit: Predictions and observations. J. Med. Èntomol. 2021, 58, 1536–1545. [Google Scholar] [CrossRef]
  95. Voyiatzaki, C.; Papailia, S.I.; Venetikou, M.S.; Pouris, J.; Tsoumani, M.E.; Papageorgiou, E.G. Climate changes exacerbate the spread of Ixodes ricinus and the occurrence of Lyme Borreliosis and tick-borne encephalitis in Europe-How climate models are used as a risk assessment approach for tick-borne diseases. Int. J. Environ. Res. Public Health 2022, 19, 6516. [Google Scholar] [CrossRef]
  96. Bisanzio, D.; Amore, G.; Ragagli, C.; Tomassone, L.; Bertolotti, L.; Mannelli, A. Temporal variations in the usefulness of normalized difference vegetation index as a predictor for Ixodes ricinus (Acari: Ixodidae) in a Borrelia lusitaniae focus in Tuscany, Central Italy. J. Med. Èntomol. 2008, 45, 547–555. [Google Scholar] [CrossRef]
  97. Bouchard, C.; Dibernardo, A.; Koffi, J.; Wood, H.; Leighton, P.A.; Lindsay, L.R. Increased risk of tick-borne diseases with climate and environmental changes. Can. Commun. Dis. Rep. 2019, 45, 83–89. [Google Scholar] [CrossRef] [PubMed]
  98. Tsoumani, M.E.; Papailia, S.I.; Papageorgiou, E.G.; Voyiatzaki, C. Climate change impacts on the prevalence of tick-borne diseases in Europe. Environ. Sci. Proc. 2023, 26, 18. [Google Scholar] [CrossRef]
Climate 12 00184 g001
Figure 1. Jenks natural breaks maps of habitat suitability (HS) for Ixodes ricinus under current bioclimatic conditions. (a) Ukraine: red, green, and blue, respectively, for HS of high (0.52–0.71), medium (0.29–0.52), and low (0.08–0.29) value (oblasts: 1—Zakarpatska, 2—Lvivska, 3—Volynska, 4—Rivnenska, 5—Ternopilska, 6—Ivano-Frankivska, 7—Zhytomyrska, 8—Khmelnytska, 9—Chernivetska, 10—Vinnytska, 11—Crimea). (b) Latvia: red, green, and blue, respectively, for HS of high (0.62–0.77), medium (0.46–0.62), and low (0.32–0.46) values (regions: 1—Kurzeme, 2—Pieriga, 3—Vidzeme). Administrative units according to OECD, 2022.
Figure 1. Jenks natural breaks maps of habitat suitability (HS) for Ixodes ricinus under current bioclimatic conditions. (a) Ukraine: red, green, and blue, respectively, for HS of high (0.52–0.71), medium (0.29–0.52), and low (0.08–0.29) value (oblasts: 1—Zakarpatska, 2—Lvivska, 3—Volynska, 4—Rivnenska, 5—Ternopilska, 6—Ivano-Frankivska, 7—Zhytomyrska, 8—Khmelnytska, 9—Chernivetska, 10—Vinnytska, 11—Crimea). (b) Latvia: red, green, and blue, respectively, for HS of high (0.62–0.77), medium (0.46–0.62), and low (0.32–0.46) values (regions: 1—Kurzeme, 2—Pieriga, 3—Vidzeme). Administrative units according to OECD, 2022.
Climate 12 00184 g001
Climate 12 00184 g002
Figure 2. Jenks natural breaks maps of habitat suitability (HS) for I. ricinus as predicted by our models by 2030. (a) Ukraine: red, green, and blue, respectively, denote areas with high (0.61–0.77), medium (0.36–0.61), and low (0.1–0.36) HS. Oblasts: 1—Zakarpatska, 2—Lvivska, 3—Ternopilska, 4—Ivano-Frankivska, 5—Chernivetska. (b) Latvia: red, green, and blue, respectively, denote areas with high (0.62–0.76), medium (0.48–0.62), and low (0.34–0.48) HS. Regions: 1—Kurzeme, 2—Pieriga, 3—Vidzeme.
Figure 2. Jenks natural breaks maps of habitat suitability (HS) for I. ricinus as predicted by our models by 2030. (a) Ukraine: red, green, and blue, respectively, denote areas with high (0.61–0.77), medium (0.36–0.61), and low (0.1–0.36) HS. Oblasts: 1—Zakarpatska, 2—Lvivska, 3—Ternopilska, 4—Ivano-Frankivska, 5—Chernivetska. (b) Latvia: red, green, and blue, respectively, denote areas with high (0.62–0.76), medium (0.48–0.62), and low (0.34–0.48) HS. Regions: 1—Kurzeme, 2—Pieriga, 3—Vidzeme.
Climate 12 00184 g002
Climate 12 00184 g003
Figure 3. Jenks natural breaks maps of habitat suitability (HS) for I. ricinus as predicted by our models by 2050. (a) Ukraine: red, green, and navy blue, respectively, denote areas with high (0.63–0.73), medium (0.45–0.63), and low (0.33–0.45) HS. Oblasts: 1—Lvivska, 2—Ivano-Frankivska, 3—Odeska. (b) Latvia: red, green, and navy blue, respectively, denote areas with high (0.69–0.8), medium (0.51–0.69), and low (0.33–0.51) HS. Regions: 1—Kurzeme, 2—Pieriga, 3—Vidzeme, 4—Zemgale.
Figure 3. Jenks natural breaks maps of habitat suitability (HS) for I. ricinus as predicted by our models by 2050. (a) Ukraine: red, green, and navy blue, respectively, denote areas with high (0.63–0.73), medium (0.45–0.63), and low (0.33–0.45) HS. Oblasts: 1—Lvivska, 2—Ivano-Frankivska, 3—Odeska. (b) Latvia: red, green, and navy blue, respectively, denote areas with high (0.69–0.8), medium (0.51–0.69), and low (0.33–0.51) HS. Regions: 1—Kurzeme, 2—Pieriga, 3—Vidzeme, 4—Zemgale.
Climate 12 00184 g003
Climate 12 00184 g004
Figure 4. Results of the analysis of SDM Maxent for I. ricinus in Europe (current, CliMond). The average training AUC for the replicate runs is 0.849, and the standard deviation is 0.001.
Figure 4. Results of the analysis of SDM Maxent for I. ricinus in Europe (current, CliMond). The average training AUC for the replicate runs is 0.849, and the standard deviation is 0.001.
Climate 12 00184 g004
Climate 12 00184 g005
Figure 5. Absolute summary plots of the datasets for Ukraine (a) and Latvia (b), where the average absolute value of the SHAP values for each variable is taken in order to obtain a bar chart as a function of the contribution of each variable to the prediction of the CliMond current model. The variables are ordered from most (top) to least (bottom) important. The y-axis represents the variables used in the study, which refer to Bio6—minimum temperature of coldest week; bio7—annual temperature range; bio9—mean temperature of the driest quarter; bio15—precipitation seasonality; bio18—precipitation of the warmest quarter; bio24—radiation of the wettest quarter; bio29—highest weekly moisture index. The x-axis represents the corresponding SHAP value.
Figure 5. Absolute summary plots of the datasets for Ukraine (a) and Latvia (b), where the average absolute value of the SHAP values for each variable is taken in order to obtain a bar chart as a function of the contribution of each variable to the prediction of the CliMond current model. The variables are ordered from most (top) to least (bottom) important. The y-axis represents the variables used in the study, which refer to Bio6—minimum temperature of coldest week; bio7—annual temperature range; bio9—mean temperature of the driest quarter; bio15—precipitation seasonality; bio18—precipitation of the warmest quarter; bio24—radiation of the wettest quarter; bio29—highest weekly moisture index. The x-axis represents the corresponding SHAP value.
Climate 12 00184 g005
Table 1. Results of the selection of uncorrelated predictors (Spearman’s r < 0.8).
Table 1. Results of the selection of uncorrelated predictors (Spearman’s r < 0.8).
Predictor DatasetSelected Subsets of Predictors
WorldClim v.2Annual mean temperature (bio1), Mean diurnal range (bio2), Temperature seasonality (bio4), Maximum temperature of warmest month (bio5), Annual precipitation (bio12), Precipitation seasonality (bio15), Precipitation of driest quarter (bio17)
ENVIREMAnnual PET *, Continentality, Emberger’s pluviothermic quotient, mean monthly PET of coldest quarter, mean monthly PET of driest quarter, PET seasonality, mean monthly PET of wettest quarter, terrain roughness index
CliMond v.1.2Mean diurnal temperature range (bio2), Minimum temperature of coldest week (bio6), Annual temperature range (bio7), Mean temperature of driest quarter (bio9), Annual precipitation (bio12), Precipitation of driest week (bio14), Precipitation seasonality (bio15), Precipitation of warmest quarter (bio18), Radiation of wettest quarter (bio24), Radiation of warmest quarter (bio26), Radiation of coldest quarter (bio27), Highest weekly moisture index (bio29)
* PET = potential evapotranspiration.
Table 2. Evaluation metrics for SDMs built using various selected subsets of predictors.
Table 2. Evaluation metrics for SDMs built using various selected subsets of predictors.
Predictor Dataset/SubsetEvaluation Metrics
Area Under Curve (AUC)SD * of AUCTrue Skill Statistic (TSS)SD of TSSContinuous Boyce Index (BOYCE)SD of BOYCE
WorldClim v.20.840.010.590.020.910.07
ENVIREM0.840.010.560.020.920.06
CliMond v.1.2 current0.790.060.500.120.910.07
CliMond v.1.2 A1B scenario for 20300.760.060.460.130.870.18
CliMond v.1.2 A1B scenario for 20500.740.080.410.120.690.16
* SD = standard deviation.
Table 3. Analysis of variable contributions in the SDM Maxent model for I. ricinus in Europe (CliMond Dataset, see Table S1).
Table 3. Analysis of variable contributions in the SDM Maxent model for I. ricinus in Europe (CliMond Dataset, see Table S1).
VariablePercent Contribution %Permutation Importance %
bio7 Annual temperature range (Bio05-Bio06) (°C)34.621.8
bio14 Precipitation of driest week (mm)15.62.1
bio24 Radiation of wettest quarter (W m−2)9.18.6
bio1 Annual mean temperature (°C)7.71.7
Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Share and Cite

MDPI and ACS Style

Tytar, V.; Kozynenko, I.; Pupins, M.; Škute, A.; Čeirāns, A.; Georges, J.-Y.; Nekrasova, O. Species Distribution Modeling of Ixodes ricinus (Linnaeus, 1758) Under Current and Future Climates, with a Special Focus on Latvia and Ukraine. Climate 2024, 12, 184. https://doi.org/10.3390/cli12110184

AMA Style

Tytar V, Kozynenko I, Pupins M, Škute A, Čeirāns A, Georges J-Y, Nekrasova O. Species Distribution Modeling of Ixodes ricinus (Linnaeus, 1758) Under Current and Future Climates, with a Special Focus on Latvia and Ukraine. Climate. 2024; 12(11):184. https://doi.org/10.3390/cli12110184

Chicago/Turabian Style

Tytar, Volodymyr, Iryna Kozynenko, Mihails Pupins, Arturs Škute, Andris Čeirāns, Jean-Yves Georges, and Oksana Nekrasova. 2024. "Species Distribution Modeling of Ixodes ricinus (Linnaeus, 1758) Under Current and Future Climates, with a Special Focus on Latvia and Ukraine" Climate 12, no. 11: 184. https://doi.org/10.3390/cli12110184

APA Style

Tytar, V., Kozynenko, I., Pupins, M., Škute, A., Čeirāns, A., Georges, J. -Y., & Nekrasova, O. (2024). Species Distribution Modeling of Ixodes ricinus (Linnaeus, 1758) Under Current and Future Climates, with a Special Focus on Latvia and Ukraine. Climate, 12(11), 184. https://doi.org/10.3390/cli12110184

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop