Introduction

The use of species distribution models (SDMs, using occurrence data) and species abundance models (SAMs, using abundance data) has grown widely in importance as a quantitative tool in climate change research and conservation planning (Evans et al. 2015). These models combine information on species occurrence or abundance with relevant predictor variables (such as habitat or climate) to infer a species-environment relationship to be used to predict a species’ distribution or abundance across space and time (Engler et al. 2017). SDMs and SAMs provide spatially explicit estimates of habitat suitability, which may be useful in prioritizing ‘hotspot’ areas to conserve, manage, or restore for focal species, assess the future impacts of environmental changes, including climate change and habitat loss, and in assessing and designing monitoring schemes for focal species (Engler et al. 2017). Depending on how SDMs and SAMs are developed, they can be used to optimize spatial predictions, investigate habitat associations, and identify thresholds of effect. As such, they can be employed for conservation planning on a landscape-scale to maximize outcomes of inherently limited funds and effort when planning occurs over large geographical areas.

An example of a conservation and management planning that occurs over a large geographical extent are the Habitat Joint Ventures of the North American Waterfowl Management Plan. These collaborative science-based partnerships were created to address local, regional, and continental goals for waterfowl and wetlands, and deliver conservation programs at large regional levels through their implementation plans. The largest of these, the Eastern Habitat Joint Venture (EHJV), covers 315 million hectares in Canada’s six eastern provinces, with more than a third of Canada’s total wetlands within its boundaries (Eastern Habitat Joint Venture 2017). The EHJV provides a collaborative forum for governmental and non-governmental partners to establish and deliver effective conservation management for migratory and non-migratory birds in eastern Canada. Its main objective is to conserve priority bird habitat, specifically wetlands and associated uplands, to benefit biodiversity while also providing ecological, social, and economic benefits to society (Eastern Habitat Joint Venture 2017). The conservation and management of eastern Canada’s wetlands, and associated habitats, are critical for the breeding, staging, and wintering of millions of waterfowl each year. As wetlands and coastal habitats may be particularly vulnerable to habitat loss or degradation and climate change (D’Arcy et al. 2005; Desgranges et al. 2006; Browne and Dell 2007; Browne and Humburg 2010; Lepage et al. 2011; Singer et al. 2016, 2020), a robust, empirical understanding of the drivers of waterfowl abundance is necessary for informed conservation planning. Habitat management for waterfowl in eastern Canada is challenging due to both the complexity of the habitat (e.g., high diversity of habitat features, competing land use practices) and the diverse relationships that each breeding waterfowl species has with their habitat due to their unique behavior and life history strategies.

The ease of conservation and management planning may vary widely depending on the ecological system being considered. For example, in the Prairie Pothole Region (PPR) of North America, dabbling duck breeding abundance and breeding success is closely tied to the number and size of wetlands as well as land cover and agricultural practices (Emery et al. 2005; Drever et al. 2007; Bloom et al. 2012; Devries et al. 2023). Activities such as wetland restoration and agriculture policy work to protect upland nests is the focus of the Canadian Prairie Habitat Joint Venture, and these activities are extremely effective at increasing landscape-level waterfowl abundance in this region (Prairie Habitat Joint Venture 2014; Devries et al. 2023). In contrast, what constitutes effective waterfowl habitat management and conservation in the EHJV is complex since there is no single habitat feature that drives waterfowl abundance analogous to prairie agriculture. Instead, waterfowl are associated with a wide variety of habitat features including lake size and configuration, wetland abundance, forestry, agriculture, fish and predator abundance, and industrial activity (Adde et al. 2020b, 2021). Additionally, the EHJV encompasses habitat at both ends of the waterfowl density spectrum, from the highest density in the south to the lower density in the north, which presents distinct habitat management challenges.

The most effective conservation and management activities may differ regionally due to the ecological features and threats inherent to the system. Across the EHJV, the southern part of the boreal forest is subject to intensive forest harvest management practices, while the northern part of the EHJV area contains some of the highest proportions of intact, primary forest and still provides good opportunities for large-scale conservation effort (Venier et al. 2018; Wells et al. 2020; Girona et al. 2023). In the southern mixed and deciduous forests biomes, habitat loss and degradation are a key challenge with 2.1 million hectares of wetlands lost in southern Ontario and Quebec along with logging, hydroelectric development, and land use changes, impacting the abundance and species composition of breeding waterfowl (Merendino and Ankney 1994; Maisonneuve et al. 2006; Eastern Habitat Joint Venture 2017). Most of the habitat management and restoration activities undertaken by stakeholders focus on the southern part of the EHJV (Eastern Habitat Joint Venture 2017). As habitat restoration and/or enhancement is costly, understanding the drivers of waterfowl abundance at the landscape-level, in addition to the local level, is important to prioritize the most cost-effective restoration projects and policy solutions. Demographic studies suggest that wetland restoration and creation in southern Ontario is more beneficial to Mallard (Anas platyrhynchos) than upland nesting habitat conservation or agricultural policy (Hoekman et al. 2006a, b). These relationships are context specific, and may exhibit threshold effects, whereby different habitat features are limited in different areas of the EHJV (e.g., areas of intense agriculture versus forestry). Given this complexity, there is a lack of published waterfowl abundance models for eastern Canada (but see Barker et al. 2014; Adde et al. 2020a).

To date, there are three abundance models for waterfowl that cover the geographic extent of the EHJV, all of which were designed to optimize spatial density predictions rather than identify habitat associations. Barker et al. (2014) modeled waterfowl abundance nationally using abundance data collected via fixed-winged/airplane aerial transect surveys and a machine learning approach with 78 predictor variables. Adde et al. (2020a) revised these models to use a hierarchical Bayesian approach and include updated covariates. Focusing on eastern Canada, Adde et al. (2020b) modeled the distribution and abundance of 12 waterfowl species using a machine learning approach and data collected via helicopter plot surveys to model the averaged-abundances over a 22 year period. Together, these abundance maps highlight areas of high abundance for waterfowl across Canada, providing information for mapping areas of interest for potential landscape enhancement that would benefit waterfowl, which will be crucial in identifying key areas for waterfowl conservation.

While the SAMs produced by Barker et al. (2014) and Adde et al. (2020b) provide valuable information on waterfowl distribution ecology and conservation in Canada, their modelling approaches make it difficult to infer habitat associations, which ultimately direct habitat management and conservation actions. While machine learning techniques have shown to predict spatial patterns of abundance better than traditional statistical models, these models can sometimes overfit the data, making the models not as useful for messy and noisy data (Thessen 2016), and they are less than ideal for making biological inferences (Barker et al. 2014). Second, limitations of the data collection must also be considered when drawing conclusions from abundance models. Barker et al. (2014) used fixed-wing surveys which are ideal for broad-scale studies but have lower species detection rates and breeding density estimates compared to helicopter plot survey data (Zimmerman et al. 2012; Cox et al. 2022). While Adde et al. (2020b) used helicopter data for their models and predictions, they averaged abundances over a 22 year period which may mask important inter-annual variation in distribution and abundances and prevents the use of annually dynamic variables.

We used a generalized additive mixed effects model approach (Wood 2017) which enables a clearer inferring of habitat associations with species’ presence/abundance, allows for a quantification of non-linear thresholds of species’ response to predictor variables, and provides additional information for applied conservation and management decisions. We combined data from four, multi-year, breeding waterfowl helicopter and ground plot surveys conducted by Environment and Climate Change Canada’s Canadian Wildlife Service (CWS) from 2001 to 2015 across the EHJV region. Using a hierarchical generalized additive mixed effects model to account for variation among survey type, years, and survey plots size, we developed predictive SAMs for 16 species of waterfowl and Common Loon (Gavia immer; Table 1) using a suite of static and dynamic environmental, climatic, and anthropogenic predictors. This work aimed to identify drivers of waterfowl species densities to forecast the impact of environmental change within the EHJV, supporting key habitat management decisions for waterfowl species breeding in the survey area. We aim to build on previous studies by developing more tailored habitat associations and abundance maps for waterfowl species across eastern Canada, including data from additional areas of the EHJV, and making predictions for species which have previously not been modeled.

Table 1 Model adjusted R2 values and predictor variables used for each species
Full size table

Methods

Study area

The EHJV is North America’s largest Joint Venture (~ 3 million km2) and includes six Canadian provinces (Ontario, Quebec, New Brunswick, Nova Scotia, Prince Edward Island and Newfoundland and Labrador). It is an area that covers both some of the most densely populated areas of the country, as well as large expanses of remote wilderness. Over 30% of the EHJV region is composed of wetlands, including coastal bays and salt marshes, lakeshore marshes, floodplain wetlands, remnant wetlands within agricultural matrix and boreal forests wetlands (Eastern Habitat Joint Venture 2017). It is also a region of historic wetland loss, occurring primarily along the sea and lake coasts (Maritimes and Great Lakes), along major river systems (e.g., St. Lawrence River), and within the productive southern agricultural areas (Eastern Habitat Joint Venture 2017).

Breeding waterfowl surveys

We used data from four different breeding waterfowl surveys conducted within the EHJV region of interest: the Eastern Waterfowl Survey (EWS; Zimmerman et al. 2012; Bateman et al. 2017), the Southern Ontario Waterfowl and Wetlands Plot Survey (SOWPS; Dennis 1974; Merendino et al. 1995), the Southern Quebec Lowland Survey (SQLS; Lepage and Bordage 2013) and the Maritime Agricultural Plot Survey (MAP; Lieske et al. 2012; Fig. 1). All surveys occurred in the spring (April–May) with similar protocols and were coordinated by the CWS. The EWS started in 1990 and integrated into the Waterfowl Breeding Population and Habitat Survey eastern survey area framework in 2004 (Zimmerman et al. 2012; Bateman et al. 2017). During the time frame of interest (2001–2015), the EWS consists of a network of 5 km × 5 km (25 km2) survey plots covering portions of Ontario, Quebec, New Brunswick, Nova Scotia, and Newfoundland and Labrador (Fig. 1) and each survey plot was surveyed two out of four years (Zimmerman et al. 2012). The SOWPS, initiated in 1971, consists of 351 survey plots (0.8 km × 0.8 km) in southern Ontario. The SOWPS is primarily ground-based (93%), with some remote survey plots surveyed via helicopter (1%) or fixed-wing (6%; Dennis 1974; Merendino et al. 1995). The SQLS, initiated in 2004, consists of 200 survey plots (2 km × 2 km) surveyed via helicopter, systematically distributed through the St. Lawrence Valley of Quebec, the Abitibi-Témiscamingue lowlands and the Lac Saint-Jean plain (Lepage and Bordage 2013). The MAP survey was a supplementary helicopter survey that took place in 2008–2012 in New Brunswick and Nova Scotia, and in 2009–2013 in Prince Edward Island, and included a total of 74 survey plots (2 km × 2 km; Lieske et al. 2018). Due to a rotating schedule of plots surveyed each year and operational changes in which plots are surveyed from the EWS and SOWPS, our analysis included 943 unique survey plots of varying sizes for all surveys from 2001 to 2015. For all surveys, species abundance is summarized as the number of groups of birds observed in each survey plot and the composition of those groups. For the analysis we used the indicated breeding pair (IBP) metric based on the number of male, female, and unknown sex individuals counted in each group of birds observed (Bordage et al. 2017). We then summed the number of IBP within each survey plot, which is used as the response variable for our abundance models.

Fig. 1
figure 1

Map of survey plots across study area, including the maritime agricultural plot survey (MAP) in red, the eastern waterfowl survey (EWS) in green, the Southern Quebec Lowland Survey (SQLS) in blue, and the Southern Ontario Waterfowl and wetlands plot survey (SOWPS) in purple. Provinces in white are part of the Eastern Habitat Joint Venture

Full size image

Environmental, climatic, and anthropogenic data

To build the SAMs, we quantified a suite of 15 predictor variables that may influence the breeding habitat use, survival, and fecundity of 16 species of waterfowl and Common Loon (Table 1). To account for varying plot sizes across surveys, we extracted covariates from a 25 km2 sample plot centered around the centroid of each survey plot; extracting covariates at the same spatial scale, regardless of the survey plot size, to standardized their interpretation. We could not use different sample plot sizes for the extraction, since the proportion of a given landscape variable at 4 km2 does not have the same biological interpretation of the proportion of the same variable at 25 km2. Similarly, the total amount of river length would be difficult to compare across survey plots of different sizes. In our model, we account for the incomplete area surveyed in the sample plots with an offset (i.e., the log of the size of the surveyplot). By including the area of the survey plots as an offset in the model we are relating the density of IBP observed with the predictor extracted at the scale of the sample plot. A main assumption of this approach is that the density of IBP is stationary within the sample plot (i.e., is the same inside or outside of the survey plot). However, given that our modelling approach already assume that the relationship between IBP density and environmental, climatic, and anthropogenic predictors is stationary at the scale of the survey area we were comfortable with this assumption.

We calculated the proportion of each landcover class within each sample plot using Ducks Unlimited Canada’s Hybrid Wetland Layer (N. Jones, unpublished report) with a 30 m resolution, similar to other recent waterfowl landscape scale/habitat analyses in Canada (Barker et al. 2014; Adde et al. 2020a, b; Singer et al. 2020). We re-classified the original 29 land cover classes into 12 broad categories, of which we quantified the proportion of six land cover classes, namely Broadleaf Forest, Coniferous Forest, Mixedwood Forest, Upland, Water and Wetland, for each survey sample plot. We used a similar approach to quantifying the mean proportion of agricultural and pasture lands from the annual crop inventory (ACI; 30 m resolution) in 2012 through 2015, the first years the ACI became available for eastern Canada (Fisette et al. 2013). We categorized all land designated as pasture/forage, fallow, or too wet to plant as “forage” and other crops as “agriculture”. Areas outside of the ACI survey coverage were assumed to have no agriculture or pasture, given that the ACI survey coverage is designed to include all agricultural areas in Canada. Each habitat spatial layer was treated independently, and we did not have any expectations for proportions of different land cover classes to add up to 100%. Areas could be classified as both wetlands and agricultural areas, or both agricultural and forested areas (for example). This is both a realistic representation the diversity of habitats that can be found within an agricultural landscape (and allowed for all the covariates to be included in the model without introducing concurvity issues.

We used the National Hydrographic Network (NHN) geodatabase (1:50,000) to quantify the length of lake shoreline, length of rivers, number of lakes and total lake area within each 25 km2square (Natural Resources Canada 2016). We extracted elevation from the Canadian Digital Elevation Model at 50 m resolution (CDEM; NRCAN 2015) by quantifying the mean elevation within each square. For the analysis, the land cover, hydrological, and elevation variables were temporally static (i.e., the same variable was used for each individual site across the 15 year period of the study).

We used the Night Lights of 2012 map from the NASA Earth Observatory (Román et al. 2018), a map of the human-made light illuminating the earth at nighttime (750 m resolution), as a proxy for urbanization at each survey site (McLaren et al. 2018). As the map is only available for a single time period (2012), the same value was used for every year of the survey for each sample plot. While there was a land cover variable available from the ACI for anthropogenic land use, this layer was not available at the same extent to the northern part of the study area, and also did not provide the same measure of density that the Night Lights spatial layer provides.

In addition to the static covariates, we extracted temporally varying covariates. We extracted climatic variables from interpolated weather station data (10 km resolution) provided by Natural Resources Canada (McKenney et al. 2011). We summarized monthly maximum temperatures (°C) and precipitation (mm) seasonally for each year, including spring (March–May), summer (June–August), fall (September–November) and winter (December–February). We used the spring precipitation from the same year as the survey year for each data point, and used previous year’s summer, fall, and winter mean precipitation as a lagged effect of precipitation. Since snowmelt affects waterfowl breeding phenology (i.e., when the birds start nesting), we extracted mean snowmelt day of year (DOY) for each sample plot for each year using snowmelt timing maps (500 m resolution) derived from MODIS (O’Leary et al. 2017). Snowmelt was further selected to determine how well the surveys are timed in relation to the arrival of the birds on their breeding grounds and timing of nesting (Petersen and Savard 2015; Bianchini et al. 2023).

Species abundance models

We built SAMs for 16 species of waterfowl and one waterbird (Table 1) using hierarchical generalized additive mixed effect models (GAMMs) with Tweedie distributions to generate expectations of how abundance was related to environmental, climatic, and anthropogenic variables across eastern Canada, using the mgcv package (Wood et al. 2016) in Program R (R Core Team 2021). We used GAMMs as opposed to generalized linear models as we had no a priori expectations that the relationships between species abundance would be linear. Along with a smooth term for each habitat covariate, we included: a tensor product of longitude and latitude as nuisance spatial variables, and survey year as random variables. We did not impose additional restrictions on the basis functions (k) and automatically selected optimal k via model estimation (Wood 2011). There were too few observations of Barrow’s Goldeneye (Bucephala islandica), Long-tailed Duck (Clangula hyemalis), and Surf Scoter (Melanitta perspicillata) to model these species so they were excluded from subsequent analyses.

We checked all model variables for concurvity using the mgcv package (Wood et al. 2016) package in Program R version 4.1.0 (2021-05-18) (R Core Team 2021). We considered covariate combinations acceptable if concurvity was < 0.7. Ultimately, we included 18 habitat covariates in our final models. Temperature covariates had high concurvity with snowmelt, and shoreline covariates were associated with lake and river covariates (> 0.7) and were excluded from the final models. Since the size of the survey plots varied between surveys, we also included the log-transformed area of the survey plot as an offset in the model. The inclusion of the offset thus allowed us to explicitly address the fact that the crew did not always survey the entire 25 km2 sample plot from which the habitat covariates were extracted. All habitat covariates were mean centered and scaled to a standard deviation of 1 after transformation and prior to modeling. We assessed model fit using the adjusted coefficient of determination (R2), which represents the proportion of variance in IBP explained by the model.

Case study: mapping wetland restoration potential

To map wetland enhancement potential for increasing waterfowl density on a landscape scale, we extracted all habitat covariates at a 25 km2 grid across the entire EHJV; dynamic variables were extracted for 2015 to represent the most contemporary period. For each of those grid cells, we predicted IBP density for each species from the GAMMs SAMs. To date, most habitat management work in the EHJV has focused on wetland habitat restoration, conservation, and policy work (Eastern Habitat Joint Venture 2017). To identify priority areas for increasing the amount of wetland to benefit current EHJV priority species (i.e., American Green-winged Teal (Ana carolinensis), American Black Duck (A. rubripes), Canada Goose (Branta canadensis), Mallard (A. platyrhynchos), and Ring-necked Duck (Aythya collaris); Eastern Habitat Joint Venture 2017), we extracted the first derivative (i.e., the slope) of the predicted smooth at the value of the proportion wetland cover for a given grid cell. As our models include a log link, this derivative represents the relative increase in the IBP density with respect to incremental increases in total wetland area for each grid cell. We multiplied the derivative by the predicted IBP density to determine the absolute increase in IBP density that correspond to an increase of 1% in the proportion of wetland area (i.e., 25 hectares) while keeping all the other covariates constant in the cell. This enabled us to evaluate the impact of conservation action focused only on increasing the amount of wetland for each cell of the grid while maintaining all the other explanatory variables in the landscape as fixed. Given that the effects of wetland area on breeding waterfowl are non-linear and show thresholds for each species (see results below), it also enabled us to assess which cells on the landscape would benefit the most from an increase in wetland area.

Results

The EWS spanned the broadest geographic range of all the surveys and had the highest species richness with 26 waterfowl species recorded on the survey, while SOWPS had the lowest with nine waterfowl species (Table 2). While the five most common species by abundance differed between surveys, the most abundant species by IBP in the survey area included Mallard, American Black Duck (hereafter Black Duck), Canada Goose, Ring-necked Duck, and Green-winged Teal (Table 2).

Table 2 Summary statistics for the number of survey years, plots surveyed and species observed for each of the four eastern plot surveys from 2001 to 2015
Full size table

SAM’s predicted abundance

Spatial predictions varied for all 17 species; however, a few similar patterns emerged across species groups (predictive maps for all 17 species can be found in online Appendix). First, Blue-winged Teal (Spatula discors), Mallard and Wood Duck (Aix sponsa) displayed similar patterns of high predicted IBP abundance in the southern portion of the EHJV, including central Ontario for Blue-winged Teal and Mallard (Fig. 2). Black Duck and Northern Pintail (A. acuta) were widespread outside of southern Ontario, though Black Duck were predicted at much higher density. Red-breasted Merganser (Mergus serrator; Fig. 3) and Northern Pintail were predicted to have their highest IBP abundances in northern Quebec and Labrador, with lower predicted IBP abundance in southern Ontario. Green-winged Teal, Canada Goose, Common Goldeneye (Bucephala clangula), Ring-necked Duck, Common Loon, and Common Merganser (M. merganser; Fig. 3) were predicted to be widespread across eastern Canada, rather than displaying high abundances in more concentrated areas such as southern Ontario or northern Quebec. Bufflehead (B. albeola) and Hooded Merganser (Lophodytes cucullatus; Fig. 3) were predicted to be widespread in forested areas, outside of Newfoundland.

Fig. 2
figure 2

Predicted density of Indicated Breeding Pairs (IBPs) for a American Black Duck, and b Mallard in Eastern Canada. Panels c–q shows the fitted relationships between the individual environmental variables and abundance for American Black Duck (pink) and Mallard (blue) based on the generalized additive model across the Eastern Habitat Joint Venture region from 2001 to 2015. Statistically significant relationships shown with darker color and solid lines and non-significant with lighter colors and dotted lines

Full size image
Fig. 3
figure 3

Predicted density of Indicated Breeding Pairs (IBP) per km2, for the three merganser species: a Common Merganser, b Red-breasted Merganser, and c Hooded Merganser

Full size image

Habitat associations

The adjusted R2 of our models ranged from a low accuracy of 0.206 for Bufflehead to a high of 0.855 for Gadwall (Mareca strepera) but averaged 0.528 (Table 1). The covariates which predicted IBP abundance varied across all species (Table 1). Figure 3 depicts habitat associations for Mallard and Black Duck; figures for all species can be found in the Online Appendix.

Water and wetland habitat associations

Many species exhibited critical thresholds for water and wetland habitat. For example, most of the dabbling ducks (Black Duck, Blue-winged Teal, Green-winged Teal, American Wigeon (M. americana), Gadwall, Mallard, and Wood Duck) and Ring-necked Duck occurred in higher density (hereafter “increased”) in sample plots (hereafter “plot”) with higher proportion of wetland habitat within the 25 km2 around the centroid of the plot. The threshold where densities peaked were species-specific and ranged from 5.7% wetland for Wood Duck to a maximum of 14.0% for Green-winged Teal. However, Common Merganser showed a threshold effect in the opposite direction, decreasing in abundance when the percentage of wetland habitat exceeded 11.6%. Abundances of Canada Goose and Common Loon also increased with an increasing proportion of wetlands, but unlike the dabbling ducks, their densities declined after a peak. Peak Canada Goose densities occurred between 13.7 and 18.5% of wetland cover and Common Loon between 4.6 and 6.9% wetland cover. Many species also increased in abundance with total lake area with species-specific thresholds. For example, above approximately 1.3 km2 of lake area per plot, the small divers (Bufflehead, Common Goldeneye, and Ring-necked Duck) as well as American Wigeon plateaued in abundance, but Common Merganser, American Wigeon and Black Duck increased in abundance until lake area reach 3.0–3.5 km2 per plot. In contrast, Common Loon and Canada Goose abundance peaked at lake areas of 3.7–5.8 and 1.5–3.8 km2, respectively. Greater Scaup abundance was positively associated with lake area with no thresholds within a 25 km2 area. The total number of lakes also had a positive effect, with a high threshold, on the IBP abundance of all species. Wood Duck was positively associated with river length, with diminishing returns after a total of 11.8 km of river habitat per plot. Other species’ relationships with river length were more complex than Wood Duck, with difficult to interpret patterns (Black Duck and Canada Goose), while other species showed a complex or non-significant pattern (e.g., Bufflehead, Common Goldeneye, Common Merganser).

Forest habitat associations

Abundances of most of the waterfowl species were significantly associated with forest features, including significant effects for coniferous (16/17 species), broadleaf (14/17 species) and mixedwood forest cover (13/17 species) on the landscape scale (Table 1). In general, coniferous forest cover was associated with declines in waterfowl abundance, while the effects of broadleaf and mixedwood forests varied by species. Mallard and Green-winged Teal tolerated some level of coniferous forest on the landscape, but decreased sharply when the proportion of coniferous forests was greater than 16.1%, and 18.7%, respectively. Blue-winged Teal immediately decreased with the proportion of conifer, with effects plateauing after 10.2%. However, Northern Pintail and Common Merganser consistently decreased with the proportion of coniferous forest. Hooded Merganser, Wood Duck and Ring-necked Duck all increased in abundance with broadleaf forest cover, until they reached a threshold of approximately 19.7%, 8.8%, and 18.5% respectively. Hooded Merganser and Wood Duck, both cavity nesting species, were also positively associated with mixedwood forest while Common Loon and Northern Pintail were negatively associated with mixedwood forest. Mallard abundance peaked when the proportion of mixedwood forest was between 38.0 and 54.0%. Black Duck was significantly associated with broadleaf, mixedwood, and coniferous forest, but the associations were complex and non-linear.

Agricultural habitat associations

The proportion of agriculture and forage lands on the landscape significantly explained the abundance of seven and 11 species respectively in eastern Canada. Green-winged Teal was most abundant when the proportion of agriculture in the plot was 4.0–14.4%, while Black Duck abundance was highest from 39.5 to 54.4%. American Wigeon and Canada Goose abundance increased with agriculture up to thresholds of approximately 6.5% and 34.9%, respectively. Black Duck, Wood Duck, Green-winged Teal, and Northern Pintail increased with forage land, while Canada Goose and Common Goldeneye decreased.

Urbanization habitat associations

We found 11 of 17 species to have a significant association with urbanization, as quantified by the amount of anthropogenic light across the landscape. Breeding abundance of Green-winged Teal, Common Goldeneye, Common Loon and Wood Duck were negatively associated with urbanization, as was Black Duck, though with a non-significant slope. Abundances of Canada Goose, Mallard and Northern Pintail were positively associated with urban areas, though for Mallard there was an inflection point after which high levels of urban habitat led to a decrease in their abundance. The four remaining species showed a more complex relationship with anthropogenic light across the landscape.

Climate habitat associations

In the eastern Canada study area, species abundance was minimally associated with precipitation (3/17 species for same year effects, 9/17 species for lagged effects) and snowmelt date (5/17 species). American Wigeon abundance on the landscape increased with high same-year spring precipitation. Common Merganser, Green-winged Teal, and Canada Goose decreased with precipitation in the previous winter and Greater Scaup (A. marila), Common Merganser and Canada Goose increased with increasing precipitation in the previous summer. Common Goldeneye abundance increased with later snowmelt dates.

Case study: benefits of increasing wetland density for five EHJV priority species

We applied our predictive models for five current EHJV priority species to assess the potential to increase their abundance through increasing wetland density. Based on our models, increasing wetland densities in the south, particularly in southern Ontario, would be predicted to primarily increase the breeding population of Mallard (Fig. 4) and Canada Goose (Fig. 6.4, Online Appendix) while increasing wetland density in the north would primarily benefit Black Ducks (Fig. 4) and Ring-necked Duck (Fig. 16.6, Online Appendix). Increasing wetland habitat along the St. Lawrence River would likely benefit all five species and specifically Green-wing Teal given it is in their core breeding region (Fig. 2.4, Online Appendix). Under the best circumstance, the predicted increases in breeding density were modest. For example, adding 25 ha of wetland habitat will increase breeding densities of Black Duck by ~ 0.27 IBP/km2 and Mallard 0.42 IBP/km2. However, those increases are ecologically relevant given that the highest IBP densities predicted for Black Ducks and Mallard were 5 and 9 IBP/km2 respectively, in the survey area, and in most cases much lower.

Fig. 4
figure 4

Relative benefits of increasing wetland habitat at the 25 km2 scale as reported in indicated breeding pairs (IBP) per km2of wetland, calculated by multiplying the first derivative (i.e., predicted smooth) at the value of the proportion wetland cover for a given grid cell by the predicted IBP for each cell for the a American Black Duck and b Mallard. Tan background represents breeding range distributions. For simplicity, only grid cells for which wetland restoration is predicted to be beneficial are plotted (i.e. grid cells where the first derivative is positive and the 95% CI does not overlap zero) and the minimum increase in IBP density has been set to 0.001 IBPs per km square

Full size image

Discussion

We found that associations between habitat features and the abundance of waterfowl on the landscape-scale across a large geographic extent of eastern Canada are often non-linear, with many of these associations showing important species-specific thresholds. This suggests that changes in the habitat features via conservation and/or management actions may result in significant changes in target species abundances in some cases, while in other cases may result in little to no additional related change in abundance. We believe that identifying those thresholds is the main strength of our modeling approach, and by leveraging our understanding of the species-specific habitat associations and thresholds, we can identify and discriminate between areas that have the highest benefits for habitat conservation from those with less beneficial effect to help prioritize decision making. We illustrated this approach by deriving a predictive layer for wetland enhancement to serve as a visual tool to identify where increasing wetland cover will likely lead to the most significant increase in breeding abundance for target waterfowl species across the study area (Fig. 4 and Online Appendix).

Predicted waterfowl abundances in eastern Canada

Two patterns emerge from mapping predicted waterfowl abundance in eastern Canada, with clear differences between northern and southern breeding areas. Not surprisingly given their respective core breeding range, Black Duck, Northern Pintail, Greater Scaup, and Red-breasted Merganser were found primarily in the north, while Blue-winged Teal, Mallard, and Wood Duck were found primarily in the south. As such, model performance greatly depends on the overlap between survey data and the core breeding range of a species. Integrating multiple surveys greatly improved the geographical extent of the sample area, and provided additional information on habitat associations that would otherwise not be included, particularly important for southern breeding species. The EWS covers a large area in eastern Canada; however, southern Ontario, the southern Québec lowlands and portions of the Maritimes (e.g., Prince Edward Island) are not sampled by the EWS and include most of the agricultural lands in eastern Canada (Atlas of Canada; https://atlas.gc.ca/). Southern Ontario wetlands and Carolinian forests were historically replaced by agricultural lands, which currently make up half of southern Ontario. Agricultural lands have in turn been declining in total area since 1920 due to urbanization (Eimers et al. 2020), impacting species habitat and abundance if thresholds are crossed. Additionally, the Quebec-based SQLS and Atlantic-based MAP contain broadleaf forests that are not well represented in the EWS, a predictor variable we found to be significant for 15 of 17 species. Broadleaf forest is expected to continue to expand northward into the conifer dominated boreal forest due to climate change (Wang et al. 2020), which our models suggest may benefit Hooded Merganser, Ring-necked Duck and Wood Duck but reduce the abundance of Black Duck, Common Loon and Common Merganser. As such, the inclusion of survey data from all four sources provided more species abundance data from a wider range of habitat conditions for our SAMs, which greatly improved our habitat associations and predictive models for waterfowl species across eastern Canada.

Our predictive maps differ in some ways from those previously created for eastern Canada. Barker et al. (2014) previously modeled waterfowl distribution and abundance across Canada using data collected from fixed-wing aircraft. While the fixed-wing data covers a large portion of the breeding area for many waterfowl (Barker et al. 2014; Adde et al. 2020a), species density is underestimated compared to helicopter plot surveys (Zimmerman et al. 2012; Cox et al. 2022). Moreover, some species, such as mergansers, can only be identified to a generic group level (Zimmerman et al. 2012; Cox et al. 2022) which results in the highly abundant Common Merganser overwhelming results for the less abundant Hooded Merganser and Red-breasted Merganser. For example, Barker et al.’s (2014) abundance map for mergansers closely resembles that of our Common Merganser predictive map while our Hooded Merganser and Red-breasted Merganser predictive maps display considerably different spatial abundance patterns than the Common Merganser. Conversely the integrated EWS and eBird models presented in Adde et al. (2020a) shows variation among the merganser spp similar to that of our models.

Habitat associations

In eastern Canada, Black Duck and Mallard are two of the most common waterfowl species harvested by hunters (Smith et al. 2022), so understanding what influences their abundance patterns will inform habitat conservation planning and harvest management decisions. Of particular interest, both species were influenced by artificial night light, which was used as a proxy for urbanization. Mallards regularly use urban habitats, such as stormwater and city ponds, and are more tolerant of human activity than Black Ducks (Longcore et al. 2020); correspondingly we found Mallard density was positively associated with artificial night light. Black Ducks appear to be less tolerant of human activity and were negatively associated with greater artificial night light, though the slope was non-significant. Both species were associated with water features (i.e., length of lake, number of lakes, and total lake area) but differed in their associations with the forest composition variables. Black Ducks tend to be found in hardwood forests and wooded uplands (Longcore et al. 2020), and we found fluctuating associations with the presence of conifer, mixedwood and broadleaf forests, suggesting Black Duck may respond to more fine scale forest features (e.g., tree composition). Mallards tend to be found in more open habitats (Drilling et al. 2020), and were negatively associated with a higher proportion of conifers and mixedwood forests in the plots.

Uniquely, our models were able to assess differences between habitat associations of the three merganser species, which have traditionally been grouped together in previous models. The Common Merganser prefers open waters, breeding in freshwater lakes and in the upper part of rivers (Pearce et al. 2020), which explains their strong, positive association with total lake area in our model. The Hooded Merganser had a positive association with the presence of conifers, mixedwood and broadleaf forests in our models, congruent with their use of tree cavities for breeding and preference for forested habitats with small ponds, marshes, or slow-moving streams (Dugger et al. 2020). Finally, the Red-breasted Merganser have the most southern winter range of all three Merganser species, wintering in protected bays, estuaries and on the Great Lakes, breeding in boreal lakes where they feed on small fish (Craik et al. 2020), congruent with a positive association with number of lakes and lake size and area in our model. We also found that the Common and Hooded Merganser abundances were associated with artificial night light and were likely more tolerant of human habitation than the Red-breasted Merganser, whose more northern breeding range was not significantly associated with night light. These habitat associations can provide key species-specific information to inform conservation management decisions for each of these species in eastern Canada. This may be particularly important, since each of the species is predicted to be impacted by climate change differently (Adde et al. 2020b). Given their divergent habitat associations, grouping them together, as is common in fixed-wing surveys, is likely to obscure and confuse any patterns in abundance, distribution, and the impact of landscape alterations.

Our models are limited by two factors. First, the resolution at which the breeding observations was compiled may not be (and in fact likely is not) the optimal scale at which the different species respond to the habitat variables. Given that the observations were compiled at the survey plot level, we modelled the habitat association at the level of the largest surveys plot (25 km2). Digitizing the physical locations of the groups observed during the surveys would allow for a more thorough investigation of the spatial scale at which species respond to the landscape during the breeding season (Jackson and Fahrig 2015; Miguet et al. 2016). Second, our models assume a stationary density within the sample plots and a stationary relationship between breeding density and habitat (i.e., a constant relationship throughout the EHJV). While the assumption of stationary density within the sample plots probably created some uncertainty in the precision of the model parameters, we were not concerned that it would bias them. However, relationships between species density and landscape variable may not be stationary over a large geographic areas is (Osborne and Suárez-Seoane 2002; Miller 2012), and drivers of distribution and population dynamics of many species of birds change over large geographical (Whittingham et al. 2007; Roy et al. 2016; Crosby et al. 2019; Doser et al. 2024). Spatial non-stationarity can be challenging to address in ecological models but there are options such as geographically-weighted regression (GWR) and spatially-varying coefficients (SVC) for linear models (Fotheringham et al. 2002; Gelfand et al. 2003). The use of SVCs is still relatively sparse for SDMs and SAMs but these models show a lot of promise (Thorson et al. 2023; Doser et al. 2024). However, there are currently no easy options to address non-stationarity for non-linear models such as GAMs. Recent developments regarding hierarchical generalized additive models might offer some options soon (Pedersen et al. 2019).

Spatially-explicit habitat management and restoration implications

Planning for effective habitat conservation management, particularly on a landscape-scale, requires an understanding of the likely outcomes for focal species weighed against the immediate and long-term costs. Some obvious habitat features, such as elevation, are difficult or cost-prohibitive to manipulate by land managers but still need to be considered when examining habitat features on the landscape. For example, most waterfowl species were associated with areas of low elevation (i.e., areas where wetlands and lakes naturally form) in our models. Similarly, managers must also balance investment in conserving or increasing habitat in areas where the abundance of a priority species is already high versus trying to manipulate the habitat feature(s) that may directly increase the abundance of a priority species in areas of low abundance. However, models, such as the ones we developed, can provide some insight by identifying areas where habitat management efforts would have large impacts on the absolute densities of specific target species.

In the EHJV, wetlands are the key habitat feature in the south targeted for restoration and management to benefit waterfowl. Given that the majority of current priority waterfowl species for the region (Black Duck, Green-winged Teal, Mallard, Ring-necked Duck and Canada Goose; Eastern Habitat Joint Venture 2017) responded positively to wetland area in our model, continuing to prioritize wetland habitat conservation is logical. Our findings provide insight into target amounts of wetland area in a given area and, more importantly, where increasing wetland area should be prioritized in order to be most effective for increasing the abundance of the current priority EHJV species. While all species benefit from increased wetland habitat on the landscape, there is no increase in waterfowl IBPs observed when the proportion of wetland area increases above 20% on the landscape at a 25 km2 scale. Furthermore, the effect of increasing the amount of wetland on waterfowl abundance is mitigated by other surrounding species-specific habitat features. For example, breeding Black Duck, Ring-necked Duck and Green-winged Teal generally showed negative associations with urbanization and intense agriculture. As such, any planning for wetland habitat restoration decisions in heavily agricultural and industrial areas in Ontario and Quebec need to consider that wetland habitat restoration in those areas is unlikely to benefit breeding Black Ducks, Ring-necked Ducks or Green-winged Teals as much as wetland habitat restoration in areas that are surrounded by undeveloped land. Instead, wetland restoration in urban and intense agricultural areas such as southern Ontario or the St. Lawrence River lowlands is more likely to benefit breeding Mallards and Canada Geese. Notably, southern Ontario was historically wetland rich, with 24.8% of the landscape classified as wetlands c.1800 but by 2002 the proportion of wetlands had declined to only 6.8% (Ducks Unlimited Canada 2010). At more northern latitudes, while breeding Black Ducks and Ring-necked Ducks would benefit from increased wetland area in the boreal forest, northern wetlands are largely still intact, and it is probably not economical to attempt to increase wetland area in this region. Conservation measures in the northern boreal could instead target preservation of widespread quality habitat for priority species through land use policy and the creation of protected areas. Based on our model predictions, increasing the area of wetland along the St. Lawrence River Valley and into agricultural mixed-wood landscapes of New Brunswick, Nova Scotia and Prince Edward Island would have the largest impact on breeding waterfowl densities due to the joint potential for wetland restoration and high breeding waterfowl densities of current EHJV priority species. Should managers wish to target other landscape variables that were included in our modeling exercise, we suggest implementing a similar model prediction-based mapping exercise. More broadly, in a world where investment in carbon compensations and land restoration are being considered as regulatory alternatives to compensate for habitat loss due to the cumulative effects of industrial, our results highlight that landscape context matters such that habitat restoration for carbon compensationpurposes must be carefully planned to achieve the expected conservation gains.

To inform EHJV policy decisions, we found positive associations between breeding Black Ducks and Green-winged Teals and pasture lands, suggesting that these species may make use of hayfields and pasture for nesting or prefer low intensity agriculture with natural features such as wetlands and forests left on the landscape (Maisonneuve et al. 2006; Lieske et al. 2012, 2018; Johnson et al. 2020; Longcore et al. 2020). Delaying mowing until after the nesting season is likely beneficial, as nest survival of many bird species declines sharply when haying begins (Hoekman et al. 2006a). Note that for Black Ducks, this relationship is primarily driven by habitat associations in the Maritimes (Lieske et al. 2012); low abundance of Black Ducks in southern Ontario is primarily driven by the latitude/longitude tensor which accounts for geographic variation that is not explained by our suite of habitat covariates, or alternatively outcompeted by the latitude/longitude tensor. Approximately 70% of the St. Lawrence River Valley and southern Ontario wetlands have been converted for agriculture and other purposes (Ducks Unlimited Canada 2010; Pellerin and Poulin 2013). Accordingly, much of the EHJV’s habitat restoration work has focused on restoring wetlands, with limited upland restoration focused on converting lands to perennial agriculture or restoring riparian areas (Eastern Habitat Joint Venture 2017). We found limited associations between breeding waterfowl densities and other types of agriculture, suggesting that many of the declines in relation to agricultural activity are driven primarily by wetland loss in eastern Canada.

Data input into our models impose limitations to interpretation. Though wetlands vary widely in quality, which in turn affects breeding productivity, we were only able to categorize areas as open water or wetland in our model. As wetland habitat maps improve at the continental level, re-analysis with more specific categories of wetlands, and the inclusion of wetland size and number (instead of or in addition to wetland density as included here) will improve our ability to identify areas with high value to waterfowl and high restoration potential (Hanson 2001). Though human light is an excellent index of urbanization, it may not capture low density human development such as cottage development very well; low density development is hypothesized to affect breeding distribution of some species particularly in the boreal forest (Ouellet d’Amours 2009; Lepage et al. 2011). We were also limited to static landcover and wetland classes; dynamic annual layers would better track habitat change. In this model we did not have any expectations for proportions of different land cover classes to add up to one, and all variables used were checked for concurvity to avoid strong interactions between land cover classes. Despite these efforts, we acknowledge that there are multiple ecological mechanisms and relationships at play that were far from fully captured by the variables used and that many of the relationships described may be caused (at least in part) with unmodelled variables.

This study is focused on breeding abundance based on an IBP metric at the 25 km2 scale, not the finer-scale attributes of a specific wetland at the scale of the observation. The presence of a breeding pair may not translate into an increase in future abundance if that pair does not successfully incubate a nest and fledge a brood. While our study focused on developing predictive models for breeding waterfowl and loons in eastern Canada, we did not examine the effect of brooding, staging or wintering habitat on species-specific populations or cross-seasonal effects on female body condition and subsequent breeding success. Future efforts to identify priority migratory staging areas for other focal species (e.g., Wood Duck) will help further refine habitat conservation.

Conclusions

Using hierarchical, non-linear modeling techniques, we were able to infer species-specific habitat associations and quantifiable thresholds of effect for 16 waterfowl and one waterbird species across a large expanse of eastern Canada, which contains approximately one third of Canada’s wetlands. Building on previous models, our work helps identify new areas that are important for breeding waterfowl, either hotspots for all species or species-specific areas. Using these models, we can further improve waterfowl abundance prediction as forest composition and availability of wetlands are altered with climate-driven change eliciting drier and warmer conditions (Gauthier et al. 2015; Adde et al. 2020b). We also demonstrated how we can go beyond static abundance maps by using the habitat associations and derivatives to help inform spatial mapping for habitat conservation, creation, enhancement, and restoration for waterfowl in Eastern Canada. As the quality of the SDMs and SAMs keep evolving, such derivative work of the model predictions should become more common in ecology as they will enable land stewards and decision makers to make more informed decisions about the investment of their limited resources. Applications like our work could be particularly useful in situations where decision makers have to choose between restoration and conservation or when they need to consider multiple species that have dissimilar habitat needs.