DOI: 10.5061/dryad.cc2fqz678
title: Earthworm abundance and availability does not influence the reproductive decisions of black-tailed godwits in an agricultural grassland

description: Maintaining the biodiversity of agricultural ecosystems has become a
 global imperative. Across Europe, species that occupy agricultural
 grasslands, such as Black-tailed Godwits (Limosa limosa limosa), have
 undergone steep population declines. In this context, there is a
 significant need to both determine the root causes of these declines and
 identify actions that will promote biodiversity while supporting the
 livelihoods of farmers. Food availability, and specifically earthworm
 abundance (Lumbricidae), during the pre-breeding period has often been
 suggested as a potential driver of godwit population declines. Previous
 studies have recommended increasing the application of nitrogen to
 agricultural grasslands to enhance earthworm populations and aid
 agricultural production. Here we test whether food availability during the
 pre-breeding period affects when and where godwits breed. Using
 large-scale surveys of food availability, a long-term mark-recapture
 study, focal observations of foraging female godwits, and tracking devices
 that monitored godwit movements, we found little evidence of a
 relationship between earthworm abundance and the timing of godwit
 reproductive efforts or the density of breeding godwits. Furthermore, we
 found that the soils of intensively managed agricultural grasslands may
 frequently be too dry for godwits to forage for those earthworms that are
 present. The increased application of nitrogen to agricultural grasslands
 will therefore likely have no positive effects on godwit populations.
 Instead, management efforts should focus on increasing the botanical
 diversity of agricultural grasslands, facilitating conditions that prevent
 hardening soils, and reducing the populations of generalist predators.
description: <b>Study
 Area</b><br> We studied the relationship between the
 reproductive decisions of godwits and earthworm abundance in the Haanmeer
 (52.9222°, 5.4353°), a 200-ha polder in southwest Friesland, The
 Netherlands. The Haanmeer polder is one of the last remaining areas with
 high densities of breeding godwits (~0.5 pairs/ha) in southwest Friesland
 — which represents the core of the godwit breeding distribution (Kentie et
 al., 2016) — and is part of a long-term landscape-scale study of godwit
 demography and breeding biology (see Groen et al., 2012; Senner et al.,
 2015 for more details). The Haanmeer consists of fields that are
 maintained under two different management schemes: (1) high-intensity
 grasslands (70 ha) and (2) low-intensity grasslands (130 ha). The two
 management types differ in a number of ways, including water levels,
 botanical richness, the number of cattle grazed, amount and type of
 nitrogen application, and mowing regime (Table 1). In particular,
 the application of cow manure — the main source of nitrogen input —
 averages ~ 70 tons/ha on the high-intensity grasslands, but only ranges
 from 0 to 25 tons/ha on the low-intensity grasslands (A. Stokman, S.
 Venema, and D. Postma pers. comm.). In addition, the water level in the
 ditches surrounding the fields is kept consistently lower in the
 high-intensity grasslands (unpubl. data Wetterskip Fryslân) resulting in a
 ~70 cm difference in spring groundwater level (unpubl. data MAV &
 NRS). These two management types and the differences between them are
 typical for our study area (see Groen et al., 2012) and the godwit
 breeding range in the Netherlands (Teunissen et al., 2012).<br>
 <br> <b>Data collection on godwit behaviour</b><i>
 </i><br> During 2008–2012, we captured 87 adult and 392
 juvenile godwits in the Haanmeer and surrounding polders. Upon capture,
 adults and chicks were individually marked, bled for molecular sexing, and
 their biometrics were measured (see Loonstra et al., 2019). The latter
 included the measurement of their exposed culmen to the nearest millimetre
 (hereafter ‘bill length’). The individual markings consisted of a unique
 colour-ring combination for adults and appropriately sized chicks, and a
 single flag with a unique alphanumeric code for smaller chicks.
             In 2013 and
 2014, we searched daily for these individually-marked godwits as they
 returned to the Haanmeer. We started with the first arrival of adults and
 continued until the first egg was found (2013: 8 Mar–10 Apr; 2014: 8 Mar–9
 Apr). When an individually-marked godwit was encountered, its colour code,
 location, and behaviour — foraging, preening, displaying, etc. — were
 recorded. During the same period, we made daily foraging observations on a
 subset of individually-marked females. For these females, we recorded the
 number of probes and successes during consecutive 3-min periods, where a
 probe was defined as a downward movement of the bill into the soil and a
 success as a swallowing motion (Senner and Coddington 2011). When
 possible, we also recorded the prey type consumed. In 2013, we observed a
 foraging individual for up to 10 3-min periods for a total daily foraging
 observation length of 30 min. In 2014, we did this for up to three, 3-min
 periods for a total daily foraging observation length of 9 min.
             In 2014, we
 also tracked the locations of three females throughout their pre-breeding
 period using 7.5 g solar-powered UvA-Bits GPS trackers (Bouten et al.,
 2013). We attached these trackers in 2013 with a leg-loop harness using 2
 mm nylon rope (see Senner, Stager, Verhoeven et al., 2018 for more
 details). These trackers stored location estimates once every 5 min when
 the battery was fully charged and once every 15–30 min in all other
 instances. One of these females was still carrying its transmitter in 2015
 and this provided us with an additional year of spatial distribution data
 during the pre-breeding period.             After the first egg
 was found, we stopped our foraging observations and started intensively
 searching for godwit nests. We floated and measured the dimensions of eggs
 from nests with completed clutches to estimate their incubation stage and
 egg volume (Liebezeit et al., 2007, Schroeder et al., 2009). We monitored
 nests every 2–3 days to determine their precise hatching date. We made a
 particular effort to associate individually-marked birds with nests. For
 nests found in the laying phase, we determined their lay date based on the
 assumption that godwits lay one egg per day (Cramp & Simmons
 1983). For nests that hatched but were not found in the laying phase, we
 back-calculated the lay date by subtracting 26 days (the combined average
 laying and incubation period; Verhoeven et al., 2020). In cases where we
 found a nest with a complete clutch, but did not know the hatching date
 due to predation or abandonment, we estimated the lay date by subtracting
 the incubation stage (as derived from the flotation method) from the date
 the nest was found.<br> <br> <b>Data collection on food
 availability</b><br> On four separate occasions, we sampled
 the entire Haanmeer for below-ground invertebrates along a 100 × 150 m
 grid (Fig. 1) as local soil invertebrate populations can vary considerably
 (Timmerman et al., 2006). In 2013, we collected 130 samples from 1–3 March
 — before the arrival of the first godwits — and 136 samples from 4–6 May,
 just after the last of our individually-marked females had begun
 incubation. We sampled six additional fields in May that had not been
 sampled in March. In 2014, we collected 135 samples from 8–10 March and
 136 samples from 1–4 May. We inadvertently skipped one sample at the end
 of a transect in March 2014.             We collected
 invertebrates in 20 × 20 × 20 cm soil samples. To estimate the depth at
 which invertebrates were present, we split each sample into four pieces
 along its horizontal axis immediately after removal from the ground,
 resulting in four equally sized slices 5 cm in height. Each slice was
 sealed separately in a plastic bag and stored between −5°C and +5°C. We
 processed the samples within two weeks by hand-sorting each slice, which
 is considered the most reliable method (Nordström and Rundgren 1972;
 Edwards and Lofty 1977). We then cleaned all invertebrates, weighed them
 to the nearest 0.01 g, and measured their (relaxed) length to the nearest
 0.01 cm.             Soil penetrability
 affects earthworm availability for godwits by influencing the ease with
 which godwits can probe the soil for invertebrates (Kleijn et al., 2011).
 We thus also measured the penetration resistance of the soil in the same
 fields in which we sampled soil invertebrates with a penetrometer with a
 1-cm<sup>2</sup> cone (Eijkelkamp, Giesbeek, penetrometer
 06.01.14). From late March to mid-May in 2013 we measured the penetration
 resistance along four 60-m long transects daily. Two transects each were
 located in high-intensity and low-intensity grasslands. We sampled every 4
 m along each transect, resulting in 15 measurements per transect. To
 ensure that the results from our original two transects adequately
 captured variation in penetrability across the study area, we added eight
 additional transects in low-intensity grasslands in 2014 for a total of
 ten transects in low-intensity grasslands and two in high-intensity
 grasslands. For these transects, we started in the middle of March and
 sampled every 2 m along 10-m long transects (six measurements per
 transect).<br> <br> <b>Statistical analysis of annual
 and seasonal differences in invertebrate biomass</b><br> The
 recorded invertebrate biomass was zero-inflated due to the number of
 samples without earthworms (Table 2) and included only positive values in
 those cases when invertebrates were recorded. As a result, we used a
 two-part binomial-gamma hurdle model. The initial binomial model predicted
 the probability of recording earthworms in a sample, using a generalized
 linear mixed model (GLMM) with binomial error structure and a logit link
 function. The gamma model then predicted the invertebrate biomass when
 invertebrates were recorded using a GLMM with a gamma error structure and
 log link function. Both models had year, month and grassland management
 intensity as categorical predictor variables, and sampling location as a
 random intercept to account for pseudo‐replication. We obtained
 chi-squared values for the significance of the predictor variables from
 likelihood ratio tests of nested models with and without the variable of
 interest. In the final model, we removed grassland management intensity as
 a predictor variable because it did not significantly improve the fit in
 either the binomial or gamma model (see Results). By multiplying the
 predictions of the two parts of the hurdle model (i.e., weighting
 predicted biomass by the probability of observing an invertebrate), we
 obtained predictions of invertebrate biomass at both sampling time points
 in each year.<br>             To additionally test for potentially
 confounding variables that might influence underlying worm abundances and,
 therefore, godwit breeding biology, we also ran the same binomial and
 gamma models with soil type — (1) more clay, (2) more sand, or (3) more
 peat (Fig. 1) — replacing grassland management intensity as a categorical
 predictor. We excluded three earthworm samples for this analysis, either
 because the sampling area had previously been filled with an unknown soil
 type (<i>n</i> = 2) or because the sample was the sole
 representative of a soil type category (more loam; <i>n</i> =
 1; Fig. 1). Soil type, however, was not a relevant predictor variable of
 worm abundance in either the binomial (χ<sup>2</sup> = 3.07,
 df = 2, <i>p</i> = 0.215) or gamma models
 (χ<sup>2</sup> = 1.76, df = 2, <i>p</i> = 0.414),
 and thus was not included in our final models (Table S1). <br> <b>Relating
 godwit behaviour to earthworm biomass</b><br> Because
 grassland management intensity and soil type were not related to earthworm
 biomass (see Results), and because prevailing management directives
 (Brandsma 1999; van der Weijden and Guldemond 2006; Kleijn et al., 2009a)
 posit a direct link between earthworm abundance and godwit behaviour and
 abundance — proposing that earthworms explain observed variation in godwit
 reproductive phenology and habitat selection — we focused our analyses on
 testing the relationship between earthworm biomass and godwit behaviour
 and abundance. To test for an
 effect of earthworm abundance on godwit foraging rates, we used only the
 biomass of earthworms in our analyses because earthworms comprised ~90% of
 prey mass found in March (see Results). Furthermore, we used only the
 earthworm abundance in the top 10 cm because those earthworms found deeper
 in the soil are not within reach of female godwits (see Results). We
 plotted the prey intake rate of female godwits during the pre-breeding
 periods of 2013 and 2014 against the average earthworm biomass in the top
 10 cm in March of the meadow in which each godwit was foraging (Fig. 2).
 We used only females that probed more than 220 times
 (25<sup>th</sup> percentile) and less than 370 times
 (75<sup>th</sup> percentile) to control for differences in
 motivation between individuals, as this is known to cause differences in
 observed intake rates (Duijns et al., 2015). We then used Holling’s type
 II functional response equation to explore the relationship between intake
 rate and earthworm density (Holling 1959). We had no empirical data for
 handling time (Th) or instantaneous area of discovery (a), and used the
 least squares method to calculate the values for both these parameters
 that yielded the best fit to our data.             We next tested for an
 effect of earthworm abundance on where godwits foraged during the
 pre-breeding period. In a generalized linear mixed model with a Poisson
 error structure and log link function, we related the number of unique
 individuals seen on a meadow during the pre-breeding period to the average
 earthworm biomass in the top 10 cm of that meadow in March, while
 accounting for the size of the meadow (continuous covariate), the year
 (two-level factor), and the fact that samples from the same meadow are not
 independent (random intercept). We ran the same model to relate the number
 of nests found in a meadow to the average earthworm biomass in the top 10
 cm of that meadow in March. For the latter analysis, we included only
 nests laid before 1 May in order to exclude replacement clutches as much
 as possible (see Verhoeven et al., 2020).             Third, we related
 earthworm abundance to godwit breeding phenology. We used linear mixed
 models to relate the (1) arrival dates, (2) lay dates, and (3) duration of
 the arrival–laying interval of females to the average earthworm biomass in
 March in the top 10 cm within 200 m of their nests. We did the same,
 separately, for the biomass within 300 m of their nests. We included year
 as a two-level factor and individual as a random intercept. We used the
 earthworm biomass from within these distances because 70% of the foraging
 locations of the GPS-tagged females were within 200 m and 88% within 300 m
 (Fig. S1a), which is similar to the spatial distribution exhibited by the
 rest of the population (Fig. S1b). In these analyses we again included
 only nests laid before 1 May and used chi-squared values to assess the
 significance of each predictor variable. Finally, we
 assessed whether earthworms influenced godwit reproductive investment. For
 this, we used linear models to relate the lay dates and average egg volume
 of all nests (including those of unknown females) laid before 1 May to the
 average earthworm biomass in the top 10 cm within 200 m and 300 m of these
 nests in March. We calculated egg volume using the formula: length ×
 breadth<sup>2</sup> × 0.52 (Schroeder et al., 2009). Female
 size and other covariates were not included in these models, as we have
 previously shown that they have little effect on egg volume (Verhoeven et
 al., 2019). Additionally, we
 tested the robustness of our worm sampling scheme and our choice to
 average worm measures across samples within 200 m and 300 m of a female’s
 nest or within the meadow in which they nested (i.e., within 500 m or
 less). Using both Moran’s I and Mantel tests, we found evidence of
 positive spatial autocorrelation or ‘clustering’ of earthworm biomass at
 these distances in March 2013 (Moran’s I for 200 m = 0.27,
 <i>p</i> = &lt;0.001; 300 m = 0.15, <i>p</i> =
 &lt;0.001; 500 m = 0.1 <i>p</i> = &lt;0.001) and March
 2014 (Moran’s I for 200 m = 0.20, <i>p</i> = 0.004; 300 m =
 0.09, <i>p</i> = 0.019; 500 m = 0.05. <i>p</i> =
 0.021). Our Mantel tests also indicated that sampling locations nearer to
 each other had more similar worm abundances, i.e., positive spatial
 autocorrelation (r<sub>March2013</sub> = 0.13,
 <i>p</i> = &lt;0.001; r<sub>March2014</sub> =
 0.11, <i>p</i> = &lt;0.001). We obtained
 F-values for the significance of all predictor variables from F-tests of
 nested models with and without the variable of interest. We ran all models
 in the R programming environment (R Core Team 2018, version 3.5.1) using
 the package “lme4” (Bates et al., 2015).
description: All files are numbered, a
 new number for each analysis/topic. Related files are given the same
 number but a different letter, i.e. files 3a, 3b, 3c, 3d are related to
 each other.  <b>1a.March2013_soil
 invertebrate sampling data</b><br> <i>This file contains
 the results from our invertebrate sampling in March 2013. It has a tab
 containing all the worms and a tab containing all the non-worms. The
 section refers to the depth of the sample; 1=0-5cm, 2=5-10cm, 3=10-15cm,
 and 4=14-20cm.</i><br> <br> <b>1b.May2013_soil
 invertebrate sampling data</b><br> <i>This file contains
 the results from our invertebrate sampling in May 2013. It has a tab
 containing all the worms and a tab containing all the non-worms. The
 section refers to the depth of the sample; 1=0-5cm, 2=5-10cm, 3=10-15cm,
 and 4=14-20cm.</i><br> <br> <b>1c.March2014_soil
 invertebrate sampling data</b><br> <i>This file contains
 the results from our invertebrate sampling in March 2014. It has a tab
 containing all the worms and a tab containing all the non-worms. The
 section refers to the depth of the sample; 1=0-5cm, 2=5-10cm, 3=10-15cm,
 and 4=14-20cm.</i><br> <br> <b>1d.May2014_soil
 invertebrate sampling data</b><br> <i>This file contains
 the results from our invertebrate sampling in May 2014. It has a tab
 containing all the worms and a tab containing all the non-worms. The
 sections refer to the depth of the sample; 1=0-5cm, 2=5-10cm, 3=10-15cm,
 and 4=14-20cm.</i><br> <br> <b>2.Bill length of
 Female godwits caught in Haanmeer</b><br> <i>This file
 contains the measured bill lengths of females in the
 Haanmeer</i><br> <br> <b>3a.2013+2014_Worms_top
 10cm_with management intensity</b><br> <i>This file
 contains the average worm biomass of worms in the top 10cm in March and
 May and has a column that identifies the management
 intensity.</i><br> <br> <b>3b2013+2014_.Worms_top
 10cm_with soil type</b><br> <i>This file contains the
 average worm biomass of worms in the top 10cm in March and May and has a
 column that identifies the soil type.</i><br> <br>
 <b>4a.2013 soil penetration data</b><br> <i>This
 file contains the soil penetration resistance data collected in
 2013.</i><br> <br> <b>4b.2014 soil penetration
 data</b><br> <i>This file contains the soil penetration
 resistance data collected in 2014.</i><br> <br>
 <b>5.2013+2014_foraging observation data</b><br>
 <i>This file contains the foraging observations made in 2013 and
 2014 and the average number of worms in the field of foraging in order to
 fit intake rates.</i><br> <br>
 <b>6.2013+2014_sightings in field and worm
 biomass</b><br> <i>This file contains the number of
 sightings in a field in 2013 and 2014, the size of the field and the
 field’s worm biomass in order to look for a
 relationship.</i><br> <br> <b>7.2013+2014_nesting
 density and worm biomass</b><br> <i>This file contains
 the number of nests in a field in 2013 and 2014, the size of the field and
 the field’s worm biomass in order to look for a
 relationship.</i><br> <br>
 <b>8.2013+2014_Arrival-Laying Date of marked females and worm
 biomass</b><br> <i>This file contains the arrival date,
 laying date and the arrival-laying interval of females as well as the worm
 biomass within 200 and 300 meters of the nest.</i><br>
 <br> <b>9a.2013+2014_all laydates and worm biomass within
 200m</b><br> <i>This file contains all of the known lay
 dates of godwits in 2013 and 2014 and the worm biomass within 200 meters
 of those nests in order to look for a relationship.</i><br>
 <br> <b>9b.2013+2014_all laydates and worm biomass within
 300m</b><br> <i>This file contains all of the known lay
 dates of godwits in 2013 and 2014 and the worm biomass within 300 meters
 of those nests in order to look for a relationship</i><br>
 <br> <b>10a.2013+2014_all egg volumes and worm biomass within
 200m</b><br> <i>This file contains all the average egg
 volumes of nests from 2013 and 2014 in which all 4 eggs were measured and
 includes the worm biomass within 200 meters of those nests in order to
 look for a relationship.</i><br> <br>
 <b>10b.2013+2014_all egg volumes and worm biomass within
 300m</b><br> <i>This file contains all the average egg
 volumes of nests from 2013 and 2014 in which all 4 eggs were measured and
 includes the worm biomass within 300 meters of those nests in order to
 look for a relationship.</i><br> <br>
 <b>11a.2013+2014+2015_Distance of tracking location from nesting
 location for FigS1</b><br> This file contains the distance
 from the nest of three females that were tracked via GPS-transmitters in
 2014 and 2015.<br> <br> <b>11b.2013+2014_Locations of
 observation and field centers for FigS1</b><br> This file
 contains the distance from the nest to the center of the field of
 observation of 42 individually-marked females that were resighted during
 the 2014 pre-breeding period.<br> <br>
 <b>12a.2013_March_sampling location and worm
 biomass</b><br> This file contains the worm biomass in the top
 10cm in March 2013 for each sampling point and the location (latitude +
 longitude) of that sampling location, used to look for the presence of
 spatial autocorrelation of worm biomass.<br> <br>
 <b>12b.2014_March_sampling location and worm
 biomass</b><br> This file contains the worm biomass in the top
 10cm in March for each sampling point and the location (latitude +
 longitude) of that sampling location, used to look for the presence of
 spatial autocorrelation of worm biomass.