Linking geomorphological processes and wildlife microhabitat selection: nesting birds select refuges generated by permafrost degradation in the Arctic

To gain better insight into the cascading impact of warming-induced changes in the physical landscape on biodiversity, it is crucial to better understand links between abiotic and ecological processes governing species distribution. Abiotic processes shaping the physical characteristics of the environment could significantly influence predator movements in the landscape and ultimately affect biodiversity through interspecific interactions. In the Arctic tundra, the main terrestrial predator (Arctic fox) avoids patches of wetlands composed of ponds with islets that can act as refuges for prey. Little is known about the geomorphological processes generating islets selected by prey species. Our study aimed to identify (i) the physical characteristics of islets selected by Arctic-nesting birds and (ii) the geomorphological processes generating islets available in the landscape. Over two breeding seasons, we determined the occurrence of nesting birds (cackling goose, glaucous gull, and red-throated loon) on islets (N=396) found over a 165 km² area on Bylot Island (Nunavut, Canada). Occupied islets were located further away from the shore (10.6 m ± 7.3 SD vs. 7.4 m ± 6.8 SD) and surrounded by deeper water (33.6 cm ± 10.6 SD vs. 28.1 cm ± 11.5 SD) than unoccupied islets. As expected, all three bird species selected islets less accessible to Arctic foxes, with nesting occurrence increasing with distance to shore and water depth around islets. Based on high-resolution satellite imagery and field observations, we found that ice-wedge polygon degradation generated the majority of islets (71 %) found in the landscape. Those islets were on average farther from the shore and surrounded by deeper water than those generated by other processes. As polygon degradation is projected to accelerate in response to warming, new refuges will likely emerge in the Arctic landscape, but current refuges could also disappear. Changes in the rate of polygon degradation may thus affect Arctic tundra biodiversity by altering predator–prey interactions.

1 Introduction

The heterogeneity of the Earth's abiotic surface, referred to as geodiversity, is increasingly gaining recognition as a pivotal force shaping the diversity of biological communities (Schrodt et al., 2019; Vernham et al., 2023). Geodiversity, defined as the natural range of geological, geomorphological, and soil features (Gray, 2004), characterizes the available physical environments and can shape species distribution (Burnett et al., 1998; Lawler et al., 2015). Therefore, inclusion of geodiversity in biodiversity research can improve our understanding of biodiversity patterns and our ability to anticipate the impact of climate changes on wildlife (Alahuhta et al., 2020; Brazier et al., 2012; Tukiainen et al., 2022). In this context, it is imperative that we establish robust connections between key abiotic processes affecting the physical landscape and ecological dynamics governing species interactions and distribution.

It is well established that the physical landscape of the Arctic tundra is strongly affected by global warming through permafrost-related changes (Farquharson et al., 2019; Jorgenson et al., 2010; Liljedahl et al., 2016). Climate change is causing deeper active layer development and thaw of permafrost in many Arctic regions (Bonnaventure and Lamoureux, 2013). Climate projections predict higher air temperatures and increased precipitation, and model results indicate that the active layer will likely deepen, and permafrost loss will continue (Farquharson et al., 2019; Lawrence et al., 2008; Shur and Jorgenson, 2007). These changes can affect the surface stability, as well as surface drainage and ponding (Lantz and Kokelj, 2008; Liljedahl et al., 2016), leading to potential alterations of habitats used by wildlife (Berteaux et al., 2017).

Predation is one of the key biotic interactions that can shape species distribution at various spatiotemporal scales (Lima, 1998, p. 199; Menge and Sutherland, 1976; Wisz et al., 2013). Physical characteristics of the environment can hinder predator movements in a landscape (Caro, 2005) and create habitat patches with reduced predation risk, which can be used by prey species to avoid predation (i.e., prey refuges, Gauthier et al., 2015; Sih, 1987). The presence of refuges in the landscape can contribute to the persistence of species vulnerable to predation and partly drive spatial distribution patterns of both predators and prey (Berryman et al., 2006; Holt, 1987; Lima, 1998).

In Arctic vertebrate communities, prey refuges can promote species occurrence and coexistence (Clermont et al., 2021; Duchesne et al., 2021; Léandri-Breton and Bêty, 2020). For example, terrestrial predators like Arctic foxes tend to avoid patches of wetlands composed of ponds with islets that can limit their movements compared to surrounding dryer habitats (Grenier-Potvin et al., 2021). Islets can thus act as important refuges commonly used by tundra prey, such as Arctic-nesting birds (Gauthier et al., 2015; Lecomte et al., 2008; Sittler et al., 2000). Birds using islets can have higher hatching success than those nesting on the shore (Gauthier et al., 2015), and species using islets can be less affected by spatial and annual variation in predation pressure (Duchesne et al., 2021). The presence of refuges like islets can thus modulate species interactions and distribution in the landscape (Clermont et al., 2021; Duchesne et al., 2021). Anticipating the impact of warming on the availability of refuges in the Arctic tundra is challenging due to our limited understanding of the abiotic processes that create the refuges preferred by various Arctic prey species.

In this study, we use a qualitative approach to investigate Arctic geodiversity–biodiversity relationships by assessing how certain geomorphological features may be linked to nest selection by Arctic birds. We aim to (i) identify the physical characteristics (distance to shore and water depth) of islets selected as refuges by Arctic-nesting birds and (ii) identify the main geomorphological processes responsible for forming islets in a high-Arctic tundra landscape (see Fig. 1). We first mapped and characterized the islets found on the southwest plain of Bylot Island, located north of Baffin Island in the Canadian Arctic. We then examined whether islet characteristics affect selection by three tundra bird species known to nest mostly on islets (cackling goose, Branta hutchinsii; glaucous gull, Larus hyperboreus; and red-throated loon, Gavia stellata). We hypothesized that birds would select islets less easily accessible to Arctic foxes (i.e., those farther from the shore and surrounded by deeper water). Using satellite imagery and field observations, we further associated each islet with a specific geomorphological or biotic process underlying its presence in the landscape. As surface hydrology, microtopography, and permafrost dynamics strongly interact in the Arctic (Khani et al., 2023; Liljedahl et al., 2016; Nitzbon et al., 2019; Woo and Young, 2006), we expected that permafrost-related geomorphological processes would generate a large proportion of the islets available in flat lowlands and upland plateaus throughout the study area.

Figure 1 Schematic representation of the link between geomorphological processes, physical characteristics of the landscape, and predator–prey interactions in the Arctic tundra. Birds are expected to select islets less easily accessible to the main predator, the Arctic fox (those farther from the shore and surrounded by deeper water), because it may reduce nest predation risk. Hence, geomorphological processes that generate physical characteristics that hinder Arctic fox movements could influence Arctic birds' distribution pattern and abundance.

2 Methods

2.1 Study area

We conducted summer fieldwork over 2 years (2018–2019) on the southwest plain of Bylot Island, a vast migratory bird sanctuary in Sirmilik National Park, Nunavut, Canada (72°54^′ N, 79°54^′ W). The study area, approximately 165 km², consists of flat lowlands and upland plateaus incised by valleys with glacial rivers. It is characterized by extensive continuous permafrost, with an active layer ranging between 30 and 100 cm deep (Fortier and Allard, 2004) and underlying permafrost reaching depths of up to 400 m (Heginbottom, 1995). Most of the area is covered with mesic tundra in the uplands and an assemblage of mesic tundra and wetlands in the lowlands (Gauthier et al., 2013, 2024). Lakes, ponds, and polygonal wetland complexes are scattered across the study area.

More than 40 bird species nest in the study area or adjacent areas, including waterfowl, shorebirds, seabirds, raptors, and passerines (Gauthier et al., 2024; Lepage et al., 1998). The majority of individuals from three of these species nest on small islets in waterbodies (Fig. 2): the cackling goose, the glaucous gull, and the red-throated loon. The Arctic fox (Vulpes lagopus), a generalist predator, is the main nest predator in our study system (Bêty et al., 2001; Giroux et al., 2012; McKinnon and Bêty, 2009). Predation is the main cause of nest failure for most avian species, and the high abundance of foxes in the study area leads to significant predation pressure in mesic tundra (Beardsell et al., 2022; Clermont et al., 2021; Dulude-de Broin et al., 2023). Avian predators, such as glaucous gulls, ravens (Corvus corax), and jaegers (Stercorarius sp.), do not represent the main cause of nest failure but can nonetheless prey upon eggs of various bird species (Gauthier et al., 1996; McKinnon and Bêty, 2009).

Figure 2 Pictures illustrating a typical waterbody with few islets (a) and a cackling goose nesting on an islet (b) on Bylot Island (Nunavut, Canada). Photo credits: Jeanne Clermont (a) and Yannick Seyer (b).

2.2 Islet characteristics and selection

We georeferenced islets in the study area using a combination of satellite image analyses and intensive field surveys conducted during the bird incubation period (between late June and mid-July). We are confident that we found the vast majority of the islets in the study area. Each year, we visited known islets, stepping on each while taking measurements. We described islets using the following characteristics: (1) distance to shore (hereafter DISTANCE: shortest distance in meters between the shore and the islet; measured on foot, ± 1 m) and (2) water depth (hereafter DEPTH: maximum water depth in centimeters recorded on the shallowest, generally the shortest, route between the shore and the islet; measured using a graduated walking stick, ± 5 cm). These two variables are the ones we aim to focus on, as we hypothesize that these characteristics can impede Arctic fox movement (Strang, 1976). The islet area (hereafter IsletArea; exposed surface of the islet in square meters) and waterbody area (hereafter LakeArea; waterbody surface entirely covered by water in square meters) were also estimated by outlining lakes and islet contours (polygons) on a satellite image (WorldView-3, color and near-infrared, 0.3 m resolution; 2 July 2020) using QGIS software (version 3.16, QGIS Development Team, 2021). Variables used to describe islet characteristics were weakly correlated (Spearman correlation coefficients varied between 0.09 and 0.44; all p<0.10; R package corrplot, version 0.92, Wei and Simko, 2021).

The occurrence of nesting birds on islets was determined annually (summer 2018 and summer 2019). We systematically visited all known islets in the study area between late June and mid-July, when most birds were incubating. When the islet was occupied by an active nest, we identified the nesting species by direct observation of incubating individuals or with egg and/or nest characteristics. In most years, nest predation is low for birds nesting on islets on Bylot Island (Gauthier et al., 2015). Occurrence was also assessed by the presence of fresh nest material and eggshells found in empty nest cups. Although we may have missed a few early-failed nests during our visits, we are confident that the vast majority of unoccupied islets (i.e., no nest was found over the 2-year study period) were not used by nesting birds.

Table 1 Brief description of the main geomorphological processes (1 to 5) and biotic processes (6) that can generate islets on Bylot Island (Nunavut, Canada). Criteria used to assign an islet to a specific process are listed in Appendix A.

2.3 Processes generating islets

We listed all potential abiotic and biotic processes that could generate islets in our study area using the high-resolution WorldView-3 satellite image and visual field observations. Based on extensive knowledge and prior research on the surface landforms in the Arctic tundra of Bylot Island (Ellis and Rochefort, 2004; Fortier et al., 2007; Fortier and Allard, 2004), we listed five main geomorphological processes that can generate islets in our study area (see Table 1). Two main biotic processes could also lead to islets used by Arctic birds (Table 1): plant succession may occur in wet plains, and red-throated loons are known to accumulate vegetation on small, submerged shoals to build nests (Bundy, 1976; Douglas and Reimchen, 1988); both processes can lead to surface aggradation and islet creation. To associate each islet with a specific process, we used a combination of criteria (see Appendix A). Criteria were mainly based on the shape of the waterbody, the nature of the surrounding terrain (e.g., littoral zone, complex polygon wetland, and watershed orientation), and the physical characteristics of the islet.

2.4 Statistical analyses

We used a permutation test to assess whether the mean characteristics of the islets selected by the nesting birds were different from those expected at random. For each characteristic, we compared the mean value obtained for the islets that were occupied at least once to the distribution of the mean value for 1000 random samples with replacement of all known islets (see details in Appendix B; R package stats, version 4.0.3, R Core Team, 2020).

We used logit-link logistic models with a binomial distribution to evaluate the influence of specific islet characteristics on the probability of occurrence of a bird species on the islet (R package lme4, version 1.1-27.1; Bates et al., 2015). An islet was considered occupied (1) when a nest was found on it at least once during the 2-year study period. Otherwise, it was considered unoccupied (0). All predictor variables were rescaled by their standard deviation. To account for a potential nonlinear effect of DISTANCE and DEPTH, we used distance-weighted functions (Miguet et al., 2017). Following Carpenter et al. (2010), we first selected the best fitting decay distance function to transform the distance to shore and water depth according to their declining effect (see the full description in Appendix C). For each of the three bird species, we then used the Akaike information criterion corrected for small sample size (AICc) to select the best supported model among a set of models describing the probability of nest occurrence according to islet characteristics (R package MuMIn, version 1.43.17, Bartoń, 2020). To account for spatial correlation between nest occurrence on islets, we incorporated the geographic coordinates of all islets into our models. The selected models were then tested with and without the coordinates. Parameter estimates were similar with and without spatial variables. All models were compared to a null model in which the probability of occurrence had no association with the variables of interest. We considered models with an AICc less than or equal to 2 to be competitive. Coefficients of the best supported model were used to visualize the results.

The presence of a bird species on an islet may influence the likelihood of finding another species on the same islet. We did not consider such interspecific interactions in our study, and we assumed that it did not affect our ability to investigate the effect of islets characteristics (DISTANCE and DEPTH) on the probability of nest occurrence. This assumption is likely valid because (i) the proportion of occupied islets is relatively low in the study area (24 %) and (ii) the study species can be found on the same islet and can nest very close to each other (minimum distances between nests: 1 m between loons and gulls, 1 m between loons and geese, and 9 m between gulls and geese). Although some very small islets could not be occupied by more than one (or two) species, the high availability of unoccupied islets in the landscape likely allowed most birds to use islets with their preferred characteristics.

3 Results

3.1 Islet characteristics and selection

We found 396 islets spread across 124 waterbodies (lakes, ponds, and wetlands) in the study area. We were able to visit and determine the distance to shore (DISTANCE) and water depth (DEPTH) for most of them (N=350 islets used in the subsequent statistical analyses). Islets were scattered throughout the entire study area, and their characteristics varied substantially (DISTANCE and DEPTH ranging from 1 to 54 m and from 3 to >41 cm, respectively; Fig. 3b). A total of 84 islets out of 350 (24 %) were occupied by a nesting bird (cackling goose, glaucous gull, or red-throated loon) at least once during the study period. The occupied islets were on average located further away from the shore and surrounded by deeper water than all available islets in the landscape (DISTANCE: occupied = 10.6 m ± 7.3 SD; available = 7.4 m ± 6.8 SD; p_DISTANCE= 0.002; DEPTH: occupied = 33.6 cm ± 10.6 SD; available = 28.1 cm ± 11.5 SD; p_DEPTH= 0.002; see Appendix B).

Figure 3 Maps showing (a) the study area (hatched area ∼ 165 km²) on Bylot Island, Nunavut, Canada; (b) the spatial distribution of islets with known characteristics (DISTANCE and DEPTH, N=350); and (c) the geomorphological or biotic processes that generated these islets based on visual field observations and analysis of a high-resolution satellite image (see also Table 1). The islets located in dense clusters were jittered in concentric circles around their centroid to reduce overlap. Geomorphological processes: polygon degradation of low-centered polygon degrading in ridge-like islet (LCP deg.); polygon degradation of flat-centered or high-centered polygon degrading in center-like islet (FCP/HCP deg.); polygon degradation with unknown shape (unknown polygon degradation); raised beach crest degradation and wetland plain degradation or glacial boulders (other process); unclassified (unknown process). Biotic processes: vegetation succession (succ.) or aggradation (aggr.). See Table 1 for more detailed descriptions of the processes. Boundaries and hydrological features were retrieved from the Government of Canada (Province and Territory, 2016; National Hydro Network, 2022). Maps are in WGS 84 and UTM zone 17N and were created by the authors.

Table 2 Generalized linear model selection of the effects of distance to shore (DISTANCE) and water depth (DEPTH) on bird nest occurrence probability on islets. Left panel reports null and competitive models (ΔAICc ≤ 2), with the number of parameters (K), change in AICc from the best supported model (ΔAICc), and Akaike weights (W). An asterisk (^*) indicates that a distance-weighted function was used for a given variable. The right panel report estimated coefficients of the model with the smallest AICc with their 95 % CI. Full model selection is presented in Appendix D.

The probability of nest occurrence on islets was best explained by the distance to shore and water depth around the islets for all three bird species (for each species, the best supported model included both DISTANCE and DEPTH; Table 2). All species selected islets less easily accessible to Arctic foxes, with nesting occurrence increasing with DISTANCE and/or DEPTH (Fig. 4). The presence of at least one weighted function in all selected models suggests that the nest occurrence probability for all species increased nonlinearly with distance to shore and/or water depth. For instance, nest occurrence probability increased sharply with distance after the first few meters and gradually stabilized after ∼ 7 m in gulls (Fig. 4a1). Out of 350 islets for which we also had DISTANCE and DEPTH estimates, we were able to estimate the IsletArea and LakeArea for 315 islets. Rerunning the analyses using this subsample did not change our main results (see full model selection in Appendix D).

Figure 4 Available islets and probability of nest occurrence on islets as a function of distance to shore and water depth (a1 – cackling goose; a2 – glaucous gull; a3 – red-throated loon). The islets used by each nesting bird species are shown using filled dark circles (legend at the bottom right; filled dark circles are the next maximum highest value, never exceeding 12 occurrences). The average characteristics (mean distance to shore and mean water depth) of islets assigned to a specific islet formation process are shown in panel (b) (error bars show 95 % confidence intervals). The number of islets associated with each formation process is indicated in the color legend below panel (b). Probabilities were derived from selected models (see also Appendix D).

3.2 Processes generating islets

Most of the islets found in the study area (328 out of 396, 83 %) were associated with a specific geomorphological or biotic process using visual criteria (see Appendix E). The vast majority (N=281, 71 %) of these islets were generated by polygon degradation (see Appendix F), with almost half (N=177, 45 %) associated specifically with low-centered polygon degradation (Fig. 3c). The same pattern was observed among the islets with known DISTANCE and DEPTH (72 % were generated by polygon degradation; see Appendix G for details and classification of all known islets). In 68 cases, we could not attribute a specific process because some islets were not clearly visible on satellite images, and field observations lacked the detail needed for a single process assignment.

Islet characteristics (DISTANCE and DEPTH) were not homogenous for islets generated by different processes (Fig. 4b). Islets derived from polygonal degradation were on average surrounded by deeper water and farther from the shore than all islets derived from other processes (all grouped together; Wilcoxon signed rank test; DEPTH: polygonal degradation = 28.7 cm ± 11.7 SD; other processes = 24.1 cm ± 11.1 SD; p_DEPTH= 0.013; DISTANCE: polygonal degradation = 8.1 m ± 7.5 SD; other processes = 5.7 m ± 3.9 SD; p_DISTANCE= 0.041. See Appendix G3 and G4 for two-by-two comparisons between all processes).

4 Discussion

The presence of refuges in the landscape can be critical for species vulnerable to predation (Berryman et al., 2006) and is known to promote species occurrence and prey coexistence in Arctic vertebrate communities (Clermont et al., 2021; Duchesne et al., 2021; Léandri-Breton and Bêty, 2020). Many Arctic-nesting birds use islets located in patches of wetlands as refuges (Dahlén and Eriksson, 2002; Mickelson, 1975; Stickney et al., 2002), but little was known about the processes that generate islets with features preferred by birds. In the present study conducted in the Canadian High Arctic, we found that islet characteristics affect nest site selection by three tundra bird species (cackling goose, glaucous gull, and red-throated loon). As expected, birds selected islets located farther from the shore and surrounded by deeper water, which are less accessible to the main nest predator (the Arctic fox). A large proportion (71 %) of islets in the landscape were generated by ice-wedge polygon degradation, which also generated islets on average farther from the shore and surrounded by deeper water compared to those generated by other geomorphological or biotic processes. Few attempts have been made to fully integrate geomorphological attributes or processes in bird nest site selection research (e.g., Eveillard-Buchoux et al., 2019, linking nest site geomorphology to cliff-nesting species preference). To our knowledge, our study is the first conducted in the Arctic that outlines the key role of polygon degradation in the origin of refuges preferred by some prey species.

4.1 Physical characteristics and nest site selection

Our results showing the effects of the water depth and the distance to shore on nest site selection are consistent with previous studies conducted on waterfowl (Giroux, 1981; Hammond and Mann, 1956; Lokemoen and Woodward, 1992) and loons (Eberl, 1993) across North America. However, very few studies were conducted on Arctic-nesting birds and at the microhabitat scale like ours (Dahlén and Eriksson, 2002; Weiser and Gilchrist, 2020). Nest site selection can be influenced by several factors that were not considered in our study. For example, site selection by red-throated loons can depend on lake or pond characteristics (e.g., bottom topography, looseness of pond floor, and distance to the ocean; Douglas and Reimchen, 1988; Eberl, 1993). Adding such variables in our analyses would likely improve our ability to explain the probability of nest occurrence on islets. Furthermore, we could not account for interannual variations in water levels for studied waterbodies. This variation may affect islet availability and characteristics between years and therefore their probability of use.

4.2 Physical characteristics and predator encounter probability

Fine-scale habitat selection is often related to predator avoidance (third scale focusing on broader habitat patch selection and fourth scale focusing on finer-scale microhabitat selection within those patches, Eichholz and Elmberg, 2014; Johnson, 1980), and our results support the hypothesis that birds select nesting sites according to physical characteristics that reduce the probability of encountering their main nest predator. Tundra bird species using islets as microhabitat refuges can partly escape predator-mediated indirect effects generated by changes in the abundance of other prey species and increase their persistence in a landscape characterized by high predation risk (Clermont et al., 2021; Duchesne et al., 2021). The use of islets or islands as refuges can increase nesting success likely due to a reduced access (Kellett et al., 2003). The quality of islets in terms of their capacity to reduce predator access should therefore be based on their physical characteristics that can impede predator movements. Several studies show a decrease in the probability of encountering terrestrial mesopredators (such as foxes, skunks, coyotes, and badgers) with an increased distance to shore, as well as increased water depth (Lokemoen and Woodward, 1992; Strang, 1976; Zoellick et al., 2004). These physical parameters likely reduce the accessibility of nests on islets because mammalian predators must swim to reach them (Mickelson, 1975).

As observed with other mesopredators, Arctic foxes are generally reluctant to swim, either in open water (Petersen, 1990) or through water channels (Zoellick et al., 2004). Moreover, walking in mud seems to be a deterrent for Arctic foxes (Sandra Lai, personal observation, 2017). A complete immersion in cold water or mud, followed by a drying or cleaning process (Dickerson et al., 2012), likely generates significant energetic costs in canids. By reducing their speed, the presence of water may also reduce their ability to successfully attack bird species that are able to protect their nests and fight back. Here we hypothesize that the maximum jumping range and leg length of foxes are likely the two main biomechanical constraints limiting their ability to reach an islet without swimming. For instance, Strang (1976) reported that most of the unsuccessful glaucous gull nests on islets were within fox jumping distance from the shore. If the islet is beyond the maximum jumping distance and the water depth exceeds the leg length, the predator is forced to fully immerse itself to reach the target islet (Zoellick et al., 2004). The nonlinear increase in the probability of occurrence of a gull nest after the very first few meters to shore likely reflects these predator biomechanical constraints and potentially outlines a mechanism explaining fine-scale islet selection based on physical characteristics. Enhancing our understanding of Arctic fox movement within wetland areas and the effects of different biomechanical limitations on their ability to access islets would enhance our capacity to assess the quality of refuges within the landscape.

4.3 Role of biotic and abiotic processes in generating potential refuges

In the Arctic, cryoturbation and frost cracking are the dominant geomorphological processes that shape the ground surface. These processes lead to pronounced microtopographic relief in the form of polygonal networks (Jorgenson et al., 2015; Minke et al., 2009). Polygons typically form due to the repeated freezing and thawing of water that is unable to drain through ground surface in polar regions. Degradation of polygons is a cyclical process typically occurring over decades, driven by the freeze–thaw cycle (French, 2017), resulting in the partial inundation of the landform. Our study highlighted the role of this main geomorphological process, ice-wedge polygon degradation, in the origin of islets as refuges selected by tundra nesting birds. This is likely a result of the study area's inherent structure, which seems to be representative of wet lowlands throughout the Canadian Arctic. Indeed, the low elevation and the predominant arrangement of plateaus, flat lowlands, and depressions throughout the southwest plain of Bylot Island have allowed for the formation of multiple polygon complexes, created by the growth of ice wedges, with a significant water supply over time (Gauthier et al., 2013).

Biotic processes such as vegetation aggradation or succession are the second most common processes that generated islets in the study area (about 10 % of those that could be classified). We may have slightly underestimated the number of islets associated with this category, as they are generally smaller and perhaps harder to interpret in the field or to classify using satellite images. Some of them may have fallen into the category of islets generated by an unknown process. Since plant succession is triggered by minor variation in water levels with the presence of colonizing plants surrounding the waterbody (Magnússon et al., 2020), islets derived from plant succession are less likely to be found in deep water or far from the shore (hence less likely to be selected by birds; see above). These islets can likely be generated over a few decades in the study area. On the other hand, islets derived from vegetation aggradation can be generated within a given year but require a biotic activity, here realized by red-throated loons. These birds are known for building up their nest by gathering mud and decaying vegetation on a shoal in shallow ponds or on emergent grasses and sedges in wet grassy shallow waters, building up “loon-made islands” (Bundy, 1976; Davis, 1972). Water depth surrounding the islets formed by such processes thus usually remained relatively shallow.

The low-lying southwest plain of Bylot Island is mainly the result of marine, fluvio-glacial, and aeolian sediment deposition over tertiary sedimentary rocks, mostly sandstone and shale (Jackson and Davidson, 1975; Klassen, 1993). Therefore, the occurrence of glacial boulders from glacial drift is rather uncommon in the landscape, which likely explains why few islets were due to the presence of such boulders in our study area. Finally, isostatic uplifting, a slow process operating over centuries and still ongoing in a part of Bylot Island, generated a succession of narrow coastal ridges from raised beaches, between which shallow wetlands were formed (Woo and Young, 2003). The degradation of coastal ridges generated few islets in the landscape, and their close parallel organization is more likely to generate islets close to shore, which are less selected by birds.

4.4 Climate change and availability of islets

Considering that ice-wedge polygon degradation can generate a high proportion of islets in the high-Arctic tundra, climate change will likely affect the availability or quality of such refuges through alterations in surface hydrology or shifts in permafrost structure. However, predicting future shifts in islet availability poses a formidable challenge due to the complex interplay among factors affecting ice-rich soil dynamics, coupled with the various temporal and spatial scales over which these changes occur (Bouchard et al., 2020; Francis et al., 2009; French, 2017; Nitzbon et al., 2019). Nonetheless, the ongoing warming trend is accompanied by a rise in extreme seasonal temperature fluctuations and hydrological flux variations, which could potentially exacerbate the natural rate of degradation of ice wedges and the underlying permafrost in the Arctic tundra (Liljedahl et al., 2016). This could rapidly lead not only to the degradation of polygonal complexes in shallow thermokarst ponds, but also to positive feedback amplifying the rate at which these changes occur (Bouchard et al., 2020; Jorgenson et al., 2010).

If we acknowledge that the degradation of ice-wedge polygons by thermokarst in ice-rich soils is a natural and inevitable process on the short and long term (French, 2017), a warming-induced increase in the rate of degradation could further influence the availability of islets, depending on the current extent of degradation observed in wetland areas. In a polygonal environment at an early degradation state, the increase in soil degradation by thermokarst processes could generate a greater number of islets by isolating polygonal emerged structures during the formation of thermokarst ponds, which can progressively coalesce (Hopkins, 1949). The opposite scenario could occur in already well-degraded environments where similar processes could accelerate ground subsidence, leading to the destruction of existing islets by the coalescence of ponds into thermokarst lakes (Bouchard et al., 2020). In this situation, the overall number of islets could eventually decrease, which would thus represent a loss of habitat structures heavily used by some bird species. The lowlands, wetlands, and complex polygonal systems generally exhibit various levels of degradation in the Arctic landscape, as observed within our study area. Additional research, such as field-based or remotely sensed islet monitoring, could provide information on past and ongoing degradation trends, as well as temporal changes in islet availability, allowing for a thorough understanding of historical patterns and dynamics. Incorporating heterogeneity in levels of degradation is necessary for a more comprehensive assessment of how warming impacts the fate of ice-wedge polygons and the availability of islets in the high-Arctic landscape.

By linking geomorphological processes and wildlife microhabitat selection, our study provides fine-grained maps of physical structures that capture ecologically relevant information and improves our knowledge of geodiversity–biodiversity patterns in the Arctic. Making such bridges between abiotic and biotic realms should ultimately improve our understanding of Arctic ecosystem trait diversity (Vernham et al., 2023). The persistence of vulnerable prey can be strongly affected by predation in the Arctic tundra (Beardsell et al., 2023), and changes in the availability of refuges could affect community trait diversity. Due to their relatively high body and egg size, birds such as cackling geese, glaucous gulls, and red-throated loons are likely easy to detect by predators like foxes (Beardsell et al., 2021). However, they do not have the defensive capabilities of larger tundra nesting species, such as greater snow geese and snowy owls (Duchesne et al., 2021). They are mainly found nesting on refuges such as islets and may not be able to persist in a landscape without islets. Climate change is irrevocably altering Arctic ecosystems through multiple mechanisms. Its effects on ice-wedge degradation and their relationship with nest site selection by birds had been little studied before. Given its influence on refuge availability through ice-wedge polygon degradation, islet formation, and changes in islets topography over time, we can reasonably conclude that global warming is likely to alter predator–prey interactions, species occurrence, and distribution in the Arctic landscape.

Appendix A: Criteria for islet classification

Table A1 List of geomorphological processes generating islets on Bylot Island and visual criteria/characteristics used to assign a given islet to a specific process. Example of islets (orange stars) identified on high-resolution satellite images are also shown. Satellite base maps and images © 2020 Maxar Technologies.

Appendix B: Permutation tests for islet characteristics

Figure B1 Permutation test comparing mean DISTANCE to the nearest shore of islets (dashed blue line) occupied by nesting birds to 1000 random samples (gray bars) of 97 out of 396 known islets (N=97; mean = 10.6 m; p_DISTANCE= 0.002${}^{**}$).

Figure B2 Permutation test comparing mean maximum DEPTH on the shortest distance to shore of islets (dashed blue line) occupied by nesting birds (N=97; mean = 33.6 cm) to 1000 random samples (gray bars) of 97 out of 396 known islets (p_DEPTH= 0.002${}^{**}$).

Figure B3 Permutation test comparing mean satellite-derived IsletArea of islets (dashed blue line) occupied by nesting birds (N=97; mean = 19.8 m²) to 1000 random samples (gray bars) of 97 out of 396 known islets (p_IsletArea= 0.102).

Figure B4 Permutation test comparing mean satellite-derived LakeArea of islets (dashed blue line) occupied by nesting birds (N=97; mean = 8230.1 m²) to 1000 random samples (gray bars) of 97 out of 396 known islets (p_LakeArea= 0.078).

Appendix C: Statistical analyses

Distance-weighted functions such as a negative exponential function paired with a distance function enable the consideration of the continuously declining effect of the surrounding landscape on an ecological response with increasing distance from the point where the response is measured (Miguet et al., 2017). Based on our hypotheses about the effects of distance from the shore and islet depth on site selection, distance-weighted functions provided an appropriate model framework for our data structure. We first selected the best fitting decay distance function to transform the distance to shore and water depth according to their declining effect as seen in Carpenter et al. (2010). We transformed each variable according to the equation ${\exp }^{-\mathit{\alpha }/\mathrm{DISTANCE}\mathrm{or}\mathrm{DEPTH}}$, with α  ranging between the minimum and the maximum rescaled distance or rescaled depth value (α_DISTANCE min = 0.04, max = 8; α_DEPTH min = 0.25, max = 4). The resulting values ranged from 0 to 1, with the highest value representing the effect of the variable at a long distance to shore or great water depth.

For each species, we then created a whole set of complete models with geographic coordinates, surface measures, and various DISTANCE and DEPTH decay distance functions. The corrected Akaike information criterion (AICc) was used to determine the best supported models. Decay distance functions in complete models presenting the lowest AICc were considered the most competitive and were then used in all competing global models (DISTANCE = ${\exp }^{-\mathrm{0.04}/\mathrm{DISTANCE}}$ for all species, DEPTH = ${\exp }^{-\mathrm{1.75}/\mathrm{DEPTH}}$ for cackling goose, and DEPTH = ${\exp }^{-\mathrm{4}/\mathrm{DEPTH}}$ for red-throated loon and glaucous gull; R package MuMIn, Bartoń, 2020).

For each species, 37 models including surface variables were built using every combination of chosen decay distance function with untransformed islet characteristics. They have been compared to 10 simplified global models excluding surface variables (see Appendix D). Final models with a AICc ≤ 2 were considered as competitive.

Appendix D: Full model selection

Table D1 Generalized linear model selection of the effects of distance to shore (DISTANCE) and water depth (DEPTH) on cackling geese nest occurrence probability on islets during year 2018 and 2019 on the Bylot Island southwest plain (N=350). An asterisk (^*) indicates that a distance-weighted function was used for a given variable. All candidate models are presented with their coefficient estimates (Int. – intercept, Long. – longitude, Lat. – latitude, NA – not applicable), number of parameters (K), change in AICc from best supported model (ΔAICc), and Akaike weights (W). Models with a ΔAICc ≤ 2 are presented in bold.

n/a: not applicable.

Table D2 Generalized linear model selection of the effects of distance to shore (DISTANCE) and water depth (DEPTH) on glaucous gulls' nest occurrence probability on islets during year 2018 and 2019 on the Bylot Island southwest plain (N=350). An asterisk (^*) indicates that a distance-weighted function was used for a given variable. All candidate models are presented with their coefficient estimates, number of parameters (K), change in AICc from best supported model (ΔAICc), and Akaike weights (W). Models with a ΔAICc ≤ 2 are presented in bold.

n/a: not applicable.

Table D3 Generalized linear model selection of the effects of distance to shore (DISTANCE) and water depth (DEPTH) on red-throated loons' nest occurrence probability on islets during year 2018 and 2019 on the Bylot Island southwest plain (N=350). An asterisk (^*) indicates that a distance-weighted function was used for a given variable. All candidate models are presented with their coefficient estimates, number of parameters (K), change in AICc from best supported model (ΔAICc), and Akaike weights (W). Models with a ΔAICc ≤ 2 are presented in bold.

n/a: not applicable.

Table D4 Generalized linear model selection of the effects of distance to shore (DISTANCE) and water depth (DEPTH), as well as islet surface (IsletArea) and lake surface (LakeArea), on cackling geese nest occurrence probability on islets during year 2018 and 2019 on the Bylot Island southwest plain (N=315). An asterisk (^*) indicates that a distance-weighted function was used for a given variable. All candidate models are presented with their coefficient estimates, number of parameters (K), change in AICc from best supported model (ΔAICc), and Akaike weights (W). Models with a ΔAICc ≤ 2 are presented in bold.

n/a: not applicable.

Table D5 Estimated coefficients and 95 %CI for the best supported model of the effects of distance to shore (DISTANCE) and water depth (DEPTH), as well as islet surface (IsletArea) and lake surface (LakeArea), on cackling geese nest occurrence probability on islets during year 2018 and 2019 on the Bylot Island southwest plain (N=315).

Table D6 Generalized linear model selection of the effects of distance to shore (DISTANCE) and water depth (DEPTH), as well as islet surface (IsletArea) and lake surface (LakeArea), on glaucous gulls' nest occurrence probability on islets during year 2018 and 2019 on the Bylot Island southwest plain (N=315). An asterisk (^*) indicates that a distance-weighted function was used for a given variable. All candidate models are presented with their coefficient estimates, number of parameters (K), change in AICc from best supported model (ΔAICc), and Akaike weights (W). Models with a ΔAICc ≤ 2 are presented in bold.

n/a: not applicable.

Table D7 Estimated coefficients and 95 %CI for the best supported model of the effects of distance to shore (DISTANCE) and water depth (DEPTH), as well as islet surface (IsletArea) and lake surface (LakeArea), on glaucous gulls' nest occurrence probability on islets during year 2018 and 2019 on the Bylot Island southwest plain (N=315). An asterisk (^*) indicates that a distance-weighted function was used for a given variable.

Table D8 Generalized linear model selection of the effects of distance to shore (DISTANCE) and water depth (DEPTH), as well as islet surface (IsletArea) and lake surface (LakeArea), on red-throated loons' nest occurrence probability on islets during year 2018 and 2019 on the Bylot Island southwest plain (N=315). An asterisk (^*) indicates that a distance-weighted function was used for a given variable. All candidate models are presented with their coefficient estimates, number of parameters (K), change in AICc from best supported model (ΔAICc), and Akaike weights (W). Models with a ΔAICc ≤ 2 are presented in bold.

n/a: not applicable.

Table D9 Estimated coefficients and 95 %CI for the best supported model of the effects of distance to shore (DISTANCE) and water depth (DEPTH), as well as islet surface (IsletArea) and lake surface (LakeArea), on red-throated loons' nest occurrence probability on islets during year 2018 and 2019 on the Bylot Island southwest plain (N=315). An asterisk (^*) indicates that a distance-weighted function was used for a given variable.

Appendix E: Processes for all known islets

Figure E1 Map showing (a) the study area (hatched area ∼ 165 km²) on Bylot Island, Nunavut, Canada, and (b) the geomorphological or biotic processes that generated these islets (N=396) based on visual field observations and analysis of a high-resolution satellite image and field observations (see also Table 1). The islets located in dense clusters were jittered in concentric circles around their centroid to reduce overlap. Geomorphological processes: polygon degradation of low-centered polygon degrading in ridge-like islet (LCP deg.); polygon degradation of flat-centered or high-centered polygon degrading in center-like islet (FCP/HCP deg.); polygon degradation with unknown shape (unknown polygon degradation); raised beach crest degradation and wetland plain degradation or boulders (other process); unclassified (unknown process). Biotic processes: vegetation succession (succ.) or aggradation (aggr.). See Table 1 for more detailed descriptions of the processes. Boundaries and hydrological features were retrieved from the Government of Canada (Provincial and territorial limits (2016) and National hydro network (2022)). Maps are in WGS 84 and UTM zone 17N and were created by the authors.

Appendix F: Polygon degradation

Figure F1 Schematic representation of the polygon degradation process in wetlands (a–b), generating islets in panel (c), pinpointed by the blue arrows. (a, b) Initiation of ice-wedge degradation can occur over a very warm summer, over two to three consecutive warm summers, or over a warming trend (5–10 years), whereas (b, c) ice-wedge melting and ridge collapse, associated with the development of ponds, occur over decades to centuries. These islets have a wide range of physical characteristics (e.g., distance to shore or water depth).

Appendix G: Variation of islet characteristics

Table G1 All known (N=396) islet characteristics and proportion of islets occupied by a nesting bird (N.occ; % occ.) at least once over 2 years for each islet formation process.

Table G2 Characteristics for all islets with known DISTANCE (in meters) and DEPTH (in centimeters; N=350) with number of islets occupied by a nesting bird (N.occ) and proportion of islets occupied by a nesting bird (% occ) at least once over 2 years for each islet formation process. Mean (meanD; meanP) ± standard deviation (sdD; sdP) and range (minD; minP – maxD; maxP) are provided for each characteristic.

Table G3 Wilcoxon signed rank test comparing mean distances to shore (DISTANCE) between categories of processes generating islets (N=350). Sample size (NP1, NP2), Wilcoxon signed rank test statistic (W_i), p values (p), and their significance (p.signif) are shown (ns – no evidence of difference in mean distances; ^* moderate evidence; ${}^{**}$ evidence; ${}^{***}$ strong evidence).

Table G4 Wilcoxon signed rank test comparing mean water depths (DEPTH) between categories of processes generating islets (N=350). Sample size (NP1, NP2), Wilcoxon signed rank test statistic (W_i), p values (p), and their significance (p.signif) are shown.