A new approach to estimate fecundity rate from inter‐birth intervals

Funded by Department of Energy and Climate Change (UK), BES, ASAB, Greenpeace, Environmental Trust, Scottish Natural Heritage, Scottish Government, Whale and Dolphin Conservation, Talisman Energy (UK) Ltd., DECC, Chevron, Natural Environment Research Council Acknowledgments Funding for this work was provided by the Department of Energy and Climate Change (UK). Photo-identification data were collected during a series of grants and contracts from the BES, ASAB, Greenpeace Environmental Trust, Scottish Natural Heritage, Scottish Government, Whale and Dolphin Conservation, Talisman Energy (UK) Ltd., DECC, Chevron, and the Natural Environment Research Council. All survey work was carried out under Scottish Natural Heritage Animal Scientific Licences. The authors have no conflict of interest to declare. We thank Mark Bravington for his helpful advice at the early stages of this work and two anonymous reviewers for their useful comments on the manuscript.


INTRODUCTION
Information on life-history parameters is essential to understand a population's dynamics in an ecological context (Baker et al. 2010). Such understanding enables population status and viability to be evaluated (Beissinger and McCullough 2002), thus providing the scientific basis to inform conservation and management actions (e.g., Runge et al. 2004, Currey et al. 2011. Noninvasive photographic identification methods have expanded the opportunities to obtain individual longitudinal data from wild populations with variable natural markings, allowing the estimation of vital rates in both terrestrial and marine species (Moss 2001, Holmberg et al. 2008, Cordes and Thompson 2013. In cetaceans, estimates of survival rate exist for several species (e.g., Smith et al. 2013, Ramp et al. 2014, Carroll et al. 2016), but information on reproductive parameters is available only for a relatively small number of well-studied populations. Fecundity is commonly estimated as the annual proportion of adult females giving birth (e.g., Wells and Scott 1990, Moss 2001, Clarke et al. 2003 or reported as the mean observed birth interval (e.g., Koelsch 2001, Henderson et al. 2014). However, these simplistic estimates of reproductive parameters tend to be biased. Estimates of birth intervals tend to be positively biased because births are missed when animals are not seen every year but negatively biased when the length of the study does not allow the longest birth intervals to be observed (Barlow and Clapham 1997).
Attempts to model birth intervals to produce unbiased estimates have been made for right whales (Payne et al. 1990, Cooke et al. 2001) and humpback whales (Barlow 1990, Barlow and Clapham 1997). These approaches utilized maximum-likelihood methods to estimate the conditional probability of giving birth at time t after a prior birth, by finding the parameter values that maximize the observed birth intervals. The approaches make different assumptions about the probability of observing births, which determine the type of data each approach can use. Barlow and Clapham (1997) put a limit of five years on the observed birth intervals for computational simplicity and because the probability of a longer birth interval for humpback whales was very small. While the above-cited approaches are appropriate ways to estimate birth intervals, the modeling frameworks do not account for individual or temporal heterogeneity in birth probabilities.
Here, we present a new approach to estimate fecundity rate also based on conditional probability of birth, from which an expected inter-birth interval (IBI) can be estimated. We model conditional birth probabilities using a generalized linear mixed-model (GLMM) framework, which can combine both categorical and binomial data and include both fixed and random effects to account for temporal and individual variability. This flexibility is advantageous when analyzing data with repeated measures of individuals in multiple years (e.g., Pomeroy et al. 1999, Ward et al. 2009). If the effects of individual and temporal variation are not of particular interest, they can be incorporated into the modeling framework as random effects, thus reducing the number of parameters to be estimated (Bowen et al. 2006, Brough et al. 2016. We illustrate the approach using individual reproductive histories from the Scottish east coast bottlenose dolphin (Tursiops truncatus) population (Wilson et al. 1999, Cheney et al. 2013 collected over a period of 24 yr. Simulations are used to validate the analytical method and investigate the most common data-driven biases related to the estimation of IBIs: the inclusion of incomplete birth intervals in the data, the effect of different probabilities of capture of a studied population, and the effect of different study lengths.

Estimation of expected inter-birth interval and fecundity rate
We define conditional probability of birth as the probability that a female will give birth t years after a previous birth, under the condition that the female survives and has not calved since the previous birth. Generalized linear mixed models were used to model conditional probability of birth as a function of the number of years since a female previously gave birth, regardless of the fate of the animal most recently born to that female. The full model included, as fixed effects, the number of years since the previous birth (YSPB), its quadratic form (YSPB 2 , to account for possible non-linearity in the relationship), and the number of calves previously born to each female (to allow experience to influence birth probability). Female identity and year were included as random effects to allow for individual and temporal variation while avoiding over-parameterization of the model.
The models were fitted assuming a binomial error distribution and with the logit-link function (Bolker et al. 2009), and model selection was based on the Akaike Information Criterion (Akaike 1973). All models were fitted using the lme4 package (Bates et al. 2015) in R (R Core Team 2015). Fitted values of the model coefficients from the most supported model were used to back-transform the probabilities of giving birth for each YSPB included in the data (i.e., probability of birth after 1, 2, . . ., t years since the previous birth occurred). We also fitted a GLM in which year and female identity could be fixed effects rather than random effects in a GLMM to investigate potential bias from not incorporating temporal and individual heterogeneity in this way.
The expected IBI for the population was estimated in two steps. First, we estimated the probabilities of each IBI in the population as: . . .
where P(IBI = 1), P(IBI = 2), . . ., P(IBI = t) are the probabilities of an IBI of 1, 2, . . ., t years in the population; and p 1 , p 2 , . . ., p t are the conditional probabilities of birth after 1, 2, . . ., t years since the previous birth occurred, backtransformed from the fitted model coefficients. The sum of the probabilities of each IBI in the population must equal one: The expected IBI was estimated as the sum of each of the IBI probabilities multiplied by the number of years in each interval: The fecundity rate, defined here as the annual probability of a mature female having a calf, was estimated as the reciprocal of the expected IBI: To provide a 95% confidence interval (CI) around the estimated expected IBI and fecundity rate, we used a parametric bootstrap based on the fitted model. In each bootstrap replicate, the variance-covariance matrix from the most supported model was used to obtain a random set of model coefficient values, using the mvrnorm function from the MASS package (Venables and Ripley 2002). Each random set of coefficient values was used to back-transform the probabilities of giving birth for each YSPB included in the data and obtain expected IBI and fecundity rate estimates as described in Eqs. 1-4 above. The 2.5% and 97.5% percentiles of 10,000 expected IBI and fecundity rate estimates generated by the bootstrap procedure defined the estimated 95% CIs.

Study population
To test the approach, we used reproduction histories of female bottlenose dolphins from the east coast of Scotland recorded between 1989 and 2012. Boat-based trips were conducted annually to collect photo-identification data (i.e., photographs of the dorsal fin of individual dolphins) following standardized protocols (Cheney et al. 2014). Most surveys occurred during the summer months from May to September; in some but not all years, surveys were also conducted during the winter months, especially in the first half of the study period (1989)(1990)(1991)(1992)(1993)(1994)(1995)(1996)(1997)(1998). Only high-quality photographs (Wilson et al. 1999) were considered to individually identify dolphins from the unique natural markings on the dorsal fin (W€ ursig and Jefferson 1990), matching them to an existing catalogue of known individuals from the population.
Mother-calf pairs were identified based on repeated observations of an adult individual with a calf in two or more trips. Calves seen as young of the year (i.e., during the first year of life) were distinguished from older calves by their smaller size, pale skin, the presence of prominent fetal folds, and nearly constant contact with the mother. When calves where first seen as older calves, their year of birth could still be determined based on their relative size and prominence of fetal lines, which in this population remain visible at least during the first two years of life. However, we only extrapolated the year of birth as far as two years as studies suggest that a calf is likely to become independent of its mother around its third year , Grellier et al. 2003.
Several problems may occur when using photoidentification data to describe individual calving histories: Births may be missed when females are not sighted in a year in which they give birth; births may occur outside the sampling season; a calf may die before the female is encountered following a birth; for a particular female, the calf of a different female may wrongly be assigned to that female (false positive) or her calf may be assigned to the wrong female or not assigned at all (false negative). In all these situations, there is a risk that observed birth intervals may be inaccurate leading to bias in estimates of IBI and fecundity rate.
Births occurring outside the sampling months or missed because the female was not seen in a year in which she gave birth could still be detected during the following sampling period if the new calf survived and was correctly assigned to the mother. However, some births may have been missed if a new calf died before it was observed or when the association with its mother could not be established. Females in this population are capable of reproducing on a two-year cycle after the death of a new born calf (Grellier 2000) but there is no record since the start of the study in 1989 of a female giving birth in two consecutive years under any conditions. Birth intervals of <2 yr have only been reported for Indo-Pacific bottlenose dolphins (Tursiops aduncus; e.g., Kogi et al. 2004, Steiner andBossley 2008). Thus, it was assumed that a female could give birth neither in the year immediately before nor in the year immediately after an observed birth year, based on a gestation period of 12 months (Schroeder 1990). Under this assumption, females that were not sighted in years before and after a known birth event were assumed to be alive and without a young of the year, as long as the female was seen alive in subsequent years. Given these assumptions, data were selected for each reproductive female to model the conditional probabilities of birth, starting with the first known calf to each female (Appendix S1).

Simulations
We used simulations to validate the proposed analytical approach and to investigate the effects of various sources of potential bias in the estimates of expected IBI and fecundity rate. Data samples were generated in two steps: (A) a biological process to generate a population of female bottlenose dolphins and their reproductive histories to represent the dynamics of the population and (B) a sampling process to generate the sighting histories of mother-calf pairs to represent the photo-identification survey effort. The parameters chosen to generate the simulated data (i.e., observed mean fecundity rate, probability of survival, and probability of capture) had been previously estimated for the study population (see Arso Civil 2015 for details). The R code used to generate the simulated data samples is given in the supporting information (Data S1).
Biological process. -We initially considered a scenario of a 50-yr study period and a maximum of 1000 females. Each year, new adult females entering the population were generated from a binomial distribution with a recruitment probability of 0.04, equal to the annual mortality rate estimated for this population (calculated as 1 À probability of female apparent survival of 0.96). The probability of a female surviving between years was drawn from a normal distribution with mean 0.96 and standard deviation (SD) 0.01. Each female could give birth to a first calf according to a binomial distribution with probability 0.25, based on the observed mean annual fecundity rate in the study population. Conditional on having given birth to a first calf, reproductive females could give birth again in a subsequent year based on the conditional probabilities of birth modeled for the study population (see Estimation of expected inter-birth interval and fecundity rate in Methods).
Sampling process. -Variation in survey effort, weather conditions, or quality of the photographs taken may affect the probability of (photographically) capturing individual females and calves. To mimic that variation, the simulated reproductive histories of females were generated in two steps: First, females were captured annually with probability drawn from a normal distribution with mean 0.74 and standard deviation 0.19, as estimated for the study population (Arso Civil 2015). Second, individual heterogeneity in capture probabilities was included through an individual coefficient of heterogeneity generated for each simulated female from a normal distribution with a mean 1 and standard deviation 0.1. These were multiplied together to give a final probability of capture for a given female in a given year. Final capture probabilities were bounded between 0 and 1.
If a female was captured in a year in which she also gave birth, we assumed the new calf was correctly recorded. This assumption is difficult to sustain for any photo-identification study of cetaceans because a calf may die before the female is encountered following the birth. However, this was not included in the simulations because it was not possible to quantify the number of births that may have been missed in the study population.
To validate the analytical method, IBI and fecundity rate estimated from the 100 simulated data samples were compared with those estimated from the study population (i.e., based on real data). We then modified the sampling process (second step above) to create different scenarios and investigate the effects of the following sources of potential bias: (1) including gaps in the data (i.e., including years when a female was not seen and thus assuming that female was alive and did not give birth to a calf in that particular year); (2) assuming different capture probabilities; and (3) varying the length of the study.
For scenario (1), results were compared between two treatments of simulated data: In the first treatment, we selected the reproductive histories excluding years in which a female was not captured, except the years immediately before or after ❖ www.esajournals.org a known birth; in the second treatment, all years were included after the first known birth, regardless of whether or not a female was captured.
For scenario (2), we sampled each simulated sample of mother-calf pairs 10 times, using a range of capture probabilities from 0.1 to 0.99 (i.e., effectively capturing all the individuals and used as a reference) with standard deviation chosen to reflect the same variation as in the study population. Estimates from simulations with capture probabilities of 0.1-0.9 were compared to those with the reference capture probability of 0.99.
For scenario (3), each simulated data sample was subsampled to select the 5, 10, 15, 20, and 25 most recent years and resulting estimates were compared to those from the complete 50-yr simulated data samples; selecting the most recent years avoided having small numbers of reproductive females at the start of each data series.

RESULTS
Between 1989 and 2012, 161 calves were born that could be associated with 84 known females. Six females that only calved in 2012 were excluded from the analysis. Year of birth was estimated for four calves first seen in 1989, but born in 1987 or 1988. Overall reproductive histories of 78 females that had 159 calves between 1987 and 2011 were selected. Each selected female produced between 1 and 6 calves (median = 2 calves, mean = 1.92 calves, SD = 1.11; Appendix S1). The number of calves associated with known females in any one year varied from 0 to 20, partially caused by changes in sampling effort (i.e., years with less effort offer fewer opportunities to detect births); from 2004 onward, all observed new calves could be associated with females. Observed IBIs for females sighted annually ranged between 2 and 9 yr (mean = 3.70, standard error [SE] = 0.20, n = 51).
The most supported model of IBI included linear and quadratic terms of YSPB as fixed effects, and female identity and year as random effects (Appendix S2), but models with either female identity or year as random effects were also well supported. Models including the number of previous calves were less well supported. Models without both YSPB and YSPB 2 were poorly supported. Estimated conditional probabilities of birth increased from YSPB = 1-6 yr and then decreased until YSPB = 9 yr (Table 1). Estimated IBI probability peaked at 4-5 yr ( Table 1). The expected IBI for the study population was estimated at 4.49 yr (95% CI = 3.94-4.93 yr), and the fecundity rate at 0.222 (95% CI = 0.218-0.253). Fitting a GLM instead of a GLMM resulted in a most parsimonious model with only YSPB and YSPB 2 as covariates. Expected IBI and fecundity rate were estimated as 3.99 yr (95% CI = 3.52-4.42 yr) and 0.250 (95% CI = 0.242-0.284), demonstrating the bias resulting from not incorporating temporal and individual heterogeneity in a GLMM framework.
Data samples simulated from the study population generated a mean expected IBI of 4.37 yr (n = 100, range 4.23-4.50 yr), only 2.8% smaller than the expected IBI estimated from the study population. All estimates from the simulations fell within the 95% CI estimated for the study population. The fecundity rate estimated from simulated data samples varied between 0.222 and 0.236, with an overall mean of 0.229, and was thus 2.8% larger than that estimated for the study population. Given that the proposed approach appears unbiased when applied to the study population, the observed IBI for this population (3.70 yr) would be negatively biased by 15-18% compared to the expected IBI estimated from simulated data and from the study population, respectively.

Data-driven biases
Including gaps in the reproduction histories, which increased the number of observed birth intervals included in the analysis by a factor of 3.7 on average, consistently resulted in overestimated expected IBI (by 36%, on average) and underestimated fecundity rates (by 26%, on average), and estimates of both parameters were generally less precise. Assuming low capture probabilities (0.1 -0.3) resulted in expected IBIs and fecundity rates that were biased, on average, by À18.6% and 24.4%, respectively, compared to assuming a capture probability of 0.99, and produced highly imprecise estimates for capture probabilities of 0.1 and 0.2. For capture probabilities between 0.4 and 0.6, the bias was reduced to À6.2% and 7.5%, respectively. High capture probabilities between 0.7 and 0.9 generated an average bias in expected IBI and fecundity rate of only À1.6% and 1.6%, respectively. Limiting the study length to 5 yr generated expected IBIs and fecundity rates that were biased by À9.5% and 10.6%, respectively, compared to the 50-yr study results. These estimated parameters were highly imprecise, and even though the probabilities of each IBI could be estimated, the coefficients of the fixed-effect parameters from the best model were not significant (P > 0.05). The bias was much reduced in 10-yr studies but estimates of expected IBIs and fecundity rates were still imprecise (bias 95% CI = À9.2% to 4.75% and À4.5% to 10.1%, respectively). Study lengths of 15 yr and greater generated unbiased and precise estimates of IBI and fecundity rate. Details on the results from all data-driven biases can be found in Appendix S3: Figs. S1-S6.

Inter-birth interval modeling to estimate fecundity rate
We developed a new approach to estimate IBI, and consequently fecundity rate, from observed birth intervals and applied it to data from the Scottish east coast bottlenose dolphin population and to simulated data. Results demonstrated that the approach produced accurate estimates of expected IBI and fecundity rate for the study population, using 24 yr of data characterized by high annual capture probability (>0.7). They also showed that observed estimates of IBI that do not incorporate missing data are negatively biased by an amount that decreases with length of study and increases with magnitude of capture probability; these results are consistent with earlier work by Barlow and Clapham (1997).
The approach could be extended to incorporate other relevant information into the modeling framework. Shorter observed birth intervals have been reported in cetaceans when pre-weaned calves die (e.g., Steiner andBossley 2008, Henderson et al. 2014). The fate of the previous calf could be included to investigate how this may influence estimated probability of birth. Fecundity has been observed to decrease with age in some odontocete species (Martin andRothery 1993, Matkin et al. 2014); female age (or time since first identification as a proxy) could be included to investigate age-related effects on probability of birth and fecundity.
Two main data-driven biases can affect estimated birth intervals: short study periods that do not allow for the longer birth intervals to be observed (Barlow and Clapham 1997) and missing births in the reproductive histories of individual females (Baker 1987), which is a function of the probability of capturing individuals. Results from our simulations showed how short study periods, missing data years, and low capture probabilities can lead to large bias in estimates of birth intervals and fecundity rate. We recommend studies of at least 10 yr long with capture probabilities of at least 0.3 to estimate birth intervals and fecundity rate using this method. We propose that the direction and size of bias obtained in each of the modeled data scenarios could be used as a reference to inform other studies regarding the likely robustness of the method to generate useful information from current data and to help plan effective future data collection regarding length of study and capture probabilities. This is clearly relevant for bottlenose dolphins but should also be informative for studies of other species with similar patterns of reproduction.

Fecundity rate in bottlenose dolphins
Our estimated birth interval of 4.49 yr (95% CI = 3.94-4.93 yr) for the eastern Scotland bottlenose dolphin population is larger than the observed average birth intervals reported for bottlenose dolphins in North Carolina, USA (2.9 yr, SD = 1.19 yr, range: 2-7 yr; Thayer 2007), and in Cardigan Bay, Wales (3.3 yr, range: 2-6 yr; Feingold and Evans 2013), but similar to the average reported in Shark Bay, Western Australia (4.55 AE 1 SE, range: 3-6.2 yr; Mann et al. 2000), and falls within the birth interval range observed in Sarasota, Florida, USA (3-6 yr; Wells and Scott 1990). These differences among populations could be caused by the above-mentioned datadriven biases. None of the studies included gaps to estimate birth intervals. However, Thayer (2007), Feingold and Evans (2013),  studies are 11-12 yr long, from which we could expect a small bias on the results based on our simulations. If the probability of capturing individuals was low, then some births could have been missed, introducing bias on the estimated birth intervals and fecundity. Differences could also genuinely reflect differences in the reproductive strategies of different populations, or it could reflect a combination of both.
Bottlenose dolphin females invest heavily in their calves, as do other marine and terrestrial long-lived mammals (Van Lawick-Goodall 1968, Lee 1987. In the study population, offspring stay strongly associated with the mothers for at least 3 and up to 8 yr (Grellier et al. 2003) but weaning time varies among and within populations . Being at the northern extreme of the coastal range of the species, individuals in the study population reach larger body sizes (3.5 m long; Scottish Marine Animal Stranding Scheme, personal communication) than in other populations. Data from stranded bottlenose dolphins in Scotland show blubber layer thickness up to 44 mm (Scottish Marine Animal Stranding Scheme, personal communication) compared to layers <22 mm from bottlenose dolphins in Sarasota biopsied during the winter (Montie et al. 2008). These extreme variations in morphology among different populations lead to variation in the energetic demands of individuals; larger females will have higher metabolic demands for growth, maintenance, and lactation. This could explain the longer IBIs and consequently lower fecundity rates in Scottish east coast bottlenose dolphins, reflecting longer lactation and resting periods to meet metabolic demands while successfully weaning calves.
Reproductive parameters estimated from the long-running study of bottlenose dolphins at Sarasota (Wells and Scott 1990) have frequently been used as representative of other populations in the absence of population-specific information on life-history parameters, including the Scottish east coast population (e.g., Thompson et al. 2000, Hall et al. 2006. The fecundity rate estimated here for the Scottish east coast population (0.222, 95% CI = 0.218-0.253) differs by 35% from that estimated for the Sarasota population (0.144, SD = 0.244; Wells and Scott 1990). This illustrates well that life-history parameters from a different population may not be a representative reflection of the dynamics of a particular study population and it may be potentially misleading when considering the viability of populations.

CONCLUSIONS
Accurate estimates of fecundity rate are key to assess the status of populations and effectively direct conservation efforts. Mean observed birth intervals, however, tend to be underestimated and consequently produce overestimated fecundity rates, which could result in overly optimistic assessments of the status of threatened populations. Our simulations showed how short study periods, missing data years, and low capture probabilities can lead to large bias in estimates of birth intervals and fecundity rate. For species with similar reproductive histories to bottlenose dolphins, we recommend studies of at least 10-yr duration with capture probabilities of at least 0.3 to estimate birth intervals and fecundity rate using this method. As datasets with individual reproductive histories are increasingly available for terrestrial and marine mammal populations, this new approach offers an opportunity to obtain unbiased estimates of reproductive parameters from populations with comprehensive data on birth intervals.

ACKNOWLEDGMENTS
Funding for this work was provided by the Department of Energy and Climate Change (UK). Photoidentification data were collected during a series of grants and contracts from the BES, ASAB, Greenpeace Environmental Trust, Scottish Natural Heritage, Scottish Government, Whale and Dolphin Conservation, Talisman Energy (UK) Ltd., DECC, Chevron, and the Natural Environment Research Council. All survey work was carried out under Scottish Natural Heritage Animal Scientific Licences. The authors have no conflict of interest to declare. We thank Mark Bravington for his helpful advice at the early stages of this work and two anonymous reviewers for their useful comments on the manuscript.