1 Introduction

Climate plays a key role in the geographic distribution of species, population dynamics, and ecosystem functionality (Easterling et al. 2000; Thomas 2010; Walther et al. 2002; Calosi et al. 2010; Sunday et al. 2010; O’Sullivan et al. 2017). Over recent decades, a global change in the climate at a large scale with slow and progressive changes in the main climatic parameters (e.g., humidity, temperature) has been observed (Easterling et al. 2000; Carey and Alexander 2003; IPCC 2014). Since the beginning of the twentieth century, the global temperature has undergone an average increase of 0.85 °C (1.2 °C in Belgium from 1900 until now) (Easterling et al. 2000; Bale et al. 2002; Vandenbohede et al. 2008; Schweiger et al. 2010). This current climate change could change ecosystems by modifying the species composition (Parmesan et al. 2000; Fitter and Fitter 2002; Hance et al. 2007; IPCC 2014).

Climate change is also related to an increase in frequency and intensity of extreme and localized weather events such as heat waves (Easterling et al. 2000; Meehl and Tebaldi 2004; Lemos and Rood 2010; Russo et al. 2015; Ragone et al. 2018). For the US National Weather Service, heat waves correspond to a minimum period of 2 consecutive days during which the temperature is higher than 40 °C. However, this definition is valid for a defined geographical area. In Belgium, the “Institut Royal de Météorologie” defines a heat wave as a period of 5 consecutive days with maximum temperatures above 25 °C with at least 3 days above 30 °C (IRM 2019). In this country, the longest heat wave of the last century has been recorded in 2018. The daily temperature exceeded the threshold of 30 °C for 26 consecutive days. The maximum temperature was recorded in July 2019 with 41.8 °C (IRM 2019). Future models show that the USA and Western Europe (especially northern France, southern England, Belgium, the Netherlands, and western Germany) should undergo more intense and frequent heat waves (Meehl and Tebaldi 2004).

These heat waves represent hyperthermic stresses for living species and are associated to direct behavioral (e.g., stupor or coma) (Hance et al. 2007; Martinet et al. 2015a) and/or physiological perturbations (e.g., ontogenic development, water balance, fertility, and immunity decrease) that could lead to an increase of mortality in several groups of organisms (Neven 2000; Parmesan 2006; Rasmont and Iserbyt 2012; Kingsolver et al. 2013). The effects of these extreme events could be more severe than those of a gradual increase in average temperature (Bale et al. 2002; Hance et al. 2007). The effects of heat waves are different depending on the physiological sensitivity of the species (Huang et al. 2006; Deutsch et al. 2008; Chen et al. 2011; Martinet et al. 2015a). Concerning pollinators, climate change is known to impact species composition and relative abundance particularly in bumblebee communities (Rasmont and Iserbyt 2012; Kerr et al. 2015; Rasmont et al. 2015, but see Herrera (2019) showing long-term stability of bumblebees in the face of steeply increasing mean ambient temperatures).

We assume that climate change will increase heat-induced mortality (Deutsch et al. 2008), testing the impacts of these effects is a crucial issue for biodiversity conservation policies.

Bumblebees are cold-adapted insects with their species diversity hotspots localized in the mountains and Arctic, subarctic, and Boreal regions (Michener 2007; Rasmont et al. 2015). They are robust and hairy heterothermic hymenopterans with strong endothermic capability. Individuals undergo thermal distress at 35 °C and die at 44 °C (Heinrich 2004), and even their endo-heterothermy system cannot protect them against such heat stress (Heinrich 1976, 2004; Martinet et al. 2015a; Oyen et al. 2016; Oyen and Dillon 2018). Previous studies suggest that hyperthermic stress could very quickly lead to fatal consequences for bumblebees (Martinet et al. 2015a; Oyen et al. 2016; Oyen and Dillon 2018). In lab conditions, Martinet et al. (2015a) showed that after successive behavioral phases (linked to physiological thresholds, Hazell and Bale 2011), hyperthermic stress leads to a heat stupor state characterized by a critical decrease of motor function leading to its death (Hutchison 1979).

Recently, several authors showed that climate change could play a key role in bumblebee decline (Kerr et al. 2015; Rasmont et al. 2015; Soroye et al. 2020). In Belgium, some bumblebee species declined sharply, especially long-tongued and specialized species, during the last century (Rasmont and Mersch 1988; Maebe et al. 2016; Vray 2018). While the abundance of some species remains stable (e.g., B. pascuorum, B. lapidarius, B. pratorum), others are declining such as B. muscorum, B. ruderarius, B. veteranus. Several species have increased their relative abundance and distribution (e.g., B. bohemicus, B. hypnorum, B. lucorum, and B. terrestris).

We aim to assess the heat stress sensitivity of ten species of bumblebee community in Belgium showing different population trends (i.e., stable, regressive, and expanding) to better understand these population trends and design specific and local conservation measures. We hypothesize that declining species have a lower resistance to heat stress than stable species and species living in a restricted number of habitat types are more sensitive to hyperthermic stress than widespread species.

2 Material and Methods

2.1 Sampling sites and specimen collection

We collected bumblebees for our experiment in 2016 from six locations in Belgium and one in the southwestern Netherlands selected for their species diversity (Vray 2018) along a climate gradient: (a) Kalmthout (51° 22′ N, 04° 28′ E 21m alt.); (b) Maasmechelen (50° 57′ N, 05° 41′ E 39m alt.); (c) Malchamps (50° 27′ N, 05° 54′ E 570m alt.); (d) Mons (50° 27′ N, 03° 57′ E 56m alt.); (e) Pecq (50° 41′ N, 03° 20′ E 23m alt.); (f) Petit Lanaye (50° 48′ N, 05° 41′ E 120m alt.); and (g) Korendijk (51° 48′ N, 04° 19′ E 0m alt.) in the Netherlands (Figure 1).

Figure 1.
figure 1

Sampling map including the six locations selected in Belgium (Kalmthout, Maasmechelen, Malchamps, Mons, Pecq, and Petit Lanaye) and one in Netherland (Korendijk).

Full size image

Bumblebees were collected by hand net and transported individually to the lab in a cooling system. Only males were used in our experiments for the following reasons: (i) they display simple and constant behavior in contrast to females which have different hormonal status (Heinrich 2004); (ii) they normally do not take shelter in thermo-regulated underground nests as females do (Heinrich 2004); (iii) males are at risk to be more strongly impacted by high air temperatures as they generally perform their nuptial behavior in open sunlight areas. Moreover, bumblebee males emerge in the mid-summer and are the caste to most likely experience the hottest conditions (Martinet et al. 2015a).

We collected 298 males belonging to ten different species including (i) widespread and ubiquitous species: (a) Bombus (Psithyrus) bohemicus (n = 11) from Maasmechelen, Malchamps, and Petit Lanaye; (b) Bombus (Pyrobombus) hypnorum (n = 30) from Pecq and Petit Lanaye; (c) Bombus (Melanobombus) lapidarius (n = 30) from Korendijk, Maasmechelen, Malchamps, Mons, and Petit Lanaye; (d) Bombus (Bombus s.s.) lucorum (n = 21) from Maasmechelen, Malchamps, Mons, and Petit Lanaye; (e) Bombus (Thoracobombus) pascuorum (n = 45) from Korendijk, Maasmechelen, Malchamps, Mons, and Petit Lanaye; (f) Bombus (Psithyrus) sylvestris (n = 24) from Maasmechelen, Malchamps, and Petit Lanaye; (g) Bombus (Bombus s.s.) terrestris (n = 37) from all sampled localities and (ii) localized and rare species (Vray 2018); (h) Bombus (Pyrobombus) jonellus (n = 56) from Malchamps, Kalmthout, and Maasmechelen; (i) Bombus (B. sensu stricto) magnus (n = 35) from Kalmthout; and (j) Bombus (Thoracobombus) muscorum (n = 9) from Korendijk (S1 Dataset). Specimens belonging to the same species were collected in different sites, but we did not expect intraspecific variation as populations recorded in Belgium should be considered as a consistent metapopulation without biogeographical structure (Maebe et al. 2016).

2.2 Protocol of experimentation

The following experiments were carried out according to the protocol established by Martinet et al. (2015a). To standardize experimental conditions, specimens were placed at 8 °C (standby temperature, Heinrich 1975) overnight in the dark with sugar syrup (Biogluc, Biobest, Westerlo, Belgium). Then, specimens were placed in controlled and constant conditions in an incubator (Herp Nursery II) at 40 °C with a humidity rate of 50–60%. This temperature was chosen to simulate the temperature of a heat wave (hyperthermic stress) (Martinet et al. 2015a).

Individuals were observed through a window in the incubator. The time before heat stupor (THS) was measured for each individual with a chronometer (± 1 min) as a measure of heat stress resistance. This THS corresponds to the time from insertion into the incubator until the heat stupor. An insect in “heat stupor” (Bodenheimer and Klein 1930; Uvarov 1931) falls upside down and loses its normal reflexes (Huang et al. 2006). The extremities of the legs are then shaken by muscle spasms that appear just before “heat coma” (Lutterschmidt and Hutchison 1997; Martinet et al. 2015a).

2.3 Body size measurement

We used intertegular distance (i.e., distance in millimeters between the two insertion points of the wings; ITD) as a proxy of the body size (Cane 1987). We measured ITD using a Facom 150-mm digital caliper (Morangis, France). To analyze heat stress resistance variation in relation to body size, we assessed the 298 bumblebee males used in this study (Tab. S1).

2.4 Statistical analyses

We performed statistical comparative analyses of the heat stress resistance using the R environment (R Development Core Team 2016). Differences in bumblebee heat stress resistance were visually assessed using the “boxplot” function with R (R Development Core Team 2018) (“graphics” R-package (Murrell 2005), “stats” R-package). The Kruskal-Wallis analyses and multiple comparison test were used to compare THS among species (“pgirmess” R-package, “funct kruskalmc,” Siegel and Castellan 1988). When a species was tested from different sampling sites, THS was compared using the Kruskal-Wallis comparison test.

To visualize the full variability of THS dataset with the two independent variables (biplot with “Site” and “Species”), we performed a principal component analysis (PCA) based on correlation distance matrix (library “FactoMineR,” “factoextra”) (Le et al. 2008). We also created a generalized linear model with Bonferroni’s adjustment (with analyses of variance) to assess the variation of heat stress resistance (THS) (i.e., dependent variable) in relation to body size (ITD as an independent variable) (library “stats”) for all specimens with “Species” and “Site” effect considered as a random factor in the model. We calculated linear models using the restricted maximum likelihood method and ensured that the normality of the residuals and independence requirements were satisfied (R Development Core Team 2017).

3 Results

Based on 298 tested males, our results show that heat stress resistance varied significantly between species (Kruskal-Wallis chi2 = 117.38, df = 9, p value < 0.001). Post-hoc multiple comparison Kruskal-Wallis test resulted in four groups (p value < 0.05): (i) B. terrestris (median 395 min a), (ii) B. lucorum (median 257 min ab), (iii) B. bohemicus (median 47 min bc) which has a higher resistance than (iv) B. hypnorum (median 81 min c), B. sylvestris (median 76 min c), B. muscorum (median 60 min c), B. lapidarius (median 58 min c), B. magnus (median 58 min c), B. jonellus (median 48 min c) and B. pascuorum (median 48 min c) (Figure 2, Tab. S1). No significant differences in THS were found between species collected from different sampling sites (multiple Kruskal-Wallis test, p values > 0.05) which confirms our hypothesis that in our case, populations should be considered as consistent metapopulation without biogeographical structure (Maebe et al. 2016). Exact p values (multiple comparison Kruskal-Wallis test) to compare each species are presented in Tab. S2.

Figure 2.
figure 2

Boxplots of hyperthermic stress resistance of ten bumblebee species (B. terrestris, B. lucorum, B. hypnorum, B. sylvestris, B. muscorum, B. magnus, B. lapidarius, B. bohemicus, B. pascuorum, and B. jonellus) recorded in Belgium. The line inside the boxes corresponds to the median, the lower-limit of boxes corresponds to the first quartile, and the upper-limit to the third quartile. The whisker extends up to the highest value of the upper limit and down to the lowest value of the lower limit. Circles correspond to outlier values.

Full size image

Based on linear models, there was no significant relation between the time before heat stupor (THS) and body size (ITD) throughout all specimens (r2 = 0.04, df = 297, p = 0.86, Figure 3). Unlike the random factor “Site” (F = 1.14, p = 0.97), “Species” factor had a significant effect on heat resistance (F = 43.73, p < 0.001). All these results were gathered and plotted on the ACP in Figure 3.

Figure 3.
figure 3

Principal component analysis (PCA, biplot) of “time before heat stupor or THS” and “inter tegulae distance or ITD” variation within ten bumblebee species. PC1 and PC2 are the first and the second principal component axes. Symbols represent the different bumblebee species. In blue, mean of THS per species; in green, the different sampling sites. The black arrows represent the contribution of variables (THS and ITD) to define principal component axes.

Full size image

4 Discussion

Our results show that males of different bumblebee species are not affected by heat stress in the same way. While a few studies have investigated the thermic resistance of bumblebees on a single species (Owen et al. 2013; Martinet et al. 2015a; Oyen et al. 2016; Oyen and Dillon 2018), our method allows an easy assessment of hyperthermic stress resistance in laboratory conditions for multiple species. Our results show that males of the widespread and common species, B. terrestris and B. lucorum, are less sensitive to hyperthermic stress. Their resistance time is up to 10 times longer than the least resistant species (Figure 2). This high interspecific variability in heat resistance in a homogeneous group like bumblebees (Figure 2 and Oyen et al. 2016) is also observed in other organisms (plants: Soliman and Sugiyama 2012; Drozdov et al. 1984; gastropods: Dong et al. 2008; frogs: Pintanel et al. 2019; review Tomanek 2010).

4.1 Tolerant species

Males of Bombus terrestris have the highest heat resistance recorded with a THS median of 395 min. While the physiological mechanisms involved in the hyperthermic resistance of bumblebees are still unknown, this could involve a more efficient protein system in Bombus terrestris (e.g., the rapid production of heat shock proteins as observed in Apis mellifera, Elekonich 2009). The high resistance of Bombus terrestris could be an explanation of its very large geographical distribution, as it occurs from the limits of the Sahara Desert in the south and crosses into the Arctic Circle in the north (Martinet et al. 2015b). According to Bale et al. (2002), species with a very wide distribution are exposed to very variable climatic parameters. The large climatic envelope in which Bombus terrestris is distributed as well as greater physiological plasticity could explain the high resistance of this bumblebee to hyperthermic stress (Rasmont et al. 2015). The positive population trends of this species nearly everywhere in Europe and its extension toward the North (Martinet et al. 2015b) could be linked to the especially high heat stress resistance of the species. Other specific adaptations of B. terrestris could reduce the negative impact of heat stress on its physiology (e.g., early emergence of queens, estivation and modification of the phenology depending on the seasonal conditions; Corbet et al. 1991; Rasmont et al. 2008).

Males of Bombus lucorum which belongs to the same subgenus as B. terrestris (i.e., Bombus sensu stricto) are the second most resistant species to hyperthermic stress. The high resistance of B. lucorum could be therefore explained in the same way as B. terrestris. However, the high resistance of B. terrestris and B. lucorum is probably not a subgeneric characteristic as B. magnus from the same subgenus, which shows a much lower resistance.

4.2 Sensitive species

The resistance time of the eight other assessed bumblebee species is short, not exceeding 80 min (Figure 2). The two sub-boreal species (i.e., B. jonellus and B. magnus) display the shortest heat stress resistance (median THS 48 min for B. jonellus and 58 min for B. magnus). Moreover, the inter-individual heat resistance variation in males of B. jonellus is lower (Figure 2). In Belgium, the distribution of these two species is limited to the most continental regions of the country (High Fens and Campine) (Rasmont and Iserbyt 2014; Rasmont et al. 2015). The low resistance of B. jonellus and B. magnus to heatwave temperature could be correlated to their boreal distribution. Once an organism is constrained to live in the boundary limits of its distribution or refuge zone (limits of its eco-climatic envelope), it is then exposed to conditions which make its physiological tolerance more challenging (Parmesan et al. 2000; Hance et al. 2007).

Rarity and the decline of species like B. jonellus, B. magnus, and B. muscorum could be explained by the low THS of males. With regard to the resistance of sensitive but common bumblebee species (B. hortorum, B. hypnorum, B. lapidarius, and B. pascuorum), their resistance time could allow them to survive to a very short extreme temperature event (i.e., several hours). These species are stable and common in Belgium (Vray 2018) with a very wide geographical area (from Mediterranean shores to the south and reaching the Barent Sea to the North). Males of these species nevertheless display a relatively low heat stress resistance.

4.3 Bumblebee conservation

Even if it has been observed that bumblebee males exposed to hyperthermic stress take refuge in the shadows to avoid excessive exposure (Heinrich 1999), heat waves could be problematic for bumblebee populations. A decrease in the male population could induce a bottleneck by reducing genetic variability (e.g., review Chakraborty and Kimmel 2001). Several authors observed in Europe a drastic regression of the montane bumblebee populations (e.g., B. monticola, B. mucidus, B. magnus, B. veteranus) during the years following an extreme climatic event (Manino et al. 2007; Rasmont and Iserbyt 2012). The regression of B. muscorum and B. magnus recorded in western Europe may be partially induced by heat waves. Indeed, the colony development and the late mating of this species occur during months where heat waves are more frequent (Rasmont and Iserbyt 2012).

It should be noted that our assessment of heat resistance is based on only one caste (even if particularly crucial for the sustainability of species) and one subpopulation of ten species (potentially thermo-differentiated). Since there are differences between males and other sex/castes in body size and physiology, further studies on queens and workers should be conducted to compare thermal sensitivity between castes in social bees. Local adaptation at a large spatial scale is also possible like it was demonstrated in other organisms (Schoville et al. 2012; Yu et al. 2018; Collier et al. 2019). However, as suggested by our results, there is no significative relation between body size of males (ITD), collecting sites, and the corresponding THS recorded for the specimen.

As bumblebees display variable thermal sensitivity depending on species, conservation plans and protection status should pay special attention to the most sensitive species (e.g., B. jonellus and B. magnus) that are already in decline. In order to preserve these species, habitats (e.g., heathlands) and nutritional resources (suitable pollen quality) of sensitive and declining bumblebee species should be conserved. To minimize the impact of heat stress on localized populations, we support the establishment of multiple and connected high-quality habitat to improve bumblebee conservation (Vray 2018). It has been shown that high-quality resources can mitigate the effect of heat stress (Vanderplanck et al. 2019). Our study also raises a question about the use of a single bumblebee species in studies focusing on environmental stressors (i.e., nutritional stress, pesticides, parasites). Bombus terrestris is the most widely used species in several bumblebee studies, but it is also the most resistant to hyperthermic stress. It is therefore important, in future research, to integrate other species in our framework for a better understanding of the impact of main decline factors on wild bee populations.