Prediction of Sites with a High Probability of Wild Mammal Roadkill Using a Favourability Function

: Roads are one of the main causes of loss of biodiversity, with roadkill one of the main causes of mortality. The aim of this research was to identify sites with a high probability of roadkill of medium and large mammals, and the environmental variables that would explain it. We used the favourability function (F) to build the predictive models. There were 57 explanatory variables, and we collected 685 records of 10 species of medium and large native wild mammals from the ECOBIO Uruguay databases. They were grouped into native forest and grassland species, according to the main habitat. Two models were developed, one with all the variables and one with the anthropogenic variables. For both groups, the model obtained with all the variables was the most signiﬁcant according to the evaluation indices used. This made it possible to identify the hot spots of roadkill (F > 0.6) for each of the groups. The anthropic variables were the ones that best explained these hot spots. This allowed the identiﬁcation of sites where the probability of roadkill is high and requires a monitoring plan to implement mitigation measures in the future.


Introduction
Road systems generate changes in the landscape and in the territory, causing direct and indirect effects on the structure and composition of ecosystems. Some of these effects include, but are not limited to, animal roadkill, barrier effect, edge effect, dispersion of exotic species, changes in microclimates, contamination of water resources, and fragmentation of ecosystems [1][2][3][4][5][6][7]. As a consequence of the above, the loss of biodiversity has significantly increased, with vertebrate species being the most affected [6][7][8][9][10].
Of the vertebrates, large mammals are among the taxa most vulnerable to the impact of roads, since they generally present low values in both abundance and reproductive rate, and they also have large territories of occupation, travelling long distances during their life cycle [2].
Research related to the impacts of roads on biodiversity began in developed countries in the Northern Hemisphere (e.g., [9][10][11][12][13]). In Latin America, these studies are very few in number, and recently, Brazil was the first country where research on this subject was conducted [2,4,[14][15][16]. In Uruguay, in 2015, the Civil Association ECOBIO Uruguay, composed of researchers, began to research this problem, generating baseline information on the presence of roadkill [17,18]. A study conducted by the NGO Vida Silvestre in 1999, in which the presence of roadkill on a main road in Uruguay was surveyed [19,20], provided a precedent in this area.
To reduce the impact of road network infrastructure on wildlife, several mitigation measures have been implemented worldwide, including wildlife crossings and signage [6,7,10,16].
In order to implement such measures, it is important to accurately identify the sites where mitigation measures should be carried out. However, it is difficult to determine these precisely, since it is necessary to collect field information, which carries a high cost. One way to overcome this difficulty is to use mathematical models to identify those sites where conditions for wildlife roadkill are present [21,22].
In recent decades, there has been an increase in the use of mathematical algorithms to model various spatial processes [18,[21][22][23][24], including the favourability function [25][26][27]. For this, explanatory variables that can be grouped into climatic, anthropic, and geographical factors [28] related to the phenomenon to be modelled, such as roadkill, are used. In the case of the favourability function, it had not yet been used to identify and predict the sites where the conditions for roadkill occur.
Based on the aforementioned points, in this work we used favourability as a tool to identify the sites with favourable conditions for roadkill of medium and large mammals in the Eastern region of Uruguay. In addition, predictive models of roadkill were generated and the most important environmental variables that explain roadkill were determined.

Study Area
The present work was carried out in the Eastern region of Uruguay, covering the ecoregions of Sierras del Este and Graben de la Laguna Merín (departments of Maldonado, Rocha, Treinta y Tres, Lavalleja, and Cerro Largo, Uruguay) reaching an area of about 3,672,457 ha (9,074,839 acres) ( Figure 1) [29]. This region presents a wide heterogeneity of ecosystems (wetlands, natural forests, coastal zone, natural grasslands, and diverse crops), a great diversity of species, and concentrates the largest number of national protected areas and international reserves such as the Ramsar sites [30].
With respect to biodiversity, the Graben de la Laguna Merín ecoregion features a richness of 57 species of mammals, of which 8 are classified as indicator species of that ecoregion [29], whereas the Sierras del Este ecoregion features a richness of 56 species of mammals, 3 of which are indicators of the same [29].
Diversity 2021, 13, x FOR PEER REVIEW 3 of 14 Figure 1. Study area in Uruguay, which includes the Sierras del Este (dark grey) and Graben de la Laguna Merín (light grey) ecoregions. We marked the roads used for modelling.

Collection of Roadkill Data
The surveys were completed on the main roads with high vehicular traffic (Annual average daily traffic of 134,487 vehicles) (Figure 1), where several routes connecting the main tourist areas of the country converge.

Collection of Roadkill Data
The surveys were completed on the main roads with high vehicular traffic (Annual average daily traffic of 134,487 vehicles) (Figure 1), where several routes connecting the main tourist areas of the country converge.
The study area covered an extension of ca. 1000 km in the Eastern region. Sampling campaigns were carried out every two months during one year (from July 2015 to June 2016). Each road was travelled in both directions by vehicle at a speed of 60 km/h (37 miles/h) [31] where all roadkill were georeferenced and photographically recorded. Once the animal was registered, it was removed from the road so as not to be recorded again in subsequent trips. Should the animal fail to be removed, it was marked with a coloured spray.
The data collected in the field were digitized and then spatialized using the open access Geographic Information System (GIS) QGIS 3.10 [32].
Then, the records obtained in the period 2017-2019 were used to validate the models generated. For this, records obtained through citizen science and the NGO ECOBIO Uruguay during that period were compiled and the abundance of roadkill per square kilometre was calculated.

Model Preparation
In order to prepare the models, we used a grid mesh made up of 177,708 cells measuring 1 km 2 , which covered the entire national territory. It was created by the Laboratory of Sustainable Development and Environmental Management of the Territory (Laboratorio de Desarrollo Sustentable y Gestión Ambiental del Territorio, LDSGAT, Montevideo, Uruguay) of the School of Sciences and included 55 environmental variables (Supplementary Materials Tables S1 and S2). Before modelling, this grid was cropped for the study area using the QGIS 3.10 software.

Explanatory Variables
The explanatory variables used in the models comprise the 55 environmental variables elaborated by LDSGAT and two road-related variables: sites where watercourses are crossed by the roadway (Curs-cam) and Annual Average Daily Traffic (AADT) for 2017 obtained from the Ministry of Transportation and Public Works of Uruguay. The variables were grouped into climatic, geographic, and anthropic factors [33,34]. In total, 57 variables were used: 32 climatic, 14 geographic, and 11 anthropic (Supplementary Materials  Table S1).

Response Variable
The roadkill data obtained were used as a response variable. The data were spatialized with the QGIS program and classified as "presence" (1: if we have at least one record) and "absence" (0: where no record was observed) on the grid. Only grids that overlapped with the roads to be modelled were considered in the analysis (number of grids = 1372). Subsequently, the species were classified into two groups: grassland and native forest species. These groups were based on the most important habitats for each species according to the literature [20,35]. A total of six species were assigned to the "native forest" group (margay (Leopardus wiedii), Geoffroy's cat (Leopardus geoffroyi), neotropical otter (Lontra longicaudis), South American raccoon (Procyon cancrivorus), nine-banded armadillo (Dasyus novemcinctus), and forest fox (Cerdocyon thous)) and four "grassland" species (Andes skunk (Conepatus chinga), lesser grison (Galictis cuja), southern long-nosed armadillo (Dasypus hybridus) and grey pampean fox (Lycalopex gymnocercus)).
Once the groups were established, we then began to unify the records with the grid mesh. For this purpose, QGIS was used, where a union by location was made between both layers. After the union, the attribute table was exported to develop the predictive models.

Predictive Models
The favourability function (FF) was used to elaborate the predictive models [25]: where y is a linear regression equation (y = α + β 1 x 1 + β 2 x 2 +···+ β n x n ), n 1 is the number of observed presences, and n 0 is the number of absences. This function is a generalized linear model (GLM) with binomial distribution and logit link function.
To build the predictive models, a script was developed using the fuzzySim package [36]. The analysis was performed sequentially using the multGLM function. First, the false discovery rate or [37] was applied to determine the proportion of falsely rejected null hypotheses [38,39]. A second filter was applied to the resulting variables, which consisted of analysing the Pearson correlation between variables with the corSelect function, considering a threshold value of 0.8 to discard those variables that were correlated [36]. The script used is presented in the supplementary material.
With the variables that passed the two filters, a forward and backward stepwise selection method was applied. After this procedure, the modelTrim function was applied to eliminate any remaining nonsignificant variables [40,41].
A total of four predictive models were constructed, two models for each of the groups ("native forest" and "grassland"): one called the complete model, in which all the explanatory variables of the three factors were used (57 variables), and another, called the anthropic model, in which only the variables related to the anthropic factor were used (11 variables). Finally, the survey package was [42] used to calculate the Wald Index [43] and obtain the relative importance of each variable in the final model.
The QGIS program (QGIS, version 3.10) was used to represent the models obtained and to create the final maps. The favourability (F) values were classified into ten categories from 0 to 1, in intervals of 0.1, for each geographic unit (cells).

Model Evaluation
Each predictive model was evaluated for performance using various indices such as Area Under the Curve (AUC), Specificity, Sensitivity, Over Prediction (OPR) and Under Prediction (UPR), and True Skill Statistic (TSS) using the package {modEvA} [44].

Model Validation
To validate what was obtained in the models, we used the roadkill records obtained in the period 2017-2019 and calculated the density of roadkill per square kilometre. For this purpose, the roadkill abundances were overlapped with the favourability values, and the spatial relationship between both layers was analysed. For this, the spatial coincidence between roadkill abundance and favourability was analysed using QGIS.
All analyses were performed in the statistical program R (version 3.5.3) [45].

Records of Roadkill Presence
A total of 685 roadkill records were obtained for the most affected medium and large mammal species, based on the number of roadkill records and those that present some degree of threat at a national and international level, as in the case of the margay or the mulita (Table 1). Of the total number of records, 186 corresponded to forest species and 499 to grassland species. Among the forest species was the margay (Leopardus wiedii), which requires conservation at a national and international level. The grassland species group included the Southern long-nosed armadillo (Dasypus hybridus), which is listed as internationally threatened (Table 1 and Figure 2). degree of threat at a national and international level, as in the case of the margay or the mulita (Table 1). Of the total number of records, 186 corresponded to forest species and 499 to grassland species. Among the forest species was the margay (Leopardus wiedii), which requires conservation at a national and international level. The grassland species group included the Southern long-nosed armadillo (Dasypus hybridus), which is listed as internationally threatened (Table 1 and Figure 2).

Predictive Models
The models obtained, both complete and anthropic, for both groups (grassland and native forest species) showed that Route 9 was the one with the largest surface area (85% in both models for both groups) with high favourability values (F ≥ 0.5). This was followed by Route 8, where for grassland species, 60% of the surface area presented values of F ≥ 0.5 in both models. However, for native forest species, the surface area of the anthropic model was 43%, and for the complete model 37% (Figures 3 and 4). On the other hand, it was observed that Route 15 presented favourability values greater or equal to 0.5 in the 30 km section that connects the city of Rocha with La Paloma (Figures 3 and 4). native forest species) showed that Route 9 was the one with the largest surface area (85% in both models for both groups) with high favourability values (F ≥ 0.5). This was followed by Route 8, where for grassland species, 60% of the surface area presented values of F ≥ 0.5 in both models. However, for native forest species, the surface area of the anthropic model was 43%, and for the complete model 37% (Figures 3 and 4). On the other hand, it was observed that Route 15 presented favourability values greater or equal to 0.5 in the 30 km section that connects the city of Rocha with La Paloma (Figures 3 and 4).

Most Important Environmental Variables
The variable Annual Average Daily Traffic (AADT) was the most important explanatory variable for both the full and anthropic models for native forest species and for the anthropic model for grassland species. The full model for the grassland species

Most Important Environmental Variables
The variable Annual Average Daily Traffic (AADT) was the most important explanatory variable for both the full and anthropic models for native forest species and for the anthropic model for grassland species. The full model for the grassland species showed that the most important variable was Precipitation in the Wettest Quarter (BIO 16) followed by the AADT variable (Tables 2 and 3). However, AADT presented a Wald Index close to the first variable and thus presented a significant importance (Table 3). Analysing the two routes with the highest roadkill favourability for both groups, Route 9 had four sections with AADT values ranging from 1440 to 2146 vehicles per day. Route 8 showed three sections with values ranging from 1440 to 2851 vehicles per day ( Figure 5).

Forestation
− 5.28 Population density − 3.93 Analysing the two routes with the highest roadkill favourability for both groups, Route 9 had four sections with AADT values ranging from 1440 to 2146 vehicles per day. Route 8 showed three sections with values ranging from 1440 to 2851 vehicles per day ( Figure 5).

Model Evaluation
The prepared models of roadkill performed well based on the indices: AUC, CCR, TSS, UPR, UPR, OPR, and Sensitivity and Specificity (Table 4). In the case of native forest species, the complete model obtained higher AUC and TSS values than the anthropic model. However, the CCR and specificity and sensitivity values were higher in the anthropic model, whereas in the case of grassland species, most of the indices (AUC, CCR, Specificity, TSS, and UPR) showed higher values in the complete model (Table 4).

Model Validation
With the roadkill mammal data obtained between 2017 and 2019, a map of roadkill abundance per square kilometre was produced ( Figure 6). The section of Route 15 that connects the city of Rocha with La Paloma was the one with the highest mortality rates during this period. This section has an AADT ranging from 2146 to 2851 ( Figure 6).
When the cells of roadkill abundances for the period 2017-2019 were overlapped with the models that were elaborated, it was observed that for native forest species, 64% (complete model) and 52% (anthropic model) matched the sites with F-favourability values of ≥0.5, whereas for grassland species, 66% (complete model) and 68% (anthropic model) matched with F-favourability values of ≥0.5. This confirmed what was obtained in the models, since the spatial overlap of the cells was higher than 60% of the surface with F-favourability values of ≥0.5. When the cells of roadkill abundances for the period 2017-2019 were overlapped with the models that were elaborated, it was observed that for native forest species, 64% (complete model) and 52% (anthropic model) matched the sites with F-favourability values of ≥0.5, whereas for grassland species, 66% (complete model) and 68% (anthropic model) matched with F-favourability values of ≥0.5. This confirmed what was obtained in the models, since the spatial overlap of the cells was higher than 60% of the surface with F-favourability values of ≥0.5.

Discussion
In this paper, the favourability function was used to model and predict the spatiotemporal patterns of the roadkill effect for six species of medium and large mammals in native forests and four in natural grasslands, and to identify (favourable) hotspots for roadkill. This represents a key contribution to the subject of road ecology in Uruguay, since it allowed us to identify the places with the highest probability of wildlife roadkill and to determine which environmental variables influence roadkill. On the other hand, the results obtained are a significant input at the national level to advance the implementation of mitigation measures in the identified hotspots.
Using the favourability function, models were obtained that provided high predictive and explanatory power according to the applied evaluations. There are research studies that used various mathematical algorithms (e.g., generalized linear models, Random Forest, Support Vector Machine, Habitat Suitability) to explain and predict the effect of roadkill on wildlife [18,21,24,[46][47][48][49][50][51][52]. However, there is no research where the favourability function was applied to predict and explain the roadkill effect; moreover, it has usually been used as a tool for the identification of favourable habitats for biodiversity [23,27,50].
For both modelled groups (grassland and forest species) the (AADT) variable (Average Daily Annual Daily Transit) was the most important for three of the four models. It was observed that the higher the vehicle traffic, the higher the probability of roadkill. This is in agreement with other authors who also identified that vehicle density is directly related to roadkill, increasing wildlife mortality with increased vehicle density [22,53,54]. In this regard, it is important to take into account the variable of vehicular traffic when implementing a mitigation measure to reduce mammal roadkill.
Roads and urban centres negatively affect the species studied through road traffic and habitat fragmentation. [6,7,10,18,49,[55][56][57][58]. This study further confirms these findings, since the most important environmental variables for both groups of mammals studied were related to anthropic factors, mainly those related to infrastructures such as the AADT, bridges, and roads.
Human infrastructures directly relate to one another. Urban areas have a higher density of road and transit infrastructures, thus increasing the probability for wildlife roadkill [6,7,10,[56][57][58]. In contrast, as we move away from urban areas the probability decreases [55,58]. However, in the present work, with the models carried out, the urbanization variable decreased the probability of roadkill. This may be due to the effects of urbanization on species distribution [57].
In relation to crossings between watercourses and road infrastructure (bridges and sewers), this variable was found to be the second most important for native forest species in the complete model. This indicates that the greater the number of bridges, the greater the probability of roadkill. Several studies highlight the importance of bridges and sewers as potential wildlife passageways [59][60][61][62]. These structures maintain landscape connectivity and strengthen permeability by allowing animal circulation. However, it is essential that such structures have particular characteristics such as the presence of dry margins underneath that allow animal circulation [59][60][61][62]. The results obtained in this study suggest that most of the bridges or sewers found in the study area did not have a dry section that allowed fauna to cross underneath. For those that did have a dry section, the section was flooded most of the year, forcing wildlife to cross over it, increasing the likelihood of roadkill. Some mitigation measures implemented in these cases include the construction of concrete platforms to allow wildlife to cross. Such platforms usually have a height such that, when the watercourse rises, the level does not affect the crossing of fauna, nor the presence of an access slope to these platforms [60][61][62]. It is worth mentioning that bridges and/or sewers are important for native forest species, since in Uruguay they are associated with native forest watercourses. Therefore, the importance of these structures is paramount to allow the mobility of mammal species. However, it is necessary to conduct studies that analyse how the fauna use these areas and if they are indeed a problem for their displacement [59,61,62].
With regard to the roads that showed the greatest surface area favourability for roadkill, Route 9 is particularly noteworthy. This highway is characterized by high AADT values, as it is one of the international corridors that connects Uruguay with Brazil [63]. In turn, it also presents other characteristics that could explain the high favourability for roadkill, such as the high flow of tourists mainly in spring-summer as it connects the main tourist resorts in the department of Rocha. In addition, it crosses the Bañados del Este Biosphere Reserve, and it is located near several protected areas included in the National System of Protected Areas [29]. Baudi et al. (2017) found that roadkill increases within protected areas because they have greater species richness and abundance. This supports the models obtained in this study since, as mentioned above, Route 9 crosses several national and international conservation areas [31].
Another important aspect of Route 9 is that the shoulders of this route have, to a large extent, natural vegetation typical of the ecosystems of the area. In this regard, in recent years, researchers and government agencies related to transportation management, public works, and environmental protection, have put the focus on roadsides as an important reservoir of public land that, if managed properly, can become a shelter for biodiversity [64,65]. Thus, the promotion of suitable vegetation on these verges can contribute, as a nature-based solution, to minimizing the impacts associated with the roads themselves [65,66], whereas at the same time, it provides a suitable habitat for other species such as micromammals, birds, amphibians, reptiles, and various groups of invertebrates that may be of interest for conservation and for which road verges can become small biological corridors [65].
Finally, it is important to further strengthen the results obtained from the models and the GIS analysis; therefore, it is advisable to monitor the locations where there is a threat of roadkill and analyse the territorial conditions of each one in order to then implement mitigation measures tailored to each site. It is also important to carry out the analyses of this study on each species separately so as to identify areas of high mortality and generate mitigation measures tailored to each species.

Conclusions
This study is a key contribution to road ecology in Uruguay, which allowed us to understand how roads affect the populations of medium and large mammals, using the favourability function as a tool to generate predictive models. Models were generated to predict and identify sites where there is a negative interaction between roads, natural habitats, and the species of medium and large mammals under study. The main environmental variables that affect roadkill are related to anthropic factors (traffic, roads, and bridges). This will make it possible to monitor these sites in order to implement mitigation measures to reduce the loss of medium and large mammals, as well as to identify sites of value for landscape connectivity. However, despite being an important input, it is necessary to survey these sites to validate the models.
Supplementary Materials: The following are available online at https://www.mdpi.com/article/10 .3390/d13110585/s1, Table S1: List of explanatory variables used to model favourability grouped according to climatic, geographical, and anthropic factors, Table S2: Records of roadkill obtained in the period 2015-2016 for the native bush and grassland species selected for this study. R modeling script used.