Summer distribution and habitat preference of beluga whale social groups in the Eastern Beaufort Sea

Two independent aerial surveys were conducted in 2019 in the southeastern Beaufort Sea, within the Inuvialuit Settlement Region, Northwest Territories, Canada. The surveys primarily covered the Beaufort Sea shelf and Amundsen Gulf, but with additional transect lines in Prince of Wales Strait and west of Banks Island (Longitude = 115°–140° W; Latitude = 68°–73° N). The southeastern Beaufort Sea consists of a continental shelf that extends 50–150 km north of the Mackenzie Delta, with the Mackenzie Canyon to the west and Amundsen Gulf to the east (Hill et al. 1991; Carmack et al. 2004) (Fig. 1). At approximately isobaths 60–100 m, the slope changes rapidly (2°–6°) leading to the deep Canadian Basin (>1000 m) (Weber 1989; Williams and Carmack 2008; Osborne and Forest 2016).

For both surveys, group size was estimated by observers on the plane. Young belugas (i.e., newborn, calf, or juvenile) were categorized as “calf” and were identified based on the relative size compared to a close adult and skin colour, as newborn belugas have dark grey skin that gradually turns white as they grow. Data from the July survey were standardized according to the time of the audio recordings between the primary and secondary observers; duplicated whale sightings were identified and removed based on the time of the sighting (maximum 15 s difference), latitude and longitude, and group size.

All analyses were conducted using R version 4.1.2. (R Core Team, Vienna, Austria). Geospatial data were transformed into an equidistant conic projection for analysis (see Supplementary Material). Transect lines were separated into 8 km segments (a trade-off between the satellite image grid cell resolution of the environmental data (see below) and the number of absences generated), merging the last two segments if the final segment was less than 8 km long. A 1 km buffer with a flat end on each side of the transect was computed, as well as the centroid of each buffered segment. For each of those segments, presence (1) or absence (0) of the following three social group types was indicated:

The SST data layers were acquired by the remote sensor MODIS onboard the Aqua satellite and available on NASA’s OceanColor Web (https://oceancolor.gsfc.nasa.gov/). SST (11 µ daytime) level-3 data were downloaded for 8 day periods from 20 July 2019 to 28 August 2019, with a 4 km cell resolution (0.0417°×0.0417°). The bathymetric (BATHY) grid was retrieved from the General Bathymetric Chart of the Oceans (GEBCO) 2020 (https://www.gebco.net/). The grid resolution was resampled from 0.00417° × 0.00417° to 0.0417° × 0.0417° so that the cell size matched the cell size of the SST grid. Bathymetric slope (SLOPE) was calculated as the degree of change in bathymetry using the GEBCO raster data set. Eight neighbours were used to calculate the slope with the function “raster::terrain” (Hijmans et al. 2021). The shortest distance to the coastline (DIST) was measured with the function “geosphere::dist2Line” (Hijmans et al. 2019). For each centroid of the buffered segment, the shortest geodesic distance to the polygon shapefile representing land boundaries was computed (Administrative boundaries in Canada—CanVec 1 M, Open Government Licence Canada).

The SST, bathymetry, and slope values were extracted at the centroid of each buffered effort segment from the raster layers with the function “raster::extract” (Hijmans et al. 2021) using the “simple” method, returning the cell value where the point falls. If the extraction returned a null value, the “bilinear” method was used to extract a value interpolated from the values of the four nearest cells. If after extracting from both methods, segments still had at least one null covariate value, they were excluded from the analysis. SST was the only dynamic variable and was extracted from the 8 day layer encompassing the day each segment was flown.

Each environmental covariate was scaled and centred around its mean and standard deviation, computed across all segments. After scaling, environmental variables were tested for multicollinearity. Because the variables were not normally distributed, Spearman’s rank-order correlation coefficient, r_S, was computed as a non-parametric measure with the function “GGally::ggcorr” (Schloerke et al. 2021). If two covariates were highly correlated, only one of the two was considered in one model. The value at which a correlation coefficient is considered high may vary according to the purpose of a study (Asuero et al. 2006). In this case, given the model design and study area, variables were considered highly correlated if the coefficient value reached ± 0.68 (Taylor 1990).

Following Pedersen et al. (2019) and McCabe et al. (2021), we used a hierarchical generalized additive model (HGAM) based on a resource selection function design. HGAMs are flexible nonlinear models that allow variation among groups in the shape of the response function (Pedersen et al. 2019). Due to the differences in survey protocol and effort, we fit independent models for each month (i.e., July was modelled separately from August) with a binomial function and complementary log-log link using the “mgcv” package (Wood 2021). In addition to the four environmental covariates, beluga group types were included as a 3-level factor variable. We fit models that account for a global smoother with a shared penalty, with restricted maximum likelihood (REML). This type of model allows a global functional response that accounts for the shared information between all groups (global smoother), in addition to group-specific functional responses that use a shared penalty to produce the same amount of smoothness for each group-level smooth (Pedersen et al. 2019). Details of the model arguments can be found in the Supplementary Materials (Table S1). We tested models with one term corresponding to a univariate smooth (i.e., a single environmental covariate, or global smoother) and a second term with the environmental covariate and the group type factor (i.e., a group-level smoother) (see eq. 1 in Supplementary Materials). We also computed bivariate smooths to test the interactive effect of two environmental variables on the response. In this case, the global smoother is an interaction term between two environmental covariates in addition to the group-level smoother including the group type factor (see eq. 2 in Supplementary Materials).

The best model was determined by the significance of each term and the lowest REML value. The deviance explained and Akaike Information Criterion (AIC) were examined to assess the model performance. The fitted values from the best model for each month were mapped on a hexagonal grid of the study area. The predicted values used the average SST of the two 8-day periods with survey effort in July (20 July 2019–4 August 2019) and the three 8-day periods for August (5 August 2019–28 August 2019) extracted from the raster layers at each point. The environmental variables used for prediction were scaled by the respective mean and standard deviation from the data used to build the models.

The best model for July included the interaction between SST and slope (REML = 441.83) (Table 2). The next best model included the interaction between SST and bathymetry (REML = 443.35). For August, the best model included the interaction between SST and bathymetry (REML = 709.71) (Table 3). The second-best model for August included bathymetry only (REML = 731.94). The model, including SST and slope, did not perform as well for August (REML = 783.79) as for July. Each model showed a global response for the observed population, that is, a similar functional response to the covariate (Pedersen et al. 2019). In the July model, preference for SST peaked around 6–10 °C for all group types (Fig. 3). In August, all beluga groups showed a high preference for areas with a temperature around 2–4 °C and bathymetry of 300–500 m (Fig. 4). All groups also had a low preference for shallow waters (<250 m) (Fig. 4). Details on the functional response curves for each environmental variable can be found in the Supplementary Material (Figs. S3, S4, and S5). The best model for each month included group type as a significant factor, also resulting in a unique response for each social group type for both environmental variables.

We examined the group-level response for each type of social group for each month and the relative preference for specific environmental conditions. In July, slope generated the most variation between group types. Individual adults and groups with calves had a higher preference for high slope (>2°), while groups of adults preferred areas of low slope (<1°). Individual belugas were found mostly at the continental slope, in waters of 6–8 °C. In comparison to the other group types, individual belugas had higher preferences for colder waters (<2 °C) (Fig. 3A). Groups of adults used habitat on the continental shelf with a wider range of preferred temperatures compared to the other two groups, between 6 and 12 °C (Fig. 3B). As with individual belugas, groups with calves also had a higher preference for high slopes and a narrower range of preferred temperatures between 6 and 8 °C (Fig. 3C).

We also investigated results from the second-best model for July, because it included SST and bathymetry, the same significant covariates as the August best model (see below) (Figs. S4 and S6). A similar global response with a high preference for SST of 7–10 °C combined with a depth of < 200 m in all group types is observed. For individual adults and groups with calves, the model also showed a high preference for a bathymetry of 400–600 m with temperatures around 4–7 °C.

In August, individual belugas had high preferences for very deep and cold waters (>1500 m, 1–2 °C), which, considering the survey coverage, is associated with the continental slope and the start of the Canadian Basin (see Fig. S1E and S1F). As mentioned above, individual adults were also found near the 500 m isobath, where temperatures were a little warmer (2–4 °C), but habitat preference was lower than for the two other groups (Fig. 4A). Groups of adults showed a high preference for the habitats preferred by individual belugas, but the former had a higher preference for the 500 m isobath (Fig. 4B). Groups with calves had a higher preference for the area around the 500 m isobath and SST of 2–4 °C (Fig. 4C). Except for one sighting, groups with calves were not found in the Canadian Basin, past the continental slope (Fig. 4F).

To interpret variability in habitat use by group type, we need to understand the composition of these groups. In this study, observations were classified into three group types: individual belugas, groups of adults only, and groups with at least one calf. Specific age or sex of the individuals cannot be confirmed with the survey methods; however, based on observations of beluga groups in other populations, we can infer possible compositions of these group types. For instance, individual belugas are most likely to be large adult males. Male belugas tend to leave mixed-age groups with females once they become independent and remain solitary or associate with other males (Krasnova et al. 2006; Colbeck et al. 2013; O'Corry-Crowe et al. 2018; O'Corry-Crowe et al. 2020). The groups composed of adults included white adults and/or subadult whales and were most likely formed by younger males (Smith et al. 1994; Suydam 2009; Colbeck et al. 2013; Michaud 2014). Groupings of males of similar size are observed in many beluga populations, exhibiting behaviours such as coordinated movement, socializing, rubbing, and aggression (Andrianov et al. 2009; Colbeck et al. 2013; McGuire et al. 2020; O'Corry-Crowe et al. 2020). However, females are not necessarily excluded from groups of adults. The last group type, groups with calves, included mother-calf dyads or groups of mixed-age whales. If calves are present in a group, the group generally includes females but does not necessarily exclude males (Colbeck et al. 2013; Krasnova et al. 2014; McGuire et al. 2020; O'Corry-Crowe et al. 2020). Interactions between mixed-age belugas can lead to the allocare of calves by other adults or juveniles (Krasnova et al. 2009; Krasnova et al. 2014; Aubin et al. 2021). The nursing period of beluga calves is considered an important social learning phase from kin or non-kin (Krasnova et al. 2009; Colbeck et al. 2013; O'Corry-Crowe et al. 2020). These three social types of associations are formed because of the benefits gained from that particular group size and composition and can lead to different behaviours (Lemieux Lefebvre et al. 2018).

In July, SST and slope explain the best distribution of group types. The preferred SST in July 2019 was between 6 and 10 °C, which is consistent with the most selected temperature range of 6–10 °C in habitat models of EBS belugas in August 2007 and 2008 (Hornby et al. 2017). This range of temperature corresponds mainly to the Beaufort Shelf and the eastern Amundsen Gulf when the survey was conducted (Fig. S1A). The SST was much colder (<2 °C) at the northern limit of the survey lines, northeast of the Beaufort Sea and in the Prince of Wales Strait (Fig. S1A). The presence of belugas corresponds to warmer temperatures on the Beaufort shelf, but not in Amundsen Gulf where only a few belugas were observed (Fig. 3D–3F). The lack of beluga detections in areas farther north could be associated with colder water temperatures or environmental variables not measured in this study (e.g., sea ice). The model also showed that compared to the two other group types, individual belugas associate more with colder temperatures. This explains the distribution of single adults in Prince of Wales Strait, the only group type observed in that area in July (Fig. 3D). Belugas travel northeast through Prince of Wales Strait to reach Viscount Melville Sound and M'Clure Strait (Richard et al. 2001). This area of the Arctic Archipelagos is heavily ice-concentrated and to date, only male belugas have been tracked in this area (Richard et al. 2001; Loseto et al. 2006; Hauser et al. 2017; Storrie et al. 2022).

Bathymetric slope was also associated with beluga distribution patterns. Individuals and groups with calves preferred the areas with a higher slope, whereas groups of adults only preferred low slope or the continental shelf. Within the July survey limits, a high slope is associated with the edge of the Canadian Basin (Fig. S1C). A steep slope or abrupt change in bathymetry, such as at the Beaufort Sea shelf-break, enhances the formation of upwelling (Williams and Carmack 2008; Kirillov et al. 2016). Upwelling brings cold nutrient-rich water toward the surface of the water column and increases local primary production (Pickart et al. 2013). Slope has been demonstrated to influence EBS beluga distribution in areas such as Barrow Canyon (western Beaufort Sea) and along the Beaufort slope margin (Asselin et al. 2011; Hauser 2016; Hauser et al. 2017). As bathymetry and slope were strongly correlated, the ecological significance of the model, including SST and bathymetry, is also interesting to examine, as the two variables do not represent the same functions. For the context of this study, slope relates to the biological productivity of upwelling topography (Pickart et al. 2013), while bathymetry links to water masses and rather describes fish communities (Majewski et al. 2017). Two bathymetric ranges had relatively high preference depending on the group type: areas of < 200 m and SST between 7 and 10 °C, and areas 400–600 m deep with water temperature between 4 and 7 °C (see Figs. S6 and S7).

Although surveys in August were conducted across the EBS, belugas were mostly observed in Amundsen Gulf. The best model for beluga distribution included SST and bathymetry. The overall average water temperature in August was colder than the average in July (Fig. S1A and S1D), and accordingly, the preferred range of SST in August was lower than in July. This is similar to Hornby et al. (2017) who report that belugas selected a colder water temperature range (2–6 °C) in August when warmer temperatures were less available. Similar to the July model, individual adults in August had a higher preference for relatively colder temperatures (<2 °C) compared to the other two group types.

In August, bathymetry showed the strongest variability among group types. All three groups showed a preference for a depth of 300–500 m, but single adults had a higher preference for areas in deep waters (>1500 m). Within the limits of the August survey, the Amundsen Gulf and the edge of the Canadian Basin (northwestern Beaufort Sea) had the combination of the two preferred conditions (300–500 m, 2–4 °C) (Figs. 4D–4F). The deepest and coldest areas were the small section above the Canadian Basin. From other habitat models, depths of 200–500 m were also selected more by EBS belugas from aerial surveys in June (Asselin et al. 2011) but differed from the preferred depth range (0–50 m) of belugas observed in August of other years (Hornby et al. 2017). The distribution of belugas in the August 2019 study area is comparable to the distribution and habitat use of tagged belugas observed in August over multiple decades (1990s–2010s). Females were observed using Amundsen Gulf in higher proportion than males, who were using deeper areas such as Viscount Melville Sound and the Canadian Basin (Richard et al. 2001; Hauser et al. 2014; Hauser et al. 2017).

Foraging opportunities potentially explain the distribution of EBS belugas in the summer of 2019. Arctic cod (Boreogadus saida (Lepechin, 1774)) is one of the most abundant fish in the Beaufort Sea and the main prey of EBS belugas (Loseto et al. 2009; Rand and Logerwell 2011; Quakenbush et al. 2015; Choy et al. 2020). Depth is a key predictor of fish assemblages in the Beaufort Sea, and Arctic cod has been primarily found at depth of 350–500 m (Majewski et al. 2016b; Majewski et al. 2017). This matches the preferred depth between 300 and 600 m for all beluga group types in July (second model) and August. The size of Arctic cod also varies with bathymetry due to the availability of their prey in the water column. Majewski et al. (2016b) found that smaller Arctic cod fed on small Calanus prey in shallow waters (<200 m) while larger cod fed on large Themisto prey in deeper waters (>200 m), particularly in the upper slope zone. In Amundsen Gulf, aggregations of Arctic cod were detected at depths of approximately 200 m in December, associated with sea ice cover, but then migrated to greater depths in the following months, up to 500–600 m in May (Geoffroy et al. 2011). Other lipid-rich EBS beluga prey species include capelin (Mallotus villosus (Müller, 1776)), found in high proportions in the diet of female and juvenile belugas, and Greenland halibut (Reinhardtius hippoglossoides (Walbaum, 1792)), a benthic fish found in high proportion in adult males diet (Choy et al. 2020). In a recent survey, capelin was only detected in the Amundsen Gulf, particularly in the Darnley Bay area, while Greenland halibut was absent from this area, found instead on the slope of the Canadian shelf (Majewski et al. 2016a), matching the distribution of the beluga group types found in this study.

In past studies of EBS beluga habitat use, bathymetric features (distance to coast, depth, and slope) were among the main predictors of the July and August distribution, supporting foraging opportunities as a main determinant of beluga habitat use (Hauser et al. 2017; Hornby et al. 2017). However, these studies did not examine the effect of group composition on the distribution. Beluga whales exhibit size dimorphism and size-related energy requirements, which can be reflected in differences in diet between size and sex (Loseto et al. 2008; Marcoux et al. 2012; Choy et al. 2020). Investigation of dive behaviour in beluga males has demonstrated deep benthic and pelagic dives at the continental slope (Storrie et al. 2022). This type of dives can be highly energy-intensive (longer post-dive surface intervals) but are most likely compensated by the abundance and/or superior quality of food (Storrie et al. 2022). This observation aligns with this study’s finding that single adult males were distributed in deeper areas and preferred the continental slope. Males can prey on higher-quality fish by accessing deeper and farther areas compared to the two other group types. Sex differences in foraging behaviour are observed in other animals, especially ungulates (Du Toit 2005). Elephant males, for instance, have more energy to devote to foraging than females, who are limited by reproductive and social needs (Shannon et al. 2006).

In July, groups with calves also preferred the continental slope. Females with calves could be looking for high-quality food as females lactate for about 2 years (Matthews and Ferguson 2015). Many female mammals have demonstrated an increased foraging activity or better food quality intake to compensate for lactation (e.g., bats (Barclay 1989), squirrels (Michener 1998), monkeys (McCabe and Fedigan 2007; Dias et al. 2011; Ruivo et al. 2017), zebras (Neuhaus and Ruckstuhl 2002), and goats (Hamel and Côté 2009)). In a mixed-age group, mothers can benefit from the presence of other group members to care for their calves while they feed (Krasnova et al. 2014).

Groups of adults had a large range of habitat preferences between July and August, from the continental shelf to Amundsen Gulf, to deeper and colder zones. These groups were composed of subadults and/or adults, but the lack of more precise size or age discrimination within groups could have the effect of hiding more specific patterns of habitat preference. As mentioned above, size and age can cause different energetic requirements for an individual and, depending on the age composition of the group, different foraging behaviours could emerge. Diet varies with beluga size, as larger belugas feed on bigger, more fatty fish and farther offshore compared to smaller belugas (Loseto et al. 2008; Loseto et al. 2009; Marcoux et al. 2012; Choy et al. 2017; Choy et al. 2020). In previous studies in the EBS, habitat use of beluga males also depended on their size, sometimes comparable to that of females with calves or sometimes to that of mixed-age groups (Loseto et al. 2006). There may be more spatial segregation by size within groups of adults, but we cannot conclude specific patterns of habitat use by this group type with the classification used in this study.

Predation risk is unlikely to explain the variability in spatial distribution among individual belugas, groups of adults, and groups with calves. As observed in other species (Heithaus and Dill 2002; Mumma et al. 2017; Iranzo et al. 2018), we expected the more vulnerable groups to avoid areas of high predation risk. The main predators of belugas are killer whales (Orcinus orca (Linnaeus, 1758)), polar bears (Ursus maritimus Phipps, 1774), and humans (O'Corry-Crowe 2018). Although a strong increase in killer whale presence has been observed in the eastern Canadian Arctic (Ferguson et al. 2010), no evidence indicates an increase in the predator presence in the western Canadian Arctic and the Beaufort Sea (Higdon and Ferguson 2009; Higdon et al. 2013). Killer whales were not detected by recent acoustic surveys of marine mammals in the Canadian Beaufort Sea either (Pyć et al. 2016). Belugas are present in polar bears’ diet in the Beaufort Sea (Boucher et al. 2019; Florko et al. 2020), although in lower proportion than seals. Polar bear predation on belugas is usually higher in spring than in summer, correlated with heavy sea ice conditions, as polar bears rely on sea ice and possibly ice entrapment events to hunt belugas (McKinney et al. 2017; Florko et al. 2020; Florko et al. 2021). The 2019 aerial surveys were not conducted in areas with high sea ice concentrations, due to the difficulty in detecting whales in sea ice habitats, or were conducted at a time when sea ice had melted. Thus, we cannot infer a predation effect on belugas in proximity to sea ice or cover within the surveyed area. The results of the models in this current study do not show specific patterns of avoiding areas to minimize predation risk, especially by groups with calves. Without drawing any definite conclusions, it seems that the risk of predation is relatively low for the EBS beluga population during the summer, and therefore would not be a major driver of the group types’ distribution.

The results of this study revealed variations and similarities in the habitat preference of the three social groups of belugas studied during July and August 2019. Our results hint at habitat preferences based on foraging and size-related energy requirements, which could create social segregation based on energetics as a by-product (Conradt 2005). Energy requirements vary with group size, sex, and age, which were not taken into account in this study. More explicit patterns of segregation could have been masked within our general group classification. For instance, within groups of adults, different groups could segregate depending on whether they are composed of young subadult males or large older males. Or, groups with calves could select different habitats if it is a mother–calf dyad or a large mixed-age group. Even more complex, there are indications that beluga populations follow a fission-fusion society framework, where groups join or separate, and group composition can change with time (Smith et al. 1994; Whitehead and Van Parijs 2010; Michaud 2014; O'Corry-Crowe et al. 2020). The timing and causes of joining or leaving other groups are not fully understood and probably relate to the benefits or costs of aggregating. Group foraging behaviour has multiple dimensions, where peers benefit from cooperation in access to food, allocare of young ones, and social learning, but also risk resource competition (Whitehead 1996; Benoit-Bird and Au 2009; Daura-Jorge et al. 2012; Vaughn-Hirshorn 2019). More detail on group composition, including sex and age, and group behaviour would improve our understanding of foraging strategies and differences in habitat use of EBS belugas.

We acknowledge that this study has limitations in the data collection and analysis. Aerial surveys capture only a snapshot in time and space of a highly mobile animal in a three-dimensional environment. Both surveys were flown only in perfect conditions (i.e., small waves, no cloud ceiling, no heavy sea ice concentration, etc.) and whale detections and correct age identification were dependent on the expertise of the observers in the aircraft. Calves could have been hidden under larger belugas or individual whales could have been part of a group that was diving at the time. However, aerial surveys have the advantage of covering extremely large areas and giving an overview of the population distribution within an extensive home range. Compared to using telemetry data (a popular method in species distribution models), aerial surveys allowed us to investigate the distribution of different group types, particularly calves who are usually not a target age class for tagging. Even though our group classification was broad and sightings could have been miscategorized, we showed for the first time social segregation based on group types instead of sex of adult belugas like in previous studies (Loseto et al. 2006; Hauser et al. 2017).

Modelling habitat use is a powerful tool but also entails limitations. Because the models used in this study are based on the used and available habitat (in our case, the segments with a group type and the total segments available), the survey effort and area covered define our available habitat. A different survey area, sampling resolution, or analytical resolution changes the available habitat and thus the estimated coefficients of habitat preference (Beyer et al. 2010; Paton and Matthiopoulos 2016). Additionally, the occurrence of animals in a habitat is not necessarily equivalent to an active selection, and a single sighting does not represent the proportion of time spent in one area (Beyer et al. 2010). Therefore, we did not consider the absolute value of preference, but the relative preference within each group. Although not perfect, we believe our models were well suited to the available data and used relevant ecological predictors, presenting useful insights into the species’ ecology. The July and August aerial surveys were conducted with different protocols (i.e., group definition, different observers, altitude flown, and area surveyed) which prevent us from comparing directly the habitat preference between the two months. For this reason, we treated the two months as independent, only using the same modelling methodology. The group classification is dependent on the observer's expertise, but these are methods currently used in visual aerial surveys of belugas to assess population abundance (e.g., Harwood and Norton 1996; Matthews et al. 2017), which was the primary objective of the aerial surveys in summer 2019. The choice of HGAM combined with geographic information systems (GIS) allowed more flexibility in the response compared to other species distribution models and allowed to model a global and group-specific response at the same time (Elith and Leathwick 2009; Pedersen et al. 2019). The results of this study still highlighted similar trends in habitat preferences in two independent models and provided important insights about the influence of environmental conditions on the EBS beluga population.

Collaboration between different institutions has allowed us to work with more information and explore different months in the summer, over different environmental conditions. This study used two different datasets collected independently but revealed similar patterns in distribution and significant environmental conditions for the EBS beluga population. Habitat selection is dynamic, especially in the marine environment and the distribution and resource selection of one population can change over years (Asselin et al. 2011; Hornby et al. 2017). Beluga societies are also dynamic and group composition can change within a season (Smith and Martin 1994; Krasnova et al. 2012; McGuire et al. 2020). Studying and understanding the many variables that influence the behaviour and movement of belugas is complex, but accessible and shared data can improve research as well as support conservation and management challenges (Reichman et al. 2011). Consultations with different knowledge holders can also provide unique perspectives on behaviour and habitat use. In the Arctic, many Indigenous communities have a long history of harvesting beluga whales, which is reflected in Indigenous knowledge (IK) (Breton-Honeyman et al. 2016a; Breton-Honeyman et al. 2021). IK observations of group behaviour, associations, and use of the environment have been recorded in Alaska, Russia, and Canada (Kilabuk 1998; Huntington et al. 1999; Mymrin et al. 1999; Huntington 2000; Lewis et al. 2009; Breton-Honeyman et al. 2016b; Waugh et al. 2018). Bridging both western science and local expertise in habitat modelling is possible (e.g., Brook and McLachlan 2009; Polfus et al. 2014) and can reveal important patterns for the assessment of climate change impacts and implementation of wildlife conservation and management policies (Robinson et al. 2017; Skroblin et al. 2021).