Proportion of killed or injured animals that fall into the search area

The proportion of killed or injured animals that fall into the search area (and do not move away) is determined by the size and distribution of the area that can be searched, the spatial distribution of the carcasses, and the proportion of injured animals that move away from the study area before they die. If you can assume that all animals that have been killed have fallen into the area that has been searched (and that no additional carcasses dead from another cause lay in the search area), you can skip this chapter and use a = 1 later. However, often only parts of the area where carcasses fall into can be searched due to inaccessibility, time and funding constraints, or because of dangers such as traffic. In many studies, the proportion of carcasses that have fallen into the search area is far lower than the proportion of carcasses removed by scavengers prior to the search or the proportion of carcasses overlooked by an observer (Warren-Hicks et al. 2013, Huso and Dalthorp 2014). Therefore, a is often among the most important parameters in terms of its influence on the resulting estimate of numbers of animals killed.

To obtain a the spatial distribution of carcasses has to be combined with the spatial distribution of the search area. In many studies, carcass density decreases with distance to the obstacle (the cause of the hazard). Carcass density may, therefore, be approximated by a function of distance from the obstacle. The density function of distance to the obstacle can be different in different directions depending on the shape of the obstacle, its height, prevailing wind directions and other factors. Carcass does not provide a general method that allows estimating a. Here, we only present an example of bat collisions at wind turbines to illustrate the principle of calculating a using information about the spatial distribution of carcasses from literature.

The data frame batdist contains examples of empirical and theoretical carcass distributions of bats at wind turbines (Table 1). The expected proportions of bats in each 10 m distance ring are given. Hull and Muir (2010) use a ballistic model to describe the spatial distribution below three different sizes of wind turbines for different carcass types and sizes, but the model does not seem to fit empirical data well (Arnett et al. 2005, Niermann et al. 2011, Huso and Dalthorp 2014). In addition, variation among sites and years may be substantial (Table 1). An empirical estimate of carcass distribution may be preferred (Dalthorp and Huso 2013).

To obtain the proportion of killed or injured animals that fall into the search area (a) , we can, e.g. multiply the proportion that was searched in each of the 10 m distance rings with the proportion of carcasses found in the corresponding ring and sum over all rings. As an example, we use the bat distribution averaged over the first three studies in batdist and multiply these proportions with the proportions of the search area per ring. The latter proportions may be obtained using a geographic information system (GIS). When sites are heterogeneous regarding visibility class, it is preferable to use the proportions of carcasses falling into the different sub-units (e.g. distance rings) as weights when averaging detection probability p or any other parameter, such as f or s, over the whole site.

If the study species may move away from the wind turbine after collision, but before dying, we have to multiply the proportion of carcasses lying in the search area (obtained as described above) by 1 minus the proportion of animals that move out of the search area. For red kites Milvus milvus the proportion of individuals that move away from the turbine after collision may be up to 15 % (Bellebaum et al. 2013). We are not aware of estimates for other taxa.

Table 1.

Estimated spatial distribution of bat carcasses below wind turbines in different studies (1 = Niermann et al. 2011,2a,b = Arnett et al. 2005, 3 = Hull and Muir 2010). Rotor = rotor diameter in m, Nacelle = nacelle height in m.

Carcass persistence time

Carcass persistence time is the duration a carcass stays on the ground before it is removed by scavengers or before it has decomposed (to a state that makes it essentially undetectable by the observers). Carcass persistence time depends primarily on the density and activity of scavengers, temperature and humidity. Thus, persistence time likely varies between sites. carcass allows the estimation of carcass persistence probabilities for multiple study sites simultaneously.

To estimate carcass persistence time, individual carcasses have to be monitored over time to record the time until they disappear. Such time-to-event data are commonly analysed using survival models to estimate (daily) persistence probability.

Survival models differ in their assumption about how removal probability changes over time. Several studies found that persistence probability of carcasses increases over time, that is, fresh carcasses are removed with a higher probability than old carcasses (Van Pelt and Piatt 1995, Erickson et al. 2004, Ward et al. 2006, Warren-Hicks et al. 2013). Bispo et al. (2013) explore the fit of different parametric survival models to carcass removal data.

At present, two different models to estimate carcass persistence time are implemented in carcass:

But other survival models may be preferred (Bispo et al. 2013). From any survival analyses, either daily persistence probabilities or the proportion of carcasses that persist until every day after death (‘survivor function’) can later be used to estimate carcass detection probability.

Here, we first estimate daily carcass persistence probability s assuming a constant persistence probability. Second, we estimate separate daily persistence probabilities for each day after death. For both analyses, we will use the data set ‘persistence’. It contains persistence times in days for 636 brown mice and a few bat carcasses from 30 wind turbines. For each carcass, we know the persistence time in days (variable ‘perstime’) if the variable ‘status’ is 1, or we know the minimal persistence time if ‘status’ is 0 for carcasses that have not disappeared until the last search (i.e. right-censored scavenging events). For example, the first brown mouse tested at turbine 5 (see first line in example data) was present at day 1 and day 2 after its placement (at day 0) but disappeared between day 2 and 3, i.e. the disappearance was observed (status = 1) at day 3. Another mouse, on line 7, persisted until the last control on day 14. Therefore, its persistence time is 14 and the status is 0 (‘disappearance not observed’). Note that the persistence times are not known exactly because a scavenging event can hardly be observed directly, but it can be recorded that the carcass has disappeared between the last and present control. The models presently implemented in carcass assume that status = 1 means that the disappearance is exactly known. Therefore, when the time between two controls is large, other methods (e.g. allowing for interval-censoring) may be preferred.

When we assume a constant persistence probability, we fit an exponential survival model and obtain one estimate for persistence probability per turbine (or site). The argument pers.const = TRUE is used to tell R that we assume a constant persistence probability.

The function returns a data frame with one row per site, the estimated daily persistence probability with the lower and upper limit of the 95% confidence interval, as well as the mean persistence time in days.

When we cannot assume a constant persistence probability, we may use the Cox proportional hazard model by setting the argument pers.const = FALSE (i.e. the default).

In this case, the function returns a list with four elements. The first element (‘time’) contains the time points for which the Cox-regression has estimated the proportion of remaining carcasses. Each of the other three elements is a matrix where rows correspond to time (one row for each time point given in the ‘time’ element) and columns correspond to sites. The values are the estimates, the lower and the upper limits of the 95% confidence intervals for the probability that a carcass persists up to the time of the corresponding row.

Both models fit reasonably well to the data (Fig. 1), and the difference between the two models is small (obviously, this may not always be the case).

When data from several sites are available, the Cox proportional hazard model estimates the shape of the disappearance probability over time (hazard function) based on the pooled data set, that is, by assuming that the proportions of disappearance probabilities ( = 1 — persistence probability) between the different sites is constant over time. Therefore, the results differ depending on whether separate analyses for each site or one analysis including all sites are done. If data for only one site is available the result from the Cox proportional hazard model equals the data (Kaplan— Meier curves in Fig. 1). Such a result may be reliable only if sample size is very large. Therefore, if data for only one site is available, or sample size is small, we recommend a parametric model such as an exponential or a Weibull model. The latter is not yet implemented in carcass. Santos et al. (2011) list carcass persistence times for different taxonomic groups. However, be aware that persistence time can strongly vary between locations and date.

The persistence probabilities estimated above are probabilities that a carcass persists for a given time interval, e.g. 24 hours. However, to estimate carcass detection probability, the relevant parameter is the proportion of persisting carcasses from all carcasses that arrived during a specific time interval. Thus, the estimates above can be used only if all carcasses arrive at the beginning of each time interval. Such an assumption may be realistic for bat carcasses since often the majority of the bat carcasses arrive in the hour after sunset (Niermann et al. 2014). However, in other situations, carcasses may arrive continuously. If we assume a constant arrival rate, the persistence probabilities can be integrated over the last time interval to obtain the proportion of carcasses persisting until the end of the time interval.

Figure 1.

Proportion of persisting carcasses over time at nine different wind turbines (nr. 11 — 19). Black = data (Kaplan—Meier curve) with 95% confidence interval (dotted lines), dark grey = estimate from exponential model assuming constant persistence probability, light grey = estimate from Cox proportional hazard model allowing for a time dependent persistence probability.

Note that this integration is done automatically in the function when the argument arrival' is set to ‘uniform’.

Searcher efficiency

Searcher efficiency is the probability that an observer finds a carcass that is in the search area at the time of the search. Searcher efficiency depends on experience of the observer, characteristics of the ground (particularly vegetation density influencing visibility), and type of object that is searched for. Furthermore, the average searcher efficiency often decreases with the number of searches the objects have been exposed to, because decaying carcasses are harder to find and because highly detectable objects will be found earlier than hidden objects. For example objects that have fallen into a hole will remain longer on the search area (Huso 2010).

Usually, experimental trials are performed to assess searcher efficiency In such trials, carcasses or carcass-like objects are distributed by a second person and the observer whose efficiency is being estimated must search for them. Such trials yield data similar to the data in the data frame searches:

The example data was collected by Niermann et al. (2011). It contains the number of detected and undetected items for each of nine observers (variable ‘person’) and three visibility classes of the ground. The function estimates searcher efficiency f for each person and visibility class, as well as for each visibility class averaged over all persons. The average may be useful, if trial data on searcher efficiency are missing for some persons in the study. The function uses the function from the package lme4 (Bates 2005) to fit a generalized linear mixed model with binomial error distribution and logit link function to model the number of detected and not detected items. From this model the fitted values (search efficiency) are extracted for each person and visibility class as well as for every visibility class averaged over all the persons. Visibility is used as a fixed effect and person as a random effect. When the data are limited to only one or two observers, observer is also treated as a fixed effect, and a generalised linear model is fitted by the function from the stats-package. Treating person as a random factor in the model has the advantage that for persons with a low number of experimental data (i.e. with weak information about their searcher efficiency) the estimate is shrunk towards the overall mean. In this way we avoid large over- or underestimation of person-specific searcher efficiencies that are due to the stochasticity of the experimental data (Gelman and Pardoe 2006).

To estimate searcher efficiency using one must provide the data (i.e. observer id, visibility for each trial, number of detected and not detected items) either bundled in a data frame or as separate vectors. If the function is applied to a data frame, this data frame must contain the following four variables with the exact names: ‘person’, ‘visibility’, ‘detected’, ‘notdetected’ (spelling and capitalization is important!). If data for only one visibility class is available (or visibility class has not been recorded), the variable ‘visibility’ has to contain the same value on each line. The function allows two (and only two) categorical explanatory variables (visibility class and observer id), thus it is not possible to include other explanatory variables or covariates.

The precision of the estimated search efficiency is obtained using Bayesian methods, because these methods are more reliable in generalised linear mixed models than frequentist methods (Bolker et al. 2009). The function sim from the package arm is used to sample from the joint posterior distribution of the model parameters (Gelman and Hill 2007). The 2.5% and 97.5% quantiles of the random draws from the posterior distributions are used as lower and upper limits of the 95% credible intervals. The number of random draws can be specified in the argument nsim. Default is 1000, but higher numbers are recommended for reliable credible intervals.

The models used by assume that the number of detected items is binomially distributed within each observer and visibility class. If the data are aggregated so that there is one row for each combination of observer and visibility class and visibility class has the same effect on searcher efficiency for each observer (no interaction between observer and visibility class), there is little danger of violating this assumption. However, if the observers have been tested in different trials so that there are two or more rows per observer and visibility class, the variance in the data may be higher than assumed by the binomial model (overdispersion). If overdispersion is likely to be present, the models have to be fitted ‘by hand’ using the functions or and checked for overdispersion as, e.g. described in Gelman and Hill (2007).

In real data, different sub-units (areas) will have different visibility classes with different corresponding searcher efficiencies and, often, also different persistence probabilities (Warren-Hicks et al. 2013). In such cases, it may be best to calculate the detection probability p^*_k for every sub-unit (1,…,K) separately. The weighted average of the detection probabilities p^*_k with weights equal to the proportion of carcasses falling into the different sub-units (a_k) can serve as an overall detection probability p. Dalthorp et al. (2014) explain how to calculate the variance for this weighted average. For simplicity, we here assume that only one visibility class was present.

Combining the proportion of carcasses in the search area a, persistence probability s and searcher efficiency f to get the detection probability p

The probability that an animal that has died is eventually found by an observer (p) is related to the proportion of carcasses that have fallen into the search area (a), the persistence time (s) and the searcher efficiency (f). These three parameters have been combined in various different ways to get p. Table 2 gives an overview over the commonly used methods.

First, some methods assume that the searches have been done regularly (constant search interval). If the searches have been done at irregular time intervals, an average search interval can be used. However, if the search interval is large compared to the carcass persistence time, mortality will be underestimated because during longer than average search intervals some carcasses will disappear without having the chance of being detected during any search. In such cases, it is prudent to use a method that allows for irregular search intervals.

Second, many studies on carcasses found that fresh carcasses are removed by scavengers with higher probability than old carcasses (Van Pelt and Piatt 1995, Ward et al. 2006, Bispo et al. 2013, Warren-Hicks et al. 2013). Therefore, it is better to use methods that allow for age-dependent persistence probabilities, if sample size is sufficient to estimate age-dependent persistence probabilities. It may not be important to account for changing persistence times with carcass age when search intervals are short compared to persistence times and searcher efficiency is high. In some studies change in persistence time with calendar date may be more important than with carcass age, e.g. when temperature or humidity substantially changed or scavengers arrived. Then the formula developed by Etterson (2013) will be a good choice.

Third, some methods do not take into account that carcasses that have been overlooked during one search can be found during subsequent searches (‘bleed-through’). Therefore, for data sets with long persistence times (relative to the search interval) and low searcher efficiency (where a high proportion of carcasses can be expected to be overlooked during the first search and found during subsequent searches) it is recommended to use a method that accounts for this situation (which may often be the case, thus bleed-through should mostly be accounted for) such as Kerns et al. (2005), Korner-Nievergelt et al. (2011) or Etterson (2013).

Fourth, if searcher efficiency is likely to change with the number of searches a carcass was exposed to or with its age the choice should be one of the methods that allow for variable searcher efficiency. This may be particularly important on heterogeneous terrain, where some carcasses fall into holes and others not, or when some animals are injured and seek cover prior to dying. In such cases, the average detectability of a carcass in a cohort decreases with the number of searches they were exposed to because easy to detect carcasses are found first. Carcass detectability may also decrease due to decay of the carcasses, i.e. dependent on carcass age, or change with calendar date, e.g. due to vegetation growth.

At present, the package carcass provides the following different methods to calculate p: Huso's estimator (Huso 2010), the formula we developed earlier (Korner-Nievergelt et al. 2011) including some extensions, and the recently published formulae by Etterson (2013). The methods currently provided by carcass are summarized in Table 3. The formula presented by Huso (2010) and implemented in the software Huso et al. (2012) require a mean persistence time estimates from a parametric model such as an exponential or Weibull model. The formula Korner-Nievergelt et al. (2011) assumes constant persistence probability and constant or decreasing searcher efficiency. In the package carcass we extend the formulas published in Korner-Nievergelt et al. (2011) to allow for time variation (based on parametric or non-parametric models) of persistence probabilities (Supplementary material Appendix 1). The formulas by Etterson (2013) are based on the same model as the one in Korner-Nievergelt et al. (2011). A mathematical proof that Eq. 15 of Etterson (2013) is equivalent to the formula Korner-Nievergelt et al. (2013) is given in Supplementary material Appendix 2. Etterson (2013) provided a generalization of this formula that allows for irregular search intervals. He further provides two different variants of time dependence in searcher efficiency and persistence probability (on carcass age or calendar date).

Table 2.

Overview over the most commonly used methods to estimate mortality from carcass searches (using unmarked carcasses). For every method, assumptions regarding the regularity of the search interval, the carcass persistence probability, and searcher efficiency are given. The 5th column indicates whether the method accounts for carcasses that are not found during one search but found in subsequent searches (termed ‘bleed-through’ by Warren-Hicks 2013). The last column names software with the specific method implemented. The first part of the table (“Formulas”) contains methods that first estimate the probability that a killed animal is found and, subsequently, use the total number of carcasses found together with the estimated detection probability to get a mortality estimate, e.g. a Horvitz-Thompson estimate or a posterior distribution using Bayes' theorem. The second part of the tables (‘Models’) contains statistical models that are fitted to the raw search data (number of carcasses found for every single search) taking detection probability into account. The table is adapted from Bernardino et al. (2013).

Table 3.

Overview of the functions currently provided by carcass for calculating the carcass detection probability p from persistence probability s, searcher efficiency f, and optional additional information such as search interval d or a vector of lengths of search intervals of different lengths J. All functions, except phuso account for bleed-through (Table 2).

Detection probability p is a function f(s,f,d,n), where s is the persistence probability, ƒ the searcher efficiency, a the proportion of carcasses that have fallen into the search area, d the search interval (usually in days), and n the number of searches. In case of irregular search intervals, a vector with the length of the search intervals J is used instead of d and n. The factor a affects p in a simple multiplicative way, therefore, we can write p^* = f(s,f,d,n), or p^*=f(s,f,J) respectively, and p = ap*.

Calculation of p^* assuming constant persistence and constant searcher efficiency

Let’s assume that searcher efficiency ƒ and persistence probability s are constant. Then we use either of the functions , or for regular search intervals or ettersonEq 14 for irregular search intervals.

The function is based on the average persistence time () which is 1/-log(s) in the case of a constant persistence probability. We can either provide s or . The search interval d has to be provided in addition to s and f. For we have to provide the total number of searches n in addition (dn is then the length of the study period), and ettersonEq14 requires a vector with the lengths of the search intervals.

The function provides 95% confidence intervals of the detection probability if the argument CI is set to TRUE and lower and upper limits of s and f are given. These confidence intervals are obtained by Monte Carlo simulations.

For the result of a Monte Carlo confidence interval can be calculated by the function CI etterson.

The function can handle irregular search intervals.

Note that when defining a constant s, it is assumed that carcasses arrive only at the beginning of each time interval. To account for continuous carcass arrival, the survivor function (proportion of persisting carcasses for each day after death) has to be integrated over each time interval and the formulas allowing for time variable persistence probabilities used (use ‘sca’ instead of ‘vars’ below).

Calculation of p^* assuming constant persistence probability and time-variable searcher efficiency

If the detectability differs between the carcasses, e.g. because some carcasses have fallen into a hole whereas others are exposed, the average detectability of a cohort of carcasses will decrease with the number of searches they were exposed to. This is because easily detectable carcasses will be found during earlier searches whereas hard to detect carcasses remain on the ground. Huso (2010) suggests that the decrease in average detectability with each search may be as high as 75%. The function permits modelling such a decrease by setting the argument .constant = FALSE. The argument k allows one to specify the factor by which the searcher efficiency f (of a cohort of carcasses) is multiplied at each search (that the cohort is exposed to). The mathematical derivation of this formula is given in the Supplementary material Appendix 1.

If the argument CI is set to TRUE, a 95% confidence interval is calculated using Monte Carlo simulations.

The function allows specifying different searcher efficiencies for all possible number of searches a carcass has been exposed to. For example, when detectability of a carcass is high during the first search and then drops markedly, because the easy to detect carcasses are found during the first search but later the change in average detectability between two searches may vanish, this pattern can be specified by a vector that contains the searcher efficiency for every search separately. Note that persistence probability needs to be specified for every possible age of a carcass in the study.

Searcher efficiency can also decrease with age of the carcass. To account for age-dependent searcher efficiency is not yet implemented. However, in case of regular (and short) search intervals, carcass age corresponds to the number of searches a carcass has been exposed to, and or may be applied. The R-package carcassb

CMR (Péron and Hines 2014) can estimate searcher efficiency dependent on two age classes (old versus fresh). The R-package that is currently under development by R. Wolpert (Table 2) will be able to account for searcher efficiency that depends on age and the number of searches a carcass has been exposed to.

The function allows for searcher efficiency that depends on calendar date. For example, the searches started in high and dense crop where searcher efficiency was low and after half-time of the study the crop was cut so that searcher efficiency suddenly increased.

For the confidence intervals, use .

The function does not allow one to specify the factor of the decrease in searcher efficiency with each search, but it seems to be fairly robust with respect to a decrease in searcher efficiency (Huso 2010, Korner-Nievergelt et al. 2011), and the new software Huso (2012) and Dalthorp et al. (2014) also allows for time dependence in searcher efficiency.

Calculation of p^* assuming time-variable persistence and constant searcher efficiency

The functions and can deal with persistence probabilities that vary with the age of the carcass, takes the proportion of remaining carcasses for each day after death in a vector as input (i.e. the ‘survivor curve’; output of persistence.prob with nonconstant persistence), whereas takes a vector with daily persistence probabilities for every possible carcass age. Both functions cannot handle missing values, thus these need to be imputed.

If the argument CI is set to TRUE, a 95% confidence interval is calculated using Monte Carlo simulations.

To use the function , e.g. because search intervals were irregular, the daily persistence probabilities instead of the survivor curve has to be given. Here, we extract the daily persistence probabilities from the survivor curve and impute the missing values for the carcasses with age higher than 14 days with the estimate for day 14.

If daily persistence probability depends on calendar date rather than carcass age, the function is used.

Also, for the last two estimates, confidence intervals can be calculated using the function .

Calculation of p^* assuming time-variable persistence and time-variable searcher efficiency

We can allow for both decreasing searcher efficiency with the number of searches a carcass was exposed to and time-variable persistence for regular ( or v2) and irregular (v2) search intervals.

Again, it is possible to specify CI = TRUE to obtain a 95% confidence interval. And, again, we recommend using more than the 100 simulations (>2000) given in the example code.

Or, when both persistence probability and searcher efficiency depend on calendar date:

And confidence intervals for the estimates obtained by the formulas of Etterson (2013) are calculated by .

Last important step to obtain p

Finally, we have to combine the probability that a carcass that has fallen into the search area is found by an observer p* with the proportion of carcasses that have fallen into the search area a to obtain the probability that an animal that has died during the study period is found by an observer p. In the simple example used here with one visibility class only (p^* is equal for all sub-units of the search area), p is the product p = ap^*. When p^* differs between sub-units (e.g. due to different visibility classes), p is a weighted average of the sub-unit-specific p^*_k with weights equal to the proportions of carcasses within the different sub-units (a_k, see section ‘Searcher efficiency’).

Estimate the number of animals that have died (N)

Assuming that a is constant over time and that no carcass has been present at day one of the study (which is the day of the first search minus the search interval d), the probability p is used to estimate the number of animals that have died N from the number of carcasses found (c). The method implemented in carcass is using Bayes' Theorem with a uniform prior distribution (Korner-Nievergelt et al. 2011, Niermann et al. 2011). Thereby, for any x, the probability that N equals x is given by:

If we assume that the detection probability p is known (without error), we use the function . The function returns the probability for each possible value of N given c (argument nf, number found) and p, if the argument dist = TRUE. If dist = FALSE, the function returns the median of the posterior distribution (expected) and the 2.5% and 97.5% quantiles (95% credible interval). The Horvitz—Thompson (Horvitz and Thompson 1952) estimate(HT.estimate) equals the number of carcasses found divided by p.

To take into account the error in the estimate of p we can use the function . This function uses Monte Carlo simulation to propagate the uncertainty of p to the estimate of N, as presented in the worked example in Korner-Nievergelt et al. (2011).

The message tells us what assumptions are made concerning the arrival of the carcasses and the time dependence of the carcass persistence probability and how we could change this if needed. Further, it plots the posterior distribution of the number of dead animals if the argument plot is set to TRUE (Fig. 2). The argument pform is one of ‘etterson’, ‘huso’ or ‘korner’. In the argument ‘x’ it is specified which posterior probabilities P(number of animals killed > x) are given in the output.

Thus, for turbine 1 with three carcasses found by the observers, we conclude based on the estimates for the proportion of carcasses that have fallen into the search area â= 0.67, the daily persistence probability ŝ= 0.80, and the searcher efficiency f = 0.79 that the true number of fatalities must be around 7 with a 95% credible interval of 3–18, taking into account the uncertainties in the estimates ŝ and Or: there is only a small probability (of 0.01) that the true number of fatalities was larger than 20.

Figure 2.

Posterior distribution of the number of fatalities.