Species distribution modelling of Bryde’s whales, humpback whales, southern right whales, and sperm whales in the southern African region to inform their conservation in expanding economies

Study area

The study area (limited by the co-ordinates 16°S, 0°E to 38°S, 80°E) encompasses the western Indian Ocean and the southeast Atlantic Ocean (Fig. 1), covering both cold and warm water habitats of these cetaceans (Findlay et al., 1992; Best, 2007). We chose this area to include all the available information on species sightings and to determine habitats suitable for whale distribution within the area.

Presence data

We collated data for all four species from opportunistic and dedicated scientific surveys covering the period from 1913 to 2016 (Table 1, Fig. 2, Fig. S1 and Table S1). This specific period was used as it included as many presence points as possible and because data falling outside these dates were not available. The majority of the data we used originated from the Whale Unit of the Mammal Research Institute at the University of Pretoria (44%) and the Ocean Biogeographic Information System Spatial Ecological Analysis of Megavertebrate Populations (OBIS-SEAMAP) (47%) (Halpin et al., 2009). Incidental sightings (7%) (Nan Rice, CEO: Dolphin Action and Protection Group, Fish Hoek) and seismic survey marine mammal observations in South Africa (2%) (Petroleum Agency of South Africa, PASA) made up the rest of the dataset.

Environmental data

We selected environmental variables based on known ecological relationships between cetaceans and environmental variables (Table 2) (Best, 2007; Torres et al., 2013; Breen et al., 2016). Topographic variables included bathymetry (metres), distance to shore (kilometres) and slope (degrees). Oceanographic variables included sea surface temperature (SST in °C) chlorophyll-a concentration (chl-a in mg m⁻³), salinity (in psu), eastwards sea water velocity (U_o in m s⁻¹) and westwards sea water velocity (V_o in m s⁻¹) (Fig. S2). We used SST, salinity, U_o and V_o data derived from satellites. Wernand, Van de rWoerd & Gieskes (2013) have shown that there has been an increase of chl-a concentration in the Atlantic Ocean and a decrease in the Indian Ocean with no overall significant chl-a trends worldwide from 1989 to 2000. SST has risen between 1.2 °C and 1.6 °C in the high latitudes and between 0.4 °C and 1 °C in the tropics and sub-tropics over the last 100 years (Deser, Phillips & Alexander, 2010). Despite these changes, studies by Purdon et al. (2020a), Purdon et al. (2020b) and Torres, Read & Halpin (2016) used averages of oceanographic variables, suggesting that they provide sufficient information for modelling baseline cetacean distribution in conjunction with presence data that dates back to 1913. In this study we averaged chl-a seasonally over the period 2002 to 2016, and salinity, U_o, V_o and SST seasonally over the period 1991 to 2016 in their original projection (WGS 84 (EPSG 4326)) and at their origional spatial resolution (∼8 km) in ArcMap (Table 2, Esri, Redlands, CA, USA).

To ensure all variables maintained the same spatial resolution and projection, we resampled them to a 0.083° (∼8 km) grid using a nearest neighbour interpolation and re-projected them to a standard WGS84, EPSG 4326 geographic projection in ArcMap. They were then clipped to the study area extent using ArcMap. We evaluated multi-collinearity between environmental variables by testing the variable inflation factor (VIF) using the ‘car’ package (Fox & Weisberg, 2011) in R (R Core Team, 2020). We retained all variables as the VIF was less than three, indicating no multi-collinearity (O’Brien, 2007).

Species distribution modelling

We used ensemble modelling implemented in the ‘SSDM’ package (Schmitt et al., 2017) within R to predict preferred species habitat using seven algorithms (Table 1). We used the algorithms as they consisted of different modelling approaches, including three regression methods, one classification method and three machine learning methods to build ensemble models for each species seasonally (Table 1). To determine seasonal variation, we followed Findlay (1989) in creating SDMs separately for summer (October to March) and winter (April to September) months.

The dataset contained both presence and absence data from scientific surveys, and marine mammal observers and passive acoustic monitors on oil and gas exploration surveys, and presence only records from opportunistic sightings. To compensate for the lack of effort data, we used pseudo-absences selected randomly in the study area (Fig. S1). The number of pseudo-absence points selected was based on the recommendations of Barbet-Massin et al. (2012) for the different types of algorithms used (Table 1). Barbet-Massin et al. (2012) suggest that, for classification and machine learning techniques, the same number of pseudo-absence points as presence points should be selected and averaged over several runs. For regression techniques, many pseudo-absence points (2% of the study area which is 25 218 188 km² for our study) should be chosen, whereas for techniques such as multi adaptive regression splines and discriminant analysis, fewer pseudo-absence points should be chosen and both averaged over several runs. We selected the number of pseudo-absence points 10 times on a random basis using the ‘SSDM’ package default settings based on the recommendations of Barbet-Massin et al. (2012) (Table 1 and Fig. S1).

To account for spatially biased data, which are prone to spatial autocorrelation, we used spatial thinning within the ‘SSDM’ package (Schmitt et al., 2017). We evaluated the algorithms and ensemble model performance using the holdout method in the SSDM package, where the initial dataset was split into separate training (70% of the data) and evaluation (30% of the data) sets (Schmitt et al., 2017). This process was repeated 10 times. We used a standard evaluation metric, area under the receiver operating characteristic curve (AUC). Scores ranged between 0.5 and 1, with those close to 1 indicating that the model has excellent performance, whereas values around 0.5 denote models that are no better than random.

Variable importance was measured by determining how much the correlation changed between predicted values before and after permuting each variable. ${\displaystyle I_v=1-Cor\left(P_f,P_v\right)}$ where I_v is index of importance of a variable, Cor is correlation coefficient, P_f is prediction from the full model, P_v is prediction after permuting or reshuffling the variable v.

We computed the partial effects of each environmental variable by predicting the response for a specific variable while holding the other variables at their mean. For example, if bathymetry was the variable that best explained the cetacean species distribution, scores of the partial effects closer to one would indicate at what depths the cetaceans were likely to be found.

To create the final ensemble model, we calculated the average of all seven algorithms if their AUC scores were over 0.7. We projected the final ensemble models for each species and season (summer and winter) on to a predicted distribution map. Probability of occurrence ranged between zero and one, with one indicating the highest probability of species presence. We created uncertainty maps using the variance of predicted probability of occurrence among the algorithms, where areas of agreement are characterised by low variance amongst algorithms (Schmitt et al., 2017).

Ensemble model and algorithm performance

Overall, we predicted that the ensemble models would produce higher AUC scores than the individual algorithms. Analysis of the results, however, indicates that, with the exception of artificial neural network, ANN, the performance of individual algorithms was equal to or better than the ensemble model (Fig. 3 and Table S2). The RF algorithm produced the highest AUC scores for all species in both seasons except for the Bryde’s whale in winter where the algorithm GBM produced the highest AUC score (Fig. 3 and Table S2).

Species distribution modelling

The most influential environmental variables driving distribution differed among species, seasons and algorithms (Fig. 4). Figure 5 shows the most influential environmental variables for the RF algorithm, and it indicates predicted preferred ranges for each species (see Figs. S3A–S3F for the most influential environmental variables for the other algorithms). Environmental variables influenced the predicted occurrence of the four cetacean species and can be seen in ensemble model predictions in Fig. 6. Seasonally, southern right whales, sperm whales and humpback whales showed a more pronounced variation in distribution than Bryde’s whales (Fig. 6). Variance in the uncertainty maps for all species during both seasons was low, with values less than 0.16, indicating close agreement between algorithms in the ensemble models (Fig. 7). When looking at the individual algorithms and their predictive ability, the RF algorithm produced probability maps that were more similar to the ensemble model predictions than the other individual algorithms did (Figs. S4A–S4H).

Model performance and sensitivity

The quality of data used in this study varies as the data are obtained from both opportunistic platforms and dedicated scientific surveys. Associated effort, which determines absence data, was not recorded in the opportunistic datasets which, despite its constraints, can be useful for cetacean modelling (Torres et al., 2013). Associated effort provides a comprehensive set of absence data over space and time, enabling more robust models to be built. The selection of pseudo-absences (in the lack of absence data) plays a pivotal role in the outcome of the predicted species distribution. By selecting different numbers of pseudo-absence points based on the recommendations of Barbet-Massin et al. (2012) for each algorithm, this study, along with spatial thinning, was able to produce ensemble models using high performing SDMs.

Certain data in this study, due to their opportunistic nature, were spatially biased (Fig. 1). Spatial thinning was used to reduce this bias, but the predicted occurrence (Fig. 6) showed strong evidence that it had not been eliminated completely. When interpreting the results of this study, readers need to be aware of this and take biases into account when applying this study for spatial planning, management and conservation.

Despite its varied data and the biases, this study produced high performing models (Fig. 3), although the ensemble models did not perform as well as expected. Our study, and a study done by Hao et al. (2020), indicates that ensemble models do not always outperform single algorithms as they produce varied results. Hao et al. (2020) suggest that ensemble models can be outperformed by single algorithms when predicting to distant areas and when single algorithms are finely tuned. In this study, the RF algorithms, on average, produced the highest AUC scores (Table S2). RF algorithms in this study have complex trees with no spatial blocking during validation and are, therefore, prone to overfitting (Hao et al., 2020). Looking at Fig. 5, the response curves of the RF algorithms indicate that overfitting was not substantial, especially when compared to the response curves of the GBM algorithms (Fig. S3C). The response curves of the GBM algorithms are flat, indicating overfitting (Derville et al., 2018). One way to address this would be to use spatial blocking (Hao et al., 2020). In spatial blocking, data are split into different blocks used for either calibrating or validating models, allowing the testing data to be spatially distant from the training data, making it more independent (Hao et al., 2020).

Key findings

Management implications

This study shows that southern Africa provides important habitats for these four cetacean species and, according to the literature, fulfil requirements during different life stages in both coastal and offshore habitats throughout the year (Best, 2007). The importance of these areas underpins their need for protection from increasing anthropogenic activities.

Golden et al. (2017) highlight the plight of ocean decision makers in the 21st century: maintaining and expanding an ocean economy to supply the demands of the rapidly growing world population while balancing resource protection through ocean governance. South Africa is not unique in this regard and is currently trying to maintain this balance in the expansion of its ocean economy (Findlay, 2018; Van Wyk, 2015) through the recently implemented MPAs (Findlay, 2020). At a 90% probability of occurrence, the MPA network in South Africa protects less than 10% of seven out of the nine odontocete species modelled habitats (Purdon et al., 2020a; Purdon et al., 2020b). Our results indicate that the predicted occurrence for the four cetaceans studied in this paper is widespread in South Africa, with a small percentage of their habitat being protected by the current MPAs. Important marine mammal areas (IMMAs) could help in identifying more areas that could be delineated as MPAs (Corrigan et al., 2014) in the future.

IMMAs do not have the same regulatory jurisdiction as MPAs over the protection of marine fauna, but they can aid in marine spatial planning, conservation and management of the oceans (Agardy et al., 2019). Delineation of IMMAs is based on four categories, (1) species status based on the International Union for the Conservation of Nature red list, (2) marine mammal distribution and abundance, (3) critical life stage habitats for behaviours such as breeding, feeding and migrating and (4) unique characteristics such as species richness and diversity (Agardy et al., 2019). Recently these criteria were used to identify two IMMAs in South Africa. The first one stretches along the 200 m isobath up the east coast of South Africa and Mozambique. This IMMA was proclaimed to identify the calving areas of humpback whales (Marine Mammal Protected Areas Task Force, 2019a). The second IMMA covers the southern coastal and shelf waters of South Africa, identifying the Indian Ocean humpback dolphin (Sousa plumbea), Bryde’s whale (B. edeni), Indo-Pacific bottlenose dolphin (Tursiops aduncus), common dolphin (Delphinus delphis) and the Cape fur seal (Arctocephalus pusillus) habitats (Marine Mammal Protected Areas Task Force, 2019b). Within these two IMMAs are not only the species that meet the IMMA criteria but others too. For example, within the second IMMA there are the calving areas of the southern right whale on top of the other eight species it was identified for. These IMMAs, however, include mainly coastal species, reducing the chances of more offshore protection when new MPAs are delineated for offshore species like the sperm whale that occurs on the shelf edge and further offshore. With more distribution modelling and better data, including age and sex classes, critical sperm whale habitats could be identified for IMMAs and eventually protected through MPAs.