# Static landscape features predict uplift locations for soaring birds across Europe

Martina Scacco, Andrea Flack, Olivier Duriez, Martin Wikelski, Kamran Safi

► To cite this version:

Martina Scacco, Andrea Flack, Olivier Duriez, Martin Wikelski, Kamran Safi. Static landscape features predict uplift locations for soaring birds across Europe. Royal Society Open Science, 2019, 6 (1), pp.181440. 10.1098/rsos.181440. hal-02394114

HAL Id: hal-02394114

https://hal.umontpellier.fr/hal-02394114

Submitted on 4 Dec 2019

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Research

<figure><img src="image_2.png" /></figure>

<figure><img src="image_3.png" /></figure>

Check for
updates

**Cite this article:** Scacco M, Flack A, Duriez O, Wikelski M, Safi K. 2019 Static landscape features predict uplift locations for soaring birds across Europe. *R. Soc. open sci.* **5**: 181440.
http://dx.doi.org/10.1098/rsos.181440

Received: 31 August 2018

Accepted: 30 November 2018

**Subject Category:**

Biology (whole organism)

**Subject Areas:**

ecology

**Keywords:**

habitat suitability, movement ecology, random forest, species distribution model, anthropogenic infrastructure, energy landscape

**Author for correspondence:**

Martina Scacco

e-mail: mscacco@orn.mpg.de

Electronic supplementary material is available online at https://doi.org/10.6084/m9 figshare.c. 4334345.

# Static landscape features predict uplift locations for soaring birds across Europe

Martina Scacco<fnref n="1" />, Andrea Flack<fnref n="1" />, Olivier Duriez<fnref n="2" />, Martin Wikelski<fnref n="1" />,<fnref n="3" /> and Kamran Safi<fnref n="1" />

<figure><img src="image_4.png" /></figure>

MS, 0000-0003-4821-4051; AF, 0000-0002-9099-2802; OD, 0000-0003-1868-9750; KS, 0000-0002-8418-6759

Soaring flight is a remarkable adaptation to reduce movement costs by taking advantage of atmospheric uplifts. The movement pattern of soaring birds is shaped by the spatial and temporal availability and intensity of uplifts, which result from an interaction of local weather conditions with the underlying landscape structure. We used soaring flight locations and vertical speeds of an obligate soaring species, the white stork (*Ciconia ciconia*), as proxies for uplift availability and intensity. We then tested if static landscape features such as topography and land cover, instead of the commonly used weather information, could predict and map the occurrence and intensity of uplifts across Europe. We found that storks encountering fewer uplifts along their routes, as determined by static landscape features, suffered higher energy expenditures, approximated by their overall body dynamic acceleration. This result validates the use of static features as uplift predictors and suggests the existence of a direct link between energy expenditure and static landscape structure, thus far largely unquantified for flying animals. Our uplift availability map represents a computationally efficient proxy of the distribution of movement costs for soaring birds across the world's landscapes. It thus provides a base to explore the effects of changes in the landscape structure on the energy expenditure of soaring birds, identify low-cost movement corridors and ultimately inform the planning of anthropogenic developments.

## 1. Introduction

All animals interact with the surrounding environment, but for some of them the role of this environment becomes particularly

<footnote n="1">¹Department of Migration and Immuno-ecology, Max Planck Institute for Ornithology, Am Obstberg 1, 78315 Radolfzell, Germany</footnote>

<footnote n="2">²Centre d'Ecologie Fonctionnelle et Evolutive, UMR 5175 CNRS-Université de Montpellier-EPHE-Université Paul Valery, 1919 Route de Mende, 34293 Montpellier cedex 5, France</footnote>

<footnote n="3">³Department of Biology, University of Konstanz, Universitätsstr. 10, 78464 Konstanz, Germany</footnote>
relevant in constraining or supporting their movement. This especially applies to aerial or aquatic animals, whose movements actively modify and are, in turn, modified by the surrounding fluid [1–3]. Air does not provide constant support against gravity and its properties vary at different temporal and spatial scales depending on turbulence. To save energy, flying animals therefore adjust timing, routes and flight modes to this turbulence [4], maximizing the advantage of horizontal and vertical air currents [5].

Soaring birds represent an extreme example of this adaptation. These large and heavy birds are particularly constrained in the use of active flapping flight, as the energetic cost of flight proportionally increases with size and weight [6]. They therefore use passive soaring-gliding flight, which is subsidized by the vertical air currents (uplifts) and may require as little energy as resting [7,8]. Body mass, wing loading and wing aspect ratio ultimately determine the cost of flapping flight, and with it a species' degree of dependence on uplifts [6,9–11]. This dependence becomes extreme in obligate soaring birds, which, due to their large size, can only fly in good uplift conditions, minimizing the use of flapping flight [12,13].

Uplifts originate from thermal convection (thermals) and/or mechanical sources (orographic uplift) [12,14]. Thermals originate from uneven heating of the earth's surface, with rapidly heated areas producing a gradient of temperature which promotes the formation of rising columns of warm air. Orographic uplifts, by contrast, result from the deflection of horizontal wind through topographic features, such as hills or ridges [12]. Thus, the occurrence of both thermal and orographic uplifts depends on a combination of local weather conditions (gradient of air temperature, wind speed and direction) and landscape features (land cover, slope inclination and aspect, roughness of the surface). The interaction of local weather conditions and landscape features thus characterizes a complex and turbulent atmospheric layer, where the spatial and temporal availability of uplifts constrains and shapes the movement patterns of soaring birds, from local scale to migratory routes [3,15–17].

Over the past decade, different models have been developed to investigate the relationship between soaring behaviour and aerial environment [8,18–21]. In many of these studies, the availability of uplifts was indirectly inferred using several weather parameters [22–25], but in recent years, these parameters have been replaced by thermal and orographic uplift potentials, as more direct estimators to quantify the probability of soaring [14,18–20,26–29]. However, some studies highlighted the inaptitude of these newly introduced variables as uplift estimators, because of the large amount of unexplained variance remaining when predicting soaring behaviour [18,19]. In fact, thermal and orographic uplift potentials are calculated based on different weather parameters [14,19], but because the uplift events are characterized by turbulences at fine spatio-temporal scale [12,19] it is challenging to predict their occurrences due to the limited spatio-temporal resolution of the available weather products.

In contrast to weather products, publicly available satellite data provide valuable static landscape information (such as land cover and elevation) at higher spatial resolution, which could be used to predict the occurrence of uplifts. Soaring birds need to locate uplifts in order to move across the landscape. Consequently, landscape features that influence uplift generation, might serve as visual cues to these birds, as they do for hang glider or paraglider pilots. By determining landscape features that birds use to locate uplifts, we may be able to predict those uplifts that are detectable and exploited by the birds.

Static landscape features alone could therefore potentially suffice in modelling the occurrence of uplifts, providing an answer to 'where', albeit not 'when', uplifts are likely to occur. Although the literature on the topic is scarce, some studies hinted at the role of static features in affecting the flight behaviour of different soaring species [25,28,30–33].

Here, we investigate to what extent static landscape features can represent the potential for generating uplifts. We explore this in an obligate soaring bird species, the white stork *Ciconia ciconia*, across the entire continent of Europe. We used first the locations of soaring and flapping behaviours of storks as an indication of the presence or absence of uplifts, and second their vertical speed as a proxy of uplift intensity. We then used only static features of the landscape to model and predict the spatial distribution of uplifts and their intensity across Europe. We also evaluated the effectiveness of these two static models by comparing their performances with the performances of two dynamic models, which included atmospheric uplift estimators used in previous studies. Finally, we explored the cost of flight (in terms of overall dynamic body acceleration) over the considered area, only based on the static landscape features. Under the assumption that soaring/flapping behaviour and vertical speed of the birds can be used as sensors of availability and intensity of uplifts, we predicted that (i) static features of the landscape (such as topography and land cover) can be used to predict the spatial availability of uplifts and to produce a static uplift suitability map at European scale; (ii) areas detected as suitable for uplifts during the first step can be further characterized in terms of uplift intensity likely to be produced in those areas; (iii) the resulting static uplift suitability map corresponds to the spatial distribution of the energetic costs of storks flying above this landscape, and thus portrays their static energy landscape.
<figure><img src="image_5.png" /><figcaption>**Figure 1.** Spatial coverage of the white storks' migration routes, relative to the extent of the environmental layers included in the model. Black lines correspond to individual stork GPS trajectories (a). Static uplift prediction maps produced using the uplift suitability model (b) and the uplift intensity model (c), projected outside the geographical range of the training set. In (b), the colour scale corresponds to uplift suitability, as predicted by the uplift suitability model; grey indicates suitable and black unsuitable cells. White represents unclassified cells (containing missing values among the predictors). In (c), cells predicted as suitable are further characterized by the predicted uplift intensity values. Colour scale corresponds to vertical speed ranging from red (high) to light yellow (low). As in (b), in (c) black represents cells that are unsuitable for uplift and white indicates unclassified cells. Latitude values outside the range of the training set were excluded from the intensity model. The two prediction maps are available at https://dx.doi.org/10.17617/3.1u.</figcaption></figure>

## 2. Methods

### 2.1. Dataset

The GPS and tri-axial accelerometry (ACC) data used in the study were collected by the Max Planck Institute for Ornithology (see [34,35]) and are deposited in the Movebank Data Repository (http://dx.doi.org/10.5441/001/1.bj96m274 [36]).

The animals were equipped as fledglings with high-resolution, solar GSM-GPS-ACC loggers (e-obs GmbH, Munich, Germany). The dataset includes 61 juvenile white storks (*Ciconia ciconia*) during their first migration (figure 1a). Because storks are diurnal, loggers provided one GPS location every 5 min between 2.00 and 20.00 GMT. If instantaneous ground speed was greater than 2 m s⁻¹, bursts of high-resolution GPS locations (1 Hz) were recorded every 15 min for 120 or 300 s. In addition to the GPS locations, ACC was recorded every 10 min for a duration of 3.8 s at a sampling rate of 10.54 Hz (40 data points per axis). High-resolution GPS recordings were collected from August to September 2014.

### 2.2. Segmentation of the flight behaviour

#### 2.2.1. Soaring flight (from GPS)

We selected high-resolution GPS bursts with a duration of at least 120 s. For each location in the burst, we calculated vertical speed and turning angle. We applied our behavioural segmentation on track segments of 15 s (average duration of one complete soaring circle [34]). We calculated the average vertical speed
<figure><img src="image_6.png" /><figcaption>**Figure 2.** Example of behavioural segmentation based on the GPS data of one stork. The figure shows the classified three-dimensional trajectory after smoothing; the red segment was identified as soaring flight, grey corresponds to gliding flight. Data for plotting the surface are provided by the EU-DEM. The black line and the red point on the ground represent the two-dimensional projection of the trajectory and the centroid of the soaring segment, respectively.</figcaption></figure>

and the absolute cumulative turning angle in these segments, and we used the expectation maximization binary clustering (EmbC) algorithm to discern the flight behaviours, introducing these two metrics as delimiters. The algorithm, implemented in the R package *EmbC* [37], efficiently detected changes in the flight behaviour, distinguishing a high turning angle (circular soaring) from two low turning angle clusters (linear flight). Based on the average vertical speed, we further differentiated the linear flight segments into gliding (linear descending flights) and linear soaring (linear ascending flights). Each 15 s segment along the animal trajectory was individually assigned to one of the behavioural classes based on its specific parameters. We applied a smoother to avoid abrupt and unnatural behavioural changes from one segment to the next along the same trajectory. Our smoother worked as a moving window: each segment assigned to a different behaviour relative to its closest neighbours was reclassified to match the modal value of two segments before and after the considered segment.

Given the high resolution of the GPS data, we could visually inspect and confirm the results of the segmentation using three-dimensional plots (figure 2). We then investigated the different classified behaviours in terms of their flight parameters, such as ground speed or vertical speed (electronic supplementary material, S1, figure S1.1).

In the subsequent steps, we wanted to contrast the use of active versus passive flight, focusing on the dichotomy soaring/flapping. We therefore did not differentiate between circular and linear soaring (both classified as soaring), and we excluded gliding segments, as they are not considered as an alternative to soaring (like flapping) but rather as its consequence [12]. In these analyses, we considered for each individual only soaring segments with a duration longer than 30 s, and treated consecutive soaring segments as different units only when separated by at least 60 s. The location of each soaring segment was defined by its centroid (mean longitude and latitude).

### 2.2.2. Flapping flight (from tri-axial accelerometry)

We interpolated the spatial location of each ACC burst based on the closest GPS locations using the R package *move* [38]. We associated each ACC burst with the height above ground corresponding to the GPS location closest in time (less than 30 s difference). The height above ground was calculated by subtracting the ground elevation value from the height above the ellipsoid.

We used ACC values to identify bursts of active flight behaviour (flapping flight). Specifically, we used overall dynamic body acceleration (ODBA), already shown to be a good proxy for energy
expenditure in soaring birds [8,39]. We quantified ODBA and dynamic body acceleration (DBA) on the z-axis following Wilson *et al.* [40] and calculated mean, sum and standard deviation, of these two variables per burst. We then used k-means clustering to categorize the bursts into three main behavioural classes based on the amount of activity recorded: least active, intermediate active and most active (electronic supplementary material, S1, figure S1.2). Within the bursts of highest activity, we wanted to isolate only the flapping behaviour marking the absence of uplifts (and to exclude, for instance, the flapping associated with taking off); we thus applied a height threshold of 100 m above ground to select our flapping locations, assuming that above this height the birds were using flapping flight only in response to the absence of uplifts.

The two flight behaviours were classified based on data collected with different instruments running on different sampling schedules (GPS for soaring and ACC for flapping). Therefore, the amount of soaring to flapping locations is not directly related to the amount of time storks spent on each flight behaviour. We thus compared the amount of time spent soaring relative to the total duration of the classified GPS segments, and the amount of time spent flapping relative to the total duration of the classified ACC bursts.

## 2.3. Environmental variables and modelling frameworks

### 2.3.1. Static predictors

We characterized the static components of the landscape in terms of elevation (digital elevation model, EU-DEM 2013), terrain unevenness (calculated as both topographic position index and roughness), unevenness in the slope (steepness of a terrain feature), aspect (compass direction faced by a slope), aspect unevenness, land cover (normalized difference vegetation index, NDVI, obtained for 2014), land use (CORINE Land Cover; CLC 2012) and presence of anthropogenic infrastructures (Global Urban Footprint, 2011). All raster layers are publicly available (electronic supplementary material, S2, table S2.1). The lowest spatial resolution was 100 m (from the CLC 2012 layer), thus we averaged cell values of higher resolution layers to match a 100 m grid. The spatial extent of the raster layers covered the southwest European countries that enclose the distribution of the storks' dataset. All the environmental layers listed above were included as predictors in our statistical models after verifying the absence of multicollinearity.

### 2.3.2. Dynamic predictors

We chose to include thermal and orographic uplift potentials in our analysis as atmospheric uplift estimators [14,18,20,26,27]. The calculation of the thermal uplift potential is based on weather data from the European Centre for Medium-range Weather Forecast Global Atmospheric Reanalysis (ECMWF) following Bohrer *et al.* [19]. The calculation of the orographic uplift potential uses ECMWF weather data and elevation from the ASTER Global Digital Elevation Model. Both thermal and orographic uplift potential are available in Movebank with a spatial resolution of 0.75° and a temporal resolution of 6 h. We associated them to our tracking data by using the Env-DATA Track Annotation service [41].

## 2.4. Modelling framework

We organized the environmental predictors in three groups, each defining a different modelling framework:

1. Static model: including exclusively static environmental predictors;

2. Dynamic model: including exclusively thermal and the orographic uplift potentials;

3. Combined model: including both static and dynamic predictors.

We used these sets of predictors for both the uplift suitability and the uplift intensity models.

## 2.5. Uplift suitability model

We used a random forest algorithm to model the effect of the three sets of environmental predictors on the occurrence of soaring (presence of uplifts) and flapping flight (absence of uplifts), using these contrasting behaviours as binary response variable. The algorithm is implemented in the R package
randomForest [42]. We manipulated the ratio between presences and absences (prevalence) and tested its effect on the model performance (see electronic supplementary material, S3). In our analysis, we included all the available data with their original (unmanipulated) prevalence values. Using regression trees, we trained each of the three models (static, dynamic and combined) with 90% of the dataset, and tested them with the remaining randomly selected 10%. The data partitioning was repeated so that each of the three models was run ten times. To evaluate and compare the predictive performance of the three models, we considered the following accuracy measures: (i) area under the curve (AUC) of the receiver operating characteristic (ROC); (ii) sensitivity, proportion of soaring locations correctly classified; (iii) specificity, proportion of flapping locations correctly classified [43]; (iv) true skill statistics (TSS: $1 - \max(\text{sensitivity} + \text{specificity})$) [44]. The contribution of each environmental variable to the final prediction was evaluated using the decrease in accuracy (increase in mean standard error) and the increase in node purity (decrease in residual sum of squares).

Next, we produced a large-scale uplift suitability map based on the static uplift suitability model. Random forest, like other machine learning algorithms, is quite unreliable when extrapolating outside the range of the predictors' values provided for training. We thus omitted (set to null) all raster cells containing environmental values outside that range, and then used each of ten runs of the static uplift suitability model to predict the uplift suitability over the area of these manipulated raster layers. The final raster prediction was derived from the pixel average of the ten predicted layers and classified into a binary map using the threshold that maximized the TSS value [43].

The only temporally related environmental variable in our model was NDVI from the year 2014; this allowed us to produce an uplift suitability map for 2014.

## 2.6. Uplift intensity model

We explored the relationship between the three sets of predictors (static, dynamic and combined) and uplift intensity, additionally including latitude among the static predictors. We used the birds' vertical speed as a proxy of uplift intensity (vertical rate of air within a thermal), assuming a higher vertical speed to indicate stronger uplift conditions.

As in the previous analysis, we considered only high-resolution GPS bursts. The vertical speed in this dataset included both negative (gliding) and positive values (soaring). We examined only the positive values (vertical speed greater than 0), because we wanted to predict uplift intensity in areas already classified as suitable by the uplift suitability model. We associated the positive vertical speed values of all individuals with their location and averaged them in a 100 × 100 m grid to match the spatial resolution of the environmental raster layers. After averaging, each cell contained a value representing the average vertical speed of all individuals during the complete temporal range in that cell. We then removed average vertical speed values exceeding the 99.97 percentile, obtaining 76 383 observations.

We used a generalized additive model (GAM) to model uplift intensity (average vertical speed) as a function of the three sets of environmental predictors to accommodate nonlinear relationships between predictors and response variable. We square-root transformed vertical speed to meet the assumptions of a Gaussian distribution of the residuals. Among the predictors, aspect was included as cyclic cubic regression spline smooth term; NDVI, elevation (DEM), roughness and latitude were included as thin plate regression spline smooth terms, given their nonlinear relationship with the response variable. We rasterized the values of thermal and orographic uplift potentials included in the dynamic and combined models to match the 100 × 100 m grid of the response variable, and included them in the models as parametric coefficients.

We used the static uplift intensity model to produce a map of uplift intensity, and enrich the binary information provided by the uplift suitability model in areas that were predicted as suitable. Raster cells containing environmental values outside the range included in the dataset were omitted; because latitude was included as predictor, the latitudinal range of the uplift intensity map was restricted to the latitudinal range of the dataset.

The models were run in R using the package *mgcv* [45]. We compared the performances of the models based on the variance explained and the AIC (Akaike Information Criterion). The relative importance of the different predictors was evaluated comparing the AIC of models containing different combinations of these predictors, computed separately for parametric coefficients and the smooth terms.

## 2.7. Static energy landscape

We quantified the relationship between the availability of uplifts along the storks' migratory routes and the energy spent travelling along these routes. We could thus test if the static uplift suitability map
produced in the previous step could convey information regarding the energetic cost of travelling across the landscape. We calculated the daily energy expenditure considering only ACC data collected when the animals were flying (height above ground > 100 m), with GPS location and ACC burst matching in time. We then calculated the mean ODBA per day along the path of each individual. The uplift suitability map was used to extract the predicted probability of uplift at the locations of the ACC bursts. We then averaged these probability values to obtain the mean daily uplift suitability, for each individual, along its migratory path (only average values computed from at least five observations were included in the model). We fitted a linear mixed effect regression model to the mean daily energy expenditure (ODBA) as a function of predicted mean daily uplift suitability. The model tested the relationship between daily uplift suitability predictions and daily ODBA based on 823 observations of 59 individuals, accounting for individual differences which were included as random effects in the model. ODBA was square-root transformed.

The importance of the predictor in explaining the daily energy expenditure was assessed comparing the AIC of the model with the respective null model. For the analysis, we used the R package *lme4* [46].

# 3. Results

## 3.1. Segmentation of the flight behaviour

We identified the location of soaring and flapping flights as proxies to detect the presence and absence of uplifts. Based on the GPS data of all individuals, we classified over 748 h of flight, of which the storks spent 297.6 h with circular and 83.5 h with linear soaring. The proportion of time spent soaring corresponded to 0.51 of their flight time (381.1 h) (electronic supplementary material S1, figure S1.3a); this proportion was similar between the 59 individuals (0.52 ± 0.07 (mean proportion ± s.d. per individual)). From the ACC data, we classified 24.3 h of flight, of which 1.3 h was spent flapping (electronic supplementary material, S1, figure S1.3b). Among all individuals, the proportion of time spent flapping corresponded to 0.05 (0.07 ± 0.05 per individual). The final dataset consisted of a total of 16 840 observations of presences and absences of uplift (15 608 soaring events marking presences and 1232 flapping events marking absences).

## 3.2. Uplift suitability model

We used multiple environmental predictors to model and predict the spatial distribution of uplifts (presence and absence data). We organized the predictors in three different modelling frameworks (static, dynamic and combined, see Methods) that we then compared in terms of predictive accuracy.

We averaged the accuracy measures of the three uplift suitability models across 10 cross-validations. The combined model (static and dynamic features) best predicted the independent test set (AUC of 0.86 ± 0.02 (mean ± s.d.)), followed by the static (AUC of 0.85 ± 0.02) and the dynamic (AUC of 0.70 ± 0.02) models (figure 3). The overall accuracy was high in all models, but both models including static variables (the static and combined models) outperformed the model based only on dynamic predictors.

We then compared the ability of the models to discriminate presences and absences. The three models returned a similarly high proportion of correctly classified soaring locations (sensitivity). They differed, however, in terms of number of correctly classified flapping locations (specificity). Again, the combined and the static model outperformed the dynamic model. To define the value of sensitivity and specificity, we used a threshold that maximized the TSS value, corresponding to 0.9 in all models. At this threshold, the static model showed a sensitivity of 0.83 ± 0.01 and a specificity of 0.69 ± 0.05. The complete output of the three models can be found in electronic supplementary material, S3, table S3.1.

In the static model, DEM, roughness (topographic heterogeneity) and NDVI contributed most to the model prediction. In the dynamic model, including only the two atmospheric uplift estimators, the thermal uplift potential contributed to the prediction more than the orographic potential; in the combined model, elevation, roughness and thermal uplift potential contributed most to the model (for more details see electronic supplementary material, S3).

Using the static uplift suitability model, we produced a map of uplift suitability covering the extent of the environmental layers (figure 1b). We classified an area of about 193 million km², of which over 81 million km² was predicted as suitable for uplifts (42% of the total area).
<figure><img src="image_7.png" /></figure>

<figure><img src="image_8.png" /></figure>

<figure><img src="image_9.png" /></figure>

**Figure 3.** ROC curves of the three uplift suitability models: static (a), dynamic (b) and combined (c). The area under the curve (AUC) represents the accuracy of the model. The accuracy was measured on both the training set (grey solid line) and the test set (red line). The dashed line represents a model whose accuracy is comparable to random (AUC = 0.5). Sensitivity and commission rate values were averaged across the 10 runs of each models (solid dots), and the error bars show their standard deviations.

## 3.3. Uplift intensity model

We then used the uplift intensity to characterize those areas identified as suitable for uplifts by the static suitability model. We used the vertical speed of the birds while soaring as a proxy for uplift intensity, and we explored the relationship between uplift intensity and the three groups of environmental predictors (static, dynamic and combined).

All three models explained very little of the total variance in vertical speed. However, here too, the combination of static and dynamic variables provided the best predictive performance (Adj.*R*² = 0.03 and AIC = 46 575.00 for the static model; Adj.*R*² = 0.03 and AIC = 49 636.88 for the dynamic model; Adj. *R*² = 0.08 and AIC = 42 003.95 for the combined model). Although in GAMs Adj. *R*² values cannot be directly compared due to the changing degrees of freedoms caused by the use of smooth terms [47], the difference in the AIC value among the three models supports the best performances of the combined model. Among the parametric predictors, the categories ‘water bodies’, ‘dumps’, ‘urban areas’ and ‘wetlands’ negatively affected uplift intensity (‘bare soil’ served as a reference), whereas thermal and orographic uplift potentials (included in the dynamic and combined models) positively affected uplift intensity (electronic supplementary material, S4, table S4.1). Aspect, NDVI, DEM, roughness and latitude were included in the models as smooth terms, given their nonlinear relationship with the response variable. Based on AIC, all these predictors contributed to explain uplift intensity. Specifically, uplift intensity was positively affected by lower latitude values, higher elevations (DEM > 2000 m), NDVI corresponding to bare soils or sparsely vegetated areas (between 0 and 0.4) and slope orientation towards SW-W (aspect between 200° and 300°) (electronic supplementary material, S4, table S4.1 and figure S4.1). Using the static intensity model, we could further characterize our uplift suitability map by predicting uplift intensity in cells already predicted to be suitable for uplifts (based on the static suitability model) (figure 1c).
<figure><img src="image_10.png" /><figcaption>**Figure 4.** Mean daily ODBA along the storks' routes as a function of the mean daily uplift suitability along those routes, as predicted by the static uplift suitability model. Coloured lines represent the different individuals accounted for in the linear mixed model. Grey points represent all observations included in the model.</figcaption></figure>

## 3.4. Static energy landscape

Finally, we quantified the relationship between the availability of uplifts along the storks' migratory routes and the energy spent travelling along these routes, to test if the static maps produced in the previous steps could convey information regarding the energetic cost of travelling across the landscape. As all uplift intensity models performed poorly in predicting the intensity of uplifts, only the uplift suitability model was included in this step. We used a linear mixed effect regression model to evaluate the role of the static uplift suitability model in conveying information about the energy expenditure of the birds (measured as daily ODBA). A negative correlation between the daily uplift suitability and the daily ODBA indicated that the birds spent more energy when flying over areas less suitable for uplifts ($ODBA = -0.67 \pm 0.07$ (estimate $\pm$ s.e.)) (figure 4). The AIC of the models was lower compared to that of the respective null model ($\Delta$AIC uplift suitability model = -54.67) (electronic supplementary material, S5, table S5.1).

## 4. Discussion

Static features of the landscape proved to be highly effective in identifying areas suitable for uplifts. However, neither static nor dynamic variables could predict the intensity of uplifts occurring in those areas. The uplift suitability predicted along the birds' migratory route using only static features, showed a clear negative relationship with the ODBA of individuals flying over those areas, indicating that birds encountering fewer uplifts along their routes experienced higher energy expenditures. This overall result validates the reliability of our static uplift suitability model, and suggests the existence of a mechanistic relationship between static landscape and energy expenditure of flying animals. We therefore propose that the static uplift availability map produced with our model corresponds to the birds' cost of transport across the landscape and can thus be considered a representation of the static energy landscape of these birds.

The possibility to describe the cost of transport in a dynamic aerial environment, only based on static features of the landscape, supports the idea that the structure of the landscape at different spatial scales could be considered as the ultimate cause for the uplifts to occur. For instance, a specific topography could represent a necessary (but not sufficient) condition of uplift occurrence, or, in other words, the potential of the landscape to produce uplifts. By contrast, local weather conditions that interact with a
specific landscape, could be considered as the proximate cause for the occurrence of the uplift, which can
define, given a suitable landscape, the temporal scale at which the uplift will, in fact, exist.

In the case of the uplift intensity, even though models including static features performed slightly better than those including only dynamic variables, the large amount of unexplained variance in all models suggested that neither static nor dynamic environmental variables were good predictors for uplift intensity. In our models, we used the birds' vertical speed as a proxy for uplift intensity. But birds' vertical speed is not only affected by uplift intensity. Their relationship is modulated by the aerodynamic performance of the bird (including wing morphology), and also by its social interactions and motivation. The ability to adjust the vertical speed within a thermal requires experience [48,49]. The storks included in our study were all juveniles during their first migration, but even among individuals of the same age, individual differences in flight performances exist, also in relation to the role of the individual within the group (leader or follower) [35]. The vertical speed of the birds might also be affected by their internal motivation to move (foraging versus migrating). During migration, birds are expected to maximize their vertical speed and travelled distance, whereas while foraging they might adopt different strategies, for instance attempting to maintain lower altitudes [12], which could explain the negative effect of land use categories such as dumps or pastures in the uplift intensity models. All these aspects could have affected the observed vertical speed of the birds and thus caused the inconsistent relationship between uplift intensity and environmental variables in our models. The spatial and temporal scales at which the uplift phenomenon occurs might have also contributed to this inconsistency. Uplifts are a turbulent and unpredictable phenomenon and they can occur at very small scale [5] as in the extreme cases of lifts produced by lines of buildings or flared methane vents [20,50,51]. The intensity of an uplift, more than the presence of an uplift, is strongly influenced by wind speed, wind direction and temperature, and thus more subject to the temporal and spatial variation of these dynamic variables.

The inadequacy of the spatio-temporal resolution of dynamic uplift estimators is not new [18], and the coarse resolution of the atmospheric data could also explain why all models including only dynamic variables performed worse than those including static variables alone, in predicting both uplift availability and uplift intensity. Nevertheless, the effect of some of the static variables included in our uplift intensity models hinted at a dependence of the uplift intensity (as in the case of the uplift availability) on the static landscape structure. This is the case, for instance, of the negative effect of water bodies, and the positive contribution of higher elevations and NDVI values corresponding to barren soils, on the uplift intensity. Also, lower latitude values positively affected uplift intensity; this result suggests a stronger thermal activity at lower latitudes, but could also indicate that young storks improved their flight performance along the route.

The static landscape features used to produce our static uplift maps are definitely not exhaustive to describe the complex fluid medium in which flying animals move, but they could represent a sufficient and efficient proxy (in terms of computational simplification) in areas and seasons where weather conditions are rather stable. The prediction maps produced by our static models are based on data from one species collected during one migratory season, but the same models could be extended to multiple soaring species and different seasons in order to generalize predictions. Such prediction maps could be used as base layers for further movement ecology analyses, and combining them with dynamic variables could provide a more accurate description of the energy available at a specific moment.

Static energy landscapes can also direct our attention to the vulnerability of flying animals to changes happening at the ground level. Anthropogenic changes in the landscape, such as deforestation, construction of wind farms and powerlines, but also roads, lines of buildings and tree rows, irrigation, or mining, could all be affecting the atmospheric environment, at a finer scale than the available weather products could possibly detect. The tight dependence of soaring birds on uplift conditions makes them particularly sensitive to changes in the landscape [19,52], in particular to anthropogenic infrastructures [10,53,54]. Our study suggests that these small changes in the landscape could affect the energy expenditure of these animals, and potentially their cost of transport over time. The static structure of the landscape and the energetic implications of changes happening on the ground should therefore be taken into account when investigating movement at larger scale, such as migratory flyways and population connectivity, and when evaluating the impact of anthropogenic infrastructure. Future studies should also focus on the interplay between vertical speed, uplift intensity and environment. By disentangling the various factors affecting this relationship, we could not only predict, based on the animal behaviour and the landscape, the quantity (availability) but also the quality (intensity) of the uplifts, and we could provide a more accurate estimation of the energetic cost of movement across the landscape.
Data accessibility. The data used for this study are deposited in the Movebank Data Repository: (http://dx.doi.org/10.5441/001/1.bj96m274) [36].

Authors' contributions. M.S., K.S. and M.W. conceived and designed the study, A.F. provided the tracking data, M.S. and K.S. carried out the analyses and wrote the first draft of the manuscript, and all authors contributed suggestions and text to subsequent drafts. All authors gave final approval for publication.

Competing interests. The authors declare that they have no competing interests.

Funding. We acknowledge funding from the Max Planck Institute for Ornithology. M.S. was supported by the German Academic Exchange Service (DAAD) and by the International Max Planck Research School for Organismal Biology.

A.F. was supported by the German Aerospace Center (DLR) and the Christiane Nüsslein-Volhard Stiftung.

Acknowledgements. We are grateful to Dina Dechmann for her valuable comments and discussions during different phases of the study and Wolfgang Fiedler for contributing during the data collection. Environmental products courtesy of European Environmental Agency (for the EU-DEM and the CLC), the Deutsches Zentrum für Luft-und Raumfahrt (for the Global Urban Footprint) and the U.S. Geological Survey Earth Resources Observation and Science Center (for Landsat Spectral Indices).

## References

1. Dickinson MH, Farley CT, Full RJ, Koehl MAR, Kram R, Lehman S. 2000 How animals move: an integrative view. *Sci. New Ser.* **288**, 100–106. (doi:10.1126/science.288.5463.100)

2. Kokshaysky NV. 1979 Tracing the wake of a flying bird. *Nature* **279**, 146–148. (doi:10.1038/279146a0)

3. Shepard ELC, Wilson RP, Gareth Rees W, Grundy E, Lambertucci SA, Vosper SB. 2013 Energy landscapes shape animal movement ecology. *Am. Nat.* **182**, 298–312. (doi:10.1086/671257)

4. Bowlin MS, Wikelski M. 2008 Pointed wings, low wingloading and calm air reduce migratory flight costs in songbirds. *PLoS ONE* **3**, e2154. (doi:10.1371/journal.pone.0002154)

5. Shepard ELC, Ross AN, Portugal SJ. 2016 Moving in a moving medium: new perspectives on flight. *Phil. Trans. R. Soc. B* **371**, 20150382. (doi:10.1098/rstb.2015.0382)

6. Pennycuick CJ. 1972 Soaring behaviour and performance of some birds, observed from a motor-glider. *Ibis* **114**, 178–218. (doi:10.1111/j.1474-919X.1972.tb02603.x)

7. Baudinette RV, Schmidt-Nielsen K. 1974 Energy cost of gliding flight in herring gulls. *Nature* **248**, 83–84. (doi:10.1038/248083b0)

8. Duriez O, Kato A, Tromp C, Dell'Omo G, Vysotski AL, Sarrazin F, Ropert-Coudert Y. 2014 How cheap is soaring flight in raptors? A preliminary investigation in freely-flying vultures. *PLoS ONE* **9**, e84887. (doi:10.1371/journal.pone.0084887)

9. Hedenstrom A. 1993 Migration by soaring or flapping flight in birds: the relative importance of energy cost and speed. *Phil. Trans. R. Soc. B* **342**, 353–361. (doi:10.1098/rstb.1993.0164)

10. Nourani E, Yamaguchi N. 2017 The effects of atmospheric currents on migratory behavior of soaring birds: a review. *Ornithol. Sci.* **16**, 5–15. (doi:10.2326/osj.16.5)

11. Agostini N, Panuccio M, Pasquaretta C. 2015 Morphology, flight performance, and water crossing tendencies of Afro-Palearctic raptors during migration. *Curr. Zool.* **61**, 951–958. (doi:10.1093/czoolo/61.6.951)

12. Kerlinger P. 1989 *Flight strategies of migrating hawks*. Chicago, IL: University of Chicago Press.

13. Spaar R. 1997 Flight strategies of migrating raptors; a comparative study of interspecific variation in flight characteristics. *Ibis* **139**, 523–535. (doi:10.1111/j.1474-919X.1997.tb04669.x)

14. Péron G, Fleming CH, Duriez O, Fluhr J, Itty C, Lambertucci S, Safi K, Shepard ELC, Calabrese JM. 2017 The energy landscape predicts flight height and wind turbine collision hazard in three species of large soaring raptor. *J. Appl. Ecol.* **54**, 1895–1906. (doi:10.1111/1365-2664.12909)

15. Vansteelant WMG, Shamoun-Baranes J, McLaren J, van Diermen J, Bouten W. 2017 Soaring across continents: decision-making of a soaring migrant under changing atmospheric conditions along an entire flyway. *J. Avian Biol.* **48**, 887–896. (doi:10.1111/jav.01298)

16. Panuccio M, Agostini N, Premuda G. 2012 Ecological barriers promote risk minimisation and social learning in migrating short-toed snake eagles. *Ethol. Ecol. Evol.* **24**, 74–80. (doi:10.1080/03949370.2011.583692)

17. Watanabe YY. 2016 Flight mode affects allometry of migration range in birds. *Ecol. Lett.* **19**, 907–914. (doi:10.1111/ele.12627)

18. Santos CD, Hanssen F, Muñoz A-R, Onrubia A, Wikelski M, May R, Silva JP. 2017 Match between soaring modes of black kites and the fine-scale distribution of updrafts. *Sci. Rep.* **7**, 6421. (doi:10.1038/s41598-017-05319-8)

19. Bohrer G, Brandes D, Mandel JT, Bildstein KL, Miller TA, Lanzone M, Katzner T, Maisonneuve C, Tremblay JA. 2012 Estimating updraft velocity components over large spatial scales: contrasting migration strategies of golden eagles and turkey vultures. *Ecol. Lett.* **15**, 96–103. (doi:10.1111/j.1461-0248.2011.01713.x)

20. Shamoun-Baranes J, Bouten W, van Loon EE, Meijer C, Camphuysen CJ. 2016 Flap or soar? How a flight generalist responds to its aerial environment. *Phil. Trans. R. Soc. B* **371**, 20150395. (doi:10.1098/rstb.2015.0395)

21. Katzner TE, Turk PJ, Duerr AE, Miller TA, Lanzone MJ, Cooper JL, Brandes D, Tremblay JA, Lemaître J. 2015 Use of multiple modes of flight subsidy by a soaring terrestrial bird, the golden eagle *Aquila chrysaetos*, when on migration. *J. R. Soc. Interface* **12**, 2–13. (doi:10.1098/rsif.2015.0530)

22. Shamoun-Baranes J, Leshem Y, Yom-Tov Y, Liechti O. 2003 Differential use of thermal convection by soaring birds over central Israel. *Condor* **105**, 208–218. (doi:10.1650/0010-5422(2003)105[0208:DU0TCB]2.0:0;2)

23. Mandel JT, Bildstein KL, Bohrer G, Winkler DW. 2008 Movement ecology of migration in turkey vultures. *Proc. Natl Acad. Sci.* **105**, 19102–19107. (doi:10.1073/pnas.0801789105)

24. Sapir N, Horvitz N, Wikelski M, Avissar R, Mehrer Y, Nathan R. 2011 Migration by soaring or flapping: numerical atmospheric simulations reveal that turbulence kinetic energy dictates bee-eater flight mode. *Proc. R. Soc. B* **278**, 3380–3386. (doi:10.1098/rspb.2011.0358)

25. Murgatroyd M, Photopoulou T, Underhill LG, Bouten W, Amar A. 2018 Where eagles soar: fine-resolution tracking reveals the spatiotemporal use of differential soaring modes in a large raptor. *Ecol. Evol.* **8**, 6788–6799. (doi:10.1002/ece3.4189)

26. Dennhardt AJ, Duerr AE, Brandes D, Katzner TE. 2015 Modeling autumn migration of a rare soaring raptor identifies new movement corridors in central Appalachia. *Ecol. Modell.* **303**, 19–29. (doi:10.1016/j.ecolmodel.2015.02.010)

27. Hernández-Pliego J, Rodríguez C, Bustamante J. 2015 Why do kestrels soar? *PLoS ONE* **10**, e0145402. (doi:10.1371/journal.pone.0145402)

28. Dodge S et al. 2014 Environmental drivers of variability in the movement ecology of turkey vultures (*Cathartes aura*) in North and South America. *Phil. Trans. R. Soc. B* **369**, 20130195. (doi:10.1098/rstb.2013.0195)

29. Oloo F, Safi K, Aryal J. 2018 Predicting migratory corridors of white storks *Ciconia ciconia* to enhance sustainable wind energy planning: a data-driven agent-based model. *Sustainability* **10**, 1470. (doi:10.3390/su10051470)

30. Hoover SL, Morrison ML. 2005 Behavior of red-tailed hawks in a wind turbine development. *J. Wildl. Manage.* **69**, 150–159. (doi:10.2193/0022-541X(2005)069<0150:BORHIA>2.0:0;2)

31. Ferrer M, De Lucas M, Janss GFE, Casado E, Muñoz AR, Bechard MJ, Calabuq CP. 2012 Weak migration. *J. R. Soc. Interface* **12**, 2–13. (doi:10.1098/rsif.2015.0530)
relationship between risk assessment studies and recorded mortality in wind farms. *J. Appl. Ecol.* **49**, 38–46. (doi:10.1111/j.1365-2664.2011.02054.x)

32. Poessel SA, Brandt J, Miller TA, Katzner TE. 2018 Meteorological and environmental variables affect flight behaviour and decision-making of an obligate soaring bird, the California condor *Gymnogyps californianus*. *Ibis* **160**, 36–53. (doi:10.1111/ibi.12531)

33. Katzner TE, Brandes D, Miller T, Lanzone M, Tremblay JA, Mulvihill R. 2012 Topography drives migratory flight altitude of golden eagles: implications for on-shore wind energy development. *J. Appl. Ecol.* **49**, 1178–1186. (doi:10.1111/j.1365-2664.2012.02185.x)

34. Weinzierl R, Bohrer G, Kranstauber B, Fiedler W, Wikelski M, Flack A. 2016 Wind estimation based on thermal soaring of birds. *Ecol. Evol.* **6**, 8706–8718. (doi:10.1002/ece3.2585)

35. Flack A, Nagy M, Fiedler W, Couzin ID, Wikelski M. 2018 From local collective behavior to global migratory patterns in white storks. *Science* **360**, 911–914. (doi:10.1126/science.aap7781)

36. Flack A, Fiedler W, Wikelski M. 2017 Data from: Wind estimation based on thermal soaring of birds. Movebank Data Repository. (doi:10.5441/001/1.bj96m274)

37. Garriga J, Palmer JRB, Oltra A, Bartumeus F. 2016 Expectation-maximization binary clustering for behavioural annotation. *PLoS ONE* **11**, 1–26. (doi:10.1371/journal.pone.0151984)

38. Kranstauber B, Smolla M, Scharf AK. 2017 *move*: Visualizing and analyzing animal track data. *R* package version 3.0.1. See https://cran.r-project.org/package=move.

39. Nathan R, Spiegel O, Fortmann-Roe S, Harel R, Wikelski M, Getz WM. 2012 Using tri-axial acceleration data to identify behavioral modes of free-ranging animals: general concepts and tools illustrated for griffon vultures. *J. Exp. Biol.* **215**, 986–996. (doi:10.1242/jeb.058602)

40. Wilson RP, White CR, Quintana F, Halsey LG, Liebsch N, Martin GR, Butler PJ. 2006 Moving towards acceleration for estimates of activity-specific metabolic rate in free-living animals: the case of the cormorant. *J. Anim. Ecol.* **75**, 1081–1090. (doi:10.1111/j.1365-2656.2006.01127.x)

41. Dodge S *et al*. 2013 The environmental-data automated track annotation (Env-DATA) system: linking animal tracks with environmental data. *Mov. Ecol.* **1**, 3. (doi:10.1186/2051-3933-1-3)

42. Liaw A, Wiener M. 2002 Classification and regression by randomforest. *R News* **2**, 18–22.

43. Franklin J. 2009 *Mapping species distributions: spatial inference and prediction*. Cambridge, UK: Cambridge University Press.

44. Allouche O, Tsoar A, Kadmon R. 2006 Assessing the accuracy of species distribution models: prevalence, kappa and the true skill statistic (TSS). *J. Appl. Ecol.* **43**, 1223–1232. (doi:10.1111/j.1365-2664.2006.01214.x)

45. Woods SN. 2003 Thin-plate regression splines. *J. R. Stat. Soc.* **65**, 95–114.

46. Bates D, Maechler M, Bolker B, Walker S. 2014 *lme4: Linear mixed-effects models using Eigen and S4*. R package version 1.0-6. See https://cran.r-project.org/package=lme4.

47. Woods SN. 2017 *Generalized additive models: An introduction with R*, 2nd edn. New York, NY: Chapman and Hall/CRC.

48. Harel R, Horvitz N, Nathan R. 2016 Adult vultures outperform juveniles in challenging thermal soaring conditions. *Sci. Rep.* **6**, 1–8. (doi:10.1038/srep27865)

49. Rotics S *et al*. 2016 The challenges of the first migration: movement and behaviour of juvenile vs. adult white storks with insights regarding juvenile mortality. *J. Anim. Ecol.* **85**, 938–947. (doi:10.1111/1365-2656.12525)

50. Shepard ELC, Williamson C, Windsor S. 2016 Fine-scale flight strategies of gulls in urban airflows indicate risk and reward in city living. *Phil. Trans. R. Soc. B* **371**, 20150394. (doi:10.1098/rstb.2015.0394)

51. Mandel JT, Bildstein KL. 2007 Turkey vultures use anthropogenic thermals to extend their daily activity period. *Wilson J. Ornithol.* **119**, 102–105. (doi:10.1676/05-154.1)

52. Marques AT, Batalha H, Rodrigues S, Costa H, Pereira MJR, Fonseca C, Mascarenhas M, Bernardino J. 2014 Understanding bird collisions at wind farms: an updated review on the causes and possible mitigation strategies. *Biol. Conserv.* **179**, 40–52. (doi:10.1016/j.biocon.2014.08.017)

53. Barrios L, Rodriguez A. 2004 Behavioural and environmental correlates of soaring-bird mortality at on-shore wind turbines. *J. Appl. Ecol.* **41**, 72–81. (doi:10.1111/j.1365-2664.2004.00876.x)

54. De Lucas M, Janss GFE, Whitfield DP, Ferrer M. 2008 Collision fatality of raptors in wind farms does not depend on raptor abundance. *J. Appl. Ecol.* **45**, 1695–1703. (doi:10.1111/j.1365-2664.2008.01549.x)
# Static landscape features predict uplift locations for soaring birds across Europe

Martina Scacco, Andrea Flack, Olivier Duriez, Martin Wikelski, Kamran Safi

## SUPPLEMENTARY MATERIAL

### List of content

**S1 - Segmentation of flight behaviour.** Containing supporting figures (Fig. S1.1, Fig. S1.2 and Fig. S1.3).

**S2 - Environmental variables.** Containing more details about the procedure we used to handle the static environmental layers; Tab. S2.1 table with environmental variables sources; Tab. S2.2 table containing the old legend and the reclassified legend of the Corine Land Cover classes.

**S3 - Uplift suitability model.** Containing additional figures (Fig. S3.1, S3.2, S3.3) and tables (Tab. S3.1) supporting the results. It also contains the detailed procedure and results of the random forest tuning and prevalence test.

**S4 - Uplift intensity model.** Containing the output of the three GAMs (Tab. S4.1) and supporting figure (Fig. S4.1).

**S5 - Static energy landscape.** Containing the output table of the linear mixed model (Tab. S5.1).
S1 - Segmentation of the flight behaviour

Supporting figures for segmentation

<figure><img src="image_11.png" /></figure>

**Figure S1.1** - Characterization of the different flight behavioural classes according to different flight parameters: cross-country speed (A), ground speed (B), vertical speed (C) and absolute turning angle (D).
<figure><img src="image_12.png" /></figure>

**Figure S1.2** - Example of segmentation of the flight behaviour detected from the ACC data of one individual of white stork. Figure (A) shows the different empirical values of the average ODBA among the three activity levels. Figure (B) shows how the three activity categories are different in terms of dynamic body acceleration measured on the three ACC axes (plot's x axis) and on the Z axis (plot's y axis).
A

<figure><img src="image_13.png" /></figure>

B

<figure><img src="image_14.png" /></figure>

**Figure S1.3** - Time allocation among flight behaviours. The barplots show the proportion of time white storks spent performing different flight behaviours (proportions obtained by cumulating the time spent by all individuals on the different behaviours). (A) Proportion of time (classified using GPS locations) spent soaring (circular or linear soaring, red bar) or using other flight types (grey bar) calculated relative to the duration of all classified GPS segments (748 h). (B) Proportion of time (classified using ACC recordings) spent flapping (blue bar) or using other flight types (grey bar) relative to the duration of all classified ACC bursts (24.3 h). The information shown in figures (A) and (B) cannot be directly compared due to the different sampling schedules of GPS and ACC.
# S2 - Environmental variables

## Static environmental variables

For descriptive purposes we group the environmental layers in two categories - surface features and land cover. All environmental layers are publicly available (Tab. S2.1).

**Surface features.** In order to characterize the surface features we used the publicly available elevation map EU-DEM (based on SRTM and ASTER Global Digital Elevation Model) [1] and we computed slope, aspect and roughness (topographic heterogeneity) using the R package raster [2]. The native spatial granularity of the elevation map is 1 arcsec (about 25 m near the Equator) but we aggregated the raster cells to match the 100 m resolution of the land cover map. Slope and aspect were computed according to Horn [3]. The roughness was calculated as the difference between the maximum and the minimum value of a cell and its 8 surrounding cells. Unevenness in the aspect, the slope and the elevation (this last one also called Topographic Position Index) were computed as the difference between the value of a cell and the mean value of its 8 surrounding cells. Highly correlated layers were excluded from the model to avoid multicollinearity (this was the case of the slope because of the high correlation with the roughness).

**Land cover.** We characterized the land cover using the Normalized Difference Vegetation Index (NDVI), CORINE Land Cover categories (CLC) and the Global Urban Footprint (GUF). The NDVI product is available as Spectral Indices product of Landsat 7 [4] with a spatial granularity of 30 m; the raster cells were resized to 100 m. We computed a summer (from June 1st to September 30th) NDVI composite for 2014 to match the temporal resolution of our tracking data. For the composite we extracted the higher monthly NDVI value of each cell and we averaged the resulting maximum monthly values. Extracting the maximum monthly value instead of the average value allowed us to avoid low values of NDVI that could be associated with errors in the
cloud cover mask of the Landsat NDVI product. The CORINE Land Cover 2012 is available from the European Environmental Agency [5] with a native spatial resolution of 100 m. We used the level 3 categories with few modifications (Table S2.2). The Global Urban Footprint is a binary thematic map with values of 1 for built-up areas (man-made building structures) and 0 for non-built-up areas. The dataset is available with a native spatial resolution of 0.4 arc seconds (about 12 m near the Equator) [6]. We resized the raster cells to 100 m computing the mean for each 9 by 9 cells (proportion of built-up areas/100 m cell).

<table><caption>**Table S2.1** - Environmental data sources.</caption><thead><tr><td>Environmental layer</td><td>source</td></tr></thead><tbody><tr><td>Digital elevation model (DEM)</td><td>EU-DEM (based on SRTM and ASTER Global Digital Elevation Model) from European Environmental Agency [1].</td></tr><tr><td>Roughness</td><td>Derived from DEM</td></tr><tr><td>Topographic Position Index</td><td>Derived from DEM</td></tr><tr><td>Slope</td><td>Derived from DEM</td></tr><tr><td>Slope unevenness</td><td>Derived from Slope</td></tr><tr><td>Aspect</td><td>Derived from DEM</td></tr><tr><td>Aspect unevenness</td><td>Derived from Aspect</td></tr><tr><td>Normalized Difference Vegetation Index (NDVI)</td><td>Landsat 7 Spectral Indices [4]. Available from U.S. Geological Service Bulk Download [9]. Data period June-September 2014.</td></tr><tr><td>CORINE Land Cover (CLC)</td><td>CLC 2012 from European Environmental Agency [5].</td></tr><tr><td>Global urban footprint (GUF)</td><td>Produced by Deutschen Zentrums für Luft und Raumfahrt, 2011 [6].</td></tr><tr><td>Thermal uplift potential</td><td>Movebank.org - Env-Data annotation service [7] based on ECMWF weather reanalysis, computed following Bohrer et al. [8].</td></tr><tr><td>Orographic uplift potential</td><td>Movebank.org - Env-Data annotation service [7] based on ECMWF weather reanalysis and ASTER DEM, computed following Bohrer et al. [8].</td></tr></tbody></table>
<table><caption>**Table S2.2:** Original and reclassified legend of the CORINE Land Cover categories.</caption><thead><tr><td>Reclassified categories</td><td>Reclassified code</td><td>CLC 2012 Legend</td></tr></thead><tbody><tr><td>Dumps</td><td>1</td><td>8</td></tr><tr><td>Artificial vegetated areas</td><td>2</td><td>10-11</td></tr><tr><td>Arable lands</td><td>3</td><td>12-14</td></tr><tr><td>Permanent crops</td><td>4</td><td>15-17</td></tr><tr><td>Pastures</td><td>5</td><td>18</td></tr><tr><td>Heterogeneous agricultural areas</td><td>6</td><td>19-22</td></tr><tr><td>Forest</td><td>7</td><td>23-25</td></tr><tr><td>Shrubs</td><td>8</td><td>26-29</td></tr><tr><td>Glaciers, snow</td><td>9</td><td>34</td></tr><tr><td>Inland wetlands, marshes</td><td>10</td><td>35-36</td></tr><tr><td>Marine wetlands, salines, inter-tidal flats</td><td>11</td><td>37-39</td></tr><tr><td>Water courses, rivers</td><td>12</td><td>40</td></tr><tr><td>Water bodies, lakes</td><td>13</td><td>41</td></tr><tr><td>Coastal lagoons, estuaries</td><td>14</td><td>42-43</td></tr><tr><td>Sea, ocean</td><td>15</td><td>44</td></tr><tr><td>Urban areas</td><td>16</td><td>1-6,9</td></tr><tr><td>Bare soil</td><td>17</td><td>7,30-33</td></tr><tr><td>No data or unclassified</td><td>NA</td><td>48-50,255</td></tr></tbody></table>
# S3 - Uplift suitability model

## Model tuning and evaluation

Different parameter values can be modified by the user in order to tune the random forest algorithm and improve its performances (such as number of variables chosen at each split, forest size and tree depth). In addition, in the field of species distribution modelling, the proportion of data belonging to different classes in the dataset used to train the model has an important effect on different accuracy measures on the model output of different machine learning algorithms [10,11]. For these reasons, we tuned the static and the dynamic landscape models, choosing the optimal values of mtry (number of variables chosen at each split) with respect to Out-of-Bag error estimate, using the package's inbuilt function *tuneRF*. Additionally, we tried different values of ratio soaring to flapping locations (named hereafter prevalence, for similarity with the definition used in species distribution modelling).

The role of prevalence on the model performances is controversial. Different studies suggested that certain accuracy measures (such as kappa) are sensitive to prevalence and that prevalence should be taken into account when evaluating the model [12,13], but few unclear suggestions have been made about whether the number of presences and absences should be manipulated in order to maximize the model performances [10,11,13]. This lack of clarity is probably due to multiple reasons, such as the algorithm used, the difference between using real absences or pseudo-absences, the difficulty to differentiate between biases in the measures used to evaluate the performance and biases in the performances itself, and finally to the uncertain effect of data manipulation on the ecological interpretation of the result (for instance when dealing with rare and specialist species versus abundant and generalist ones, and the need to compare the predicted prevalence with the observed prevalence). For this reasons we decided not to manipulate the prevalence in our analysis but to perform a
prevalence manipulation test. The prevalence manipulation test was performed on 100
trees. Before running the model with different prevalence values, the complete
dataset was randomly partitioned in test set (20% of the data) and train set (80% of
the data), separately for soaring (presences) and flapping (absences) locations in
order to maintain the same ratio of presences to absences in both the original dataset
and the test set; from the train set we then manipulated the number of presences to
meet the values of prevalence (n. presences/n. absences) we wanted to test. The
same data partitioning procedure was repeated 10 times. The data partition used to
test the model was separated from the training set before manipulating the
prevalence in the training set; in this way we could ensure a constant size for the test
set during the evaluation of each prevalence value. During the models' evaluation we
considered the following accuracy measures: AUC or area under the ROC curve,
sensitivity and specificity; sensitivity and specificity are threshold dependent
measures and the chosen values correspond to a threshold that maximize the sum of
specificity and sensitivity (True Skill Statistics or TSS).

The results of the tuning procedure showed an increase in the model accuracy (AUC
and TSS) with increasing prevalence. Sensitivity and specificity showed variable values
at different values of prevalence, but they both showed a slight positive trend and the
difference between the two slightly decreased, with increasing prevalence (Fig. S3.1).
In contrast, the threshold that maximized the TSS increased, being around 0.5 with a
prevalence equal to 2, and around 0.9 with a prevalence of 12, which is in agreement
with the fact that the best performances in our model were associated with high
thresholds. This didn't affect our model accuracy, but its robustness. In fact, all models
including the original presence/absence ratio (12.66), the threshold at which the
number of flapping locations correctly classified matched the number of soaring
locations correctly classified was really close to 1, between 0.9 and 0.95; a threshold
of 0.9 favoured the identification of soaring (with higher commission rate) whereas a threshold of 0.95 favoured the identification of flapping (with higher omission rate). This means that this model prediction, although accurate, can be considered sensitive compared to a prediction obtained with a lower threshold, and a small change in the threshold leads to really different predicted results. For this reason, if the sample size of the least represented class allows it, we recommend the use of balanced classes in species distribution and habitat suitability modelling, by randomly subsampling the more represented class to about twice the size of the least represented class, to be able to choose a more centered threshold and increase model robustness.

<figure><img src="image_15.png" /><figcaption>**Figure S3.1** - Prevalence test. Effect of different prevalence values (ratio soaring to flapping locations) on different accuracy measures (A), and on the threshold that has to be chosen in order to maximize True Skill Statistics, calculated as (Sensitivity + Specificity) - 1 (B).</figcaption></figure>
**Table S3.1** - Random forest evaluation of the three uplift suitability models, based on the test set and averaged across the ten cross-validations.

<table><caption>Table S3.1 A - Static uplift suitability model, based only on static features (AUC±sd = 0.851±0.022).</caption><thead><tr><th>Threshold</th><th>Sensitivity</th><th>Specificity</th><th>Commission Error</th><th>Omission Error</th><th>TSS</th></tr></thead><tbody><tr><td>0</td><td>1 (0)</td><td>0.001 (0.003)</td><td>0.999 (0.003)</td><td>0 (0)</td><td>0.001 (0.003)</td></tr><tr><td>0.05</td><td>1 (0)</td><td>0.061 (0.022)</td><td>0.939 (0.022)</td><td>0 (0)</td><td>0.061 (0.022)</td></tr><tr><td>0.1</td><td>1 (0)</td><td>0.097 (0.033)</td><td>0.903 (0.033)</td><td>0 (0)</td><td>0.097 (0.033)</td></tr><tr><td>0.15</td><td>1 (0)</td><td>0.122 (0.025)</td><td>0.878 (0.025)</td><td>0 (0)</td><td>0.122 (0.025)</td></tr><tr><td>0.2</td><td>1 (0)</td><td>0.162 (0.044)</td><td>0.838 (0.044)</td><td>0 (0)</td><td>0.162 (0.044)</td></tr><tr><td>0.25</td><td>1 (0)</td><td>0.185 (0.05)</td><td>0.815 (0.05)</td><td>0 (0)</td><td>0.185 (0.05)</td></tr><tr><td>0.3</td><td>1 (0)</td><td>0.187 (0.049)</td><td>0.813 (0.049)</td><td>0 (0)</td><td>0.187 (0.049)</td></tr><tr><td>0.35</td><td>1 (0)</td><td>0.211 (0.047)</td><td>0.789 (0.047)</td><td>0 (0)</td><td>0.21 (0.047)</td></tr><tr><td>0.4</td><td>0.999 (0.001)</td><td>0.257 (0.046)</td><td>0.743 (0.046)</td><td>0.001 (0.001)</td><td>0.255 (0.046)</td></tr><tr><td>0.45</td><td>0.997 (0.002)</td><td>0.341 (0.054)</td><td>0.659 (0.054)</td><td>0.003 (0.002)</td><td>0.338 (0.054)</td></tr><tr><td>0.5</td><td>0.997 (0.002)</td><td>0.356 (0.059)</td><td>0.644 (0.059)</td><td>0.003 (0.002)</td><td>0.352 (0.059)</td></tr><tr><td>0.55</td><td>0.996 (0.002)</td><td>0.36 (0.06)</td><td>0.64 (0.06)</td><td>0.004 (0.002)</td><td>0.356 (0.06)</td></tr><tr><td>0.6</td><td>0.994 (0.002)</td><td>0.367 (0.058)</td><td>0.633 (0.058)</td><td>0.006 (0.002)</td><td>0.361 (0.058)</td></tr><tr><td>0.65</td><td>0.992 (0.002)</td><td>0.376 (0.058)</td><td>0.624 (0.058)</td><td>0.008 (0.002)</td><td>0.368 (0.058)</td></tr><tr><td>0.7</td><td>0.988 (0.003)</td><td>0.393 (0.056)</td><td>0.607 (0.056)</td><td>0.012 (0.003)</td><td>0.38 (0.056)</td></tr><tr><td>0.75</td><td>0.977 (0.005)</td><td>0.414 (0.054)</td><td>0.586 (0.054)</td><td>0.023 (0.005)</td><td>0.391 (0.054)</td></tr><tr><td>0.8</td><td>0.956 (0.003)</td><td>0.468 (0.06)</td><td>0.532 (0.06)</td><td>0.044 (0.003)</td><td>0.424 (0.06)</td></tr><tr><td>0.85</td><td>0.914 (0.006)</td><td>0.55 (0.051)</td><td>0.45 (0.051)</td><td>0.086 (0.006)</td><td>0.464 (0.051)</td></tr><tr><td>0.9</td><td>0.827 (0.009)</td><td>0.685 (0.048)</td><td>0.315 (0.048)</td><td>0.173 (0.009)</td><td>0.512 (0.045)</td></tr><tr><td>0.95</td><td>0.625 (0.014)</td><td>0.862 (0.041)</td><td>0.138 (0.041)</td><td>0.375 (0.014)</td><td>0.487 (0.036)</td></tr><tr><td>1</td><td>0.113 (0.01)</td><td>0.991 (0.009)</td><td>0.009 (0.009)</td><td>0.887 (0.01)</td><td>0.103 (0.011)</td></tr></tbody></table>
<table><caption>Table S3.1 B - Dynamic uplift suitability model, based only on thermal and orographic uplift potentials (AUC±sd = 0.695±0.024).</caption><thead><tr><th>Threshold</th><th>Sensitivity</th><th>Specificity</th><th>Commission Error</th><th>Omission Error</th><th>TSS</th></tr></thead><tbody><tr><td>0</td><td>1 (0)</td><td>0 (0)</td><td>1 (0)</td><td>0 (0)</td><td>0 (0)</td></tr><tr><td>0.05</td><td>1 (0)</td><td>0 (0)</td><td>1 (0)</td><td>0 (0)</td><td>0 (0)</td></tr><tr><td>0.1</td><td>1 (0)</td><td>0 (0)</td><td>1 (0)</td><td>0 (0)</td><td>0 (0)</td></tr><tr><td>0.15</td><td>1 (0)</td><td>0.002 (0.004)</td><td>0.998 (0.004)</td><td>0 (0)</td><td>0.002 (0.004)</td></tr><tr><td>0.2</td><td>1 (0)</td><td>0.004 (0.004)</td><td>0.996 (0.004)</td><td>0 (0)</td><td>0.004 (0.004)</td></tr><tr><td>0.25</td><td>1 (0)</td><td>0.007 (0.005)</td><td>0.993 (0.005)</td><td>0 (0)</td><td>0.007 (0.005)</td></tr><tr><td>0.3</td><td>1 (0)</td><td>0.018 (0.012)</td><td>0.982 (0.012)</td><td>0 (0)</td><td>0.018 (0.012)</td></tr><tr><td>0.35</td><td>1 (0)</td><td>0.025 (0.016)</td><td>0.975 (0.016)</td><td>0 (0)</td><td>0.025 (0.016)</td></tr><tr><td>0.4</td><td>0.999 (0)</td><td>0.03 (0.016)</td><td>0.97 (0.016)</td><td>0.001 (0)</td><td>0.03 (0.016)</td></tr><tr><td>0.45</td><td>0.999 (0.001)</td><td>0.039 (0.013)</td><td>0.961 (0.013)</td><td>0.001 (0.001)</td><td>0.039 (0.013)</td></tr><tr><td>0.5</td><td>0.998 (0)</td><td>0.053 (0.008)</td><td>0.947 (0.008)</td><td>0.002 (0)</td><td>0.051 (0.008)</td></tr><tr><td>0.55</td><td>0.996 (0.001)</td><td>0.069 (0.015)</td><td>0.931 (0.015)</td><td>0.004 (0.001)</td><td>0.065 (0.014)</td></tr><tr><td>0.6</td><td>0.994 (0.001)</td><td>0.096 (0.016)</td><td>0.904 (0.016)</td><td>0.006 (0.001)</td><td>0.09 (0.017)</td></tr><tr><td>0.65</td><td>0.99 (0.002)</td><td>0.127 (0.017)</td><td>0.873 (0.017)</td><td>0.01 (0.002)</td><td>0.117 (0.018)</td></tr><tr><td>0.7</td><td>0.984 (0.003)</td><td>0.157 (0.026)</td><td>0.843 (0.026)</td><td>0.016 (0.003)</td><td>0.141 (0.026)</td></tr><tr><td>0.75</td><td>0.973 (0.005)</td><td>0.203 (0.033)</td><td>0.797 (0.033)</td><td>0.027 (0.005)</td><td>0.176 (0.033)</td></tr><tr><td>0.8</td><td>0.954 (0.006)</td><td>0.276 (0.039)</td><td>0.724 (0.039)</td><td>0.046 (0.006)</td><td>0.23 (0.04)</td></tr><tr><td>0.85</td><td>0.919 (0.009)</td><td>0.358 (0.047)</td><td>0.642 (0.047)</td><td>0.081 (0.009)</td><td>0.277 (0.049)</td></tr><tr><td>0.9</td><td>0.844 (0.008)</td><td>0.487 (0.045)</td><td>0.513 (0.045)</td><td>0.156 (0.008)</td><td>0.331 (0.045)</td></tr><tr><td>0.95</td><td>0.459 (0.018)</td><td>0.737 (0.033)</td><td>0.263 (0.033)</td><td>0.541 (0.018)</td><td>0.196 (0.039)</td></tr><tr><td>1</td><td>0.026 (0.005)</td><td>0.985 (0.014)</td><td>0.015 (0.014)</td><td>0.974 (0.005)</td><td>0.011 (0.015)</td></tr></tbody></table>
<table><caption>Table S3.1 C - Combined uplift suitability model, based on both static and dynamic predictors (AUC±sd = 0.862±0.016).</caption><thead><tr><th>Threshold</th><th>Sensitivity</th><th>Specificity</th><th>Commission Error</th><th>Omission Error</th><th>TSS</th></tr></thead><tbody><tr><td>0</td><td>1 (0)</td><td>0.001 (0.003)</td><td>0.999 (0.003)</td><td>0 (0)</td><td>0.001 (0.003)</td></tr><tr><td>0.05</td><td>1 (0)</td><td>0.051 (0.019)</td><td>0.949 (0.019)</td><td>0 (0)</td><td>0.051 (0.019)</td></tr><tr><td>0.1</td><td>1 (0)</td><td>0.075 (0.032)</td><td>0.925 (0.032)</td><td>0 (0)</td><td>0.075 (0.032)</td></tr><tr><td>0.15</td><td>1 (0)</td><td>0.107 (0.027)</td><td>0.893 (0.027)</td><td>0 (0)</td><td>0.106 (0.027)</td></tr><tr><td>0.2</td><td>1 (0)</td><td>0.133 (0.024)</td><td>0.867 (0.024)</td><td>0 (0)</td><td>0.133 (0.024)</td></tr><tr><td>0.25</td><td>1 (0)</td><td>0.157 (0.02)</td><td>0.843 (0.02)</td><td>0 (0)</td><td>0.157 (0.02)</td></tr><tr><td>0.3</td><td>1 (0)</td><td>0.169 (0.024)</td><td>0.831 (0.024)</td><td>0 (0)</td><td>0.169 (0.024)</td></tr><tr><td>0.35</td><td>1 (0)</td><td>0.192 (0.023)</td><td>0.808 (0.023)</td><td>0 (0)</td><td>0.191 (0.023)</td></tr><tr><td>0.4</td><td>0.999 (0)</td><td>0.234 (0.03)</td><td>0.766 (0.03)</td><td>0.001 (0)</td><td>0.234 (0.03)</td></tr><tr><td>0.45</td><td>0.999 (0.001)</td><td>0.287 (0.041)</td><td>0.713 (0.041)</td><td>0.001 (0.001)</td><td>0.285 (0.04)</td></tr><tr><td>0.5</td><td>0.998 (0.001)</td><td>0.319 (0.032)</td><td>0.681 (0.032)</td><td>0.002 (0.001)</td><td>0.317 (0.032)</td></tr><tr><td>0.55</td><td>0.998 (0.001)</td><td>0.341 (0.037)</td><td>0.659 (0.037)</td><td>0.002 (0.001)</td><td>0.339 (0.036)</td></tr><tr><td>0.6</td><td>0.997 (0.001)</td><td>0.356 (0.034)</td><td>0.644 (0.034)</td><td>0.003 (0.001)</td><td>0.352 (0.033)</td></tr><tr><td>0.65</td><td>0.995 (0.001)</td><td>0.375 (0.037)</td><td>0.625 (0.037)</td><td>0.005 (0.001)</td><td>0.37 (0.036)</td></tr><tr><td>0.7</td><td>0.992 (0.002)</td><td>0.392 (0.038)</td><td>0.608 (0.038)</td><td>0.008 (0.002)</td><td>0.384 (0.038)</td></tr><tr><td>0.75</td><td>0.983 (0.003)</td><td>0.417 (0.039)</td><td>0.583 (0.039)</td><td>0.017 (0.003)</td><td>0.4 (0.038)</td></tr><tr><td>0.8</td><td>0.962 (0.005)</td><td>0.471 (0.036)</td><td>0.529 (0.036)</td><td>0.038 (0.005)</td><td>0.433 (0.033)</td></tr><tr><td>0.85</td><td>0.92 (0.008)</td><td>0.557 (0.038)</td><td>0.443 (0.038)</td><td>0.08 (0.008)</td><td>0.477 (0.035)</td></tr><tr><td>0.9</td><td>0.829 (0.014)</td><td>0.706 (0.036)</td><td>0.294 (0.036)</td><td>0.171 (0.014)</td><td>0.535 (0.036)</td></tr><tr><td>0.95</td><td>0.609 (0.014)</td><td>0.878 (0.032)</td><td>0.122 (0.032)</td><td>0.391 (0.014)</td><td>0.488 (0.022)</td></tr><tr><td>1</td><td>0.055 (0.004)</td><td>1 (0)</td><td>0 (0)</td><td>0.945 (0.004)</td><td>0.055 (0.004)</td></tr></tbody></table>
Supporting figures for the models' output

<figure><img src="image_16.png" /></figure>

**Figure S3.2** - Accuracy of the three uplift suitability models, static (A), dynamic (B) and combined (C), in terms of sensitivity (proportion of soaring locations correctly classified, in green) and specificity (proportion of flapping locations correctly classified, in red) at different thresholds values. The solid points represent the value of Sensitivity and Specificity, averaged across the ten runs of each model, at a threshold that maximize the value of the True Skill Statistics.
A - Static model

Variables

<figure><img src="image_17.png" /></figure>

Variables

<figure><img src="image_18.png" /></figure>
**B - Dynamic model**

Variables

Orographic uplift

Thermal uplift

<figure><img src="image_19.png" /></figure>

Mean decrease in Accuracy
(%IncMSE)

Variables

Orographic uplift

Thermal uplift

<figure><img src="image_20.png" /></figure>

Decrease in Node Impurity
(IncNodePurity)
C - Combined model

<figure><img src="image_21.png" /><figcaption>**Figure S3.3** - Variable importance. Contributions of the different variables to the static (A), dynamic (B) and combined (C) uplift suitability models, measured in terms of mean decrease in accuracy (left) and decrease in node impurity (right).</figcaption></figure>
S4 - Uplift intensity model

<table><caption>**Table S4.1** - Summary of the three uplift intensity models (generalized additive models).</caption><thead><tr><th></th><th colspan="3">Response: sqrt(Vertical Speed)</th></tr><tr><th></th><th>Static<br>Uplift intensity model</th><th>Dynamic<br>Uplift intensity model</th><th>Combined<br>Uplift intensity model</th></tr></thead><tbody><tr><td colspan="4"><strong>Parametric coefficients:</strong></td></tr><tr><td colspan="4">Estimate (C.I.)</td></tr><tr><td>Intercept</td><td>0.898 (0.866, 0.929)</td><td>0.765 (0.759, 0.771)</td><td>0.719 (0.688, 0.751)</td></tr><tr><td>Aspect unevenness</td><td>0.0001 (0.00000, 0.0001)</td><td></td><td>0.0001 (0.00000, 0.0001)</td></tr><tr><td>DEM unevenness</td><td>0.001 (-0.0002, 0.002)</td><td></td><td>0.001 (-0.0004, 0.002)</td></tr><tr><td>Slope unevenness</td><td>0.005 (0.002, 0.008)</td><td></td><td>0.005 (0.002, 0.008)</td></tr><tr><td>GUF</td><td>0.013 (-0.001, 0.028)</td><td></td><td>0.022 (0.008, 0.037)</td></tr><tr><td>CLC 1 - dumps</td><td>-0.155 (-0.310, -0.0004)</td><td></td><td>-0.039 (-0.190, 0.112)</td></tr><tr><td>CLC 2 - art. veg. areas</td><td>0.063 (0.023, 0.103)</td><td></td><td>0.055 (0.015, 0.094)</td></tr><tr><td>CLC 3 - arable lands</td><td>-0.018 (-0.049, 0.014)</td><td></td><td>-0.019 (-0.050, 0.012)</td></tr><tr><td>CLC 4 - perm. crops</td><td>-0.023 (-0.055, 0.009)</td><td></td><td>-0.016 (-0.047, 0.016)</td></tr><tr><td>CLC 5 - pastures</td><td>-0.034 (-0.068, 0.0001)</td><td></td><td>-0.028 (-0.062, 0.005)</td></tr><tr><td>CLC 6 - het. agr. areas</td><td>-0.017 (-0.049, 0.015)</td><td></td><td>-0.012 (-0.044, 0.019)</td></tr><tr><td>CLC 7 - forest</td><td>0.014 (-0.018, 0.046)</td><td></td><td>0.009 (-0.022, 0.041)</td></tr><tr><td>CLC 8 - Shrubs</td><td>0.024 (-0.008, 0.057)</td><td></td><td>0.022 (-0.011, 0.054)</td></tr><tr><td>CLC 10 - wet. marshes</td><td>-0.050 (-0.088, -0.012)</td><td></td><td>-0.026 (-0.064, 0.012)</td></tr><tr><td>CLC 11 - marine wet. salines</td><td>0.002 (-0.032, 0.036)</td><td></td><td>0.007 (-0.027, 0.040)</td></tr><tr><td>CLC 12 - rivers</td><td>-0.020 (-0.070, 0.029)</td><td></td><td>-0.026 (-0.075, 0.023)</td></tr><tr><td>CLC 13 - wat. bodies lakes</td><td>-0.102 (-0.152, -0.051)</td><td></td><td>-0.093 (-0.143, -0.044)</td></tr><tr><td>CLC 14 - coast. lagoons est.</td><td>0.060 (0.014, 0.105)</td><td></td><td>0.069 (0.024, 0.114)</td></tr><tr><td>CLC 15 - sea</td><td>0.006 (-0.450, 0.462)</td><td></td><td>-0.004 (-0.449, 0.441)</td></tr><tr><td>CLC 16 - urban areas</td><td>-0.013 (-0.046, 0.020)</td><td></td><td>-0.017 (-0.049, 0.016)</td></tr><tr><td>Therm. uplift pot.</td><td></td><td>0.124 (0.119, 0.129)</td><td>0.161 (0.156, 0.166)</td></tr><tr><td>Orog. uplift pot.</td><td></td><td>0.042 (0.029, 0.056)</td><td>0.063 (0.049, 0.077)</td></tr><tr><td colspan="4"><strong>Smooth terms:</strong></td></tr><tr><td colspan="4">F (Effective df)</td></tr><tr><td>s(Aspect)</td><td>4.057 (3.040)</td><td></td><td>2.080 (2.545)</td></tr><tr><td>s(DEM)</td><td>49.019 (6.978)</td><td></td><td>40.020 (6.901)</td></tr><tr><td>s(Ruggedness)</td><td>1.863 (4.443)</td><td></td><td>4.510 (4.875)</td></tr><tr><td>s(NDVI)</td><td>18.068 (6.710)</td><td></td><td>21.650 (4.277)</td></tr><tr><td>s(Latitude)</td><td>60.916 (8.881)</td><td></td><td>91.420 (8.885)</td></tr><tr><td>Observations</td><td>76383</td><td>78307</td><td>74946</td></tr><tr><td>Adjusted R2</td><td>0.032</td><td>0.032</td><td>0.083</td></tr><tr><td>AIC (df)</td><td>46575.00 (54.745)</td><td>49636.88 (4)</td><td>42003.95 (53.19181)</td></tr></tbody></table>
A

<figure><img src="image_22.png" /></figure>

B

<figure><img src="image_23.png" /></figure>

**Figure S4.1** - Non linear relationship between some of the environmental predictors included as smooth terms in the static uplift intensity model (generalized additive model) and the response variable (vertical speed). In (A) the effect of NDVI and aspect, whereas in (B) the effect of DEM and latitude on the vertical speed. In lighter colours, combination of variables' values that predict a higher vertical speed.
## S5 - Static energy landscape

<table><caption>**Table 5.1** - Output of the Linear Mixed Model analysing the mean daily ODBA as a function of the mean daily uplift suitability (as predicted by the static uplift suitability model).</caption><thead><tr><th></th><th>Response variable:<br>sqrt(Mean Daily ODBA)</th></tr></thead><tbody><tr><td><strong>Fixed effects:</strong></td><td></td></tr><tr><td>Estimate (C.I.)</td><td></td></tr><tr><td>Intercept</td><td>3.548 (0.166)</td></tr><tr><td>Mean daily uplift suitability</td><td>-2.252 (0.189)</td></tr><tr><td><strong>Random effects:</strong></td><td></td></tr><tr><td>Variance (St. Dev., Corr.)</td><td></td></tr><tr><td>Intercept</td><td>0.990 (0.995)</td></tr><tr><td>Mean daily uplift suitability</td><td>1.248 (1.117, -1.00)</td></tr><tr><td>Residual</td><td>0.055 (0.234)</td></tr><tr><td>Observations</td><td>823</td></tr><tr><td>Groups (Individuals)</td><td>59</td></tr><tr><td>AIC</td><td>46.233</td></tr><tr><td>AIC (null model)</td><td>114.254</td></tr><tr><td>ΔAIC = AIC - AICnull</td><td>-68.021</td></tr><tr><td>Marginal R<sup>2</sup></td><td>0.306</td></tr><tr><td>Conditional R<sup>2</sup></td><td>0.426</td></tr></tbody></table>
## REFERENCES - SUPPLEMENTARY MATERIAL

1. EEA. 2013 European Environment Agency - Digital Elevation Model over Europe (EU-DEM). *Accessed March 2017*. See https://www.eea.europa.eu/data-and-maps/data/eu-dem.

2. Hijmans RJ. 2016 raster: Geographic Data Analysis and Modeling. R package version 2.5-8. See https://cran.r-project.org/package=raster.

3. Horn BKP. 1981 Hill shading and the reflectance map. *Proc. IEEE 69*, 14-47.

4. Masek JG, Vermote EF, Saleous N, Wolfe R, Hall FG, Huemmrich F, Gao F, Kutler J, Lim T. 2006 A Landsat surface reflectance data set for North America, 1990-100. *IEEE Geosci. Remote Sens. Lett.* **3**, 68-72. (doi:10.1109/LGRS.2005.857030)

5. EEA. 2012 European Environment Agency - Corine Land Cover 2012 raster data. *Accessed March 2017*. See https://www.eea.europa.eu/data-and-maps/data/clc-2012-raster.

6. DLR. 2011 Earth Observation Center - Global Urban Footprint. *Accessed March 2017*. See http://www.dlr.de/eoc/en/desktopdefault.aspx/tabid-11725/20508_read-47944/.

7. Dodge S et al. 2013 The environmental-data automated track annotation (Env-DATA) system: linking animal tracks with environmental data. *Mov. Ecol.* **1**, 3. (doi:10.1186/2051-3933-1-3)

8. Bohrer G, Brandes D, Mandel JT, Bildstein KL, Miller TA, Lanzone M, Katzner T, Maisonneuve C, Tremblay JA. 2012 Estimating updraft velocity components over large spatial scales: Contrasting migration strategies of golden eagles and turkey vultures. *Ecol. Lett.* **15**, 96-103. (doi:10.1111/j.1461-0248.2011.01713.x)

9. U.S. Geological Survey. In press. USGS Earth Resources Observation and Science (EROS) Center Science Processing Architecture (ESPA) On Demand Interface. *Accessed March 2017*. See https://espa.cr.usgs.gov/.

10. Barbet-Massin M, Jiguet F, Albert CH, Thuiller W. 2012 Selecting pseudo-absences for species distribution models: how, where and how many? *Methods Ecol. Evol.* **3**, 327-338. (doi:10.1111/j.2041-210X.2011.00172.x)
11. Franklin J. 2009 *Mapping species distributions*. Spatial inference and prediction. Cambridge University Press.

12. Allouche O, Tsoar A, Kadmon R. 2006 Assessing the accuracy of species distribution models: Prevalence, kappa and the true skill statistic (TSS). *J. Appl. Ecol.* **43**, 1223-1232. (doi:10.1111/j.1365-2664.2006.01214.x)

13. Franklin J, Wejnert KE, Hathaway SA, Rochester CJ, Fisher RN. 2009 Effect of species rarity on the accuracy of species distribution models for reptiles and amphibians in southern California. *Divers. Distrib.* **15**, 167-177. (doi:10.1111/j.1472-4642.2008.00536.x)