Alaska's climate sensitive Yukon–Kuskokwim Delta supports seven million Arctic-breeding shorebirds, including the majority of six North American populations

ABSTRACT Baseline information about declining North American shorebird populations is essential to determine the effects of global warming at low-lying coastal areas of the Arctic and subarctic, where numerous taxa breed, and to assess population recovery throughout their range. We estimated population sizes on the Yukon–Kuskokwim Delta (YKD) in western Alaska on the eastern edge of the Bering Sea. We conducted ground-based surveys during 2015 and 2016 at 589 randomly selected plots from an area of 35,769 km2. We used stratified random sampling in 8 physiographic strata and corrected population estimates using detection ratios derived from double sampling on a subset of plots. We detected 11,110 breeding individuals of 21 taxa. Western Sandpiper (Calidris mauri), Red-necked Phalarope (Phalaropus lobatus), Dunlin (subspecies C. alpina pacifica), and Wilson's Snipe (Gallinago delicata) were the most abundant taxa. We estimated that ∼7 million individual shorebirds were breeding on the entire YKD in 2015 and 2016. Our surveys of this region provided robust population estimates (coefficient of variations ≤ 0.35) for 14 species. Our results indicate that the YKD supports a large proportion of North America's breeding populations of the Pacific Golden-Plover (Pluvialis fulva), the western population of a Whimbrel subspecies (Numenius phaeopus hudsonicus), a Bar-tailed Godwit subspecies (Limosa lapponica baueri), Black Turnstone (Arenaria melanocephala), a Dunlin subspecies (C. alpina pacifica), and Western Sandpiper. Our study highlights the importance of breeding shorebirds of this relatively pristine but climatically sensitive deltaic system. Estuaries and deltaic systems worldwide are rapidly being degraded by anthropogenic activities. Our population estimates can be used to refine prior North American population estimates, determine the effects of global warming, and evaluate conservation success by measuring population change over time. LAY SUMMARY Many shorebird populations are declining and estimates of population sizes and trends are essential for effective conservation action. Limited data indicated previously that the Yukon–Kuskokwim Delta (YKD) had high densities of breeding shorebirds, but the area had never been thoroughly surveyed. Our results indicate the relatively pristine but climatically sensitive YKD is important to breeding shorebirds in North America, with 21 species and ∼7 million individuals breeding in the region. This area is home to most of North America's breeding Pacific Golden-Plovers, Black Turnstones, Western Sandpipers, and subspecies of Whimbrel, Bar-tailed Godwit, and Dunlin. Habitat loss and degradation at estuaries and delta ecosystems could be driving population declines for many species. Our results contribute to refining American shorebird population estimates, determining the effects of global warming, identifying species and sites of conservation concern both on theYKD and elsewhere, and evaluating conservation success by measuring population change over time. RESUMEN Es esencial la información de base sobre el declive de las poblaciones de aves playeras en América del Norte para determinar los efectos del calentamiento global en las áreas costeras bajas del Ártico y subártico, donde numerosos taxones crían, y para evaluar la recuperación de las poblaciones a lo largo de sus rangos. Estimamos el tamaño de las poblaciones en el Delta de Yukon-Kuskokwim en el oeste de Alaska, en el borde este del Mar de Bering. Realizamos censos terrestres durante 2015 y 2016 en 589 parcelas seleccionadas al azar de un área de 35.769 km2. Utilizamos un muestreo aleatorio estratificado en ocho estratos fisiográficos y corregimos las estimaciones poblacionales utilizando tasas de detección derivadas de un muestreo doble en un subconjunto de parcelas. Detectamos 11.110 individuos reproductores de 21 taxones. Calidris mauri, Phalaropus lobatus, la subespecie C. alpina pacifica y Gallinago delicata fueron los taxones más abundantes. Estimamos que ∼7 millones de aves playeras estaban criando en todo el Delta de Yukon-Kuskokwim en 2015 y 2016. Nuestros censos en esta región proporcionaron estimaciones poblacionales robustas (CV ≤ 0.35) para 14 especies. Nuestros resultados indican que el Delta de Yukon-Kuskokwim alberga una gran proporción de las poblaciones reproductoras de América del Norte de Pluvialis fulva, la población occidental de la subespecie Numenius phaeopus hudsonicus, la subespecie Limosa lapponica baueri, Arenaria melanocephala, la subespecie C. a. pacifica y C. mauri. Nuestro estudio destaca la importancia para las aves playeras reproductoras de este sistema deltaico, relativamente prístino, pero climáticamente sensible. Los estuarios y los sistemas deltaicos en todo el mundo están siendo degradados rápidamente por actividades antropogénicas. Nuestras estimaciones poblacionales pueden utilizarse para refinar las estimaciones poblacionales originales de América del Norte, para determinar los efectos del calentamiento global y para evaluar el éxito de la conservación mediante la medición del cambio poblacional a lo largo del tiempo.


LAY SUMMARY
• Many shorebird populations are declining and estimates of population sizes and trends are essential for effective conservation action.
• Limited data indicated previously that the Yukon-Kuskokwim Delta (YKD) had high densities of breeding shorebirds, but the area had never been thoroughly surveyed.• Our results indicate the relatively pristine but climatically sensitive YKD is important to breeding shorebirds in North America, with 21 species and ~7 million individuals breeding in the region.This area is home to most of North America's breeding Pacific Golden-Plovers, Black Turnstones, Western Sandpipers, and subspecies of Whimbrel, Bar-tailed Godwit, and Dunlin.• Habitat loss and degradation at estuaries and delta ecosystems could be driving population declines for many species.

INTRODUCTION
Most North American shorebirds have declined significantly since the 1970s (Rosenberg et al. 2019, Smith et al. 2020, 2023), but we have incomplete information about the size and trend of many populations (Bart et al. 2007, Andres et al. 2012b, Smith et al. 2020).Large-scale conservation efforts-including the Atlantic, Pacific, and Midcontinent shorebird conservation initiatives in the Americas-depend on accurate population assessments to determine species and sites of high conservation priority, help identify threats, and evaluate conservation success by measuring population change over time (Atlantic Flyway Shorebird Initiative 2015, Senner et al. 2016, Angarita-Martínez 2021).Programs such as the Western Hemisphere Shorebird Reserve Network and Important Bird and Biodiversity Areas require population size estimates to identify sites that meet biogeographic population thresholds for protection (i.e., ≥1% of a geographic population).Population size, distribution, and trend data are crucial for assessing whether a species should be classified as a species of conservation concern (U.S. Fish and Wildlife Service 2021) or listed as a threatened or endangered species (Smith et al. 2018, IUCN 2021).Often, little or no conservation action is taken when the population size and trend of a species are unknown.In contrast, reliable population size and trend data often allow optimal allocation of resources to conserve species most in need (Gerber et al. 2018).
Since the inception of the Program for Regional and International Shorebird Monitoring (PRISM) in 2002 (Skagen et al. 2003), biologists working on the Arctic portion of PRISM have completed studies to estimate the abundance of breeding shorebirds in northern Alaska and across Arctic Canada (Bart and Johnston 2012, P. A. Smith, personal communication).The areas surveyed in northern Alaska include the Arctic National Wildlife Refuge (NWR; Brown et al. 2007), the Teshekpuk Lake Special Area of the National Petroleum Reserve-Alaska (Andres et al. 2012a), and other portions of the Arctic Coastal Plain (Johnson et al. 2007, Bart et al. 2013).While prior PRISM surveys also included ~50,000 km 2 of Arctic and subarctic tundra in western Alaska, <1,000 km 2 of the Yukon-Kuskokwim Delta (hereafter YKD) were included in those efforts (McCaffery et al. 2012).Much of what we know about shorebirds in the YKD region comes from autumn staging studies (Gill and Handel 1990, Handel and Gill 2010, Ruthrauff et al. 2021) and more focused speciesspecific breeding studies (Handel and Gill 1992, Ruthrauff and McCaffery 2005, Johnson et al. 2009, Jung et al. 2016).
Obtaining accurate and precise estimates of population sizes for breeding shorebirds on the YKD is important for several reasons.First, most of the suspected important breeding areas for shorebirds in Arctic and subarctic North America have been surveyed except for the YKD, which may support a large fraction of breeding populations of multiple species and subspecies, for example, Bar-tailed Godwit (Limosa lapponica baueri; McCaffery et al. 2012) and Black Turnstone (Arenaria melanocephala; Handel and Gill 1992).Second, the lack of baseline population size and distribution information makes it difficult to ascertain likely impacts from climate change (e.g., sea-level rise and storm surge), which are projected to degrade or reduce the amount of shorebird habitat in the vast low-lying areas of the Yukon, Kuskokwim, and adjacent rivers (Jorgenson andEly 2001, Terenzi et al. 2014).Impacts from other anthropogenic factors on the YKD (e.g., mineral extraction, and oil, gas, and wind development) also are hard to evaluate without accurate information on the birds breeding there (Alaska Shorebird Group 2019).Deltaic systems, such as the YKD, are among the most productive ecosystems in the world and are disproportionately important to shorebird populations (Gill andHandel 1990, Butler et al. 2001).It is essential to understand the vulnerability of these ecosystems to a changing climate and threats to global shorebird populations (Overeem andSyvitski 2009, Murray et al. 2019).
In 2015 and 2016, we conducted the first extensive, statistically robust surveys of breeding shorebirds across the YKD.The objectives of the study were to (1) obtain population estimates of shorebird taxa breeding across the YKD, (2) compare our estimates to prior North American estimates for shorebirds breeding primarily in the YKD, (3) determine the importance of the YKD to breeding shorebirds in North America, and (4) evaluate the importance of these findings considering direct and indirect anthropogenic threats to the YKD and other deltaic systems throughout the world.

Study Area
Our study area included most of the Yukon-Kuskokwim coastal lowlands physiographic province (Wahrhaftig 1965), which was formed by the Yukon and Kuskokwim Rivers and their tributaries.Habitats on the YKD are predominantly graminoid meadows interspersed with numerous tidal rivers, sloughs, and shallow water bodies (Jorgenson 2000).The western portion of the YKD is extremely flat, with only ~1-m increase in elevation within 7.5 km of the coast (Jorgenson and Ely 2001).The inland areas are drier, hilly in some places, and dominated by dwarf shrubs, mosses, and lichens (Jorgenson 2000).Several volcanic cinder cones, craters, and lava flows are scattered across the lowlands and uplands away from the immediate coast.The boundary of Yukon Delta NWR overlays most of the YKD, although large areas within the boundary of the refuge are not under federal jurisdiction (Figure 1).

General Survey Approach
We followed the Arctic PRISM approach, which uses double sampling to estimate population size (Bart andEarnst 2002, Bart andJohnston 2012).Fieldwork consisted of a combination of rapid surveys by one team (rapid surveyors) at a large sample of plots and more intensive surveys by a second team (intensive surveyors) at a subsample of those same plots.At the subsample of plots, we estimated a detection ratio as the number of breeding individuals detected by a rapid surveyor on one survey divided by the final tally of breeding individuals determined after repeated intensive surveys (usually >30 survey hours per plot).The detection ratio, which accounts for imperfect detection during rapid surveys, was used to adjust the count data for the large sample of plots that received only rapid surveys.We used stratified random sampling based on major habitats to (1) select sites for the rapid surveys, (2) generate stratum-specific density estimates, and (3) by area extrapolation, estimate population totals for the YKD.Taxonomic names for species follow the American Ornithological Society's Checklist of North and Middle American Birds (Chesser et al. 2020).Because previous estimates of shorebird population sizes were generated for full species, subspecies, and biogeographic populations, we used the scientific names of subspecies and English names of populations designated in Andres et al. (2012b) for comparisons to previous estimates.For general presentation and discussion, we use "taxa" for the species, subspecies, or biogeographic populations of shorebirds breeding in the YKD (Supplementary Material Table 1); there is no overlap in breeding populations within the YKD.

Sampling Frame and Stratification
We primarily followed the methods of Bart et al. (2012a) for small-scale stratification to sample large landscapes.To delineate sampling units, we first overlaid the Yukon Delta NWR with a grid of 390,276 16-ha plots (400 × 400 m each) in ArcGIS 10 (Environmental Systems Research Institute, Redlands, CA; Alaska Albers Equal Area Conic NAD83).Next, we used available data on habitat type (Ducks Unlimited 2011) and our knowledge of shorebird breeding habitat requirements to exclude plots unlikely to contain suitable breeding habitat.Specifically, we excluded plots containing >60% clear and turbid water, >60% mudflat/sandbar, >50% unsuitable vegetated habitats (i.e., forested, tall shrub, and alpine dwarf shrub), >5% unsuitable non-vegetated habitats (i.e., rock gravel, urban, snow/ice, burn, and streams), or >5% image obstruction (e.g., cloud and terrain shadow).It is possible that most suitable breeding habitat for alpineassociated shorebirds, and some suitable breeding habitat for boreal-associated species (e.g., Tringa species), was excluded with this approach, but it reduced our sampling frame to a manageable area and included most shorebirds.Given the challenges of obtaining permission to access lands of >50 Alaska Native Tribes and villages, we also limited our spatial sampling to public lands administered by the U.S. Fish and Wildlife Service (USFWS), which we identified using the National Wildlife Refuge Land Status Data for Alaska (U.S. Fish and Wildlife Service 2014).After removing areas unlikely to contain suitable breeding habitat (32,406 plots, 8% of the total) and suitable Alaska Native lands (134,315 plots, 37% of suitable habitat), our final sampling frame included a total of 223,555 plots (35,769 km 2 ; Table 1).
To classify sampling units (grid cells) to strata, we used a modified version of the physiographic regions and subset classifications defined for the Yukon Delta NWR by Jorgenson and Roth (2010).We identified eight physiographic regions as strata: Tidal, Coast, North Coast, Coastal Plain, Lowlands, Riverine, Uplands, and Mountains.Coast and North Coast were subsections of the Jorgenson and Roth (2010) Coast physiographic region; North Coast was delineated as a separate stratum because we expected shorebird density to be higher in this area.This approach improved the precision of our estimates, and provided separate population totals for this area as requested by the Yukon Delta NWR staff (Thorsteinson et al. 1989).The dominant physiographic region (% area) in each grid cell was assigned as the stratum classification.We classified any areas of the YKD not included by Jorgenson and Roth (2010) into one of the eight strata based on habitat types delineated by Ducks Unlimited (2011) and aerial imagery (see Figure 1 for distribution of the strata).
We divided the spatial sampling frame into two sections: (1) a central coast section that was accessed by boat in 2015, and (2) the remaining section which was accessed by helicopter in 2015 and 2016; this partitioning, and the use of boats, was necessary because helicopter use in a section of Alaska Native lands on the central coast (1,426 km 2 ) was not approved at the time of our 2015 survey.To accommodate the limited spatial extent observers could survey by boat in one day, we used a two-stage sampling strategy in the central coast and selected a spatially balanced, simple random sample at each stage using the generalized random tessellation stratification (GRTS; Stevens and Olsen 2004).A two-stage design was beneficial because it facilitated concentrated survey effort in areas accessible by boat yet maintained a representative sample of the area where helicopters were restricted.For the two-stage design, we first created a 1,600 × 1,600 m grid of primary units (8,913 cells in total) and then established a grid of sixteen 400 × 400 m plots (secondary units) within each primary unit.In the first stage, we selected a spatially balanced sample of 52 primary units across the entire central coast section.In the second stage, we selected a spatially balanced sample of three plots (plus one oversample plot to replace plots that were inaccessible [e.g., moose or bear prevented access] during fieldwork) to be surveyed within each randomly selected primary unit.These procedures identified a total of 156 randomly selected plots in the central coast section that were to be surveyed by boat in 2015.For the helicopter survey section, we used the GRTS procedure across all 8 physiographic strata to select a spatially balanced, stratified random sample of 300 and 328 plots in 2015 and 2016, respectively.Allocation of sampling effort across the eight strata was close to proportional allocation, although it included minor adjustments to reduce the variance of the estimates for as many taxa as possible based on exploratory analysis with available data on shorebird density (McCaffery et al. 2012).We also selected 120 and 180 oversample plots in 2015 and 2016, respectively.

Rapid surveys
Each plot was surveyed once by a single observer from May 15 to June 9, 2015, or May 15 to 25, 2016, when shorebirds were most actively displaying as they established territories and initiated nests (Nebel and McCaffery 2003).Each survey lasted 1.25 hr and consisted of an observer systematically traversing the plot and recording the location of all nests, probable nests, pairs, males, females, birds of unknown sex, and groups of shorebirds on a plot map.We differentiated bird activity that occurred within and outside the plot and excluded birds likely not using the plot at all (i.e., birds flying over), allowing us to estimate the number of breeding individuals (defined below) whose territory centroid was within the plot.Observers summarized their plot information immediately after completing the survey.Each member of the rapid survey team also surveyed each of the eight intensive plots in 2016, using these same rapid survey methods (see below) and with no prior knowledge of the species or number of birds on the intensive plots.
Following Andres et al. (2012a), our method to enumerate the number of breeding individuals from the rapid-survey field observations depended on whether the species was monogamous or polyandrous (we did not detect any polygynous species breeding in the study area).For monogamous species, we enumerated breeding individuals on the plot by summing and then doubling each instance of a pair, nest, probable nest, lone breeding male, lone breeding female, and lone breeding individual of unknown sex.This approach assumed all single males showing territorial behaviors, and others recorded as breeding on the plot, were already mated or would be later in the breeding season.For polyandrous species, we used the "adjusted total" metric of Andres et al. (2012a), which enumerated breeding individuals by doubling only the number of pairs and males, nests, and probable nests, but then adding in all other breeding individuals detected alone.We did not double the number of single females (i.e., the non-incubating sex) because it was uncertain whether these individuals had 0, 1, or >1 mates.Biases in adult sex ratios in shorebirds (Székely et al. 2006) support our assumption that there is not always one mate for each single bird.Finally, three species and Pectoral Sandpiper [C.melanotos]) were considered passage migrants because they have never (Red Knot, Sanderling; J. Johnson and B. McCaffery, personal observation), or only very rarely (Pectoral Sandpiper;Farmer et al. 2020), been observed breeding on the YKD.A small fraction of the Pectoral Sandpipers that we encountered may have remained to breed in the study area; however, migrating individuals also exhibit breeding behavior, and it was not possible to distinguish breeding birds from migrants (McCaffery et al. 2012).Therefore, we treated all Pectoral Sandpipers in our study as migrants.Because the large passage population of this species on the YKD was previously undocumented, we estimated the total number of Pectoral Sandpipers in the study area by summing the total number of individuals observed on each plot regardless of their breeding status.We did not estimate the number of migrant Red Knots or Sanderlings.Jorgenson and Roth's (2010) physiographic regions.Survey plots (16 ha) were randomly selected from the USFWS portion of the frame.Alaska Native lands were not surveyed but are included for extrapolation purposes."No.plots surveyed" is the number of rapid surveys conducted at randomly selected plots in each stratum and "% of plots surveyed in stratum" is the stratum sampling intensity (USFWS only).

Intensive surveys
To estimate a detection ratio using the double sampling approach, we surveyed 8 (16-ha) intensive plots in 2016.One set of intensive surveyors made repeated visits to the plots throughout the nesting season.The intensive survey information was used to estimate the actual number of breeding individuals on the plot, which was then compared to the estimated number of breeding individuals thought to be present on these same plots by a set of different observers who conducted only rapid surveys on the plot.Ideally, intensive plots would have been chosen randomly, but this was not possible given the logistics of establishing and maintaining field camps and the need to have sufficient shorebirds nesting within the plot to generate a detection ratio.We located four intensive plots near Boot Lake (62.1020°N, 164.4870°W), a float plane accessible lake thought to have moderate densities of shorebirds, and four near the Kanaryarmiut Field Station (60.3633°N, 165.1255°W),where prior studies found high densities of nesting shorebirds (Figure 1).
To estimate the actual number of breeding individuals in a plot, intensive surveys were conducted for 4 hr each day from May 14 to June 10, 2016 at Boot Lake and from May 17 to June 21, 2016 at Kanaryarmiut.Using intensive area searches and rope dragging (Saalfeld and Lanctot 2015), we found nests by flushing adults from nests or following adults back to nests.In addition, observers recorded all observations of paired birds, territorial males or females, and birds engaged in territory boundary disputes.These observations provided clues for finding nests in future surveys and to estimate the number of pairs likely nesting on a plot (i.e., those pairs whose nests were never found but were believed to be present).Collectively, our intensive surveys maximized our likelihood of enumerating all the breeding pairs on the intensive plots (Smith et al. 2009), which was multiplied by two to yield the number of breeding individuals.

Doubling sampling and detection ratio
To increase our sample size for detection ratios, we combined our 2016 data from intensive plots with identical doublesampling data from previous studies at the YKD and other PRISM sites in Alaska.We assumed that the factor with the greatest impact on the probability of detection is the behavior displayed during the nesting season and that these behaviors would be similar throughout a species' breeding range.Therefore, we combined our detection ratio data with similar data collected using the same protocol at intensive nest plots on the YKD (9 plots in 2001 to 2002; McCaffery et al. 2012), on Alaska's North Slope near Utqiaġvik (6 plots in 2014 to 2015; R. Lanctot, personal communication), the Colville River (20 plots, 1998to 2001;Bart and Johnston 2012), and the Canning River (7 plots, 2002, 2004;Brown et al. 2007).We estimated a taxon-specific detection ratio using a linear mixed model with a random effect for taxon.The response variable in these models was the log-scale ratio (R) of the estimated number of breeding adults from rapid survey to actual breeding adults from intensive surveys (R = log[(detected breeding adults + 1)/(actual breeding adults + 1)]); a constant was added to the numerator and denominator because the log of zero is undefined and to avoid division by zero.We fit a null model with no effects and five different models using species (random effect), year (random effect), and site (fixed effect) as predictor variables in one-and twofactor combinations (Supplementary Material Table 2).For example, the model with species and year effects is as follows: where β 0 is the average detection ratio among all observations; δ j and γ k are random effects of species j and year k for observation i, respectively (normally distributed with mean equal to 0 and constant variance); and ε i is normally distributed random error with mean equal to 0 and constant variance.The linear mixed models were implemented with the lmer function in the R package lme4 (Bates et al. 2015).We used Akaike's information criterion (AIC c ) to evaluate support in the data for each model and made inference from the model with minimum AIC c (Burnham and Anderson 2002).We estimated a detection ratio for nine species using a back-transformation of fitted values, R j = exp Ä β 0 + δ j ä for species j.There were 10 monogamous species that bred in the study area for which we had little (n < 30) or no data on detection ratio (Table 2); for these species, we used the overall average detection ratio from the linear mixed model above (exp [ β 0 ]).We did not adjust the number of breeding individuals of polyandrous species (i.e., Red Phalarope [Phalaropus fulicarius] and Red-necked Phalarope [P.lobatus]) with a detection ratio because we wanted our estimates for these 2 species to be comparable to the "adjusted total" metric of Andres et al. (2012a).

Estimating population size
For the sample of plots in the central coast section surveyed by boat in 2015, we estimated the population total for each taxon using estimators for a two-stage design with simple random sampling at each stage (Thompson 2002: 145).In this design, there are N total primary units (256 ha each) in the population and M i grid cells in the primary unit i, of which m i are sampled; the estimated total for the primary unit i is as follows: where ȳi = (1/m i ) mi j=1 y ij = ŷi /M i.The population total for all N primary units in the population, of which n were sampled, was T = (N/n) n i=1 ŷi, with variance where and for i = 1, …, n,, and µ 1 = 1 n n i=1 y i (Thompson 2002: 145).For the sample of plots accessed by helicopter in 2015 and 2016, we used conventional stratified random sampling to estimate the population total for the 2 years combined (Thompson 2002: 118).In this design, the population total in each stratum, T h, was T h = N h ȳh, where N h was the number of grid cells in the stratum h, and ȳh was the sample mean for the stratum h, in which n h grid cells were sampled: ȳh = (1/n h ) n h i=1 y hi .We assumed that the mean shorebird density on Alaska Native lands was the same as the mean density on lands administered by USFWS and extrapolated our densities to include these areas; N h thus included USFWS and Alaska Native grid cells across all strata (a total of 357,870 grid cells; see Table 1).The population total for all L = 8 strata was T = L h=1 N h ȳh, with variance where was the sample variance from the stratum h.
Population totals were adjusted for imperfect detection using the detection ratios (R) estimated by double sampling: T * j = T j /R j , where T j is the unadjusted estimate of the population total and R j is the detection ratio for taxon j or the average detection ratio (footnotes in Table 3).Following Goodman (1960), the variance of T * j was estimated using the coefficient of variation (CV): 2 ).Approximately 95% confidence intervals for the population total of taxon j was calculated using T * j ± 1.96 × » Var(T * j ).As the helicopter and boat surveys were independent, the combined population total for the entire area was equal to the sum of the two population totals; the variance of the com-bined population total was equal to the sum of the variances.Below we present population size and density estimates for 2015 and 2016 combined by summing the population totals for (1) the 2015 boat-based sampling area on the central coast and (2) the 2015 to 2016 helicopter-based sampling area (all rapid counts from both years).This approach allowed us to accomplish one of our objectives, i.e., provide a single estimate for our study area that incorporates variation across these 2 years to the extent possible (n = 589 total plots).To estimate overall density, we divided the population size by the area (km 2 ) of the sampling frame.We calculated the variance of the overall density estimate using the delta method (Powell 2007).Stratum-specific density estimates were created using stratum means (ȳh), adjusted with the detection ratio as described above for population totals (T * j ), and expressed as birds per kilometer square.The mean stratum density (ȳh ) was calculated using a combination of helicopter data, which were from a stratified design, and boat-based data, which were post-stratified using the same eight physiographic strata.The variance of detection-adjusted, stratum-specific density was calculated using the delta method (Powell 2007).Density estimates are reported as mean ± SE.Differences in density among strata were evaluated with a Kruskal-Wallis (KW) test.When the KW test indicated significant differences among strata, we used the Dunn test for multiple comparisons and adjusted the P-values in the pairwise comparisons with the Benjamini-Hochberg adjustment to maintain an error rate of 0.05 (dunnTest function in the package FSA in R 4.0.2;R Core Team 2020).
TABLE 2. Estimated detection ratios for shorebird taxa breeding in Alaska.Data were collected using double sampling at 17 intensive plots on the Yukon Delta NWR (2001 to 2002, 2016) and 33 plots on the North Slope of Alaska (1998Alaska ( to 2004Alaska ( , 2014Alaska ( to 2015)).The detection ratio estimates (response variable) are from a linear mixed model with random effects of species and year as predictor variables.Each detection ratio sample reflects a single rapid surveyor visiting a single intensive plot.Detection ratio estimates <1 indicate that rapid surveyors detected fewer breeding birds compared to the final estimate of breeding birds after repeated intensive surveys of the same plot."-" = insufficient data to estimate a detection ratio (sample size < 30); for these species we used the average detection ratio (Avg.)across all taxon to adjust population totals (Table 3).LCL and UCL are respectively lower and upper confidence limits.Scientific names are given in Supplementary Material Table S1.1).Density and population size estimates are for breeding individuals and were adjusted using detection ratios (see footnotes and Table 2).Number of individuals detected is the total number of detections on all plots combined (breeding and nonbreeding birds, excluding birds seen flying over the plot).Rows in bold were species that had large CV (>0.35) for population estimates; these estimates should be interpreted with caution.LCL and UCL are lower and upper confidence limits respectively.Scientific names are presented in Supplementary Material Table S1.* 95% LCL replaced with number of breeding individuals recorded in the field because lower confidence limit was less than zero.

Species
J. E. Lyons et al.
Breeding shorebird populations the Yukon-Kuskokwim Delta, Alaska 9

RESULTS
We conducted rapid surveys on 589 of the 223,555 plots available to be surveyed (0.26%, Table 1).The number of surveyed plots per stratum ranged from 21 to 180, with a sampling intensity relative to plots available between 0.09% (Uplands stratum) and 1.45% (North Coast stratum; Table 1).We recorded 11,110 breeding individuals of 21 taxa; Western Sandpiper (C.mauri), Red-necked Phalarope, and Dunlin (C.alpina pacifica) were the most abundant breeding shorebirds encountered (2,424; 3,593: and 1,890 birds detected, respectively; Table 3).No other taxon had >720 individuals detected, and 12 taxa had <100 individuals detected.Individual taxa were detected on between 0.2% and 67.4% of the plots (Table 3); taxa that were observed in higher numbers were observed in more plots (Pearson r = 0.91, t 20 = 9.9, P < 0.001).Between all surveys on the YKD and the North Slope of Alaska (50 intensive plots), we recorded a total of 851 detection ratio samples from 12 species (Table 2).Detection ratios were estimated from 34 to 166 rapid surveys per species; for most species, data were collected at both the North Slope and YKD.Our model selection results indicated an additive model with species and year as predictor variables explained the most variation in our detection ratio data, as this model had the minimum AIC c and nearly all the model weight (Supplementary Material Table S2).Using this model, we estimated the detection ratio for nine species with more than 30 samples and the overall average detection ratio for all species combined (Table 2).Detection ratios varied from 0.45 for Rock Sandpipers (C.ptilocnemis tschuktschorum) to 1.20 for Long-billed Dowitchers (Limnodromus scolopaceus), indicating rapid surveyors were likely to undercount the former and overcount the latter.
The total estimated number of shorebirds in the study area, obtained by summing the taxa-specific estimates, was ~7 million individuals (95% CI: 5,572,645 to 8,404,500; Table 3).For taxa encountered on ≥20 plots, estimated population totals ranged from 10,186 individuals (Red Phalarope) to 3,518,195 (Western Sandpiper).In general, taxa with the lowest densities were also found on the fewest plots, which resulted in high CVs (Table 3).Conversely, population estimates for taxa found on more plots were more precise; 14 taxa had a CV of ≤0.35 (Table 3).In Table 3, we included only taxa with ≥1 breeding individual detected; we excluded 2 nonbreeding Semipalmated Plovers (Charadrius semipalmatus) and a single nonbreeding (flyover) Bristle-thighed Curlew (Numenius tahitiensis).Similarly, we did not include in Table 3 the observation of species considered passage migrants: 306 Red Knots (n = 3 plots), 3 Sanderlings (n = 1 plot), and 1,038 Pectoral Sandpipers (n = 110 plots).The minimum passage population estimate for Pectoral Sandpipers using the stratified random sample estimator (without adjustment for detection ratio) was 298,572 birds (95% CI: 219,511 to 377,633).
Density of the 4 most abundant taxa ranged from ~10 birds km -2 (Wilson's Snipe Gallinago delicata) to ~62 birds km -2 (Western Sandpiper; Table 3).Dunlin (~12 birds km -2 ) and Red-necked Phalarope (~17 birds km - 2 ) were the other 2 species among the 4 highest densities.Density estimates for all other taxa were ≤3.4 birds km -2 (Table 3).Among the 21 taxa observed breeding in the study area, density varied among strata for 17 (81%) taxa (Supplementary Material Table S3).Density of all taxa combined varied significantly among strata, generally declining with distance from the coast, as the highest densities (mean ± SE) were found in the Tidal (~259 ± 47 birds km -2 ) and Coast strata (~232 ± 52 birds km -2 ), followed by the Coastal Plain (~195 ± 60 birds km -2 ; Figure 2; Supplementary Material Table 3).Dunlin and Semipalmated Sandpiper (C.pusilla) occurred in significantly higher densities in the Tidal stratum than other strata, and Black Turnstones and Rednecked Phalaropes occurred in significantly higher densities in Tidal and Coast strata.Density of Red Phalaropes was highest in the Coast stratum.Western Sandpipers were concentrated in the Coastal Plain, where their density was significantly higher (~114 ± 24 birds km -2 ) than other strata; however, this species was relatively widespread and abundant in all strata (overall density ~62 ± 12 birds km -2 ; Table 3; Supplementary Material Table 2).Black-bellied Plovers (Pluvialis squatarola squatarola) were also concentrated in the Coastal Plain, whereas Pacific Golden-Plovers (Plu.fulva) occurred in significantly higher density in the Uplands stratum (Supplementary Material Table 3).American Golden-Plovers (Plu.dominica) were relatively scarce in our sampled area, but like Pacific Golden-Plovers, were concentrated in the Uplands stratum.Density of Long-billed Dowitchers (~2 ± 0.4 birds km -2 ), while moderate compared to other taxa in this study, was significantly higher in the North Coast stratum (~11 ± 2 birds km -2 ) than other strata.Whimbrels (N.phaeopus hudsonicus) tended to have higher densities in Riverine and Uplands strata, and Lesser Yellowlegs (Tringa flavipes) were only recorded in the Riverine stratum.

DISCUSSION
Using a design-based approach to sample 8 physiographic strata, we provide the first large-scale assessment of shorebird populations on a 57,259 km 2 subarctic coastal delta in western Alaska.We estimated that ~7 million birds of 21 taxa were breeding on the YKD from 2015 to 2016.Western Sandpipers, Red-necked Phalaropes, Dunlin, and Wilson's Snipe had the largest populations (in decreasing order of abundance), and the other breeding taxa were relatively less abundant.Our effort was the first on the YKD to include most shorebird habitat types, not just the coastal areas where previous surveys found very high densities of breeding shorebirds (McCaffery et al. 2012).Our population estimates for many taxa were relatively precise for ecological data (i.e., CV < 0.35).Our taxon-and stratum-specific density estimates are useful for comparison with other areas, evaluating the relative importance of habitat types within the YKD, and informing North American population estimatesan important aspect of setting conservation priorities and measuring trends over time.
The most abundant taxa in the coastal habitats were Dunlin, Black Turnstone, Long-billed Dowitcher, Semipalmated Sandpiper, Red Phalarope, and Ruddy Turnstone (A. interpres interpres), whereas those most common in the interior habitats were Red-necked Phalarope, Least Sandpiper (C.minutilla), Pacific Golden-Plover, and Whimbrel.A few other taxa tended to be common both near the coast and inland, including Western Sandpiper, Wilson's Snipe, Black-bellied Plover, and Bar-tailed Godwit.The Hudsonian Godwit (L.haemastica), Lesser Yellowlegs, Greater Yellowlegs (Tringa melanoleuca), and Solitary Sandpiper (Tringa solitaria tended to occur in most inland regions, usually occupying the Lowlands or Riverine strata which typically include tall shrublands and forests.We also detected three species that we considered passage migrants: Pectoral Sandpiper, Red Knot, and Sanderling.Our estimate of 298,572 Pectoral Sandpipers migrating through the YKD was noteworthy because most of the population (1.6 million, Andres et al. 2012b) is thought to migrate through the mid-continent of North America (Farmer et al. 2020).Our observations, combined with prior reports of birds moving northward through British Columbia (Campbell et al. 1990) and Cook Inlet, Alaska (Farmer et al. 2020), suggest a Pacific Coast migration corridor that may connect wintering birds to their breeding locations in Alaska or Siberia.
Only one other study provided population sizes for breeding shorebirds on the YKD.McCaffery et al. (2012) investigated a relatively small area (853 km 2 ) near Hazen Bay that consisted of only 2 strata (wet and moist), compared to our larger (57,259 km 2 , 65× larger) and more diverse (8 strata) study area.Given the vast differences in the amount of area sampled, it is not meaningful to compare our population totals with those of McCaffery et al. ( 2012), but a comparison of bird density estimates may be useful.Including all shorebird taxa, McCaffery et al. (2012) estimated densities of ~300 birds km -2 in the moist stratum and ~400 birds km -2 in the wet stratum.These density estimates are greater than our most dense stratum (Tidal: ~259 birds km -2 ; Coast: ~232 birds km -2 ; and Coastal Plain: ~195 birds km -2 ).These differences could be due to declines in bird density between the early 2000s and mid-2010s or, in some cases, may be related to our use of taxon-specific detection ratios (and no detection ratio adjustment at all for phalaropes) versus McCaffery et al.'s (2012) use of an overall average detection ratio of 0.81 for all taxa.Perhaps more likely, the differences may simply reflect the inherent higher quality of the habitats around Hazen Bay compared to our larger study area.
The total number of shorebirds breeding on the YKD (~7 million) is equal to or greater than the numbers breeding on the North Slope of Alaska (~6.2 million), which includes the National Petroleum Reserve-Alaska (NPR-A; ~4.5 million), the Prudhoe Bay oil field and surrounding areas (~1.4 million), and Arctic NWR (~308,000, Bart et al. 2013).The overall size of the North Slope study area of Bart et al. (2013) was 1.32× larger than our YKD study area (57,259 km 2 vs.  1).Physiographic strata are described by Jorgenson and Roth (2010).Note differences in y-axes.Bar charts of bird density in 8 physiographic strata for Dunlin, Western Sandpiper, Whimbrel, and all species in the study combined.
73,248 km 2 ), confirming higher densities of shorebirds occur on the YKD.A closer look, however, indicates that shorebird densities on the YKD (~123 birds km -2 ) are higher than Prudhoe Bay (81 birds km -2 ; Bart et al. 2013) and the Arctic NWR (27 birds km -2 ; Brown et al. 2007), but lower than the NPR-A (151 birds km -2 ; Bart et al. 2013) or the smaller Teshekpuk Lake Special Area (126 birds km -2 ; Andres et al. 2012a).Like the coastal gradient found on the North Slope (Brown et al. 2007, Johnson et al. 2007, Andres et al. 2012a), overall shorebird density and diversity on the YKD were greatest near the coast and declined with distance from the coast.
Because much of the YKD has not been previously surveyed and no large-scale surveys have been completed elsewhere in western Alaska (McCaffery et al. 2012), our estimates for taxa that breed primarily in western Alaska are substantially greater than the previous population sizes reported by Andres et al. (2012b).For example, our estimate of Alaskabreeding shorebirds is greater than available estimates for Pacific Golden-Plover (2.5×), pacifica Dunlin (1.2×), western Whimbrel (2.3×), and Black Turnstone (1.5×).Based on our results, North America estimates for these taxa could be revised substantially upward from those presented in Andres et al. (2012b).These estimates are a testament to the importance of the YKD for breeding North American shorebirds.Our estimate for baueri Bar-tailed Godwits (113,624) is comparable to the most recent winter-ground count of 126,000 birds (Schuckard et al. 2020), after adjusting for juveniles that remain on the wintering ground and birds that breed in other parts of Alaska.In addition, our estimate for this taxon was only slightly larger than a recent post-breeding count on the YKD by Ruthrauff et al. (2021;100,926 birds).Because our current estimate for Western Sandpiper on the YKD is 98% of the previous total estimate, the total population estimate may likely increase when data from the YKD are considered.Overall, our design-based estimates should contribute to providing more robust estimates of several North American breeding shorebird populations.Some of the differences between our estimates and previously generated population estimates are due to differences in methods.One key difference is whether a sampling frame and a sampling plan were employed for estimation.We used design-based survey sampling methods (Thompson 2002) to maintain inference to an un-surveyed area, and our assumptions with these methods are described below.Other studies have delineated all known areas of occurrence, assuming complete knowledge of range boundaries and the spatial distribution of individuals and conducted a census that assumes all individuals in the population have been exposed to sampling effort (i.e., there is no un-surveyed area).This approach uses non-random site selection and assumes no site selection bias (Fournier et al. 2019, Mentges et al. 2021).Accounting for individuals present during the survey but not detected by the observer is another important consideration for population estimation (Elphick 2008, Nichols et al. 2009).Some population studies have assumed perfect detection of all individuals present.In contrast, we used double sampling (Bart and Earnst 2002) to adjust our count data for imperfect detection of monogamous species (both under-and over-detection bias).We did not use detection adjustments with phalaropes in our study due to the complexities of determining breeding status for polyandrous species.Had we applied the average detection rate, our population size estimates for phalaropes would have been larger: 23,259 vs. 10,104 breeding individuals for Red Phalarope, and ~1.6 million vs. 948,206 for Rednecked Phalarope using the monogamous and polyandrous approach, respectively.
We encountered multiple logistical constraints and relied on important assumptions of our design-based methods.One of our logistical constraints was lack of access to Alaska Native land, which encompassed 21,490 km 2 (35%) of our entire study.To estimate population totals for the entire YKD (Federal and Alaska Native lands combined), we assumed that average shorebird density was similar in the areas we could and could not survey.Alaska Native lands are spread throughout our study frame (Figure 1) and very little of the area (Federal or Alaska Native) is developed.We chose not to access Alaska Native lands for this study given the large number of jurisdictions (villages) present on the YKD and the urgency in obtaining population estimates given the declines in many shorebird taxa (Andres et al. 2012b).Another constraint was the limited amount of double-sampling data available during our study to estimate detection ratios.To increase our sample sizes, we relied on data from previous studies on the YKD (McCaffery et al. 2012) and the North Slope of Alaska (Brown et al. 2007, Andres et al. 2012a).We consider this an appropriate use of available data because our analyses based on 50 intensive plots found taxon effects (as expected) and year effects, but no evidence that the detection ratio differed between the North Slope and the YKD (Supplementary Material Table S2).In addition, we assumed, for monogamous species, that solitary birds classified as breeding on the plot were mated, which may have resulted in over-estimation in some species.Finally, when creating the sampling frame for our study area, we used a set of assumptions and decision rules to eliminate habitats we considered unsuitable for most shorebird species.These procedures reduced the inclusion probabilities for plots that were suitable for a small number of species found predominantly in forest, tall-shrub, alpine, and gravel barren habitats, resulting in unreliable population estimates for Hudsonian Godwit, Solitary Sandpiper, and Lesser and Greater Yellowlegs.
Our results provided a broad-scale baseline to understand the future impacts of direct and indirect anthropogenic change on the YKD.As climate change accelerates and sea level rise continues, low-lying coastal habitats (i.e., all our habitat strata aside from uplands and mountains) are vulnerable to increased flooding, storm surge, salinization, and sedimentation (Jorgenson and Ely 2001).Beringia, including the YKD, is most vulnerable to loss of suitable climatic conditions for endemic breeding shorebirds within the world's Arctic and subarctic regions (Wauchope et al. 2017).In the YKD, breeding shorebirds occurred in the highest densities in the low-lying Tidal and Coast strata.Black Turnstone, Dunlin, and Red-necked Phalarope all occurred in higher densities in these strata compared to strata at higher elevations.While widely distributed in the region, Western Sandpiper occurred in exceptionally high density in Tidal and Coast strata (~33 to 51 birds km -2 ).The YKD is a critical breeding area for Black Turnstone, Dunlin, and Western Sandpiper because our results indicate that the YKD supports a large portion of the North American populations of these species.Our results also provide important baseline information that could help support resource managers to prioritize areas for conservation considering proposed anthropogenic factors (e.g., mineral extraction, and oil, gas, and wind development), which are on the YKD and threaten shorebird populations in the region (Alaska Shorebird Group 2019).Finally, our population estimates could help co-manage Alaska Native subsistence harvest when used in a harvest-theoretic framework to estimate sustainable mortality limits for shorebird populations (Watts et al. 2015, Naves et al. 2019).
Deltaic systems like the YKD, especially the regions closest to the coast, are some of the most important breeding areas for shorebirds, and other aquatic birds, in the Arctic and subarctic (Bart et al. 2013).These river deltas also provide crucial, initial stopover, and staging sites for post-breeding shorebirds departing their breeding areas (Gill and Handel 1990, Taylor et al. 2010, Ruthrauff et al. 2021).Worldwide, shorebirds rely on estuarine deltas for critical food resources during migration (Butler et al. 2001, Gill et al. 2013).These important ecosystems are vulnerable to a variety of direct and indirect threats, such as land reclamation, contamination from upstream sources, sea level rise, and tidal surges (Jorgenson and Ely 2001, Galbraith et al. 2002, Yang et al. 2011), and many are in danger of collapse in the 21st century (Overeem and Syvitski 2009;MacKinnon et al. 2012).It is essential to continue to generate information on the abundance, distribution, and status of shorebirds to develop effective strategies to mitigate threats to species dependent on these important ecosystems.The lack of dedicated funding for conducting PRISM surveys on the YKD (and Alaska generally) limit surveys to every 15 to 20 years; a greater frequency could increase our ability to assess changes in population size on a fine time scale.More frequent, species-specific surveys would be beneficial when immediate threats are identified on the YKD.Maintaining shorebird populations on the YKD, particularly endemic taxa and populations with most of their breeding range on the YKD, will require a concerted global effort.

FIGURE 1 .
FIGURE 1. YKD, Alaska, study area showing the physiographic strata, the location of rapid survey study plots in 2015 and 2016, and the location of intensive study plots near Boot Lake and Kanaryarmiut that were surveyed in 2016.Coast and North Coast strata were pooled ("Coast") for clarity of presentation.The area surveyed by boat in 2015 is delineated as "Central coast." Areas not surveyed included lands of >50 Alaska Native Tribes and villages and Alaska State land.A map of the Yukon-Kuskokwim Delta, Alaska, showing the locations of surveyed plots.

Figure 2 .
Figure 2.Estimated bird density (and 95% confidence interval) in 8 physiographic strata, showing 3 selected species (with different spatial distributions) and all species combined.Dunlin were largely concentrated near the coast, Western Sandpipers were widespread and found in relatively high density in all strata, and Whimbrel had highest density in the interior.These density estimates, which were adjusted for imperfect detection using a detection ratio, were based on data from shorebird surveys at the YKD, Alaska, 2015 to 2016.Strata are arranged from low elevation (Tidal) to high elevation (Mountains).Stratum code definitions: TD = Tidal; CO = Coast; NC = North Coast; CP = Coastal Plain; LO = Lowlands; RV = Riverine; UP = Uplands; MT = Mountains (see Figure1).Physiographic strata are described byJorgenson and Roth (2010).Note differences in y-axes.Bar charts of bird density in 8 physiographic strata for Dunlin, Western Sandpiper, Whimbrel, and all species in the study combined.

TABLE 1 .
The stratified random sampling frame for shorebird surveys on the YKD, Alaska, in 2015 and 2016 included 8 strata.Strata were based on modification of Estimated density and population size of shorebirds on the YKD, Alaska, 2015 to 2016.Estimates are for USFWS and Alaska Native lands combined, a total of 57,259 km 2 (Table