Predicted effects of landscape change, sea level rise, and habitat management on the extirpation risk of the Hawaiian common gallinule (Gallinula galeata sandvicensis) on the island of O‘ahu

Study area

We studied the population of gallinules on O‘ahu, Hawai‘i (21°28′N, 157°59′W), which persists in fragmented subpopulations along the island’s coastal plain (Fig. 1). Wetland habitats supporting breeding subpopulations of gallinules include state and federal wildlife refuges that are actively managed for waterbirds, gardens, agricultural areas, and golf courses. O‘ahu is the most developed island of Hawai‘i with the highest human populations, and it has experienced rapid landscape change (Giambelluca, 1986; Giambelluca, 1996; Klasner & Mikami, 2003), resulting in a highly diverse and increasingly urbanized landscape matrix. Most wetlands support small numbers of gallinules (five–30 individuals), although two larger wetland complexes may support subpopulations of ∼100 (C van Rees, 2018, unpublished data). To our knowledge, the O‘ahu population of gallinules is an effectively closed system, with no published records of movements between O‘ahu and Kaua‘i (Banko, 1987; van Rees et al., 2018a), although biannual waterbird surveys show extremely rare occurrences of gallinules on the Big Island, Maui, and Molokai, so the possibility of movement cannot be excluded. Demographic data used for this study were taken from van Rees et al. (2018a), and were collected from fourteen wetland sites on the island from different time periods from 1979–2017.

Baseline population model

We used Vortex 10 (Lacy, 1993; Lacy & Pollak, 2017) to generate an individual-based stochastic simulation model of O‘ahu’s gallinule population. We chose this approach to explicitly model how the small size of many of our subpopulations affects extirpation risk, and the movement of individuals between subpopulations (Lacy, 2000; Walters, Crowder & Priddy, 2002). We used a time frame of 160 years, the span of ∼40 generations recommended by O’Grady et al. (2008) and Reed & Mccoy (2014), based on our estimated generation time of approximately four years as calculated in Vortex (using the Euler equation, Lacy, Miller & Traylor-Holzer, 2017). We set the duration of each year in days as 365 days. We defined a subpopulation as extirpated (or, in the software, the “extinction definition”) when only one sex remained; the same criterion was used for island-wide extirpation. We did not model population harvest or supplementation. Accordingly, the order of events in each Vortex year (as described in the program) was EV—Breed—Mortality—Age—Disperse—rCalc—Ktruncation—UpdateVars—Census. For each scenario, we ran 1000 iterations, and recorded the following metrics of extirpation risk for each subpopulation and for the total island population: (1) probability of extirpation by 160 years, (2) mean population size of extant populations (and subpopulations) at 160 years, and (3) deterministic and stochastic growth rates for that scenario. We did not include state variables as defined in the Vortex simulation model.

Subpopulations and carrying capacity

We delineated subpopulations based on genetic evidence from van Rees et al. (2018b), combining subpopulations with no evidence of genetic structure (e.g., nonsignificant or near-zero F_ST) into wetland complexes (Fig. 1). The carrying capacity of each subpopulation was estimated as the maximum recorded count of gallinules in the last 12 years, using data from both Hawai‘i’s biannual waterbird survey and playback surveys (following DesRochers, Gee & Reed, 2008) that we conducted for a separate project (C van Rees, 2018, unpublished data). The carrying capacities of wetland complexes were the sum of carrying capacities of their constituent wetlands (Table 1). The starting size (time = 0) of each subpopulation was set to 80% of K, which is similar to current estimated population sizes, and we assigned the model to start populations at a stable age distribution. We chose the default value (0.5) for correlation of environmental variance between subpopulations given an absence of quantitative data on this parameter and a general understanding that while habitats are geographically proximate due to the small size of the island, that precipitation conditions can differ dramatically based on orographic rainfall (Lau & Mink, 2006). We set the correlation between reproduction and survival to 0.1 due to repeated observations of failed nests and depredated broods with no evidence of increased adult mortality (C van Rees, pers. obs., 2018, M Silbernagle, USFWS ret., pers. comm., 2017).

For almost all wetland sites, gallinule counts from unpublished playback survey data yielded the highest numbers, and these values were consequently used. These counts likely have much lower detection error than waterbird survey counts (DesRochers, Gee & Reed, 2008), but our carrying capacities and starting population sizes nonetheless assume high detection and could be underestimates of true population size.

Survival and reproductive system

We generated baseline model parameter values using information on reproduction and survival from van Rees et al. (2018a), with support and supplementation from data on related taxa (e.g., Common moorhen, Gallinula chloropus; Table 2). Even though Hawaiian gallinules do not breed until age 2, there is no evidence of age structure in survival rate beyond year 1 (van Rees et al., 2018a). Accordingly, we modeled survival using two age classes: hatch year (hatching to age 1), and after hatch year (age 1+ or adult). For adult (after hatch year, AHY) birds, we used the less conservative of two survival estimates generated by van Rees et al. (2018a), because it better accounted for the extremely poor detection in this taxon. We calculated first-year survival by combining survival rates to fledging with adult survival rates for the remainder of year 1 (10 months). We used information on fledgling survival from van Rees et al. (2018a), which reported ∼41% mortality within the 60-day period up to fledging; similar to observations in G. g. cachinnans in North America (Greij, 1994). Assuming adult mortality rates for the remaining 10 months of year 1, the total mortality in year 1 was estimated as 67% (Table 2). We estimated mean adult mortality to be 26.8%, using pooled data from van Rees et al. (2018a), and calculated variance in adult survival using the binomial equation (yielding a value of 3.3%), because the standard deviation from van Rees et al. (2018a) appeared to be inflated by detection errors. We chose a value of 15% for the variance in juvenile to reflect the large variability in chick survival observed in the field. The sensitivity of our model to this parameter estimate and those above was evaluated to assess the degree to which our assumptions and parameter uncertainty affect the outcome of the model and its predictions (see Sensitivity Analyses).

We used the long-term monogamy mating system in Vortex, given observations of two closely monitored pairs (van Rees et al., 2018a) and general knowledge of the species (Bannor & Kiviat, 2002). We set the age of first breeding for males and females to 2 years (van Rees et al., 2018a). We used a maximum age of 10, which is the oldest observed Common gallinule (Gallinula galeata) from mark-resighting data in the United States (Clapp, Klimkiewicz & Kennard, 1982); Hawaiian gallinules have not been recorded as living this long (van Rees et al., 2018a). We assumed no reproductive senescence, since no evidence of this has been observed in Hawaiian gallinules (van Rees et al., 2018a). Because we had no information on the potential impacts of inbreeding depression in Hawaiian gallinules, and because it is considered of potentially lower importance for island bird species (Duncan & Blackburn, 2007) we did not include it in our model.

We defined offspring in the model as chicks, and used 8 as the maximum number of progeny per brood, which is the highest observed brood size on O‘ahu (van Rees et al., 2018a). We specified the distribution of brood sizes based on data from 103 broods on O‘ahu collated by van Rees et al. (2018a) (Table 2) and assumed a 1:1 sex ratio at hatch in the absence of evidence to the contrary. We set the maximum number of broods in a given year to 4, basing our estimate on field observations of two closely monitored pairs (van Rees et al., 2018a), as well as accounts of this and conspecifics under natural conditions (Benthum, 1931; Nagata, 1983). We modeled the distribution in number of broods produced per female per year using a binomial distribution with the estimate of nest success for O‘ahu (0.65; van Rees et al., 2018a) as the probability parameter, and the maximum number of broods in a year (4) as the number of trials (Table 2). Because Vortex does not accommodate additional modeling steps for nesting (and thus counts only chicks as offspring) our approach to modeling reproduction accounted for nest failure and hatch rates implicitly. Observed nest success rates were used to generate the number of successful broods per year, and our observed data on brood sizes rather than clutch sizes internalized hatch rates into reproduction in the model. We used sensitivity analyses of the distribution of brood sizes and mean brood size to examine the potential impacts of changing nest success rates as well as hatch rates on the outcomes of our model. Notably, the effects of these parameters could not be distinguished from uncertainty in clutch size, though this is one of our best-studied parameters (van Rees et al., 2018a), and we see little biological reason it should be changing across time or differing between subpopulations.

We did not make young dependent on their parents for multiple years, because Hawaiian gallinules can feed independently within a month of hatching (Chang, 2010). In the absence of evidence of polygyny, and assuming that all males of breeding age have the potential to compete for breeding vacancies, we set the percentage of breeding-age males included in the model’s breeding pool to 100%. All reproduction and survival parameters were identical between subpopulations, because we currently lack site-specific data on Hawaiian gallinule vital rates (van Rees et al., 2018a).

Density dependence

Given observations that Hawaiian gallinules are aggressively territorial (Chang, 2010; C van Rees, pers. obs., 2018), and population models in other rails that took territoriality into account (Wanless, 2002; Hockey, Wanless & Brandis, 2011), we included density dependence in our model. We did this using the density dependence function provided by Vortex, which varies the probability of an individual female breeding in a given year based on the population size with respect to carrying capacity (Lacy & Pollak, 2017): (1)${\displaystyle P\left(N\right)=\left(P\left(0\right)-\left[P\left(0\right)-P\left(K\right){\left(\frac{N}{K}\right)}^B\right]\right)\frac{N}{N+A}.}$ Here P(N) is the proportion of females that breed when the population size is N, P(K) is the proportion that breed when the population size has reached carrying capacity, and P(0) is the proportion of females breeding at low densities. We used a value of 0 for A, the Allee parameter, because we have seen no evidence for behaviors that would lead to Allee effects in this species (e.g., no dependence on group defense from predators; large habitats where finding mates would be difficult), and a value of 20 for the steepness parameter B, which made the function show few density dependent impacts until the population is above 0.8K (80% of carrying capacity). We assumed this high ceiling because the only published research on Hawaiian gallinules indicates that they are not food limited (DesRochers, SR & Reed, 2010). This study focused only on plant-based nutrition, which is the majority of the adult diet. We acknowledge, however, that the degree to which animal prey might possibly be limiting, especially at certain life stages (e.g., nesting or chick-rearing), has never been investigated. Empirical evidence for food limitation and food-driven density dependence in waterbirds appears to be restricted to dabbling ducks in boreal wetlands (e.g., Sjöberg et al., 2000; Elmberg et al., 2005), which are likely ecologically very different than our system, a highly cursorial rallid in tropical wetlands. Accordingly, given the best available knowledge of our system, we assumed that Hawaiian gallinules would likely not experience density dependent negative feedback until territorial behaviors (e.g., chasing and fighting) began negatively affecting survival and reproductive success. We set P(0), the baseline probability of a female breeding in a year at 0.90, and P(K), the probability of a female breeding under maximum density dependence, at 0.33. We estimated P(K) based on our observations that large family groups in densely populated wetlands typically had at most four non-breeding adult helpers, of which 2 were likely females assuming a 1:1 sex ratio, which implies that 1 in 3 adult females per territory would be breeding at high density (C van Rees, pers. obs., 2018).

Dispersal

We estimated dispersal using unidirectional gene flow data from van Rees et al. (2018b) calculated using software MIGRATE (Beerli & Felsenstein, 1999; Beerli & Felsenstein, 2001). We estimated the number of migrant individuals per generation (m) from one wetland to another as m − Mμ, where M is the effective emigration rate between one subpopulation and another, and µis the mutation rate. We used μ = 10⁻⁴, which is a standard mutation rate for microsatellite dinucleotide repeats (Vigouroux et al., 2002; Marriage et al., 2009). We then converted the resulting per-generation estimate to per-year using Hawaiian gallinule generation time (four years) and to a percentage, to determine the percent of individuals from each wetland moving each year to each of the other wetlands. We used gene flow estimates to estimate movement rates instead of measures of genetic differentiation (e.g., F_ST) because of their ability to provide separate values for different directions between the same subpopulation pair, and because fewer assumptions were made in converting those estimates to actual probabilities of movement (Beerli & Felsenstein, 1999; Beerli & Felsenstein, 2001; Whitlock & McCauley, 1999). We input these values as a matrix of pairwise inter-wetland percentage movement rates, included both males and females as dispersing sexes in the absence of evidence for sex-biased dispersal, and set the age range for dispersal as one to five years, based on observations from mark-resight studies (van Rees et al., 2018a). We estimated high percent survival of dispersers (95%) given our current understanding that gallinules likely disperse at night to avoid predation and fly at higher altitudes where they are unlikely to be struck by vehicles (van Rees et al., in press).

Catastrophes

Generally, wetlands on O‘ahu are subject to few catastrophic events (typically, hurricanes), which we view as unlikely to have significant effects on survival or reproduction of Hawaiian gallinules. Observations indicate that adult Hawaiian gallinules are unlikely to be killed by flood events on O‘ahu (M Silbernagle, USFWS ret., pers. comm., 2017), although flooding is more common and possibly a greater threat on Kaua‘i (Htun et al., 2016; K Uyehara, USFWS, pers. comm., 2018) Hawaiian gallinules can re-nest quickly (van Rees et al., 2018a), indicating that nest and chick losses due to storm-related flooding might be quickly compensated by repeated breeding attempts. Additionally, Hawaiian gallinules are aseasonal breeders, apparently breeding year-round (DesRochers et al., 2009), so even widespread losses at a single time point would represent reproductive loss for only a small portion of the total annual breeding window, rather than disrupting a limited breeding season. Three major hurricanes have directly hit the main Hawaiian Islands with enough proximity to affect O‘ahu in the last 68 years (Central Pacific Hurricane Center, 2017). We generated a per-year probability using this value (0.04). Because subpopulations are spread across three different coasts of the island, we estimated that a hurricane could cause total reproductive failure due to flooding to a maximum of 50% of the population. Due to their aseasonal breeding habits, however, such a catastrophe would probably eliminate only about one fourth (25%) of the year’s breeding attempts in affected subpopulations, given that the combined nesting, incubation, and fledging time of Hawaiian gallinules is around 90 days. This would result in reducing the reproductive output of 50% of the island’s subpopulations by 25%, so we estimated that a hurricane would reduce reproduction by 12.5% in the year that it struck. We set catastrophes to reduce survival by 5% to account for the possibility of a small number of individuals being killed by flooding or during dispersal from flooded areas.

Sensitivity analyses

We followed the sensitivity analysis protocols of Mortensen & Reed (2016), who conducted sensitivity testing using three approaches: (1) perturbation analysis, (2) relative sensitivity or elasticity, and (3) the logistic regression approach (Cross & Beissinger, 2001). Their perturbation approach involves systematically changing a single parameter to see how much a parameter value can be changed until a population either declines to extinction or persists throughout the study period (if it is going extinct under baseline conditions). This method explores potentially extreme parameter values in order to ensure proper model behavior (e.g., to ensure that a modeled population can go extinct, or persist, given parameter values of sufficient magnitude) and to identify the threshold value (or set of values) at which model behavior changes from extinction to persistence (or vice versa). This threshold can be used to determine how close current estimated parameter values are to a major change in population behavior, and whether management or future conditions could possibly result in a shift between extinction and persistence We used perturbation analyses of subpopulation carrying capacity and dispersal rates as a test of two competing management strategies, wetland management and connectivity conservation, to assess their relative value as strategies for decreasing island population and subpopulation extirpation risk. In particular, we were interested in whether extirpation risk could be altered dramatically with realistic shifts in these parameters, or if the changes required are beyond the capability of current management. The logistic regression approach uses logistic regression to examine the relationship between the value of a given parameter and the probability of extinction given a large number of samples of possible parameter values and a binary outcome of extinct or not extinct at the end of the study period (McCarthy, Burgman & Ferson, 1995; Cross & Beissinger, 2001). Notably, each of these methods of sensitivity analysis offers a different advantage and yields different information about the model (Mortensen & Reed, 2016; Manlik, Lacy & Sherwin, 2017). The principle advantage of elasticity analysis is its comparability to other studies as a widely-used metric, while perturbation analysis reveals the points at which model behavior changes across parameter values, and provides perspective on whether system shifts (e.g., extinction vs. persistence) are feasible with realistic parameter values. The logistic regression approach allows for the examination of each parameter’s effect on model outcome while taking the effects of all others into account (Cross & Beissinger, 2001; Mortensen & Reed, 2016). We performed sensitivity analyses on all subpopulations separately, and on the entire O‘ahu population combined.

We performed perturbation analyses, following Reed, Murphy & Brussard (1998), on mean and variance of juvenile and adult mortality, percentage of females breeding in a given year, carrying capacity, mean brood size, and dispersal rate (Table 2). We varied both survival parameters from 0-1 in increments of 0.1, and percentage of females breeding from 0-100% in increments of 10%. In order to simulate habitat management, we varied K from 5-100 in increments of 5 for all subpopulations that had a probability of subpopulation extirpation >0 in our baseline model, which included all populations other than the Windward and Kahuku complexes (Fig. 1). We did not increase the carrying capacities of subpopulations with ∼0% extirpation risk under baseline conditions, because no changes in carrying capacity could further reduce extirpation risk in these subpopulations. We varied the distribution of brood sizes (which we used as a proxy for reproductive success) based on a normal distribution with a mean that we varied between 1 and 8 chicks (standard deviation of 1), encompassing the range of observed values on O‘ahu (van Rees et al., 2018a), and changed mean brood size by 0.5 chicks at a time. We changed dispersal rate by applying a multiplier across all inter-wetland movement rates, thus keeping relative dispersal rates the same and increasing overall movement and population connectivity on O‘ahu. We varied this multiplier from 2 to 100, effectively varying annual probability of individual dispersal by two orders of magnitude, but maintaining relative rates. We used this broad range to reflect current uncertainty about the extent to which connectivity might be altered by management, and uncertainty over the true mutation rate of microsatellite markers used to estimate per-generation movement rates. Each scenario was run for 1000 iterations and 160 years.

We conducted elasticity analysis by changing each parameter by ±10% of its mean value, and dividing the difference between the stochastic lambda (λ) of positive and negative scenarios by 0.2 times the stochastic lambda of the baseline scenario, according to the equation (λ₊ − λ₋)∕(0.2∗λ₀) (Cooper, Walters & Priddy, 2002), where λ₊ and λ₋ are the positive and negative scenarios, respectively, and λ₀ is the baseline scenario. We calculated stochastic lambda using the stochastic instantaneous growth rate (r) provided by Vortex for each scenario (using λ = e^r). We assessed the relative sensitivity of our modeled populations to mean and variation of juvenile and adult mortality, the distribution of number of broods per female per year, the average brood size, the percentage of breeding-aged females breeding in a year, the carrying capacity K of each subpopulation, and population connectivity (Table 2). We varied the distribution of the number of broods per year by subtracting 10% of the proportion of brood numbers falling in each category (zero, one, two, three, and four broods per year) and adding it either to the next highest or lowest category, depending on the direction being tested. Accordingly, the distribution was shifted to higher or lower values by 10%. We called this parameter “distribution of no. broods per year”, which represents the probability distribution of the number of broods per year a simulated pair might have (i.e., the probability of having zero, one, two, three, or four successful nests) given that they attempted to breed. In our sensitivity analysis for this parameter, we changed the distribution by adjusting the frequency of higher vs. lower numbers of broods per year. We changed the distribution of brood sizes per year by modeling brood size as a normal distribution, and adding or subtracting 10% from the mean value. We altered dispersal rate by adjusting the overall dispersal multiplier by ±10%. We used Cooper, Walters & Priddy (2002)’s rule of thumb for assessing the relative sensitivity of model outcomes to changes in parameter values, whereby any parameter with a sensitivity value of >1 or <−1 was considered to have a disproportionate effect on population growth rate.

For logistic regression analysis, we used Latin Hypercube sampling in Vortex to randomly generate parameter sets selected from uniform distributions that we determined using observed and feasible values for Hawaiian gallinules. We ran 10 iterations for each parameter set, resulting in 10,000 total simulations for regression analysis. We then performed logistic regression using the ‘car’ package (Fox & Weisberg, 2011) in R 3.2.2 (R Core Team, 2015), treating extirpation (or persistence) at 160 years as the dependent variable. We conducted logistic regression using carrying capacity, mean juvenile and adult mortality, environmental variation in juvenile and adult mortality, percentage of females breeding, and connectivity as explanatory variables. Each parameter was varied according to a uniform distribution, bounded where appropriate based on knowledge of feasible values. We varied mean juvenile mortality rate from 0–1 because of our extremely poor knowledge of the parameter and its large apparent variation in the field, and adult mortality from 0.02–0.65, based on the potential range of annual adult apparent mortality estimates found in van Rees et al. (2018a). We varied the dispersal multiplier from 2 to 80, again reflecting poor knowledge of possible values (Table 2). We combined these parameters as predictor variables in a single generalized linear model with probability of extirpation (PE) as the response variable. We compared the explanatory value of different parameters using their standardized regression coefficients, calculated by dividing the regression coefficient by its standard error.

Climate scenarios

We used readily available spatial data on sea level rise and maps of the location and extent of modeled Hawaiian gallinule habitats to estimate the reduction in habitat area expected from future sea level rise on O‘ahu. Wherever existing habitat features overlapped with future areas of inundation (indicated by the distribution of future inundated areas), we counted that area as habitat loss. We assumed that reductions in habitat would correspond with a proportional reduction in carrying capacity; in other words, that gallinule population densities are uniform throughout their habitats. The spatial data that we used to approximate sea level rise on O‘ahu were created by the National Oceanic and Atmospheric Administration’s (NOAA) Office of Coastal Management, and are available through their sea level rise data portal (https://coast.noaa.gov/slrdata/). These sea level rise projections do not account for the impacts of erosion, island subsidence, wetland migration through accretion, or human management for sea level rise mitigation, but represent a baseline model of potential habitat losses due to inundation. We used these models to generate an estimate of the potential magnitude of reduction of gallinule habitat and carrying capacity on O‘ahu under projected sea level rise scenarios, and to gauge the relative threat of sea level rise compared to other potential factors affecting extirpation risk in this taxon.

We estimated changes in habitat area (carrying capacity) at two scenarios, 0.914 m (3 feet, hereafter ∼1 m) and 1.829 m (six feet, hereafter ∼2 m), given a limited number of scenarios for which data were available on the NOAA sea level rise portal. These values correspond approximately to the range of sea level rise projected by Vermeer & Rahmstorf (2009) for the year 2100, which other researchers have found to be more predictively robust than the IPCC (2007) projections (Rahmstorf, Perrette & Vermeer, 2011; Kane et al., 2015a). Notably, Vermeer and Rahmstorf’s (2009) estimates are only for a 100-year projection, so we chose the highest available sea level rise value from available NOAA maps (1.829 m) for our 160 year time frame, representing a highly conservative estimate. Our smaller value (0.914 m) corresponds to a conservative value for an 80-year projection (Vermeer & Rahmstorf, 2009).

We combined our own maps of gallinule habitats on O‘ahu with data from the U.S. Fish and Wildlife Service’s National Wetlands Inventory (USFWS, 2017) to generate outlines of the habitats supporting all subpopulations modeled in this study. We converted these polygon maps to 3m resolution rasters using the Polygon to Raster tool in ArcMap 10.4.1 (ESRI, 2012), and for each subpopulation, added this map to one of the two sea level rise maps using the raster calculator. We then calculated the area of overlap to estimate the proportion of habitat pixels that would be inundated with salt water. Using this value as an ending carrying capacity at either 80 or 160 years, we designed a power function of the form aX^b to approximate the shape of the sea level rise curve depicted in IPCC (2007) and Vermeer & Rahmstorf (2009). We applied this to the carrying capacity of affected populations in Vortex to approximate the change in carrying capacity across time according to (2): (2)${\displaystyle K_t=K_0-\left[K_0\ast L\left(\frac{a{t}^b}{Y}\right)\right],}$ where K_t is carrying capacity at time t, K ₀ is the original carrying capacity, L is the total proportion of habitat area lost at year 80 or 160 (depending on the simulation), t is the current year, Y is the last year of the simulation (80 or 160), and a and b are shape parameters used to approximate the pattern of projected global sea level rise. We determined values of a and b separately for 80 year and 160 year scenarios to maintain curve shape while passing through a different point at t = 80 or t = 160 (for 80 years, a = 0.06 and b = 1.613; for 160 years, a = 0.02 and b = 1.7385). We ran sea level rise scenarios using baseline model parameters. Carrying capacities of individual subpopulations were altered independently according to separate analyses of their potential area loss. Both scenarios were run with 1,000 iterations, with the ∼1 m scenario projecting for 80 years, and the ∼2 m scenario for 160 years (Table 3).

Model projections

Our baseline model showed that O‘ahu’s island-wide Hawaiian gallinule population has a high likelihood of persistence over the 160 year time frame (Table 4). Probability of extirpation varied strongly between subpopulations, however, with the large subpopulations having probabilities near 0, and the small populations having probabilities near 1.0. Among these small subpopulations, median times to extirpation ranged from four to seven years and mean times from five to 8.9 years. The longest times to extinction among any iterations for small populations were on the order of 50 years. Medium subpopulations with nonzero probabilities of extirpation had mean times to extirpation ranging from 55.5–87.7 years; the median time to extirpation for Waimea valley was 53 years. The extirpation of small subpopulations appears to drive the very slight decreases in mean estimated island-wide population size at year 160, though average increases in population size in medium and large populations compensate for much of this loss. The stochastic population growth rate (r) for the total island population was 0.25, indicating rapid growth, although this ranged from 0.09 to 0.24 in different subpopulations, where larger subpopulations exhibited higher rates.

Model sensitivity

Because subpopulations of similar sizes showed similar trends in stochastic population growth rate and probability of extirpation, we do not show results for all subpopulations, but instead show representative examples of each size class (Figs. 2–5). Perturbation analysis showed that large and medium subpopulations transitioned rapidly from low probabilities of extirpation to high (P =  ∼ 1.0) probabilities when juvenile mortality rose above 80% (Fig. 2). Perturbation of adult mortality showed a transition point for large subpopulations from low to high probability of extirpation at 50% mortality, although this point shifted toward lower values in medium-sized subpopulations (Fig. 3). Small subpopulations’ extirpation probabilities remained at 1.0 for all values of juvenile and adult survival. The range of parameter uncertainty for adult survival encompassed values that were meaningful for medium-sized subpopulations (i.e., there were large differences in probability of extirpation across values ±1 SD from our parameter estimate), but this was not the case for large or small subpopulations.

Large subpopulations and the overall island population transitioned from near 0 extirpation probability to near 1.0 extirpation probability when the percent of females breeding decreased below 30% (Fig. S4). This transition was more gradual for medium subpopulations, which increased steadily from about 70% of females breeding, reaching an extirpation probability near 1.0 at 30%. For mean brood size, large subpopulations and the overall island population had extirpation probabilities near 1.0 when mean brood size was zero, but decreased to an extirpation probability of ∼0.0 when mean brood sizes were one or higher (Fig. S5). Medium-sized subpopulations had a probability of extirpation near 1.0 for mean brood sizes of zero and one, but declined to ∼0.25 with a mean brood size of two, and remained near zero at all higher values. Extirpation probability of small populations remained at 1.0 for all simulated values of percent of females breeding and brood size.

Our investigation of management scenarios in which connectivity was increased showed that with large (e.g., 10 to 25 times) increases in the dispersal multiplier, the extirpation risk of small subpopulations could be reduced by 40–75% (Fig. 4), but that it had little to no effect the extirpation risk of medium and large subpopulations, or the overall island population. The probability of extirpation of the two smallest subpopulations, Klipper and Keawawa, declined rapidly with increases of two to 20 times current dispersal, and showed slower declines after that point. For medium-sized subpopulations, extirpation risk declined sharply from two to 12 times detected dispersal rates, at which point probability of extirpation was near zero. For the habitat management scenario, perturbation of the carrying capacity of small and medium subpopulations showed consistently that carrying capacities above 15 individuals led to rapid decline in extirpation probability, reaching 0 by around ∼30 individuals (Fig. 5).

Conventional, relative sensitivity analysis showed that none of the tested vital rates had a disproportionate (>1.0) effect on stochastic annual population growth rate (λ). Juvenile mortality had the largest effect on the island-wide population persistence (−0.61; Table 5), though its effects were smaller for small subpopulations (e.g., −0.46 for Klipper). Mean brood size (a proxy for reproductive success, including nest failure and hatch rates) and percent females breeding had the next largest effects, with 0.30 and 0.23, respectively, with greater sensitivity to brood size among small subpopulations, and greater sensitivity to percentage of females breeding among large subpopulations. Adult mortality had a moderate effect on population growth rate (−0.17), which was greater for smaller subpopulations (e.g., −0.30 for Keawawa and −0.54 for Lotus Farm). Dispersal rate had a negligible effect on λ for the island-wide population, but small populations showed sensitivities up to 0.14. Environmental variation in adult and juvenile mortality and carrying capacity both also had little effect on λ, with slightly stronger effects (e.g., 0.03–0.06) on small subpopulations. Our logistic regression analysis showed that mean juvenile mortality, mean adult mortality, and mean brood size accounted for the most variability in observed extirpation probability of the overall island population (Table 5). Variance in juvenile mortality accounted for a greater proportion of variability than did variance in adult mortality, and both carrying capacity and dispersal rate explained very little variation. The importance of dispersal and carrying capacity was much larger for small subpopulations, with standardized coefficients as much as seven times larger than those for the island-wide population (e.g., Keawawa, −21.73) for dispersal, and three times higher for carrying capacity. The p values for all covariates in our logistic regression model were statistically significant (p < 0.0001 in all cases).

Sea level rise scenarios

Our Hawaiian gallinule habitat maps encompassed ∼430 ha of wetland habitat on O‘ahu. In the 80 year (∼1 m sea level rise) scenario, a total of 36 ha of wetland habitat (8% of total) are predicted to be lost due to salt water inundation, with a sharp increase in the 160 year (∼2 m sea level rise) scenario, in which 239 ha (56%) are predicted to be lost. Due to their landscape context, four wetlands (Keawawa, Klipper, Lotus Farm, Waimea Valley) were unaffected by projected sea level rise in either scenario, and Turtle Bay was unaffected in the 80-year scenario. In contrast, the Windward, Kahuku, and Pearl Harbor Complexes are predicted to lose 5.6%, 17%, and 11%, respectively, of their total carrying capacity in the 80-year scenario. In the 160 year and ∼2 m sea level rise scenario, the Windward Complex is predicted to lose 54% of its carrying capacity, the Kahuku Complex 51% (Fig. 6), the Pearl Harbor Complex ∼100%, and Turtle Bay 99%.

Simulated island-wide gallinule populations on O‘ahu had 0.0 probability of extirpation under both our 80-year and 160-year sea level rise scenarios, though the mean population size at the end of the simulation was 11% lower than the starting size at 80 years, and 57% lower in 160 years (Table 4). Probabilities of subpopulation extirpations were not changed by sea level rise at 80 years, although the ending population sizes of the three largest subpopulations were smaller than in the baseline scenario. Pearl Harbor and Turtle Bay had population declines in excess of 80% over the 160-year (∼2 m rise) timeframe, while the Windward and Kahuku complexes showed declines of 40% and 37%, respectively. The probabilities of extirpation of the Windward and Kahuku complexes remained approximately 0.0 after 160 years, but they increased dramatically for both Pearl Harbor (from 0.05 to 0.77) and Turtle bay (from 0.08 to 0.80). Small populations were equally likely to be extirpated during either sea level rise scenario as in the baseline scenario. Stochastic population growth rates showed little change between climate change scenarios, both for the overall island population and individual subpopulations.

Sea level rise

We found that a large portion of current gallinule habitat is expected to disappear due to sea-level rise. Based on our population model, however, this reduced carrying capacity is not expected to cause extirpation of the overall island population. There are reasons, however, to believe that our assessment is conservative. The ∼2m sea level rise scenario that we chose for our 160 year projection is probably a substantial underestimate of sea level by 2180 (∼160 years from present). Additionally, we did not take into account habitat degradation through salinization. Coastal freshwater wetlands on O‘ahu could be salinized either by storm surge or by saltwater intrusion in underlying basal aquifers (Kane et al., 2015), reducing their habitat value to gallinules. The impacts of salinization with sea level rise could thus be quite severe and drive reductions in carrying capacity well beyond the spatial extent modeled in our projections. Finally, our assumption that gallinule densities are uniform across habitats can be reasonably challenged. We have found that many of the areas most threatened by sea level rise are also those with the highest population densities of Hawaiian gallinules (van Rees & Reed, 2016–2017, unpublished data). For example, Hamakua Marsh (part of the Windward Complex) and James Campbell National Wildlife Refuge (part of the Kahuku Complex) account for 51% and 76% of the carrying capacities of their respective wetland complexes, and will be the first parts of their respective wetland complexes lost to sea level rise. The effects of climate change on gallinule carrying capacity, and therefore on extirpation risk, therefore, is likely to be greater than that estimated in this paper.

Wetland migration due to soil accretion is a potential mitigating effect against habitat loss for some coastal wetland types, in which wetland distributions shift inland and upland with rising sea levels (e.g., Traill et al., 2011). This phenomenon may have little benefit on O‘ahu, where most areas farther inland from wetlands are either densely developed or feature dramatically sloped topography that transitions directly from the coastal plain into mountains, leaving little room for wetland migration. Additionally, soil accretion takes place after inundation, and is thus of value primarily for salt water wetland habitats (Kirwan et al., 2016). The persistence of the O‘ahu’s two major strongholds for Hawaiian gallinules, the Windward and Kahuku complexes, may depend strongly on the migration of wetland habitats to higher elevations and increased management of wetland habitats within those regions that are less vulnerable to rising sea levels. Both Kawainui marsh in the Windward Complex and the Shrimp Farms in the Kahuku Complex are positioned inland and adjacent to high-density managed sites from which gallinules might easily emigrate under sea level rise. Management (which is currently limited at both sites) to increase carrying capacity at these higher-elevation sites could create habitat capable of supporting a large portion of the gallinules currently found in their respective wetland complexes. Similar inland wetland alternatives are not available for the Pearl Harbor Complex or Turtle Bay, meaning that more intensive and perhaps economically infeasible measures like habitat creation and land acquisition might be necessary to mitigate sea level rise impacts there.

The high population growth capacity of Hawaiian gallinules implies that management strategies that increase carrying capacity of existing wetlands that may be less threatened by climate change (on O‘ahu or other islands) may compensate for anticipated losses at sites with greater risk. Such restoration, (e.g., Kawainui marsh, US Army Corp of Engineers (USACE), 2008), and potential creation of new habitats would require a solid understanding of the habitat requirements of Hawaiian gallinules. Unfortunately, no quantitative data on habitat correlates of gallinule abundance or breeding success have been published, so a better understanding of determinants of habitat quality may be a research priority for this subspecies.

Sensitivity analysis and management strategies

Both conventional and regression sensitivity analyses highlighted the importance of juvenile survival for population persistence in Hawaiian gallinules. This finding is especially relevant to management, given the general belief that exotic invasive predators have larger impacts on juvenile than adult mortality, and that predator management is one of the longest-implemented and logistically feasible management strategies for this subspecies (VanderWerf, 2012). Reproductive output and mean adult mortality were also shown to be important, and they are affected by the same management strategies. Consequently, the prevailing emphasis on predator control in current management plans for this subspecies should be considered highly justified in light of our population viability analysis. The high importance and uncertainty of our estimates of juvenile mortality make it a top priority for future field research on Hawaiian gallinules. Related to this, however, are unanswered questions pertaining to predator control to protect Hawaiian gallinules and other endemic waterbirds. Notably, which predators are having the largest impacts and at what life stage, how these effects differ across wetlands, islands, and years, and which management strategies are most effective at reducing the impacts of these predators. The management of endangered and introduced invasive species is constrained by finite budgets, so more applied research will be necessary to clarify the exact management implications of our results. It is additionally unclear if introduced predators are the main driver of juvenile mortality. Although this notion is supported by numerous anecdotal observations (C van Rees, pers. obs., 2018, A Nadig, USFWS, pers. comm., 2016) and prevailing beliefs (USFWS, 2011; VanderWerf, 2012) on O‘ahu, it is unsupported by empirical evidence. Accounts of emaciated chicks on Kaua‘i (K Uyehara, pers. obs., 2018) might imply that access to food also influences chick and juvenile survival rates. Understanding the mechanisms driving survival in these early life stages is essential to an effective management solution to low survival rates.

Our sensitivity analysis of mean brood size and the distribution of number of broods per year showed that the processes intrinsic to those parameters (nest success and hatch rates) have an important but secondary impact on the viability of Hawaiian gallinule populations on O’ahu. The importance of nest success and hatch rate for population viability is sensible, being that they represent a key stage in the gallinule’s reproductive cycle, but other idiosyncrasies of the gallinule’s reproductive biology may be reducing its importance. The readiness with which Hawaiian gallinules (and congeners) may limit the impact of individual nest failures (van Rees et al., 2018a), and a year-round nesting season indicates that multiple broods are easily possible in any given year. Accordingly, birds may continue attempting to breed until a nest is successful (as has been observed in one closely monitored pair, C van Rees, pers. obs., 2018), at which point juvenile survival becomes the final determinant of reproductive success.

The lower importance of dispersal rates on population growth rate and extirpation risk relative to juvenile and adult survival to population persistence in our study corresponds with findings in other systems and reinforces the notion that connectivity typically contributes less to population viability than demographic rates (South, 1999). In both regression-based and conventional sensitivity analyses, however, dispersal had some impact on the extirpation risk of small subpopulations, amounting to as much as half of the sensitivity value of directly altering carrying capacity for those populations. The extirpation of these subpopulations under low-dispersal scenarios caused a slow decline in overall island population throughout the study period, which was ameliorated by increasing the dispersal multiplier. Interestingly, dispersal multiplier values >10 caused smaller increases in mean island population size at 160 years, indicating that, at a certain point, gallinule emigration to isolated populations was detrimental to the overall population. This is similar to observations by Medici & Desbiez (2012), who examined metapopulation persistence in lowland tapirs (Tapirus terrestris). Our perturbation analysis also showed steady decreases in extirpation probability across a wide range of connectivity values. For small and especially medium populations, even modest (2–5×) increases over current dispersal rates can produce nontrivial changes in the probability of persistence of small Hawaiian gallinule subpopulations. These increases are within current variation in dispersal rates between populations, implying that they are within the realm of possibility for management. It is important to recognize that, according to our models, management actions that directly increase the carrying capacity of small subpopulations would be generally more effective at reducing extirpation risk. While creating new habitat for Hawaiian gallinules on the O’ahu might be prohibitive in social and economic cost due to constraints on land ownership and availability, the restoration of existing habitats would be a highly effective approach to increasing carrying capacities that would avoid such the complications of land acquisition. Increasing connectivity through dispersal might be a viable additional option for decreasing extirpation risk under future sea level rise conditions, where large subpopulations will be reduced to small and medium-sized ones, and more direct and effective methods of decreasing extirpation risk will become more challenging. Recent work by van Rees et al. suggests that stream networks and green water management infrastructure (e.g., drainage swales) proposed as solutions for climate change adaptation on O‘ahu may increase landscape permeability to gallinules, providing a feasible measure for managing connectivity for this subspecies.