1 Introduction

Annual global application of agrochemicals has amounted to 2.7 × 106 tons in recent years (OECD/FAO 2017), mostly glyphosate (Ferreira et al. 2017). The environmental implications of the increased use of this active substance have led to a growing number of glyphosate-tolerant native plants (Ferreira et al. 2017), negative consequences for many other groups of organisms, and social impact, among others (Lupi et al. 2019). In Argentina, the most widely used practice is industrial agriculture with direct sowing (SD) and agrochemical pesticides, accounting for 90% of the agricultural area (Aapresid 2018). This model relies exclusively on the application of increasing amounts of pesticides to control weeds and pests, with glyphosate, 2,4-D and atrazine being the most widely used systemic herbicide (Casafe 2014). Industrial agriculture in Argentina has reduced natural and semi-natural habitats, generating habitat fragmentation and significantly affecting floral diversity and abundance (Pengue 2005; Galetto and Torres 2020). This has resulted in a net impact on the abundance and health of wild and managed bees (Dolezal et al. 2019; St Clair et al. 2020).

Annual losses of honeybee (Apis mellifera) colonies in Argentina range from 15 to 30% (Maggi et al. 2016; Requier et al. 2018). These losses are mainly ascribed to the presence of parasites as well as the application of agrochemicals in the areas surrounding the apiaries (Maggi et al. 2016).

The main routes of bee exposure to agrochemicals occur by contact with the drift of contaminated spraying or dust, and by the ingestion of residues in vegetation and water bodies (Botías and Sanchez-Bayo 2018). Pesticides are considered one of the main abiotic stressors for pollinators (IPBES 2016). Medici et al. (2019) carried out a pesticide survey in the largest agricultural exploitation area of Argentina, and identified the provinces with the highest risk for beekeeping, due to the presence of one or more agrochemical residues in honey. The province of Buenos Aires manages the largest number of beehives in the country and reports the highest number of honey pollutants (Medici et al. 2019). Similar results were also observed by Villalba et al. (2020). Indeed, the presence of these compounds poses a significant risk not only in terms of the activity but also the health of the colonies themselves.

During foraging flights, honeybees settled in agricultural ecosystems may encounter glyphosate residues both on cultivated and surrounding flowers (native or exotic) (Farina et al. 2019). This systemic herbicide has been considered non-toxic to bees based on acute contact and oral toxicity tests conducted (> 100 µg/per bee) (Frasier and Jenkins 1993). Nevertheless, a recent study has found that, under laboratory conditions, a low concentration of ̀glyphosate (2.5 mg/L) affects the flavor’s taste, responsiveness and learning performance of managed bees (according to proboscis extension response (PER) tests) (Balbuena et al. 2015). Other studies have reported deleterious effects of glyphosate on bee gut microbiota (Motta et al. 2018) as well as a negative impact on brood development (Farina et al. 2019), and difficulty in establishing non-elemental associations that could affect resource collection by foragers bees (Herbert et al. 2014). Intensive agriculture can contribute to pollination decline, exemplified by the high annual losses of honeybee colonies in regions dominated by annual crops (Dolezal et al. 2019).

The configuration and heterogeneity of the landscape in which the bee colonies are located, along with the characteristics of agricultural management, have direct and indirect implications on the survival of bees (e.g., Dicks et al. 2016; Dolezal et al. 2019; St. Clair 2020). Dolezal et al. (2019) showed that intensively farmed areas can provide a short-term feast unable to sustain the long-term nutritional status of colonies and suggested that partial restoration of biodiversity in such landscapes may provide relief from nutritional stress. The study conducted by St. Clair et al. (2020) postulated that fruit and vegetable farms provide some benefits for honeybees; however, they do not benefit wild bee communities. These authors further suggest that incorporating natural habitats, rather than diversified farming, in these landscapes may be a better choice for wild bee conservation efforts.

The aim of this work was to analyze the relationship between the landscape configuration and the presence of pesticides in honey since the characteristics of the crop matrix provide proxies for pesticide contamination (Medici et al. 2019). In particular, we tested the relationship between the presence of glyphosate and its metabolite Aminophosphoric Acid (AMPA) residues in honey and some metrics of the landscape configuration. As the cropland proportion in the landscape grows with increasing volumes of herbicides applied per site, values of glyphosate and AMPA in honey are expected to be higher. At the same time, as the functional homogeneity or heterogeneity of the landscape dominated by croplands is higher, increasing values of glyphosate and AMPA in honey are expected because honeybees can easily explore the whole area due to less barriers are present.

2 Materials and methods

A total of thirty samples of honeys from different sites (apiaries) from the southeast of Buenos Aires province were obtained during 2019 and 2020 harvests (Table I). Honey samples were acquired from each studied site, collecting them from batch of the entire apiary: so, the honey sample do not represent a beehive but the whole apiary. Honey was collected in the same season at the end of the honey harvest (Bicudo de Almeida-Muradian et al. 2020). Apiaries were opportunistically sampled, considering the possibilities of contacting the producers and, in some cases, some sites were sampled twice either being sampled in two different years (same seasons) or when skinny hives needed a second extraction.

Table I. Transitions and chromatographic parameters used for glyphosate and AMPA detection
Full size table

Each apiary was georeferenced by GPS. Once the samples were obtained (one per apiary), they were frozen at -20 °C until their analysis. We identified a total of 22 different sites (apiaries) to test the relationship between landscape configuration and pesticide residues in honey, eight of the samples were obtained from the same apiary (during the same or different season; Figure 1).

Figure 1.
figure 1

Supervised classification of each of the 22 sites presented in Table II for the study period using Sentinel 2 image.

Full size image

2.1 Analysis of glyphosate and AMPA residues in honey

A method of analysis of glyphosate residues and its main degradation metabolite, aminomethylphosphonic acid (AMPA), in honey was developed.

Extraction: 100 gr of sample was homogenized before the analysis by stirring thoroughly. From this homogenate, 1 g of honey was weighed in a 50 mL centrifuge tube on an analytical balance at 0.1 mg, and 5 ml of water was added. Subsequently, 10 µL of concentrated hydrochloric acid and 5 ml of methanol were added. It was shaken 1 min and centrifuged 2 min at a speed greater than 1500 rpm. 1 ml of the supernatant was transferred to a glass vial and 200 µL of borate buffer pH 9 and 500 µL of 9-fluorenylmethoxycarbonyl chloride (FMOC) solution were added. Subsequently, it was placed in an oven at 50 °C for 90 min. It was then acidified with 10 µL of concentrated formic acid and filtered with a syringe filter.

Calibration curve: The necessary volumes of working solutions were added to 1 g of honey to obtain concentrations of 10, 25, 50, 75, and 100 µg/kg.

Chromatographic conditions: For the analysis, a Waters Acquity H-Class Liquid Chromatograph was used coupled to a Xevo TQ-S Spectrometer. An Acquity BEH C18 100 × 2.1 mm column was used, with a guard column, 1.7 µm particle size. The mobile phases were: A: Water: Methanol (95: 5) 10 mM Ammonium formate; B: 10 mM Methanol Ammonium formate. The oven temperature was 40 °C, and the injection volume was 50 µl. MS detector conditions: The TQS detector mode was ESI (-), MS/MS, with a day flow of 600 L/h and a temperature of 400 °C. The source temperature was 125 °C. To increase the likelihood of a positive identification, an accurate number of MRM transitions were acquired for each target compound (Table I) to generate a library searchable product ion spectrum, designed to minimize the risk of false positive or negative detection. The duration of the chromatography run was 12 min.

Once the run was finished, the result was evaluated and quantified with the chromatography software (Waters Target Lynxs v4.1.SCL945.SCN960). The results obtained were sample-weight corrected. The results were accepted if pesticide recovery in the spiked was in the range of 80–110%. The detection limit of glyphosate and AMPA obtained was 1.0 µg/kg.

2.2 Habitat configuration

Of the selected study area, honey samples were collected during autumn in sites sowed by summer crops (soy, corn, sunflower) and winter crops (oats, barley, wheat). Among summer crops, soybean represents 65% of the total planting area, while corn (with the best contribution to soil cover) only 13% (Ligier et al. 2018). Croplands could be a favorable habitat for honeybees during the soybean flowering (January to mid-February), but it would not be later, during fruit maturing. Visits of honeybees to soybean flowers have been reported for different regions of Argentina (Zelaya et al. 2018; Huais et al. 2020).

Two radii (1000 m and 2250 m) surrounding each apiary were selected, because different studies have shown that, on average, honeybees do not often forage farther than 2 km from their hives (Carr-Markell et al. 2020). The 1000 m and 2250 m radii are equivalent to a landscape of 9.87 and 22.21 km2 area, surrounding each apiary, respectively (Figure 1). Cloud-free, Sentinel-2, Level-1C products were acquired from the United States Geological Survey (USGS; http://earthexplorer.usgs.gov/) in January 2020 to assess and classify the area during the same period for most of the honey samples obtained. The spectral bands used in this study included blue (0.49 μm), green (0.56 μm), red (0.66 μm), and near infrared (0.842 μm) with a spatial resolution of 10 m. Digital numbers of the Sentinel 2 imagery were converted to top-of-atmosphere reflectance. Numerous sites selected during field reconnaissance and high-resolution images in Google Earth were used as the training sites for the supervised image classification. The classification of images was performed using the semi-automatic classification plugin (SCP) in QGIS version 3.12 (Quantum GIS Development Team 2020). To reach adequate results, different successive classifications were necessary, masking or adding areas, combining and removing training sites. Post-processing was performed by correcting classification. Land cover characteristics were assessed for circular buffer areas of increasing radius (1000 m and 2250 m) around each sampling site. The land use categories considered were: 1) Urban (impervious surfaces such as buildings, roads); 2) Open vegetation (grass, low vegetation, isolated trees); 3) Forest (continuous arboreal vegetation), 4) Agricultural, and 5) Water, as water bodies. The supervised classification of each site is shown in Figure 1.

We applied the spatial pattern analysis software FRAGSTATS 4.2 (McGarigal et al. 2012) to calculate the landscape metrics of each class type and total landscape. We chose two class-level metrics: Class area (CA, in hectare) and percentage of landscape (PLAND, in percentage) for each of the land use categories. These are landscape composition measures that calculate how much of the landscape is comprised of a particular patch type. Also, the Contagion index (CONTAG, in percentage) was calculated. This index consists of the sum, over patch types, of the product of 2 probabilities: (1) that a randomly chosen cell belongs to a patch type i (estimated by the proportional abundance of patch type i), and (2) the conditional probability that given a cell is of patch type i, one of its neighboring cells belongs to patch type j (estimated by the proportional abundance of patch type i adjacencies involving patch type j). The product of these probabilities equals the probability that 2 randomly chosen adjacent cells belong to patch type i and j. Contagion index measures both patch type interspersions (i.e., the intermixing of units of different patch types) and patch dispersion (i.e., the spatial distribution of a patch type). A landscape in which the patch types are well-interspersed will have lower contagion than a landscape in which they are poorly interspersed. Higher values of contagion may result from landscapes with a few large, contiguous patches (agricultural patches in our case, i.e. a higher homogeneity in the landscape because predominance of croplands), whereas lower values generally characterize landscapes with many small and dispersed patches. Thus, contiguous patches of agriculture will have greater contagion than a landscape in which the patch types (i.e. croplands) are fragmented into many small patches. The contagion index represents the observed level of contagion (homogeneity/heterogeneity) as a percentage of the maximum possible given the total number of patch types (McGarigal et al. 2012). Spatial correlation was performed using Moran's spatial autocorrelation index with the Crime stat V software (Levine 2015).

Lineal Mixed Models (LMMs) were used to determine the effects of the landscape configuration on each of the dependent variables from the honey samples: Glyphosate (µg/kg) and AMPA (µg/kg). Among the CA land use categories, we selected those directly related with herbicide application: the amount of the cropland area (CA AGRO) in hectares (at 1000 m and 2250 m radii) and percentage of croplands (PLAND AGRO) of this class area for 1000 m and 2250 m radii in the landscape. Due to PLAND AGRO had a significant positive correlation with CA AGRO, CA AGRO was not modeled. Thus, the independent variables included in the models were percentage of croplands and the contagion index. Year was included in the models as a random factor. We log-transformed all 30 data of the response variables (Table 2) to meet assumptions of normality and performed all statistical analyses in R v.4.0.4 A (R Core Team 2021). For each main model (one per radius), we used all the independent variables with the function lmer and from lme4 package (Bates et al. 2015), and then we eliminated one by one and compare different models with the likelihood ratio test with the function ANOVA from stats package (R Core Team 2021).

Table II. Values obtained for glyphosate and AMPA from 30 honey samples obtained from 22 sites, which metrics for the landscape configuration at 1000 m and 2250 m radii are presented. In particular, the hectares of croplands (CA AGRO 1000 and CA AGRO 2250), the percentage of agriculture in the landscape (PLAND AGRO 1000 and PLAND AGRO 2250), and the measures of homogeneity or heterogeneity of the landscape (CONTAG for both radii) are presented for each sample
Full size table

3 Results

Table II lists the values obtained for glyphosate and AMPA residues for all the sites during 2019 and 2020. A total of 30 honey samples belonging to the 2019 and 2020 harvests were analyzed in various areas of south-eastern Buenos Aires. The results obtained for each honey sample from georeferenced apiaries and their association with some metrics of the landscape configuration are presented in Table II and Figure 1. A frequency of 50% of glyphosate positivity was found. In turn, AMPA was found in 30% of the samples (Table II). The range of values found in the samples was between 2 and 27.5 µg/kg for glyphosate and 1.9 – 18.1 µg/kg for AMPA. The selected metrics for agriculture showed some variations between sites. Percentage of cropland index calculated on agriculture class at 1000 m and 2250 m radii ranged between 1.16 and 77.27% and between 16.6 and 76.13%, respectively; the amount of cropland area ranged between 4.6 and 310 ha for 1000 m radius and between 333 and 1534 ha for 2250 m radius (Table II).

The contagion index at a landscape level was greater than 50 percent in all sites, thereby indicating a high level of clustering between pixels of the same coverage (most of them with agriculture as the dominating land use category). Values for the contagion index at short scale (1000 m radius) ranged from minimum values of 51.9% at site 18 to maximum values of 92% at site 15, with most sites showing values over 60% (16 of 22 sites; Table II). Values for the contagion index at a larger scale (2.25 km radius) ranged from minimum values of 43.39% at site 3 to maximum values of 93.04% at site 22, with most sites showing values over 60% (16 of 22 sites; Table III).

Table III. Statistical results obtained for best fitted LMMs using all data
Full size table

3.1 Habitat configuration

Some sites that were sampled twice during the same or different seasons displayed a variety of trends. Some of them showed the same consistent pattern (absence of pesticide residues in site 10 for sampling 2019 and 2020; or presence in sites 8, 11 and 12 for sampling 2019 and 2020; Table II), while others showed the presence of residues or not between sampling times (sites 5 and 9; Table II).

The tests for spatial autocorrelation (Moran´s index) showed values close to 0 (see Table S1, Supplemental Material), which indicates that there is no spatial autocorrelation in the observations and that they are spatially distributed independently of each other.

The best fitted LMMs models with significant results selected only percentage of croplands to explain the variability of the response variables (Table III). LMMs results showed that the variable percentage of cropland significantly explained the amounts of glyphosate (p = 0.04272) and AMPA (p = 0.0101) honey residues in 1000 m radius of the landscape configuration and the residues of AMPA (p = 0.01211) in 2250 m radius (Table III). All relationships were negative (Figure 2A-C) meaning that the concentration of glyphosate and AMPA in honey fell 0.16 µg/kg and 0.17 µg/kg, respectively, when the percentage of cropland is increased 1% in the 1000 m radius, and AMPA concentration fell 0.17 µg/kg when percentage of cropland increased 1% in the 2250 m radius. Table III shows statistical results for all models display in this work. Because site 3 could be an outlier (very low percentage of croplands in the surroundings of the apiary), we run the LMMs models without this data (Figure 2D-E) and the tendencies were similar (Table S2, Supplemental Material).

Figure 2.
figure 2

Significant estimates of the models that best fit. With all data: A: Scale 1000 for glyphosate; B: Scale 1000 for AMPA; C: Scale 2250 m for AMPA. Without the site (3) because it showed the lowest percentage of cropland: D: Scale 1000 for glyphosate; E: Scale 1000 for AMPA; F: Scale 2250 m for AMPA.

Full size image

When sites with apiaries showing honey pesticide residues were mapped, an interesting pattern was observed at a regional level, with most sites (except one) with residues in honey are regionally concentrated (Figure 3). Nevertheless, some adjacent sites also showed no residues of pesticides in honey (Figure 3).

Figure 3.
figure 3

Map of 22 sites with apiaries showing honey pesticide residues. Source: Google satellite-QGIS.

Full size image

4 Discussion

Our results indicate that glyphosate was present in 50% of the samples analyzed while AMPA was found in a lower percentage. These results are in line with those reported by Rubio and Kamp who informed 59% positive samples for glyphosate in a survey carried out in the USA, but lower than those reported by Thompson et al. (2019) who informed 98 and 99.5% of samples positive for glyphosate and AMPA in Canada, respectively.

We expected a higher presence of glyphosate and AMPA as the proportion of croplands or their closeness (i.e., a higher homogeneity for this land use category) in the landscape are increased. However, our analysis found that glyphosate and AMPA residues (level and presence) are negatively related to some metrics of habitat configuration mostly in the areas surrounding the apiaries (1000 m radius). Our results contrast with those reported by Berg et al. (2018) were, among all the variables explored, agricultural land use showed a strong positive correlation with glyphosate concentration. During the last three decades, Argentina increased dramatically the total annual amount of glyphosate released to the environment (Galetto et al. 2021), mainly due to the expansion of the agriculture frontier boosted by the biotech soybean adoption produced largely under no-till practice (78% of the total seeded area) (SAGyP 2016). This use of glyphosate does not occur only for weed control, but also for chemical fallow. In addition, the amount released per hectare to control weeds generated tolerance or resistance of many native species (Ferreira et al. 2017). Primost et al. (2017) state that, under current practices, application rates of this herbicide are higher than dissipation rates, so they proposed that glyphosate and AMPA should be considered “pseudo-persistent” pollutants and that a revision of management procedures, monitoring programs, and ecological risk for soil and sediments should be also recommended. In view of honeybee’s broad foraging preferences and land use changes play a key role in colony survival and honey production (Smart et al. 2016; De Groot et al. 2021), colonies have the potential to be used as environmental pollution monitoring tools, reflecting the use of agrochemicals in the area (Shimshoni et al. 2019). Nevertheless, our results suggest that the increasing amount of croplands along with the intensification of industrial agriculture is not enough to explain the relationship within glyphosate and AMPA residues in honey. In particular, when we mapped those sites with apiaries showing honey pesticide residues, it is found a regional concentration of the sites with residues in honey, often independent of landscape configuration (i.e., weak negative relationships between the pesticides in honey and metrics of landscape). These trends suggested that more detailed studies in particular regions of landscape could be useful to better understand the relationship of agricultural practices and the presence of pesticides in honey. Further detailed studies should be conducted in the areas surrounding the apiaries in order to determine whether glyphosate is in fact a contributing factor to colony losses previously reported in Argentina (Maggi et al. 2016; Requier et al. 2018), based on pesticide prevalence in the environment, as well as on our findings in this regard in honey samples.

Usually, glyphosate is applied by airplanes in those croplands far away of urban areas and by terrestrial devices close to urbanizations. These methodologies imply different rates of glyphosate drifts that impact differentially on the croplands and vegetation of semi-natural habitats that could be related to the residues in the apiaries according to their placement, independently of the landscape configuration and honeybee behavior. Nevertheless, the results of this work strengthen the assumption that honey could be a good environmental sentinel of agrochemical contamination (Medici et al. 2019). All food destined for human or animal consumption in the European Union (EU) is subjected to a maximum limitation of pesticide residues to protect human and animal health (EUR-LEX). Hence, better understanding how environmental factors modify the presence of these herbicides in the colonies would provide valuable information for decision-making and planning at a national and global scale.

5 Conclusions

The presence of both glyphosate and AMPA in honey, even at very low levels, identifies an important pathway whereby pesticides migrate from site of application to the hive and into the human food supply (Damalas and Eleftherohorinos 2011; Bargańska et al. 2016; Mullin et al. 2010). The method validated in this work is simple and may be used to control the quality and safety of honey in terms of glyphosate and AMPA residues, since the MRL, for glyphosate in honey according to the European Union, is fixed at 0.050 mg kg-1.

A geospatial analysis, like the one performed in this study with some metrics of landscape configuration, can help honey producers to think about on pesticide exposure risks related to industrial agriculture. When bees are used for the commercial exploitation of honey production, hive placement on semi-natural habitats should be optimized to reduce their exposure and then bee colonies could be less impacted by pesticides (Mullin et al. 2010; Hartel and Ingolf 2014). As it has been postulated by Berg et al. (2018) linking spatial and temporal dynamics of flowering crops, agro-environmental schemes, and pesticide applications will lead to a better understanding of environmental risk assessment, management of pollination services, and biodiversity protection (Mullin et al. 2010; Donkerskly et al. 2018).