Removal modelling in ecology: A systematic review

Removal models were proposed over 80 years ago as a tool to estimate unknown population size. More recently, they are used as an effective tool for management actions for the control of non desirable species, or for the evaluation of translocation management actions. Although the models have evolved over time, in essence, the protocol for data collection has remained similar: at each sampling occasion attempts are made to capture and remove individuals from the study area. Within this paper we review the literature of removal modelling and highlight the methodological developments for the analysis of removal data, in order to provide a unified resource for ecologists wishing to implement these approaches. Models for removal data have developed to better accommodate important features of the data and we discuss the shift in the required assumptions for the implementation of the models. The relative simplicity of this type of data and associated models mean that the method remains attractive and we discuss the potential future role of this technique.


Introduction
Population density, the number of animals per defined area, is commonly used as a simple measure of how an organism interacts with local conditions. If circumstances are undesirable for the species, the density will be low, whereas if conditions are good the density will be high [1].
One of the most interesting problems for the ecologist is the estimation of the density of a particular animal species [2] and a tremendous diversity of techniques are available for estimating the size of populations [3][4][5]. Despite frequent refinements to make some of the more sophisticated techniques more realistic for particular situations (see for example, [6][7][8][9][10]) the removal method remains exceptionally popular.

Synthesised findings
The reviewed literature was published from 1939 to 2019 and interestingly there have historically been long gaps in publications on this topic. However in recent years there has been a more constant stream of published papers, suggesting a resurgence of interest in the method (Fig 2). It has not been possible to determine if this increase in publications over time is a noteworthy trend or if this is confounded with the general increase in numbers of publications more generally. If there is a genuine increase in use of this methodology, it is potentially an indication of the role of removal modelling in studies of reintroduction, especially when translocated individuals are removed from endangered populations [26], and the adaptations of model collection protocols to adopt a "removal" design for other data types such as occupancy and distance sampling, as will be discussed in the Adapting sampling schemes using removal theory Section f this paper.
A strong representation of removal model publications was observed in the fields of Statistics (29 publications), Fisheries and Ecology journals (20 papers each of the disciplines) and 13 papers published in ornithological journals.

Methodological contributions
Early model developments. The basic principle of the removal method is that a constant sampling effort will remove a constant proportion of the population present at the time of sampling. Thus, if the total population size is N and p denotes the probability of capture of an individual during a particular sampling occasion [21], the expected number of captures will be Table 1. Parameters used to categorise removal models included within the database.

Category Definition
Year of publication Year of publication as it appears in the final print

Reference
Author-year citation style Species Scientific name of the species included in the analysis Taxonomic Group Taxonomic group of the species included in the analysis (Class and Order) given by: pN, p(N − pN) and p[N − pN − p(N − pN)] for the first, second and third sampling occasions, respectively. Population estimates can be obtained either by plotting catch per unit effort of collection as a function of total previous catch (see for example [2,27,28]) or by obtaining maximum likelihood estimates [20,21,29,30]. This model can be formalised by defining the corresponding likelihood function, which we denote by L(N, p;data). Suppose x t denotes the number of individuals captured at sampling occasion t = 1, . . ., T, where T denotes the number of sampling occasions. Using a binomialformulation, if N t individuals are still in the study at sampling occasion t, we can define the probability of removing x t individuals by where N 1 = N and N t ¼ N À P tÀ 1 k¼1 x k for t � 2. Then, the probability of removing x 1 , . . ., x T individuals on occasions 1, . . ., T is given by Alternatively, we could specify that the N individuals within the population belong to one of T + 1 categories: either they are captured on one of occasion 1, . . ., T, or they could never be captured. Let π t denote the probability that an individual is captured at occasion t p t ¼ ð1 À pÞ tÀ 1 p and let π 0 denote the probability that an individual is never captured: Let n denote the number of individuals never captured, which is given by n ¼ N À S T t¼1 x t , then the likelihood can be expressed as: Early developments of removal models were methodological in nature, to overcome issues which these days are simple to deal with because of computing power.
[21] formalised the conditional binomial removal model of [20], providing an asymptotic variance of the abundance estimator. Further, they demonstrated how graphical methods can be used to estimate the parameters of capture and abundance. The likelihood of [20] is weighted with a beta prior by [31], which results in estimators with lower bias and variance. [32] proposed an improved confidence interval for abundance for small populations, whilst [33] proved that the profile loglikelihood for the removal model is unimodal and demonstrated that the likelihood-ratio confidence interval for the population size has acceptable small-sample coverage properties. Similarly, [34] proposes a profile likelihood approach for estimating confidence intervals which showed improved performance.
Validity of model assumptions. The model assumptions required for these early models were very restrictive. Specifically, capture probability, p is assumed to be constant both across all individuals and for each sampling period. [20,21,29,30] have typically included diagnostic tests of assumptions or information on necessary sample sizes in their studies. [14] included a table of percentage errors to be expected if p varies during sampling. [6] presents a large number of capture-recapture models for closed populations which can be used to estimate abundance. These models accommodated temporal, behavioural and heterogeneity effects in capture probability. However, such generalisations for removal models are not possible as each individual is only captured once. [35], noting that the assumptions of the removal model are often violated, proposed the non-parametric jackknife estimator as an alternative to the removal model.
[31] proposed a standard test that combines testing for arrivals or departures of individuals in the population with testing for equal catchability. The test entails examining trends in the number of removals over time: "When the expected third catch as determined from the first two catches is larger than the observed third catch, emigration or a decreasing probability of capture is indicated. When the opposite condition exists, immigration or an increase in the probability of capture is indicated". Both the presence and absence of trends yields ambiguous information. Significant trends in the observed data can be caused by the population being open or by unequal catchability. An absence of trends implies either that the assumptions were met or that migration balanced changing probabilities of capture over time.
[36] indicates that unequal catchability tends to be the rule in biological populations. Hence testing for equal catchability is crucial unless one adapts the model. Assuming that equal catchability exists when it does not leads to underestimation of the population size. One procedure that avoids this bias is to identify subsets of the population that are equally catchable and to obtain separate population estimates for each subset [37]. However, because such subdivision of the data greatly decreases the precision of the estimate of the total population size, one should not divide the population unnecessarily if the assumption of equal catchability of all the individuals is met. Therefore, employing a test of equal catchability is a crucial step in any population analysis, even if failure to reject the null hypothesis of equal catchability is an ambiguous result [38].
Two aspects of equal catchability are important for the removal method: equal catchability among groups and equal catchability in all sampling occasions. The first is tested analogously to the test for marked and unmarked captures. If groups with different catchability are identified, separate population size estimates are made for each such group. Individual differences in catchability unrelated to a particular group membership are still possible, but [32] showed that unless these differences are great, their effect on the population estimate is small. [39] investigated the robustness of the removal model to varying behaviours exhibited by fish using simulation.
The second assumption, that catchability remains constant in all sampling periods can be tested by the χ 2 chi-squared test given by [21] or further test given by [6] or [40]. Conclusions drawn from any of these tests will be accurate only if the population remains closed during sampling. Use of a physical barrier around the study site, if possible, can reduce arrivals and departures from the population. Alternatively, independent verification by sampling of marked animals will allow the validity of the assumption of population closure to be investigated [37].
The closure assumption of the removal model has been relaxed in [41], where a model was proposed which allows for population renewal through birth/immigration as well as for population depletion through death/emigration in addition to the removal process. The arrivals of new individuals are modelled by an unknown number of renewal groups and a reversible jump MCMC approach is used to determine the unknown number of groups. Within this paper however it is assumed that any emigration from the population is permanent, and this assumption has been relaxed in [13] which presents a robust design, multilevel structure for removal data using maximum likelihood inference. The implemented hidden Markov model framework [42], allows individuals to enter and leave the population between secondary samples. A Bayesian counterpart to this robust design model is presented in [43].
Change in ratio, index-removal and catch-effort models. Change in ratio models for population size estimation are closely related to removal models [44,45]. The model requires that the population can be sub-divided into distinct population classes and the removals will be performed in such a way that the ratio of removals of the sub-populations will be the same as the underlying ratio of the sub-classes within the population.
We can generalise the basic removal model likelihood function of Eq (2) by extending the definition of the parameters and the summary statistics. Specifically, suppose the population is sub-divided into G mutually exclusive and exhaustive groups, and let N(g) denote the unknown abundance of sub-population g = 1, . . ., G. We now record x t (g) individuals of subpopulation g being removed at occasion t. The likelihood becomes, where N 1 (g) = N(g) and N t ðgÞ ¼ NðgÞ À P tÀ 1 k¼1 x k ðgÞ for t � 2. Catch-effort models are a straight-forward extension of basic removal models which allow capture probability to be related to sampling effort. If catch per unit effort declines with time, then regressing accumulated removals by catch per unit effort allows the starting population to be estimated. This approach however strongly relies on the assumption that if more effort is put into capturing the individuals then a higher proportion of the population will be caught and if this is not satisfied estimators might be appreciably biased [46]. More generally we can extend the removal likelihood of Eq 2 to define a functional form of capture probability: where p t is the capture probability at occasion t which can be linked to a recorded covariate of survey effort, denoted by w t . Possible forms for the functional form might be logit(p t ) = α + βw t , where α and β are parameters to be estimated, or p t ¼ 1 À exp À yw t , where θ is a single parameter to be estimated has been used for fisheries applications where w t denotes the amount of time spent fishing. Alternatively if w t denotes the number of traps on occasion t and each animal is assumed to be caught in any trap with probability θ, p t ¼ 1 À ð1 À yÞ w t [5].
Indeed the logistic form of time-dependent capture probability can also be used to model time-variation in capture probability as a function of climatic conditions-see for example [47]. When sampling is with replacement and the sampling efforts are known, [48] modelled the survey sampling process as a Poisson point process where each animal is counted at random with respect to increments of sampling effort and it is assumed that the encounter probabilities for each individual are independent.
[49] propose a class of catch-effort models which allow for heterogeneous capture probabilities.
The index-removal method makes use of the decline in a measure of relative abundance due to a known removal. The relative abundance is measured in surveys before and after the removal [50].
[51] proposed an index-removal estimator which accounts for seasonal variation in detection.
Further model developments.
[52] demonstrates why auxiliary information is beneficial in removal studies and how to incorporate the extra information into the model and [15] extended the idea of incorporating sub-class level information by proposing a conditional likelihood approach for incorporating auxiliary variables. Capture probabilities are directly estimated from the conditional likelihood and then abundance estimates can be obtained using a Horvitz-Thompson-like approach.
[53] relaxed the assumption of [20] that traps are not limited in capacity by developing a model in which traps have a reduced capacity to catch once they have been filled. The model assumes that the probability that a specific animal will be caught is proportional to the number of traps that are unoccupied and is often referred to as the proportional trapping model. The theory of analysing multiple types of data in an integrated model within ecology has gained traction in recent years-see for example [58]. Early ideas of combining data types has been found in the removal literature. For example, [37] proposed combining capture-recapture and removal methods for fish removals when sampling is over a limited study period and [59] showed how the proportional trapping model can be extended to include data on non-target species.
Removal models have been presented as a class of hierarchical models, for example [16] present a hierarchical removal model where the sites are assumed to have several distinct subsites located spatially. Suppose removals occur at S sites, then records are made of x st , the number of individuals removed from site s = 1, . . ., S at occasion t = 1, . . ., T. Following the multinomial form of the likelihood of Eq 3 the probability of observing a sequence of removal counts, where N s denotes the abundance at site s, p s denotes the capture probability at site s and n s ¼ N s À P T t¼1 x st . Within the hierarchical formulation, a probabilistic formulation, defined by density function g(N|ψ), with parameter ψ, specifies the variation in abundance among the S spatially distinct sub-populations in the sample. The site-specific removal counts (Eq 6) can be combined with this model by integrating over possible values of N s : The likelihood function is then the product over the observations from the S sites, which assuming independence is defined by This model is in fact a multinomial N-mixture model [60] and it has been shown that the removal N-mixture model outperforms the standard N-mixture model using simulation [61]. In practice, a value K has been used in place of the infinite sum in Eq 8 when evaluating the likelihood, however what value of K is appropriate is subjective, and it has been shown there can be some problems with proposed values [62].
[63] has demonstrated that the multinomial N-mixture model for removal data, with negative binomial mixing distribution, has a closedform likelihood and therefore no numerical approximations are required to fit the model.
Adapting sampling schemes using removal theory. Little work has been found which investigates study design for removal surveys.
[64] explored how to optimally allocate total sampling effort for multiple removal sites by maximising the Fisher information of the constant capture probability in the classic removal model. This approach was extended by [65] to allocate effort between primary and secondary sampling occasions within the robust design removal model [13].
An adaptation of survey design for various data types have been augmented with the concept of removal studies. For example, [66] described a time-removal model that treats subsets of the survey period as independent replicates in which birds are 'captured' and mentally removed from the population during later sub-periods. This method has been implemented in many subsequent papers, see for example [67,68]. Further, the removal study design has been proposed for occupancy surveys [69], whereby once a site has been observed as occupied by a species no further surveys are required [70].
[71] developed a spatially explicit temporary emigration model permitting the estimation of population density for point count data such as removal sampling, double-observer sampling, and distance sampling.

Software
There has been a considerable amount of recent work on developing software to make complex statistical models accessible to the wider ecological community. Much software has been developed to estimate population parameters, including abundance and demographic parameters which account for imperfect detection.
Capture [6], was developed to compute estimates of capture probability and population size for closed population capture-recapture data and given the basic removal model is a special case of closed capture-recapture model M b , Capture can be used to fit the geometric removal model. RCapture [132] is an R package [133] for fitting models to capture-recapture data. As well as open, closed and robust-design versions of these models based on multinomial likelihoods it is the only software which also implements a log-linear modelling algorithm to estimate the parameters.
Mark [134] provides a wide-range of models which can be fitted to more than 65 different data types to estimate several population parameters from the encounters of marked animals. Typically, parameters are obtained by method of maximum likelihood estimation through numerical methods (Newton-Raphson by default). However, an MCMC algorithm has been added to provide estimates using a Bayesian framework [135]. [136] demonstrate how removal models can be fitted using Mark for fisheries data. RMark [137] is a software package for the R computing environment that was designed as an alternative interface that can be used in place of Mark's graphical-user-interface to describe models with a typed formula so that models do not need to be defined manually through the design matrix. At the time of writing, RMark supports fitting 97 of the 155 models available in Mark. R packages marked [138] and unmarked [139] can also be used to fit the standard removal model and multinomial N-mixture removal model, respectively-see [61]. There is also more specialised software that has arisen for specific applications. Removal Sampling v2 [140] was designed to estimate population size from removal data and [105] apply this software in order to analyse the effectiveness of stream sampling methods for capturing invasive crayfish.
In addition there is of course well-known software which accommodates the removal study design when fitting models to other data types. In particular Presence and RPresence [141] for occupancy surveys and Distance [142] for distance sampling.

Discussion
Early removal models were simplistic and did not adequately account for potential variability exhibited by the underlying population, however computational and methodological advances give the possibility of more complex models, increasing opportunities, scenarios and accuracy in population estimation. The framework we have presented here is designed to assimilate the use of removal models in order to assist future practitioners in the effective application of removal models.
In our review we have not only presented a list of papers published, differentiating application and methodological advances, but also we have explained the evolution of the model. We have shown how the model has been developed since [2] presented the first case with the evolution in likelihood function from the basic to, for example, that proposed by [49], accounting for heterogeneous capture probability, and more recently the work of [63], theoretically developing the multinomial N-mixture model for removal data.
We have shown how the model assumptions have been adapted, trying to fit the model to different scenarios such as unequal catchability [36], non closure of population [41], heterogenity accross sites [16] and temporary migration [13,43].
Software development, means that even the complex models described in this paper are accessible to ecologist, meaning that maximum utility can be obtained from removal data.
There are several advantages for non-specialists that wish to apply removal methods: there is a vast array of available models for removal data, with the possibility of selecting the approach where the model assumptions best align with their particular study; there is no restriction to frequentist or Bayesian paradigms; there are several software packages and R code accompanying publications of more recent development to investigate where model assumption might fall short.
In our research we noted in earlier papers a thorough assessment of effects if model assumptions were violated but this rigour was not found in late methodological papers, except in some cases through part of a simulation. New methods developed in this research field have been motivated by unique aspects of particular data sets, and therefore nuances of a case study should be embraced rather than avoided in order to encourage methodological advances.
There is a worldwide interest in identifying tools for effective estimation of species population size and removal models show great potential for application in a wide range of situations, such as species relocation projects and control of exotic invasive species. Even if the original aim of this method was to deplete a species, and some of the studies included in this review use the method for this purpose [2,22,107,108,113,117,119], most of the studies included in this review were focused on understanding abundance of a particular species due to general interest (for example [78, 81-84, 88, 143]). However, as a management tool it has been used not only to estimate abundance [12,144], but also, to estimate catchability [35,66,68,109,118], migration [13,43,71,85] [23]. The potential of removal models to facilitate the estimation of population size in the source population whilst also obtaining a pool of individuals to translocate/reintroduce means that such models will remain important and will likely be further developed.
Species relocations are becoming more prevalent in conservation worldwide [145][146][147]. They are performed in several countries on an extensive range of species including plants [145], amphibians and reptiles [148,149]. There are many studies of translocated species and the success of reintroductions, including settlement, survival and reproduction of translocated individuals and their effects on the viability of the reintroduced population [150][151][152][153][154][155]. However, there is less information regarding the impact of translocations on the source or donor population [156,157]. These impacts can dramatically affect community stability, which is especially important when translocated individuals are from endangered populations [25]. The main components that can affect the stability of a population are: resistance, that is the ability to maintain its current state when subjected to a perturbation [158]; amplitude, that will determine, after some alteration, if it will return to its original state [159]; elasticity is the property that will determine the rate of return to its initial configuration when the perturbation exceeds the resistance of a community, but not its amplitude [160]. Removal data and removal models may be a powerful tool in order to understand and manage these populations.
As we have shown in this review there are several methodological tools available for practitioners. Although a deep mathematical/statistical knowledge is not needed to apply these methodologies due to the prevalence of software, some assessment should be made of the appropriateness of the statistical models in order to obtain robust estimates of parameters of interest. When removal models were first proposed two core assumptions were key: populations had to be closed and capture probability had to be constant across both time and individuals. We have shown that these assumptions are no longer necessary and more general models exist. However we have to recognise that the ability to estimate all parameters from a very complex model will depend on the data you have available. Substantial numbers of papers exist introducing the concept of removal study design in other types of data collection studies however we have found little work on study design of removal studies themselves, an exception being [64]. Some papers have demonstrated the power of collecting additional information during removal studies, such as sub-group information or spatial information, and it is likely this type of adaptation to basic removal data that will facilitate the fitting of more complex models. Addressing the important aspect of removal study design more generally is an area of active research.