A habitat-based approach to determining the effects of drought on aridland bird communities - Roberts, Samuel; Thoma, David; Perkins, Dusty; Tymkiw, Elizabeth; Ladin, Zachary; Shriver, Gregory (Dryad, 2021)
Aridland breeding bird communities of the United States are among the most
vulnerable to drought, with many species showing significant population
declines associated with decreasing precipitation and increasing
temperature. Individual breeding bird species have varied responses to
drought which suggests complex responses to changes in water availability.
Here, we evaluated the influence of water deficit, an integrative metric
of drought stress, on breeding bird communities within three distinct
aridland habitat types: riparian, pinyon-juniper, and sagebrush shrubland.
We used 12 years of breeding bird survey data from 11 National Parks and
Monuments in the Northern Colorado Plateau Inventory and Monitoring
Network (NCPN). We used a multivariate community level approach to test
for the effect of annual water deficit on the bird communities in the
three habitats. We found that bird communities responded to water deficit
in all three habitat types, and 70% of the 30 species-habitat combinations
show significant relationships between density and variation in water
deficit. Our analyses revealed that the direction and magnitude of species
responses to water deficit was habitat dependent. The habitat specific
responses we observed suggest that the adaptive capacity of some species
depends on the habitat in which they occur, a pattern only elucidated with
our habitat-based approach. The direction and magnitude of the
relationships between predicted densities and annual water deficit can be
used to predict the relative sensitivity of species within habitat climate
changes. These results provide the first attempt to determine how the
indirect effect of changes in water availability might affect aridland
breeding birds in distinct habitat types. Linking breeding bird density to
annual water deficit may be valuable for predicting changes in species
persistence and distribution in response to climate change.
Study Design and Habitat Descriptions We used 12 years (2005-2015 and
2017) of breeding bird survey data to estimate bird densities in three
habitat types (riparian, pinyon-juniper, and sagebrush shrubland) in 11
National Parks and Monuments of the National Park Service’s Northern
Colorado Plateau Inventory and Monitoring Network (Figure 1). We
identified pinyon-juniper and sagebrush shrubland habitat types on
land-cover maps from the Southwest Regional Re-GAP Analysis Project (Lowry
et al. 2005) or using vegetation-association maps (Daw et al. 2017).
Riparian habitat occurs along perennial streams as narrow strips of
habitat, surrounded by dry uplands, and contains multiple layers of
canopy, which may experience different degrees of dryness depending on
ground water levels and rooting depths. The riparian habitat is dominated
by willow (Salix spp.), Tamarisk (Tamarix ramosissima), and isolated
stands of Fremont cottonwood (Populus fremontii) and boxelder (Acer
negundo). Pinyon-juniper habitat is dominated by two-needle pinyon (Pinus
edulis) and junipers (Juniperus spp.), which varied in relative abundance.
The shrub layer in pinyon-juniper habitat varies throughout the NCPN, but
is often dominated by sagebrush (Artemisia spp.), mountain mahogany
(Cercocarpus spp.), jointfir (Ephedra spp.), and cliffrose (Purshia spp.).
Sagebrush shrubland habitat is dominated by sagebrush, primarily big
sagebrush (Artemisia tridentata) and prairie sagewort (A. frigida), with
rabbitbrush (Chrysothamnus spp.), greasewood (Sarcobatus spp.), and other
shrub species interspersed. Riparian transects ranged in elevation from
1,283 – 1,901 m, pinyon-juniper transects ranged in elevation from 1,393 –
2,402 m, and sagebrush shrubland transects ranged in elevation from 1,666
– 2,447 m. Transects within each habitat type were established across a
large latitudinal gradient to allow for habitat comparisons and to guard
against variation that occurs across broad geographic gradients (Figure
1). Water Deficit Estimates Water deficit is a measure of drought stress
and is the amount of additional water vegetation would use if it was
available (Stephenson 1998). The use of water-year, defined as October 1 –
September 30 of the following year, in the estimation of water deficit is
routine in studies of vegetative responses to climate in the desert
Southwest, as it helps account for precipitation legacies that influence
plants during the growing season (Reichmann et al. 2013, Bunting et al.
2017, Thoma et al. 2019). We used a monthly water balance model to
estimate water-year water deficit using temperature and precipitation at
the center of 45 bird survey transects (see below) following the methods
of Lutz et al. (2010). Precipitation was partitioned into soil moisture,
the quantity of water stored in the top meter of soil at the end of each
month, and runoff which included overland flow plus infiltration below the
rooting zone. Runoff is the proportion of water that is not available for
plant growth. Maximum storage in the top meter of soil was defined by
water holding capacity obtained from soil surveys (Soil Survey Staff
2019). Potential evapotranspiration (PET, mm) was the amount of water that
could be evaporated or transpired with available energy if soil moisture
was unlimited. Actual evapotranspiration (AET, mm) was the estimated
monthly loss of water from soil via evaporation and transpiration, limited
by availability of soil moisture. Water deficit (mm) was calculated as the
difference between PET and AET (Stephenson 1998). We used Daymet daily
temperature and precipitation data at a 1-km grain as the climatic input
to the model (Thornton et al. 2016). These data were co-located with bird
point count transects and thus represent the local elevation and latitude
effects on precipitation and temperature at 1-km resolution. Within each
1-km grid cell of a temperature and precipitation time series the water
balance model calculated heat load due to slope and aspect obtained from
30 m digital elevation model at the center of each transect (U.S.
Geological Survey, 2017). Heat load was calculated as a scaling factor
used to adjust PET up on south aspects or down on north aspects, thus
accounting for slope and aspect interactions (Lutz et al. 2010 after
McCune and Keon 2002).   Breeding Bird Community We sampled the
breeding bird community at 675 unique sampling locations on 45 transects
during the breeding season (May 1st through July 15th), the timing of
which varied by transect location. The elevation and latitudinal position
of each transect was considered during the scheduling of each field season
such that all transects were surveyed within their peak breeding season
and after the passage of migratory birds. Transect orientation in
pinyon-juniper and sagebrush shrubland habitats were randomly determined
following protocol procedures outlined in Daw et al. (2017), the riparian
transects surveyed in this study consisted of narrow strips of habitat,
requiring transect placement to be centered along the length of the
habitat and to follow the natural orientation of the riparian zone. All
riparian habitat was narrow enough that the entire width of the riparian
zone was included in the point radius of each survey (see below) and that
dry upland habitat was also included in surveys of riparian habitat at
some locations. Each transect consisted of 15 point count locations, with
points spaced 250 m apart. Observers visited one transect per day, and
initiated sampling within 30 minutes before sunrise and completed the
survey within five hours after sunrise. Observers walked the transect,
navigated to each point using a GPS unit and completed a five-minute point
count survey. We did not conduct surveys if winds exceeded a four on the
Beaufort Scale (13 – 18 mph) or precipitation was more than a drizzle.
During each point count, observers recoded the species and distance (m) to
all individuals detected. Observers estimated the distance to each bird or
cluster of birds using a laser rangefinder (Simmons LRF 600; Simmons
Outdoor Products, Overland Park, KS, USA). From 2005 – 2013, field teams
surveyed each transect twice annually. Starting in 2014, survey effort was
reduced to one visit per year because a sufficient number of detections to
estimate density had been obtained for most bird species detected and to
reduce costs. The number of completed surveys varied by year, but
>75% of possible surveys were completed annually. For more detailed
information about survey protocols, refer to Daw et al. (2017). We first
estimated transect-level detection-adjusted densities for each bird
species using the Distance package (Miller et al. 2019) in program R
(3.6.1, http://www.r-project.org, accessed July 10, 2019). To avoid
including double-counted individuals in our analyses, we only included
detections within 125 m from the point center in our analyses. The total
area covered by each transect was 73.6 ha (i.e. 15, 125 m radius points
per transect). Next, we used the detection-adjusted density estimates in a
multivariate analysis to determine the bird community response to changes
in water deficit using the manyglm function in the R package mvabund (Wang
et al. 2012). The mvabund package provides a novel set of hypothesis
testing tools that is a flexible and powerful framework for analyzing
multivariate abundance data (Wang et al. 2012). We used the manyglm
function to simultaneously fit general linear models to each species using
water deficit as the common predictor variable (Wang et al. 2012). This
model-based approach alleviates the mean-variance relationship problem
associated with distance-based community level metrics (Warton et al.
2012). To meet the data structure requirements of mvabund, we sought to
maximize the number of transect-level density estimates available using
Distance by constraining our analyses to only include the 10 species
within each habitat type with the greatest number of transect-level
density estimates available (Table 1). Species with the greatest number of
detections allow for the greatest number of transect-level density
estimate, therefore, the species included in each habitat are not the
species with the highest densities, rather the species with the greatest
number of detections. For species-transect-year combinations where
detections were too few to estimate density, we summed the unadjusted
counts and divided by the average detection probability for that
species-year. This was done for 6% of the density estimates to satisfy the
need for a response variable for all species-transect year combinations.
Because there is likely a time lag in the bird community response to water
deficit driven changes in vegetation, we compared two time-lag models (a
one-year lag and a two-year lag effect of annual water deficit) against a
null model within each habitat and selected the model with the lowest
Akaike Information Criterion (AICc, Burnham and Anderson 2002). We did not
explore other lags because vegetation response to legacy effects of
precipitation shortfalls revert to average condition within two years in
semiarid grasslands (Thoma et al. 2016). Once we determined the most
appropriate time lag within each habitat type, we then tested for
community-wide responses to annual water deficit within each habitat using
ANOVA on the ‘manyglm’ object and assessed the resulting likelihood ratio
tests (LRT; Warton 2011) and resampled P-values (Wang et al. 2012). If we
detected a significant community level response to water deficit within a
habitat type, we then used p.uni = “adjusted” argument to determine the
effect size and direction of the response for each species. Finally, we
predicted the relationships between species’ densities and annual water
deficit using the ‘predict’ function on the best model within each habitat
type.
The code and model to produce the water deficit values are not included
here because they are in the process of being published.
Keywords: Ecology