Translator Disclaimer
1 June 2011 Evaluation of alternative management strategies of muskrat Ondatra zibethicus population control using a population model
Author Affiliations +

Muskrats Ondatra zibethicus are considered a pest species in the Netherlands, and a year-round control programme is in effect. Currently, the agency responsible for the management of muskrat populations in the Netherlands (the LCCM) is preparing for field studies to compare alternative strategies of control. In order to decide on the specific design of such field studies, a population dynamic model was built. The model compares the current management strategy with alternatives in which the effort is focused in space or in time. The model allows us to prioritise future research questions. The major gaps in knowledge at this moment are: 1) insight into the costs of harvesting at different harvest rates, and 2) the relationship between population density on the one hand and (financial damage or) safety risk on the other hand. We suggest continuing the current management, and to test our hypothesis that intensifying harvest will lead to lower numbers of animals killed in the medium term than more extensive harvest rates. The muskrat control programme offers excellent opportunities for applied biological studies of which the benefits are likely to outweigh the costs.

Muskrats Ondatra zibethicus are and have been exploited for fur in much of their native range in North America, but are treated as a pest species in many other places. In contrast to conservation, agriculture or fisheries, the objective of pest species management is to minimise the population. This is the situation for the muskrats in the Netherlands. The reason for this is a presumed relationship between financial damage or safety risks and population size of muskrats. Damage is caused by their digging activities and the biggest safety risk is the potential flooding of low lying country due to dike failure. The target of muskrat population management is to keep muskrat damage and required repair below a publicly acceptable level. Secondary objectives of the management are to reach this target population size with a minimum of resources, little animal suffering and limited impact on other species.

Whatever the objective, a good understanding of the population dynamics can help in choosing the most promising management options from the available alternatives. From the vast literature on exploitation of populations, it is clear that timing, intensity and spatial implementation of harvesting are crucial in determining the success (Clark 1986, Boyce et al. 1999). The strategy to manage muskrat populations in the Netherlands over the past four decades has been to kill/capture animals during the whole year, all over the country, and to do so without differentiation in space and time (Barends 2002). With one exception (the Oostvaardersplassen), this includes National parks and other nature areas. The intensity of harvesting expressed as a proportion of the population captured over a given period of time is uncertain, given that population numbers are unknown. In terms of resources invested, it has been quite intensive (14 man−2 on average, at a cost of 35 million euro per year in 2007; LCCM 2008). No poison is used.

Currently, the responsible body for management of muskrat populations (LCCM) in the Netherlands is reconsidering this basic strategy. It is preparing for field studies to compare alternative strategies of harvesting that focus all available effort and resources either in space or in time (van Vliet & Lengkeek 2007). In order to decide on the relevance of such field studies and their specific design, we built a population dynamic model in our study. The model simulates the developments in population sizes and total harvest under two basic alternatives. The first alternative is called ′space-differentiated harvesting′ and refers to intensive capture efforts near infrastructure that is sensitive to damage by muskrat digging, such as dikes (BCM 2006) or water treatment facilities (Kadlec et al. 2007). The idea behind this is that the acceptable level of damage away from such objects is higher, such that higher population levels of muskrats can be tolerated. The second alternative is called ′time-differentiated harvesting′ or seasonal harvesting in which no effort is made in summer and autumn, but much effort in winter and spring. Harvesting in these seasons is suggested to be more effective, given the reproductive biology of muskrats (van Vliet & Lengkeek 2007). Similarly, in North America, there are seasons for muskrat trapping, although the reason in North America is to harvest muskrat populations in a sustainable way rather than to reduce population levels effectively.

Muskrats are medium-sized rodents with a semi-aquatic lifestyle and a strong social organisation. Basic reviews are given in Boutin & Birkenholz (1987), Heidecke & Seide (1990) and Perry (1982). The history and result of muskrat introductions in Europe, as well as their dispersal rates and the impact of muskrats on biota and their habitats in northwestern Europe, are discussed in Danell (1996). However, few studies report on rigorous field experiments, and many suffer from difficulties in estimating population size. There is a clear need to design experiments that manipulate harvest proportions and strategies, while studying their effect on demographic parameters and dispersal (Boutin & Birkenholz 1987, van Vliet & Lengkeek 2007). The model study we report in this article has been performed with the objective to guide the development of such field experiments. We aim at: 1) predicting muskrat population development under the basic and the two alternative strategies of harvesting discussed above, 2) providing an evaluation of the alternative strategies based on the model predictions, and 3) clarifying the need and usefulness of field experiments dealing with harvesting strategies of muskrats, by identifying gaps in the knowledge. These aims serve our main objective of suggesting improvements of current plans for field studies to compare alternative strategies.


The basic model

We evaluated the alternative harvesting strategies for muskrats using a metapopulation simulation model. We required our model to have a spatial as well as a temporal component, to include the essential demographic processes and to be flexible in comparing the different strategies of population management. We therefore used a stage-structured stochastic dynamic metapopulation model. It was formulated in an existing model environment for reasons of robustness and transparency. Based on a review in Lindenmayer et al. (1995) comparing ALEX, VORTEX and RAMAS, we judged that Ramas Metapop 5.0 would be the most suitable model environment for our purpose, given the requirements presented above. Ramas Metapop 5.0 is well described (Akcakaya 2008), simplifying communication about calculations and technical details of the model.

The model simulates population development in a series of subpopulations and compares the outcome of different harvesting strategies. In each simulation, relevant aspects of the model were changed, such as season of harvesting, location of harvesting and intensity of harvesting (the proportion of animals harvested). The details of the general model are given in the appendix and in Akcakaya (2008). Here, we present the essential choices, such as those related to temporal and spatial scale as well as vital rates. As we are interested in seasonal harvesting, we used time steps equal to one season of three months. We distinguished seven stages to account for juveniles in three seasons and adults in four seasons. Juveniles reach adult size and weight by winter (Vincent & Quéré 1972), so we decided not to model juveniles in spring. Although in muskrats the sex ratio in adults is often biased towards males (Vincent & Quéré 1972, Parker & Maxwell 1984), we assumed a 50% ratio in adults and sex ratio was not modelled explicitly.

Because of our interest in a strategy of harvesting that is differentiated in space, we modelled a metapopulation that started with 49 subpopulations in 49 cells. The cells were arranged in a square of 3.5 × 3.5 km. The centre of each cell was 500 m from the centre of a neighbouring cell (Fig. 1). This level of scale was chosen in relation to known values of home range and population density, as well as ranges of dispersal by muskrats, and scale and grain of the Dutch landscape. Each subpopulation was defined by the same parameters (Table 1), for reasons of simplicity. These parameters are: type of density dependence (see below), survival and recruitment per stage (the stage matrix), initial abundance and, under assumptions of density-dependence, carrying capacity (K) and maximum population growth rate (Rmax). Within simulations, there can be variation between subpopulations due to demographic and environmental stochasticity (see Appendix I and Akcakaya 2008). The cells are connected by dispersal. Dispersal to neighbouring cells is defined as the proportion of individuals from a subpopulation that move to another population at a given distance. It is calculated from a distance-dispersal function defined by four parameters, i.e. three parameters determining the shape of the function (a, b & c) and maximum distance (Dmax).

The input parameters were obtained from the literature as far as possible (see Table 1). In order to deal with uncertainty, we used a range of values for each parameter involved and grouped parameter values that have positive and negative effects on population viability. We labelled these groups ′pessimistic′ and ′optimistic′, respectively (see Table 1). This short-cut, suggested by Akcakaya (2008), eases the burden of making large numbers of simulations varying each uncertain parameter individually. We took demographic parameters from three different sources (Smith & Jordan 1976, Clay & Clark 1985, Clark & Kroeker 1993). These sources allow the calculation of seasonal recruitment and survival rates, based on mark-recapture estimates over multiple years. We estimated dispersal parameters from data presented in Caley (1987) and Mallach (1971). We fitted the distance-dispersal function to their data by eye. Dispersal between cells is symmetric, and relative dispersal is assumed to be highest in juveniles in autumn and winter (see Appendix I, following Caley 1987). We estimated the capacity (K) based on Leboulengé & Leboulengé-Nguyen (1981), Doude van Troostwijk (1976) and Brooks & Dodge (1986). The estimate of maximum growth rate (Rmax) has no reference to a source. It is consistent with the preliminary results of a non-linear regression analysis on capture data for the Netherlands (Bos et al. 2009), to which its values were calibrated. In our basic model simulations, we assume that density dependence operates, affecting all stages and all vital rates. This is in line with published findings (Errington 1954, 1963, Clark & Kroeker 1993). The type of density dependence was assumed to be ′scramble competition′. Density dependence of this type refers to situations where population growth diminishes at high densities of animals due to worsening returns, i.e. all animals receive fewer resources. In order to test the model robustness for this choice, we also ran simulations assuming ′contest competition′ or ′no density dependence′.


Validation of our models has been done in a qualitative way, by comparing model output with published information (i.e. Parker & Maxwell 1984, Errington 1954) on patterns of population trajectories, comparisons of harvested vs unharvested populations, effects of seasonal harvesting and values of harvest proportions.

Analysis of the model and evaluation of harvesting strategies

Each simulation lasted 40 time steps and was run 100 times. Forty time steps is the equivalent of 10 years. Short to medium time horizons are appropriate when models are used to compare alternative management options (Akcakaya 2008). One hundred runs is considered suitable for our aims. Our model results were compared using the expected minimum abundance, a parameter used in Population Viability Analysis (PVA; Akcakaya 2008). The expected minimum abundance is an index of propensity to decline. It can be used to compare results of two models: the area between two interval extinction risk curves is equal to the difference between the two expected minimum abundance estimates (McCarthy & Thompson 2001). Note that normally ′risk′ is the probability of an unwanted event, but in this case, the event of extinction is not unwanted. Other relevant model outputs are population density and harvest results.

First, we demonstrate the effect of harvesting in different seasons, as it is one of the basic assumptions that such effects are relevant (van Vliet & Lengkeek 2007). Over the whole metapopulation, we simulated a harvest in a specific season at proportions of 10 and 25% of the numbers present, and we compared the results in terms of population viability. We considered proportions of 10 and 25% moderate and intensive, respectively. We then proceeded with the evaluation of alternative strategies defined in our introduction. These alternatives are:

  1. Year-round sustained effort → harvesting in all 49 cells (see Fig. 1; all cells);

  2. Space-differentiated harvesting → harvesting only in five adjacent cells (see Fig. 1; cells labelled ′S′, this pattern mimics a dike); and

  3. Time-differentiated harvesting → harvesting only in winter and spring.

We ran simulations for each strategy at different proportions of harvesting. The year-round sustained effort strategy, at a proportion of 0% is of course equivalent to a no-harvesting scenario (i.e. 4). If harvest is relevant for a given strategy, it takes place during a time-step, before dispersal. We illustrate the effect of harvesting proportion in each strategy on: population density, total harvest and final yearly harvest. The strategies are ranked according to their outcome in terms of resulting reduction in population viability and number of animals harvested. A reduction in population viability is the biological objective of the management. This biological objective should be weighed against other objectives, such as financial costs, social constraints and ethical or practical considerations (Akcakaya 2008). Ethical considerations relate to the fact that animals are killed per se, the way that they are killed, and the rationale behind it (Warren 2007) and to some extent also the sheer numbers. In our study, the number of animals harvested is taken as a proximate for both financial costs and ethical considerations.


General model output

Without harvesting, the metapopulation in our model quickly grew to high levels, with strong seasonal fluctuations and a seasonal maximum determined by K (carrying capacity; Fig. 2). Growth is generally positive in spring and summer but negative in autumn and winter. The range of uncertainty in the parameter values is such that the estimated average maximum levels differ very much between the optimistic and pessimistic parameter groups. The equilibrium situation is reached in approximately four years. When we mention population densities in the rest of the text, we refer to autumn population densities. They are the peak numbers and are relevant for our problem. They are also indicative for the average population over the year.

Harvesting in different seasons

Our model predicts a strong effect of harvesting in the different seasons on population viability. Harvesting in the autumn or summer (square and round symbols in Fig. 3) has less effect on the minimum expected abundance than harvesting the same percentage in winter or in spring (triangles in Fig. 3). In addition, the number of animals harvested is lower when implementing harvest in winter or in spring. So, for a given number of animals killed, the predicted reduction in population viability is greater when harvesting in winter or spring.

Year-round harvesting - effect of harvested proportion

With an increase in harvest proportion, the average autumn population density declines (Fig. 4). The harvest is highest at intermediate intensity of harvesting. Above 15% harvesting per time step, the population is overharvested from an exploitation point of view, given the parameters of growth and carrying capacity chosen in this model. For convenience, we classify proportions of harvesting of up to 15% as ′moderate′ and proportions > 15% as ′intensive′. Note that the average harvest over the entire run of 40 time steps differs from the harvest in the final year, especially at high proportions of harvest. This difference is caused by the difference in initial and equilibrium abundance of muskrats.

Alternative strategies - space or time differentiated

The resulting population density in the central population for the year-round and the alternative strategies as defined in our introduction is shown in Figure 5. The central population is a population of interest in each of the strategies, including the space-differentiated strategy. Each strategy has an impact on the population density, but the effect of harvest proportion is of overriding importance. Under time-differentiated harvesting, the autumn population density is higher than under year-round harvesting. This shows that it is not useless to harvest in summer and autumn, although the effectiveness per animal was shown to be low for these seasons (see Fig. 3). Under space-differentiated harvesting, the population in the central cell is also reduced, but not as much as with year-round harvesting at the same intensity (see Fig. 5). For the space-differentiated strategy, it is of course useful to examine muskrat density in the remainder of the metapopulation as well. This is shown in Figure 6. Where high harvest proportions under the year-round strategy result in very low population densities both in the central and the metapopulation. This is not the case in the space-differentiated strategy, where the metapopulation is predicted to have a high population density, close to what would be achieved under no harvesting.

Evaluation of strategies

In Figure 7, we ranked the major alternative strategies on the basis of 1) our biological objective, i.e. the predicted reduction in population viability, and 2) an approximation of the ′costs′ associated with each strategy, i.e. the final yearly harvest. The graph allows us to jointly consider the two criteria (Akcakaya 2008). Note that the space-differentiated scenario is not properly captured in this diagram; one also has to consider the density or viability of the population in the main area of interest, where harvesting takes place.

We need a selection criterion to be able to select the most appropriate alternative. These selection criteria can be simple, e.g. lowest viability or lowest cost. But one could also include constraints in the consideration, which can be drawn into the graph as horizontal or vertical lines. Finally, one can use the cost:benefit ratio as a selection criterion, while satisfying either the cost or the viability constraint. This is matter of public debate, but we will elaborate upon it to some extent in the discussion section (below).

As we deduced already from Figure 5, the strategy of year-round harvesting is predicted to result in the greatest reduction in average autumn population density and thus in population viability (see Fig. 7). Now that we have concluded that harvest proportion has such an overriding effect, we felt obliged to represent each strategy for intermediate proportions of harvesting and for high proportions of harvesting. Of course these are only two out of a range of possibilities. With increasing harvest proportions, the viability is further reduced in each strategy (see Fig. 7, in which the filled symbols are above their open counterparts). For the year-round strategy, the final yearly harvest also decreases with harvest proportion (see Figs. 4 and 7, and compare intermediate to high harvest intensity year-round). But in the space-differentiated strategy as well as in the time-differentiated strategy, the final yearly harvest increases when increasing the harvest proportion to 25% (see Fig. 7).

The year-round strategy at high harvest intensity (25%) leads to the greatest reduction in viability. If we agree that the number of animals harvested is a good approximation of ′costs′, the year-round strategy at high harvest intensity also has the lowest cost:benefit ratio. This is indicated by the line that connects this point to the origin (see Fig. 7; all other points lie below this line). Note that the year-round strategy at intermediate harvest intensity (10%), results in the highest cost:benefit ratio, in comparison to the other alternatives. To use financial costs as a criterion, one has to assume a relation between the financial costs of each strategy and the amount of animals harvested. This is a very critical relation that we will discuss below.


The ranking of alternative strategies is robust to different assumptions of density dependence (Table 2), and the ranking does not change when pessimistic or optimistic input parameters are used in the model.



Spatially explicit models are important tools for investigating scale-related questions in population ecology (Dunning et al. 1995). Our paper is an example of the application of such a model in the context of pest management. In the case of the muskrat, sufficient information is available (e.g. Errington 1963, Perry 1982, Boutin & Birkenholz 1987, Heidecke & Seide 1990), and the model provides us with a framework into which this information can be incorporated in a transparent way. Our objective was not to precisely predict population trajectories or abundances, and we strongly emphasise that the model results, as with any model, need to be interpreted with substantial caution. Nonetheless, as a tool to think quantitatively, to specifically formulate hypotheses and research questions, and as a means to communicate effectively, quantitative models are indispensable. The model we present here conforms to these objectives.

Robustness of results

All assumptions rank the management options in a similar order, in terms of population viability. In that sense, our model is robust. But this does not mean that all other model output, for each of these runs with different assumptions than our best estimates, is credible. For example, under the assumption of no density dependence, the metapopulation size keeps increasing at an exponential rate and so does the yearly harvest (data not shown). This pattern of unchecked population growth has not yet been observed in any other vertebrate species than humans (Allee et al. 1949); the evidence that density-dependent mechanisms affect the population dynamics of muskrats is strong (Errington 1954, Simpson & Boutin 1989, Clark & Kroeker 1993, Erb et al. 2001).

Comparisons with field examples - validity of the model

The model provides qualitatively realistic results, though of course quantitative validation is required. Population trajectories in the literature show great seasonal fluctuations (Errington 1954, Clark & Kroeker 1993) and so does our model. In the absence of harvest, the densities may become very high, but they are bounded by a maximum (Errington 1954, Clark & Kroeker 1993, Danell 1996). This maximum may vary strongly between habitat and years, by orders of magnitude. Such great variation in carrying capacity was not implemented in our model, and in that sense, our results should be regarded as somewhat idealised. On the other hand, the Dutch habitats are also more stable than the prairie potholes and river systems in the U.S., where such great annual fluctuations have been observed, as they are carefully managed for agricultural production, and the climate is moderated by the presence of the sea. Generally speaking, population models predict that harvested populations have lower average densities than unharvested populations (Boyce et al. 1999). Our output is in line with that finding, but for muskrats, no rigorous field studies are at hand to support it, except the work of Errington (1963) on fur refuges in Iowa. Errington (1963) stated that within the refuges, densities are generally higher than outside. There is one National Park in the Netherlands where a core area is left unharvested, but no data on population numbers have been collected from this area so far. Parker & Maxwell (1984) report on an experiment with controlled harvesting in different seasons. Their results, in line with our predictions, show that combined harvesting in spring and autumn leads to stronger effects upon the population than harvesting in either autumn or spring alone.

Reported harvest proportions vary between 50 and 90% of the autumn population (Smith et al. 1981, Parker & Maxwell 1984, Clay & Clark 1985). Some of these proportions are suspect, because the underlying estimate of the population suffers from difficulties (Boutin & Birkenholz 1987). Parker & Maxwell (1984) harvested about 60% either in autumn or spring, without an obvious decline in population or harvest. Harvesting at those proportions in both seasons caused the population to decline. In autumn, 74% is suggested to be the optimal sustained yield according to a model by Smith et al. (1981). Most of these values can be compared directly to the proportions we used for the autumn, because the literature estimates often refer to animals that were harvested during the autumn trapping season. Our model suggests only a moderate decline in population viability relative to an unharvested situation, at an harvest proportion of 25% in autumn alone (see Fig. 3). A year-round harvest at 35% per time step results in a strong decline in our model. During autumn and winter alone, this corresponds to a compounded harvest of 56% of the autumn population. These two observations are not inconsistent with the above mentioned proportions in the literature. For the other seasons, the comparison is problematic for one reason, because in our study the proportions refer to one time step rather than a full year, and the animals harvested in spring and summer are to some extent born in the year following the autumn that we should refer to.

Current muskrat management

The current situation in Flanders (Belgium) is probably a good representative of year-round intensive management at high harvest proportions. The catch history indicates that muskrat densities in Flanders are probably much lower than in the Netherlands (i.e. in 2008). In Flanders, muskrat management has strongly intensified over the past 15 years. The only reliable data have been collected as of the year 2000, when > 39,000 animals were harvested by the Ministerial agency, the Vlaamse Milieu Maatschappij (VMM), along 2,866 kilometers of major waterways. Other agencies, private parties and the local councils were also campaigning against muskrats at that time, killing an unrecorded amount of animals. By the year 2008 the VMM was the sole responsible party harvesting muskrats. The catches had been reduced during 2000-2008, whereas the area searched had increased dramatically. Currently, the managed populations are presumably small in Flanders. All the effort, which is still year-round and intensive, yielded 5,300 animals in 2008, which corresponds to 1.9−2 or 0.4−1 of shorelength. There are no claims of damage and more or less independent checks of population density by controlling teams of trappers also point at low densities. Interestingly, a large proportion of the catches is made in 5-km buffer zones within the boundaries of Flanders. Current labour investment for muskrat management in Flanders (Belgium) appears not to be significantly higher per unit area than for the Netherlands (12 man−2 in Flanders vs 14 ± 2.1 (S.E.) man−2 (N = 12 provinces) in the Netherlands).

The current Dutch strategy is also a year-round sustained effort. There has been quite some variation in harvest and effort, between different regions within the country over the past 20 years. In general, the current situation in the Netherlands can be considered an intensive management at high harvest proportion and consequently low population sizes. But there are areas where harvest proportion actually appears to be intermediate rather than high. On average, 8 ± 21 animals are caught per km2 or 1.1 ± 4.8−1 of shorelength (Bos et al. 2009). This is considerably more than in Flanders, but similarly to Flanders, the Dutch trend is declining.

We believe it is possible to increase harvest proportion everywhere in the Netherlands to such a level that populations decline to publicly agreed low levels, without an exceptional rise in costs at the medium or long term. This belief is supported by experience from Dutch and Flemish professional trappers (H. Post, pers. comm.), and by data series of effort and catches in the Netherlands for different provinces and habitat types (Bos et al. 2009).

Role of predation

Interestingly, muskrat numbers in Poland have strongly declined without an intensification of harvest effort. A recent analysis by Brzezinski et al. (2010) documents the decline in the hunting bags of muskrats in Poland after the invasion of American mink Neovison vison and provides an analysis that identifies mink predation as one of the most important factors affecting muskrat numbers. Note that this refers to a situation in which the muskrat is unadapted to mink predation (c.f. Soper & Payne 1997). In mainland North America, American mink are also primary predators of muskrat, but they are especially successful when muskrats are stressed for reasons of food-shortage or disease (Shier & Boyce 2009). In the Netherlands, no viable population of American mink has established yet, in spite of a more or less continual release of individuals from fur farms over the past decades (Muskens & Dekker 2010). The chances of future establishment of American mink are limited. The American mink is not protected and can legally be killed. Besides, mink farming in Germany and Belgium has come to a halt, and maybe this will happen in the Netherlands as well. The native counterpart of American mink, the European mink Mustela lutreola is extinct in the major part of its ancient range. Another potential predator of muskrat is the European otter Lutra lutra, formerly extinct in the Netherlands but now being reintroduced. At the moment, the role of predation in population regulation can only be small in comparison to the effects of trapping.

Management alternatives

In deciding upon the most promising alternative strategy, one has to consider harvest, effort, population viability and average population size (given in Figs. 5 and 7), but also practical feasibility of a strategy, the damage level that is publicly acceptable and the effects upon other natural values. These will be discussed below in relation to the different alternative strategies.

Practical feasibility and damage

Obviously, the year-round strategy is feasible, as it is practiced already. The Flemish example illustrates that it is possible to intensify the management within this year-round strategy, which means increasing the harvest proportion. It is unclear, however, to what extent harvest proportion can be further increased before yearly financial costs become an inhibiting factor. There shall be an optimal proportion of harvesting, determined by the relation between yearly financial costs on the one hand and harvest proportion or target average population size on the other. The target population size, in turn, depends on the relation between average population size and damage, in association with the publicly acceptable level of damage for a particular landscape or province. These aspects are identified as the major gaps in knowledge that hamper proper policy making at this moment.

The space-differentiated strategy is infeasible, when implemented with a low surface:edge ratio. This is the case in a situation where harvesting is limited to many small or long sensitive objects, a situation which was mimicked in our model strategy of space-differentiated harvesting. Under those circumstances, immigration will have a disproportionally large effect. However, if implemented on a much larger scale, for example the scale of whole provinces on the higher grounds, it might be promising. In that case, we predict average high population levels in the unharvested areas, and a need for extra harvest effort in a strip of several kilometers along the border to the unharvested region. On the higher grounds, the risk of flooding will only be limited, but it is poorly known what other damage will result from higher population levels in agricultural and urban areas and what the publicly acceptable level of damage is. Although such information is not properly documented, experience suggests that damage will increase with higher muskrat density. For this reason, it is more advisable to experimentally study this relationship between damage and population size under controlled circumstances first, before implementing space-differentiated harvesting over large areas.

The time-differentiated strategy is unfavourable in terms of cost:benefit ratio in comparison to the year-round strategy under the harvest proportions studied. This is because, under these proportions, the metapopulation is able to recover during the no-harvest seasons. But the notion that harvesting animals in winter and spring is predicted to affect population viability more per animal harvested, than harvesting in summer and autumn is important. The principle may be applied in practice by hiring unskilled labour to support the professional trappers in winter and spring, yielding a ′mixed′ strategy. Practically, it is infeasible for any organisation to maintain high quality staff (and trapping requires skill) for parts of the year only. Besides, all capital investments will be underused under a pure time-differentiated strategy. It is either an expensive or an infeasible strategy to embark on in its pure form.

Our study suggests that it is locally possible to increase harvest proportions of muskrats in the Netherlands, and that doing this year-round will be the most effective way to keep muskrats at a low population level. Relative to harvesting at intermediate harvest proportions, this leads to less animal killing.


Our model points at differences in effectiveness of the alternative strategies and allows us to evaluate these differences quantitatively with respect to muskrat numbers. Now, in order to also quantify cost effectiveness of each alternative, increase precision and validate the model, clearly identified field information needs to be gathered. Our model thus clarifies the need and usefulness of field experiments dealing with harvesting strategies. It also allowed us to focus such experiments with respect to spatial scale, the variables that need to be manipulated and the parameters that need to be measured, which was our main objective.

We need insight into the costs of harvesting at different proportions for different seasons, landscapes and population densities. We also need to establish the relationship between population density on the one hand and economical damage or safety risk on the other hand. A final constraint that needs to be quantified is the publicly acceptable level of damage per region of interest. These aspects are identified as the major gaps in knowledge that hamper proper policy making at the moment. The benefits that can be derived from guiding expensive control programmes like these with information derived from well-designed field studies, are likely to outweigh the costs of such research.


We wish to thank Jelmer van Belle, Paul Thomas, Niek Bosma and two anonymous reviewers for their input into our study and comments on earlier drafts of this article.



H. R. Akcakaya 2008. RAMAS Metapop Viability Analysis for Stage-structured Metapopulations, user manual for version 5. Applied Biomathematics. Setauket, New York, USA. pp.  Google Scholar


W. C. Allee, A. E. Emerson, O. Park, T. Park, and K. P. Schmidt . 1949. Principles of Animal Ecology. W.B. Saunders Company. Philadelphia and London. pp.  Google Scholar


F. Barends 2002. The Muskrat (Ondatra zibethicus): expansion and control in the Netherlands. Lutra 45:97–104. Google Scholar


BCM 2006. Gevolgen van graverij door Muskusratten en Beverratten voor de veiligheid van waterkeringen. DHV. Amersfoort, the Netherlands. pp. (In Dutch).  Google Scholar


D. Bos, J. van Belle, P. W. Goedhart, S. van Wieren, and R. C. Ydenberg . 2009. Populatie dynamica van Muskusratten. Huidige en alternatieve strategie van bestrijding in Nederland. (In Dutch with a summary of two chapters in English: Evaluation of current and alternative population management of muskrat Ondatra zibethicus in the Netherlands). Altenburg & Wymenga, ecologisch onderzoek. Veenwouden. A&W-rapport 1250,. pp.  Google Scholar


S. Boutin and D. Birkenholz . 1987. Muskrat and Round-tailed Muskrat. In: M. Novak, J. A. Baker, M. E. Obbard, and B. Malloch . (Eds.);. Wild Furbearer Management and Conservation in North America pp. 315–325. Google Scholar


M. S. Boyce, A. R. E. Sinclair, and G. C. White . 1999. Seasonal compensation of predation and harvesting. Oikos 87:419–426. Google Scholar


R. P. Brooks and W. E. Dodge . 1986. Estimation of habitat quality and summer population density for Muskrats on a watershed basis. Journal of Wildlife Management 50:269–273. Google Scholar


M. Brzezinski, J. Romanowski, M. Zmihorski, and K. Karpowicz . 2010. Muskrat (Ondatra zibethicus) decline after the expansion of American mink (Neovison vison) in Poland. European Journal of Wildlife Research 56:341–348. Google Scholar


M. J. Caley 1987. Dispersal and Inbreeding Avoidance in Muskrats. Animal Behaviour 35:1225–1233. Google Scholar


W. R. Clark 1986. Influence of Open Season and Weather on the Harvest of Muskrats. Wildlife Society Bulletin 14:376–380. Google Scholar


W. R. Clark 1994. Habitat selection by muskrats in experimental marshes undergoing succession. Canadian Journal of Zoology 72:675–680. Google Scholar


W. R. Clark and D. W. Kroeker . 1993. Population-Dynamics of Muskrats in Experimental Marshes at Delta, Manitoba. Canadian Journal of Zoology-Revue Canadienne de Zoologie 71:1620–1628. Google Scholar


R. T. Clay and W. R. Clark . 1985. Demography of muskrats on the Upper Mississippi River. Journal of Wildlife Management 49:883–890. Google Scholar


K. Danell 1996. Introductions of aquatic rodents: Lessons of the muskrat Ondatra zibethicus invasion. Wildlife Biology 2 (4):213–220. Google Scholar


W. J. Doude van Troostwijk 1976. The Muskrat (Ondatra zibethicus L.) in the Netherlands, its ecological aspects and their consequences for man. PhD thesis. Leiden University, Leiden, the Netherlands. pp.  Google Scholar


J. B. Dunning, D. J. Stewart, B. J. Danielson, B. R. Noon, T. L. Root, R. H. Lamberson, and E. E. Stevens . 1995. Spatially Explicit Population-Models - Current Forms and Future Uses. Ecological Applications 5:3–11. Google Scholar


J. Erb, M. S. Boyce, and N. C. Stenseth . 2001. Population dynamics of large and small mammals. Oikos 92:3–12. Google Scholar


P. L. Errington 1954. On the Hazards of Overemphasizing Numerical Fluctuations in Studies of Cyclic Phenomena in Muskrat Populations. Journal of Wildlife Management 18:66–90. Google Scholar


P. L. Errington 1963. Muskrat populations. Iowa State University Press. Ames Iowa, USA. pp.  Google Scholar


D. Heidecke and P. Seide . 1990. Bisamratte Ondatra zibethicus (L.). In: H. Stubbe (Ed.);. Buch der Hege, Band 1: Haarwild. VEB Deutscher Landwirtschaftsverlag. Berlin, Germany. pp. 640–666. (In German).  Google Scholar


R. H. Kadlec, J. Pries, and H. Mustard . 2007. Muskrats (Ondatra zibethicus) in treatment wetlands. Ecological Engineering 29:143–153. Google Scholar


LCCM 2008. Landelijke jaarverslag 2007 muskusratten bestrijding. Waterschap Rivierenland. Tiel, the Netherlands. pp. (In Dutch).  Google Scholar


E. Leboulengé and P. Y. Leboulengé-Nguyen . 1981. Ecological Study of A Muskrat Population. Acta Theriologica 26:47–82. Google Scholar


D. B. Lindenmayer, M. A. Burgman, H. R. Akcakaya, R. C. Lacy, and H. P. Possingham . 1995. A Review of the Generic Computer-Programs Alex, Ramas/Space and Vortex for Modeling the Viability of Wildlife Metapopulations. Ecological Modelling 82:161–174. Google Scholar


N. Mallach 1971. Markierungsversuche zur Analyse des Aktionraums und der Ortsbewegungen des Bisams (Ondatra zibethicus L.). Anzeiger für Schädlingskunde und Pflanzenschutz 9:129–136. Google Scholar


M. A. McCarthy and C. Thompson . 2001. Expected minimum population size as a measure of threat. Animal Conservation 4:351–355. Google Scholar


G. J. D. M. Muskens and J. J. A. Dekker . 2010. Amerikaanse nerts. In: N. Huisenga, R. Akkermans, J. Buys, J. van der Coelen, B. Morelissen, and L. Verheggen . (Eds);. Zoogdieren van Limburg. Verspreiding en ecologie in de periode 1980-2007. Stichting Natuurpublicaties Limburg. Maastricht, the Netherlands. pp. 347–349. (In Dutch).  Google Scholar


G. R. Parker and J. W. Maxwell . 1984. An Evaluation of Spring and Autumn Trapping Seasons for Muskrats, Ondatra zibethicus, in Eastern Canada. Canadian Field-Naturalist 98:293–304. Google Scholar


H. R. Perry Jr. 1982. Muskrats. In: J. A. Chapman and G. A. Feldhamer . (Eds.);. Wild Mammals of North America. The Johns Hopkins University Press. Baltimore, USA. pp. 282–325. Google Scholar


C. J. Shier and M. S. Boyce . 2009. Mink Prey Diversity Correlates with Mink-Muskrat Dynamics. Journal of Mammalogy 90:897–905. Google Scholar


M. R. Simpson and S. Boutin . 1989. Muskrat, Ondatra zibethicus, Population Responses to Harvest on the Old Crow Flats, Yukon-Territory. Canadian Field-Naturalist 103:420–422. Google Scholar


H. R. Smith and P. A. Jordan . 1976. An unexploited population of muskrats with unusual biomass, productivity, and body size Connecticut Department of Environmental Protection. Report of investigations No. 7,. pp.  Google Scholar


H. R. R. Smith, J. Sloan, and G. S. Walton . 1981. Some management implications between harvest rate and population resiliency of the muskrat (Ondatra zibethicus). In: J. A. Chapman and D. Pursley . (Eds.);. Worldwide furbearer conference, Frostburg, Maryland, USA pp. 425–442. Google Scholar


L. R. Soper and N. F. Payne . 1997. Relationship of introduced mink, an island race of muskrat, and marginal habitat. Annales Zoologici Fennici 34:251–258. Google Scholar


F. van Vliet and W. Lengkeek . 2007. Alternatieve strategieën voor bestrijding muskusratten Bureau Waardenburg, Culemborg. Rapportnr 07-182,. pp. (In Dutch).  Google Scholar


J-P. Vincent and J-P. Quéré . 1972. Quelques Données sur la reproduction et sur dynamique des populations du Rat Musqué Ondatra zibethica L. dans le Nort de la France. Annales de Zoologie Ecologie Animale 4:395–415. (In French).  Google Scholar


C. R. Warren 2007. Perspectives on the ′alien′ versus ′native′ species debate: a critique of concepts, language and practice. Progress in Human Geography 31:427–446. Google Scholar


Appendix I

In this appendix, we summarise the basic model that was developed in this study, under the assumption of density dependence using the best parameter group. All model variations can be understood from the tables and the text in the main document, in combination with the information presented here.

Models name:

Date: 08/04/2009 11:48:57

Program: RAMAS Metapop version 5.0

Title: mag ik effe vangen v4 best estimates

Comments: 49 cells of 500 × 500 m = 25 ha

Replications: 100

Duration: 40 time steps (120.0 months)

Stage structure

There are seven stages, three for juveniles (Juv; summer = su, autumn = aut and winter = wi) and four for adults (Ad; spring = spr).


Stage matrix



Proportion of each stage matrix element that is survival (as opposed to fecundity).



Demographic stochasticity is used;

Environmental stochasticity distribution: Lognormal;

CV for dispersal: 0.05;

Extinction threshold for metapopulation = 0;

Explosion threshold for metapopulation = 0;

When abundance is below local threshold: count in total;

Within-population correlation: All uncorrelated (F, S, K where F = fecundity, S = survival and K = carrying capacity).

Standard deviations matrix



There are no catastrophes.

Initial abundances


Spatial structure

There are 49 populations (see Fig. 1 for coordinates).


There are 2,352 migratory/dispersal connections among the 49 populations (100% of the 2,352 possible connections);

The dispersal rates range from 0.00022 to 0.02263;

All migration/dispersal rates are symmetric (same in both directions).


Populations have uncorrelated fluctuations (independent environments).



Initial abundance is 20;

Local threshold is 0.0.

Density dependence

Density dependence type is Scramble;

Density dependence is based on the abundances of all stages;

Density dependence affects all vital rates;

Max. growth rate (Rmax) is 1.3;

Carrying capacity (K) is 28;

Standard deviation of K is 9.0;

Density-dependent dispersal as a function of source population size (slope) is 0.0.

Population management


all populations;

10% of individuals from each stage;

at each step.


[1] Edited by Associate Editor: Klaus Hackländer

Figure 1.

The metapopulation consists of 49 subpopulations (at the start) in 49 cells, arranged in a square of 3.5 × 3.5 km. The centre of each cell is 500 m from a neighbouring cell. Indicated are several (groups of) subpopulations that play a particular role in the simulations; the central population (C), the five cells that together form a ′sensitive object′ (S), to be harvested in the management alternative, ′2. Space-differentiated harvesting′.


Figure 2.

Metapopulation density trajectory of an unharvested population under three groups of parameters (best, optimist and pessimistic estimates). Error bars give standard deviations. The simulations refer to the average Dutch habitat with 100 m water edge/ha. The fluctuations refer to a positive growth in spring and summer and to a negative growth in autumn and winter.


Figure 3.

Effect of season of harvesting and harvest proportion (10 and 25%) on population viability expressed as the minimum expected abundance on an area of 625 ha. Horizontal error bars give standard errors.


Figure 4.

Effect of harvest proportion on autumn population density in the final year in the central population (1st y-axis), average yearly harvest and harvest in the final year (2nd y-axis) under a year-round strategy. Error bars give standard errors. All estimates refer to the average Dutch habitat with an average of 100 m water edge/ha.


Figure 5.

Effect of harvest proportion on the autumn population abundance in the final year in the central population (C in Fig. 1) for four different strategies of harvesting. Error bars give standard errors. All estimates refer to the average Dutch habitat with an average of 100 m water edge/ha.


Figure 6.

Effect of harvest proportions on the autumn population density in the central population and the metapopulation in the final year, for the space-differentiated (squares) and the year-round strategy (circles). In the year-round strategy, harvesting takes place in the entire metapopulation (all 49 cells), whereas under space-differentiated harvesting, only in five adjacent cells is harvested (S in Fig. 1)


Figure 7.

Ranking of scenarios on the basis of population viability (% reduction in minimum expected abundance relative to the unharvested situation) and ′costs′ in terms of number of animals harvested. Time-differentiated harvesting refers to havesting in winter and spring only. From the strategies plotted, the year-round strategy at high harvest proportion has the lowest ′cost:benefit′ ratio (i.e. ′number of animals killed′: ′reduction in viability′) as indicated by the line that connects this point to the origin. It also yields the greates reduction in viability.


Table 1.

Input parameters chosen to run the model comparisons. The sources that have been used for reference are given in brackets behind each estimate and include: 1 Mallach (1971), 2 Caley (1987), 3 expert view (DB), 4 Clay & Clark (1985), 5 Clark & Kroeker (1993), 6 Smith & Jordan (1976), 7 Clark (1994), 8 Leboulengé & Leboulangé-Nguyen (1981), 9 Doude van Troostwijk (1976), 10 Brooks & Dodge (1986), 11 c.f. Bos et al. (2009). Abbreviations in the column Parameters are: Ad = adult, Juv = juvenile, aut = autumn, wi = winter, spri = spring and su = summer.


Table 2.

Reduction in population viability for the scenarios year-round sustained effort, time and space-differentiated harvest. Each scenario, except the no-harvest strategy, is represented twice, i.e. at intermediate (10%) and high (25%) proportions of harvest. The table is sorted on the basis of reduction in viability.

Daan Bos and Ron Ydenberg "Evaluation of alternative management strategies of muskrat Ondatra zibethicus population control using a population model," Wildlife Biology 17(2), 143-155, (1 June 2011).
Received: 21 December 2009; Accepted: 1 March 2011; Published: 1 June 2011

Back to Top