PollMap: a software for crop pollination mapping in agricultural landscapes

Ecosystem service mapping is an important tool for decision-making in landscape planning and natural resource management. Today, pollination service mapping is based on the Lonsdorf model (InVEST software) that determines the availability of nesting and floral resources for each land cover and estimates pollination according to the foraging range of the desired species. However, it is argued that the Lonsdorf model has significant limitations in estimating pollination in a landscape that can affect the results of this model. This paper presents a free software, named PollMap, that does not have the limitations of the Lonsdorf model. PollMap estimates the pollination service according to a modified version of the Lonsdorf model and assumes that only cells within the flight range of bees are important in the pollination mapping. This software is produced for estimating and mapping crop pollination in agricultural landscapes. The main assumption of this software is that in the agricultural landscapes, which are dominated by forest and agriculture ecosystems, forest patches serve only as a nesting habitat for wild bees and the surrounding fields provide floral resources. The present study provided new software for mapping crop pollination in agricultural landscapes that does not have the limitations of the Lonsdorf model. We showed that the use of the Lonsdorf model for pollination mapping requires attention to the limitations of this model, and by removing these limitations, we will need new software to obtain a reliable mapping of pollination in agricultural landscapes.


Introduction
Ecosystem service mapping is a powerful platform to support the articulation of this knowledge in landscape planning and natural resource management (Fernandes et al., 2020). In this regard, researchers are interested in developing methods to quantify ecosystem services and show the scope of these services on maps (Picanço et al., 2017). These maps provide a spatial representation of ecosystem services that can be used to identify ecosystem services hot spots and predict the risks related to these services (Fernandes et al., 2020). Therefore, ecosystem service mapping helps decision-makers to identify key areas of different ecosystems providing ecosystem services and consider different trade-offs in land-use planning (Cunningham et al., 2018). To map pollination service, Lonsdorf et al. (2009) presented a model that examines the spatial arrangement of nesting and foraging habitats of wild bees. In this model, for each land cover, the availability of nesting and floral resources is determined and pollination is estimated according to the foraging range of the desired species.
Bees need nesting and foraging habitats, and the suitability of the nesting habitats depends on the floral resources surrounding the nests at bee foraging ranges. Therefore, bee abundance is limited by nesting and foraging habitat (Lonsdorf et al., 2009). Crop pollination is different from forest pollination as prior is an ecosystem service and the latter is an ecosystem function (Rahimi et al., 2021d). The Lonsdorf model estimates both ecosystem function and service and shows them on a map. This model is available in the InVEST software (Integrated Valuation of Ecosystem Services and Trade-offs) (Sharp et al., 2014) and has been applied in many studies for estimating pollination. However, it is argued that the Lonsdorf model has significant limitations in estimating pollination in a landscape that can affect the results of this model. For example, one of the important limitations of the Lonsdorf model is the lack of a behavioral mechanism to describe the foraging behavior of pollinators in the landscape so that his model does not show enough efficiency in complex or heterogeneous landscapes (Olsson et al., 2015). The Lonsdorf model only considers distance in determining bee foraging patterns. It is assumed that regardless of the quality of floral patches, the visiting rate probability from nesting patches towards floral patches decreases with increasing distance, meaning that bees are equally spread to surrounding areas (Fernandes et al., 2020).
However, bees are central place foragers (Bell, 1990), and studies have shown that they distinguish between high-quality and poor resources and locate their nesting habitat accordingly (Collett et al., 2013). Therefore, there is a trade-off between proximity and quality in visiting farms that are near nesting habitats than those with high-quality but farther away (Fernandes et al., 2020). To overcome this limitation, Olsson et al. (2015) proposed a behavioral mechanism based on central place theory (CPF). The strength of the CPF model is that it weighs the cost of traveling versus the rewards obtained from far patches in the landscape. In CPF-based models, bees do not visit patches beyond the foraging distance. However, the Lonsdorf model, assumes that some bees visit distant areas that are beyond bees flying ability (Olsson et al., 2015).
The abundance of pollinators in a landscape is influenced by biotic and abiotic factors (Klein et al., 2007;Rahimi et al., 2021b). To include mentioned abiotic factors in mapping relative pollination potential, an applied methodology adopted by Zulian et al. (2013), known as the ESTIMAP model (Ecosystem services mapping at European scale), was introduced to map pollination services across Europe. In addition to the factors required to implement the Lonsdorf model in InVEST software, Zulian et al. (2013) added several other factors to the output map obtained from the Lonsdorf model. They included the effects of some components like roads, water bodies, and climate in the process of pollination mapping. Recently, the ESTIMAP model is receiving more attention and new modifications and innovations have been added to this model to reduce the limitations of the Lonsdorf model (Cunningham et al., 2018;Fernandes et al., 2020).
In the landscapes that nesting patches have a higher value than surrounding foraging habitats, the central cells of the patches always receive higher scores than marginal cells (Rahimi et al., 2021d). Therefore, in the output maps of this model, it can be seen that the central cells of the patches are most likely to be pollinators. In other words, ecotones receive lower values than interior parts of nesting and foraging patches, while ecotones support a higher diversity of pollinators (Kells and Goulson, 2003;McKechnie et al., 2017;Potts et al., 2016;Svensson et al., 2000;Westphal et al., 2003). Thus, another significant limitation of the Lonsdorf model is the overestimation of pollinators' presence in inland habitats that is inconsistent with empirical studies. Numerous studies have shown a high abundance of pollinators in ecotones, and especially near the forest edge, in agricultural landscapes. For example, as the distance from forest edges increases, pollination, and crop production decrease, indicating the effect of forest edges on crop pollination (Carvalheiro et al., 2010). To overcome this limitation of the Lonsdorf model, Zulian et al. (2013), Cunningham et al. (2018), and Fernandes et al. (2020) differentiated forest edge at distances of 25 to 200 m from forest cores and increased the score of the edges compared to the inner forest.
Generally, the limitations of the Lonsdorf model are as follows: (1) this model assumes that bees are spread to surrounding areas with equal probability without considering floral patch quality, (2) the model includes areas that are beyond the flight ranges of bees in its calculations, (3) the model does not consider environmental factors that play a key role in the presence of pollinators, and (4) the model assigns higher scores on the interior parts of patches than the margins. From our point of view, these limitations have a significant impact on the results of the Lonsdorf model, and currently, only In-VEST software is widely used to map the pollination service. Therefore, the present study aims to provide new software for mapping crop pollination in agricultural landscapes that does not have the limitations of the previous model. Based on the model presented by Fernandes et al. (2020), this software assumes that bees selectively move from nesting habitat to floral resources. We also used multi-criteria evaluation (MCE) to include environmental variables to pollination modeling.

Estimating pollination based on the modified Lonsdorf model
The Lonsdorf model first measures the desirability or quality of patches that are suitable for the bee nesting habitat according to the floral resources around these patches (Equation 1). In the assessment of floral resources around nests, near pixels weigh more than those around distant ones according to the expected foraging ranges of the bee. The result is a map showing the desirability of nests between 0 and 1. In the next step, the model predicts the relative abundance of visiting bees in agricultural fields according to the desirability of the nests (Equation 2).
where FAj is the floral resource in neighbor cell j. Dij is the Euclidean distance between cell j and i. αβ is the foraging range for species β. The numerator is the total weight of the distance between all cells of floral resources adjacent to nesting patches. The quality of these cells (FAj) is between 0 and 1. The relative pollinator abundance in each nest site i (PA iβ ) is thus derived by the map algebra between nesting suitability (NS iβ ) and (FQ iβ ) (Equation 2): To determine the bees' abundance on farms, similar to the framework of Equation 1, cells from farms that are closer to the nesting habitats receive more probability of wild bee abundance. In this equation (Equation 3), the normalized proportion scores each cell according to the distance weighted PA summed across all J cells visiting cell i. As mentioned in the introduction section, one of the limitations of the Lonsdorf model is considering the equal probability for bees to spread toward the floral patches regardless of their quality. To overcome this limitation, Fernandes et al. (2020) modified the foraging model in Equation 4 according to the distance from nests and respective floral quality.
where the normalized proportion scores each cell according to the Euclidean distance between cell i and PA source cell j and the floral availability of cell i in comparison to the overall FQ available in the neighborhood J (Fernandes et al., 2020). In this study, we used the modified Lonsdorf model for pollination mapping. Therefore, using Equation 4, one of the most important limitations of the Lonsdorf was solved. We included only cells that are within the flight range (α) of bees in the software calculations. This flight range is defined by the user. We also include a condition to our modeling as nesting patches receive their values only based on the surrounded floral resource. With this condition, because the central cells are farther from the floral resources than the marginal cells, they always receive less scores.

Including environmental factors in pollination mapping
To reach a specific objective, it is clear that several criteria need to be evaluated. Such an approach is known as Multi-Criteria Evaluation (MCE). In multicriteria evaluation, there is a concern about combining the different factors to form a single index of evaluation. In the present study, a weighted linear combination (Eastman, 2012) was used to combine the additional criteria to the Lonsdorf model output. In the weighted linear combination method (Equation 5), criteria are categorized into factors and constraints (Eastman, 2012). The factor is a measure that increases or decreases the appropriateness of an option for the intended purpose. Criteria are usually measured at different scales; therefore, factors must be standardized before combination (Eastman, 2012). There are different procedures for standardization, typically using the minimum and maximum values as scaling points. In MCE, fuzzy set membership is used in the standardization of criteria. Fuzzy sets are sets without sharp boundaries, meaning that there is a gradual transition between membership and non-membership of a location in the set (Zadeh, 1996). A fuzzy set is characterized by a fuzzy membership grade that ranges from 0 to 1, indicating a continuous increase from non-membership to complete membership. We used a linear membership function characterizing a linear relationship between pollinators and environmental factors (Fig. 1).
where, S is suitability map, W i is the weight of criterion i, and X i is the criterion i.
In PollMap software, for each factor, according to the linear membership function (Fig. 1), two inflection points are considered, which are displayed with letters a and b. If the user enters the value of a less than b, the software standardizes the desired factor as a monotonically increase. For example, if the roads negatively affect pollinators, then by increasing the distance from the road, the habitat suitability for pollinators increases. In this case, the value of a must be assigned less than b. If the low values of an environmental factor are more suitable for pollinators, the value of a must be greater than b. environmental factors have different effects on determining the pollinators' presence and their relative importance is different from each other. For this, a weight is usually assigned to each factor. For example, The Analytic Hierarchy Process (AHP) is a tool for weighting that is done through pairwise comparisons and judgments by weighting experts. The weight of the factors is from 1 (extremely important) to 9 (extremely important) (Saaty, 2008). Assignment of weight to different layers is based on the literature and the experts' opinion. Therefore, the user must assign a weight between 0 and 1 according to the relative importance of each layer. The sum of all weights must not exceed 1. According to Equation 5, the weight of the layers is multiplied by them and all the layers are added together.

Input data
PollMap is a free software for Windows (32 or 64 bit), developed using the MATLAB framework. The main assumption of PollMap software is that consistent with other studies (Mitchell et al., 2015;Rahimi et al., 2021a), in the agricultural landscapes, which are dominated by forest and agriculture ecosystems, forest patches serve only as a nesting habitat for pollinators especially aboveground nesting bees that nest in tree cavities or bee hotels (Rahimi et al., 2021c). On the other hand, the surrounding agricultural fields provide floral resources. Therefore, to run this software, two layers of nesting habitat and floral resources are needed separately. Nesting and foraging habitats might have several patches with different fitness to bees, which is determined according to the Lonsdorf model between 0 and 1. Therefore, the user must determine the fitness value of each patch for input layers in advance. It is suggested that the user first prepare two layers of nesting and foraging habitat separately in Arc GIS software as a shapefile and determine the fitness of each patch for bees and then convert the desired layers to a raster tif format based on the fitness field. Figure 2 shows the PollMap software screen. The software screen is divided into four general sections. Section A is related to the land use map layer and the user calls the desired layers, i.e., floral resources and nesting habitat by load button in this section. After calling the desired layers, the user can display the layers using the imshow button. In the box embedded for alpha (α), the user must specify the alpha value per cell unit. Therefore, to provide the pollination map, the amount of alpha and two mentioned maps are needed. The box marked with W is assigned for weights of factors and when the user wants to add environmental layers to the model. Section B is created to incorporate environmental factors into the pollination map. In this section, the user can add seven layers to the model. After calling each layer, the user must specify the values of a and b to standardize each layer and must display it with the imshow button. Section C has two buttons, Run 1 and Run 2. Run 1 is used to estimate the pollination map using two layers of floral resources and nesting habitat. For adding one or more environmental factors to the pollination map, the user must use Run 2. After obtaining the results, the user can save the pollination map in a specific path using the save button. Figure 3 shows examples of the Lonsdorf model and the modified version. Figure 3A shows the original landscape. It is assumed that each farm cell has a random value between 0 and 1 (Fig. 3B). The modified model shows well that the edges of the forest patches have a higher score than the central cells (Fig. 3C2), and farms near forest patches indicate more pollination rate than the center cells of the forest patches (Fig. 3D2). Figure 4 shows a comparison of the Lonsdorf and the modified models in considering the quality of floral resources. Figure 4A shows a real landscape in which forest patches are highlighted in blue. The yellow matrix (crop 1) indicates an agricultural field that has a moderate quality for attracting pollinators (0.5), the green patches (crop 2) show fields that have a high quality of floral resources (0.9). Figure 4B shows the nesting suitability of forest patches, in which marginal cells scored higher than central cells. Figure 4C shows the pollination map, assuming that the bees are equally distributed into their surroundings. Figure 4D shows the pollination map obtained from the modified model that takes into account the quality of the fields in bees' movement. Comparison of Figure 4C and D show that the modified model can well consider the quality of farms in providing the pollination map.

Software description
To better illustrate the differences between the Lonsdorf and the modified models, we separated the two subsets  Fig. 4C1 and D1 show a yellow forest patch with a farm in the center of the patch, which is shown in red on Fig. 4D1. In Fig. 4D1 subset, this farm has obtained a high score because it has a high quality of floral resources and is surrounded by forest patches, but in Fig. 4D1 subset, the quality of this farm is not considered and bees are assumed to distribute equally, and therefore it scored a low value. Subsets Fig. 4C2 and D2 also show how the modified model considers farm quality in pollination mapping. In Fig. 4C2 subset, farm cells received an equal service from forest patches, while the Fig. 4D2 subset shows that there are high-quality farms around forest patches, and as a result, these farms score higher than the surrounding farms.  Figure 5 shows an example of adding a road factor to a pollination map. Figure 5A shows a pollination map obtained concerning floral resources and nesting habitat. Figure 5B shows the distance from road map, which is standardized based on fuzzy membership functions between 0 and 1. This standardization includes numbers from 0 to 516 m, which is the maximum distance from the road. Figure 5C assumes that roads can affect the presence of pollinators up to 50 meters away. In Fig. 5D, this distance is considered to be 100 m. Adding these maps to the pollination map has resulted in significant changes to the final map.

Discussion
The present study introduces new software for crop pollination mapping based on the modified Lonsdorf model that does not have the limitations related to previous software (i.e., InVEST). In this study, we introduced four basic limitations related to the Lonsdorf model based on other studies (Fernandes et al., 2020;Olsson et al., 2015;Rahimi et al., 2021d;Zulian et al., 2013) and showed that changing some of the assumptions of the Lonsdorf model can provide different estimates of crop pollination.
Consistent with other studies PollMap assumes that in the forest environment, only the edges near the fields are attractive to pollinators (Schulp et al., 2014). This claim has been confirmed by other experimental studies. For example, Solitary bees build their nests at the edge of the forest, where the shade is less and the canopy is less dense (Klein et al., 2003). The richness of bumblebee species increases at the edge of agriculture and semi-natural habitats (Aguirre-Gutiérrez et al., 2015). The highest number of queen bees has also been observed at the forest edge (Svensson et al., 2000). Wild bee diversity decreases in forest environments with Fig. 5 Examples of the including road factor in the estimating pollination. In this figure, (A) shows the pollination map based on land use only, (B) shows the map showing distance from road standardized between 0 and 1, (C) adding road factor to pollination map assuming that roads up to 50 meters away negatively affects the presence of pollinators, and (D) assuming that roads up to 100 meters away negatively affects the presence of pollinators dense stand and large canopy (Hanula et al., 2015;Odanaka and Rehan, 2020;Roberts et al., 2017) because in dense forests the temperature is lower and sunlight does not reach understory plants (de Lima et al., 2020). In Fig. 3, we show how this limitation of the Lonsdorf model can lead to over estimates of pollination at the center of forest patches that are far from farms or floral resources.
PollMap take into account abiotic factors that affect the presence of pollinators in a landscape. There are many environmental variables that affects the presence of pollinators. For example, the ambient temperature affects pollinators' activity considerably (Zulian et al., 2013). Altitude affects local weather and as it increases, the population of many pollinators decreases (Devoto et al., 2005 293;Hodkinson, 2005 294;Kimball, 2008;Totland, 2001 296). Wetlands and riverside areas also play an important role in pollination by providing diverse nesting and foraging habitats for pollinators (Cole et al., 2015;McInnes, 2018). Roads also have a considerable effect on insects such as bumblebees (Kallioniemi et al., 2017;Keller and Largiader, 2003). Phillips et al. (2020) examined 141 studies related to the effects of roads on pollinating insects. They found that traffic and pollution caused by roads negatively affected pollinators. Muñoz et al. (2015) also examined the effects of roads on insects and reviewed 50 studies and found that fragmentation, pollution, accidents, and traffic caused by roads adversely affected the diversity and abundance of insects. Therefore, it is necessary to include the role of abiotic factors in the presence of pollinators. In the present study, we also showed in Fig.  5 that including roads to pollination map can significantly change the pollination rate at the landscape level.
According to the CPF model, PollMap assumes that there is a maximum distance for the central place foragers that they do not go beyond this distance for foraging. Therefore, near the nest patches, all highquality patches are visited, but at distances away from the nests, only the best patches are used. A comparison of Fig. 4C and D shows that cells that are far from forest patches have a score of zero in Fig. 4C, but in Fig. 4D, which is the result of the Lonsdorf model, these cells have a value greater than zero.

Conclusion
The present study emphasizes that the Lonsdorf model has significant limitations that challenge the accuracy of the maps obtained from this model. Therefore, we tried to provide more reliable results by providing compact software concerning these limitations. The results of the alternative model to the Lonsdorf model showed that farms near forest patches show higher pollination rates than those far away. In general, the present study showed that the use of the Lonsdorf model for pollination mapping requires attention to the limitations of this model, and by removing these limitations, we will need new software to obtain a reliable mapping of pollination in agricultural landscapes.

Software availability
Software and sample data are available at https://github. com/ehsanrahimi666/PollMap.git