Formal Model Article Format
Print
Formal Model Article Format
The Formal Model for the butterfly Pieris napi (Lepidoptera, Pieridae) agent-based model in the Animal Landscape and Man Simulation System (ALMaSS)
expand article infoChristopher John Topping, Xiaodong Duan
‡ Aarhus University, Aarhus, Denmark
Open Access

Abstract

We present a formal model for Pieris napi, the green-veined white. This model is intended for inclusion in the Animal Landscape and Man Simulation System as the basis for regulatory assessment of the impacts of pesticides on butterfly pollinators. We propose implementing the model using an individual-based format, with the added complication of a dynamically coupled individual-based model for parasitoids and hyperparasitoids. The model’s main drivers are weather, temperature and the distribution of larval food plants and nectar forage resources in space and time. A prototype model description is presented, describing the full model ready for implementation. The model considers individuals at all life stages, from eggs to adults and utilises thermal performance models to represent development. Movement is modelled in detail, integrating dispersal, egg laying and foraging. Mortality sources include parasitoids, background mortality, slow development and pesticide and farm management mortality. A simple toxicological model is described as a basis for future expansion.

Key words:

Agent-based model, ALMaSS, dispersal, foraging, parasitoids, Pieris napi

Introduction

This paper follows the formal model format (Topping et al. 2022) and describes a prototype model prior to implementation for the green-veined white. The green-veined white (Pieris napi) is widespread across Europe and Asia, including the Indian subcontinent, Japan, the Maghreb and North America. It is found in meadows, hedgerows and woodland glades. However, unlike its close relatives, the large and small whites (Pieris brassicae and Pieris rapae) rarely lay eggs on the garden or crop varieties of crucifer. It is not considered a pest species (Balachowsky and Mesnil 1936; Ebert and Rennwald 1993). Host plant choice is likely related to the plant volatiles produced and P. napi has a different preference from P. rapae (Okamura et al. 2019). Hence, although morphologically similar, some responses of the two species will be quite different.

Pieris napi is a relatively well-studied species, which allows the creation of a relatively detailed and well-founded formal model. The model is designed for inclusion in the ALMaSS (Animal Landscape and Man Simulation System (Topping et al. 2003; Topping 2022). Therefore, we rely on this system’s features and adapt to the functional limitations imposed by the framework (see Topping and Duan (2024a, 2024b)).

Aims and purpose

The P. napi model is designed to represent day-flying lepidopteran pollinators in European agricultural landscapes. The initial application of the model is to provide a representative of this group for use in a systems-based approach to regulators risk amendment for pollinators impacted directly or indirectly by agrochemical use, primarily by pesticides. As such, we have a specific section in the formal model for the implementation of exposure and effects of pesticides. The model will be developed in the ALMaSS framework for this purpose at spatial scales of typically 100, but up to 2500 km2, with a resolution of 1 m2 and at daily time steps.

Theoretical framework and modelling approach

The model aims to represent the reproduction, mortality, development and movement of P. napi using an individual-based approach within the ALMaSS framework in a detailed manner that is suitable for potential use in pesticide regulatory environmental risk assessment for pollinators.

As with other published pollinator models in ALMaSS (Duan et al. 2022; Ziółkowska et al. 2023), we use an agent-based modelling approach (Grimm and Railsback 2005). Each butterfly is represented as a separate code object. At any point in time, each object has a particular behavioural state and characteristics such as age, size, location and pesticide loading. These characteristics are important, as they can also influence the behavioural states, such as size being linked to reproductive capacity.

The ALMaSS modelling environment also provides a detailed spatio-temporal representation of the landscape from which individual butterflies can obtain the information necessary to simulate their behaviour to a high degree of spatial resolution. The details of the landscape representation in ALMaSS and the interface to agent-based models are described in Topping and Duan (2024a, 2024b). This representation describes spatial landscape heterogeneity through a detailed raster land-cover map with a spatial resolution of 1 m2 and a polygon map of habitats. Farms are represented as the collection of fields managed by that farm. These units are classified into farm types (e.g. cattle or arable farms), which determines which crops are grown and which affects pesticide use and other crop management characteristics. Each crop management is modelled throughout the year as an individually tailored management model. The cropping system is represented as a pluri-annual crop rotation, based on the proportion of each crop type grown by a farm classification type, arranged in an agronomically sensible order.

Crop and non-crop vegetation grows daily depending on management, weather and soil conditions. A pollen and nectar model is also implemented for flowering vegetation (Ziolkowska et al., in prep). This vegetation model provides the daily pollen and nectar availability per unit area, which can be foraged by the butterflies for nectar.

All life stages of the butterfly will be represented as individuals. Their internal development and interactions with each other and the environment are covered in each section below.

Framing the model

This model aims to describe population changes in time and space and does not include all possible known types of behaviour unless these significantly influence these processes. For example, population genetics and behaviour related to genetic makeup are not considered even though they may affect population processes in the long term.

Other features of the ecology and behaviour are excluded because of a lack of clear scientific support for the consequences. These include the fact that, when mating for the first time, the males transfer a large ejaculate that represents, on average, 15% of their body mass. This ejaculate contains nutrients essential for egg production and somatic maintenance and this ‘gift’ may decrease the number of matings for females (Kaitala and Wiklund 1994). However, the consequence of this behaviour seems unclear at the population level, despite genetic consequences and, thus, we did not include this.

Other known types of behaviour that were excluded are the behaviour of the larvae on the food-plant and the location of eggs when laid. Both have implications for the individual organism, but occur at a spatial scale (sub-plant) that is beyond the resolution of the model. We also excluded long-distance dispersal events, which may be rare, but are likely, because of the landscape scale at which the model is expected to be used.

Although development is included in terms of thermal performance curve functions, the effect of temperature and food plant quality on growth is not included, although a density effect is. This omission could potentially be important because it could alter the size of the butterfly and, therefore, the toxicological impact of pesticides. However, as noted, we cannot include food-plant quality in a mechanistic way, nor is there enough evidence to link temperature and size. This feature could easily be addressed if new information were to be found.

The only mechanistic predation/disease factor included in the current version is parasitoid-butterfly interactions. Therefore, we ignore deterministic predator mortality, such as insects feeding on eggs and bacterial, protozoan and viral diseases. All of these drivers can be important in the population control of butterflies (Dempster 1983; Courtney 1986), but there is no clear evidence that these are reliable density-dependent drivers in P. napi. Thus, the current model includes mortality from these factors as a stochastic background mortality. This assumption can be altered if new information suggests a predictable mechanistic relationship. This is in agreement with Dempster (1983), who concluded that, for lepidopteran populations, regulation was primarily through resource limitation and not intraspecific competition, but the natural enemies have the potential to act in density regulation.

In this formal model, we have defined a simple toxicological model. This is intended as a starting point only for future carefully developed toxicological models for regulatory use. This may include expansion of the current model to include sub-lethal effects and more detailed internal toxicology, for example, including toxicokinetic/toxicodynamic models.

A further, but important, simplification is that we currently assume that females can be mated and do not model males individually. This is simply for ease of implementation and if, for instance, genetics were to be included, males would need to be added to the adult class.

Overview of processes

In Europe, P. napi overwinters as a pupa, with the first butterflies emerging during the spring between mid-April and late May. They are not restricted to specific habitat types and form continuous populations in large areas of heterogeneous habitat (Hanski and Kuussaari 1995). Males are first to emerge (Meyrick 1927), followed by females 4–7 days later; both sexes fly by day. After mating, eggs are laid singly on crucifers and develop according to temperature. Larvae hatch, feeding on the leaves of the host plant. When host plants are abundant, the females prefer to lay eggs on smaller hosts, which provide a more suitable microclimate for faster larval development and higher survival during development (Forsberg 1987). Each of the five larval instars lasts 4 to 5 days according to temperature (Buckler 1885). The number of broods per year varies with location, but two or three broods in central Europe are normal. In colder climates, a single brood is possible.

Growth and development

Earlier studies reported that the threshold temperature for larval development is 8 °C, 398 degree-days are required to complete development and that, at temperatures between 15 and 30 °C, the rate of development is directly proportional to temperature (Stepanova 1962). However, it is well-known that insects grow faster under fluctuating temperatures compared with constant temperatures (e.g. Sanderson (1910)). Therefore, since these fluctuating temperatures will occur in natural conditions, they are a better basis for modelling insect responses than artificial constant temperature experiments. When considering development, Ratte (1985) describes the “non-linear temperature-velocity relationship”, meaning that fluctuating treatments should be “normal”, whereas constant temperature insect development studies were essentially conducted under “abnormal” conditions. In addition, a key trait of insects is that exothermic development is non-linear and asymmetrical (Angilletta Jr. 2009). As a result of this inequality and environmental fluctuations, the resulting dynamics are complex, but mathematically explicable. In temperature variations, animals can reach critical points like the critical thermal minimum or maximum, where processes such as protein denaturation or water freezing occur. Small temperature increases can push insects beyond their upper limits, while changes at lower temperatures are slower. The disparity in temperature effects is explained by Jensen’s inequality (Ruel and Ayres 1999), with fluctuating temperatures having a more significant impact, particularly with larger amplitudes.

In the P. napi model, we will use a thermal performance curve model to describe growth and temperature relationships. As originally proposed, the model has the form of a modified Gaussian curve (Angilletta 2006). This describes a Gaussian increase in performance up to a thermal optimum (To) and a quadratic decline towards a critical maximum (CTmax), at which the growth rate is zero (Equation 1).

r(T)=e-T-To2σ2,ifTTo1-T-ToTo-CTmax2,ifT>ToandTCTmax0,ifT>CTmax Eq. 1

This function requires the estimation of three parameters, To, the thermal optimum at which development is fastest, CTmax, the critical maximum temperature at which development is zero and σ, which describes the slope of the Gaussian component of the function.

However, for P. napi, a subsequent model was proposed describing the same shape of curve, but using four, rather than three parameters. The Lobry–Rosso–Flandrois (LRF) model named by Ratkowsky and Reddy (2017) was proposed to represent microbial growth by Lobry et al. (1991) and Rosso et al. (1993). The four parameters represent three cardinal temperatures (the maximum temperature at which development can occur Tmax, the minimum development temperature Tmin and the optimum development temperature Topt) and the specific growth rate at the optimum (Uopt). Uopt is independent of the three cardinal temperatures and T is the temperature (Equation 2).

r(T)=uopt(T-Tmax)(T-Tmin)2(Topt-Tmin)[(Topt-Tmin)(T-Topt)-(Topt-Tmax)(Topt+Tmin-2T)] Eq. 2

This model was used by von Schmalensee et al. (2021) to describe the thermal performance curves for eggs, larvae and pupae of Pieris napi.

The thermal performance curve function was used to fit performance curves to a set of Swedish populations of P. napi, yielding a table of parameter values for each life stage (Greiser et al. 2022) (Table 1) (sexes were not specified separately).

Table 1.

Parameters describing thermal performance models as used by Greiser et al. (2022) for Pieris napi. Model uncertainty in average percentage prediction error is presented using 90% highest density intervals (HDIs).

LIFE STAGE PARAMETER VALUE LOWER HDI UPPER HDI
Egg Tmin 1.936 0.755 3.004
Egg Topt 30.450 29.520 31.799
Egg Tmax 35.988 32.759 39.618
Egg uopt 0.354 0.337 0.373
Larva Tmin 1.542 0.407 2.824
Larva Topt 29.854 29.021 30.829
Larva Tmax 35.162 32.843 38.033
Larva uopt 0.092 0.084 0.101
Pupa Tmin -0.215 -1.654 1.195
Pupa Topt 29.309 28.657 30.086
Pupa Tmax 34.969 32.635 37.811
Pupa uopt 0.180 0.157 0.205

Eggs

Physically, the eggs are elongated, tapering to a blunt point, approximately 1.25 mm long, 0.45 mm wide and less than 0.2 mm wide across the top. There are 13 or 14 strong ridges longitudinally and many fine cross ribs (Capinera 2020).

Eggs are laid singly on the host food-plant and are assumed to develop following the thermal performance curve using parameters from Table 1 and Equation 2. The size of the hatched larvae is related to the size of the egg in butterflies (Fischer et al. 2002), although size changes might be in specific body parts as an adaptation to food-plant specifics (Ohata et al. 2011).

Fischer et al. (2002) provide a very strong linear relationship between egg mass (megg) and hatchling mass (mhatch) (p < 0.0001) for Bicyclus anynana and, reading from the graph, we can recreate the relationship (Equation 3).

mhatch=-0.1+1.046154megg Eq. 3

It seems likely that other butterflies follow a similar linear relationship, although we did not find similar data for P. napi. In the absence of further information, our assumption is that egg size in P. napi follows the same relationship.

Implementation in the model

We assume the thermal performance relationships from Equation 2 and parameters from Table 1. We assume no effect of egg size on development time, but that egg size will be proportional to larval size on hatching using Equation 3.

Larvae

In P. rapae, larvae reared at low temperatures produce larger pupae and adults than those reared at high temperatures (Jones et al. 1982). In P. napi, we did not find evidence to suggest whether the same trend occurs. Larvae reared at higher temperatures are recorded as having a reduced immuno-competence (Bauerfeind and Fischer 2014a). Although larvae follow the expected temperature-related growth function, it has been noted that extreme temperature events have a disproportionate impact on juvenile survival, body size and longevity (Bauerfeind and Fischer 2014b).

Larval density alters larval survival and final pupal mass, although it has no effect on pupation time; high density reduces pupal size by approximately 10% (Kivelä and Välimäki 2008).

Implementation in the model

We assume the thermal performance relationships from Equation 2 and parameters from Table 1. However, we need to include a density-dependent function for growth rate. This requires a measure of larval density per unit area and a rate per unit density decrease in growth rate. Neither parameter is available in literature; hence, these must be estimated by reverse parameter fitting to result in pupal weights matching the weights and range of weights recorded by Kivelä and Välimäki (2008). We can express this as Equation 4, where G is the amount of mass added during a day, D is a density measure per unit area, d is a fitting parameter representing the change in growth addition per unit density, P is the proportion of growth for the larval stages for that day and M is the maximum pupal mass:

G=D×d×P×M Eq. 4

Equation 4 can be fitted by assuming D is a minimum threshold (Dl), at which point the sum of G over the development period should be equal to M. At the other extreme, with D above an upper threshold (DU), the sum of G would be equal to 0.9M, if we assume a 10% reduction for high densities (Kivelä and Välimäki 2008). Values outside the range Dl to DU would be assumed to be constrained to be within this range, resulting in a maximum reduction in the daily mass added of d and a minimum of zero.

Pupae

The initiation of diapause in P. napi is primarily triggered by changes in day length (Friberg et al. 2011; Lehmann et al. 2016), while its cessation is predominantly influenced by temperature (Posledovich et al. 2015; Lehmann et al. 2016; Lehmann et al. 2017). Consequently, pupae maintained at higher temperatures (≥ 15 °C) will remain in diapause, unable to terminate it, whereas those exposed to colder temperatures (≤ 10 °C) will undergo diapause termination, a process taking approximately three months (Lehmann et al. 2017). Since diapause can commence quite early in the summer, the necessity for lower temperatures probably serves as a safeguard mechanism to prevent diapause termination during late summer or early autumn (Lehmann et al. 2018). In addition to this pattern, in a study of butterflies in Vermont and Massachusetts (Van Driesche et al. 2004), there is an indication that the quality of the food-plant can also trigger diapause. Thus, food-plant deterioration is identified early by the larvae, triggering a higher chance of early diapause. However, this factor and day length must co-vary and are, therefore, difficult to separate. Recently, von Schmalensee et al. (2024) provided a model for P. napi linking diapause termination and post-diapause development as directly sequential, continuous and temperature-dependent rate processes.

The conditions needed to switch between direct development and diapause in pupae are experienced in the larval stages, with the likelihood of switching to diapause inversely related to the age of the larvae experiencing a short photoperiod (Friberg et al. 2011) (Table 2).

Table 2.

The percentage of P. napi pupae continuing with direct development when exposed to short photoperiods as a larvae, from Friberg et al. (2011).

LARVAL STAGE % ENTERING DIAPAUSE
I 3.3
II–III 6.3
IV 31.4
V 87.5

Implementation in the model

Like eggs and larvae, we assume the thermal performance relationships from Equation 2 and parameters from Table 1, which describes development. On hatching, we assume that 50% are male and these are then removed from the model. However, there is one key difference with respect to diapause. Diapause onset is based on the day length and PhoD describes a threshold below which any developing larvae will potentially enter diapause on pupation. If the day length is larger than PhoD, then, at this point, each larva will take a development bifurcation test using the probabilities given in Table 2.

As described by Lehmann et al. (2016), we assume a bifurcation of developmental pathways with an alternative developmental pathway induced by short photoperiods (PhoD) that results in temperature-independent metabolic suppression and maintenance. During this phase, the pupae will simply suspend development, regardless of temperature. If a period of low temperature is experienced, development still occurs, but at a very low rate. High temperatures increase the development rate until the development follows the same processes as directly developing pupae. This means that, as described by Petersen (1949), the first generation adults can be created from a combination of pupae from 2nd and 3rd generations if a partial third generation can be obtained before day length forces pupal diapause. PhoD will need to be fitted during the pattern-orientated modelling phase to ensure the correct long-term and seasonal dynamics under different light regimes.

We will implement the Thermal Performance Curve model (Equation 2) provided by von Schmalensee et al. (2024). The study provides two curves, one for the probability of breaking diapause with temperature, using a mirrored TPC, the other for post-diapause development (Table 1).

Foraging

The adults of this first generation may restrict themselves to taking nectar from one host species (Lees and Archer 1974). Dover (1989) reported that, in southern England, the first generation fed exclusively from Sinapis arvensis, although later in the season, adults of the second (summer) generation feed on up to 11 different host plants). Flower consistency of this type may reflect the best available source of food since this behaviour seems to be plastic, even if the response time to change is slower than seen in other species, such as honey bees (Goulson and Cory 1993). We make the assumption that suitable nectar sources are based on the general availability of nectar per unit area and that the local continuity of resources can explain flower consistency.

To our knowledge, there is no detailed study on the energetics of P. napi. However, a simulation study indicates mechanical limitations on the range of nectar sugar concentrations and nectar extraction times available to butterflies (Kingsolver and Daniel 1979). This model was adapted for Pieris rapae, predicting an overall optimal range of nectar concentration of 31–39% sucrose for butterflies, which was in agreement with previously reported laboratory values (Daniel et al. 1989). In Speyeria mormonia, an experimental approach considered body mass, resting metabolic rate and lifespan, flight metabolic rate, egg number and composition and food intake across the adult lifespan (Niitepold and Boggs 2015). Different flight treatments did not affect body mass or lifespan. Still, they did alter food intake, showing that, when food resources were abundant, butterflies living in a continuous meadow landscape resisted flight-induced stress, exhibiting no evidence of a flight-fecundity or flight-longevity trade-off.

A study, also using Speyeria mormonia, measured the volume of nectar consumed per day for males and females (Boggs and Ross 1993). On average, females consumed approximately 25 µl per day, with males closer to 20, both estimated from the graphs presented. Body mass of S. mormoria is approximately 120 mg (Niitepold et al. 2014). P. napi wingspans are typically in the range of between 40 and 52 mm wide (Emmet and Heath 1990), with a forewing length of 23.8 mm in summer females (Shkurikhin and Oslina 2016). This compares to S. mormoria with a wing length of 26.6 mm in females (Boggs 1987). Thus, we might scale the expected daily intake by 23.8/26.6, giving us an estimate of approximately 22.4 µl per day nectar intake.

The dependence of butterfly flight on weather has been known for a long time (Clench 1966) and this provides limits to when the butterflies can be actively flying.

Implementation in the model

Since the foraging movements are intertwined with reproduction and dispersal, we consider foraging movement under the heading ‘Dispersal’ below.

For food intake rate, the sucrose percentage of nectar can be used to limit the choice of forage for the butterfly. Nectar is available from ALMaSS, as is the sucrose concentration. The studies indicate that the rate of intake in butterflies is limited. Based on P. rapae, the sucrose/sugar concentration between 31–39% will be foraged by P. napi. Further, multiple foraging bouts need to be considered during the day and, in reality, the butterflies will combine foraging and reproduction and integrate the activities together. From the modelling perspective, this integration is not easy, so we suggest viewing this as several feeding bouts where the maximum intake of nectar is limited. This is an artificial construct, but ideally, it will represent a realistic pattern of movement. Thus, we assume that Nf represents the number of bouts and the nectar intake is 22.4/Nf, which is assumed as the volume of the stomach. If we further assume that nectar in the feeding range is between 31–39%, then we can assume that 22.4 µl is the daily required nectar volume. When the stomach’s glucose level drops below zero, it triggers the foraging behaviour for nectar.

Reproduction

Egg production

Females may mate once or multiple times, but there is no indication that this has an impact on their longevity. It is strange that no increase in mating frequency with age occurs since Karlsson (1998) showed that a virgin male of P. napi transfers a large and nutritious ejaculate containing 14% nitrogen, equivalent to that in approximately 70 eggs. This nuptial gift can be seen as a nutritious resource (Wiklund et al. 1993; Wiklund et al. 1998), but, as argued by Bergström and Wiklund (2005), can equally be seen as a parental investment, leading to there being no clear effect on fecundity. Nevertheless, it is suggested that larval competition could be strong enough to affect female mating strategy and drive the divergence of strategies (Kivelä and Välimäki 2008). As a consequence, Kivelä and Välimäki (2008) suggest that the offspring of monandrous females start to develop in relatively low densities and they are relatively large when the offspring of highly polyandrous females start to hatch. However, there is much variation within populations and the mechanisms are still not completely clear. Hence, we assume that mating frequency is not important in egg production and that it need not be simulated as part of the model; rather we will represent this as stochastic variation in the timing of egg production.

In butterflies, the egg size itself appears to be related to temperature and size and age of the mother. The size of eggs in Pieris rapae is inversely correlated with both the age and the size of the mother (Jones et al. 1982). The relationship appeared to be linear and Jones (loc cit) fitted a regression of egg weight in milligrams (megg) on maternal pupal weight in milligrams (mpupa) and maternal age in degree days (a) had a linear relationship described by:

megg=0.125-0.0012mpupa-0.0044a Eq. 5

The regression was significant, but also included a large degree of variability for any age-size combination. However, in P. napi, Karlsson and Johansson (2008) found that temperature affects egg size with the size decreasing with increasing temperature, but found no significant effect of female size (p = 0.13). Unfortunately, no data on how size was affected by temperature are available. Additionally, in this later study, the age of the females was not considered. However, Bauerfeind and Fischer (2008) also investigated P. napi and concluded that “the role of maternal body size as a constraint on egg size has been previously over-emphasised and other factors are likely to be more important”. Hence, the situation in P. napi is unclear.

In Pieris rapae, total egg production is linearly related to pupal weight (Jones et al. 1982). Maternal size and lifetime fecundity were also significantly positively correlated in P. napi in the study by Bauerfeind and Fischer (2008). In the same study, P. napi egg size contributed significantly to the variation in egg number, showing a positive correlation.

Whether the female mates once or multiple times appears to affect total fecundity. Wiklund et al. 01993) found that monandrous first-generation butterflies had a total fecundity of 247 +/- 47 compared to polyandrous females with 440 +/- 49. Two second generation experiments yielded 280 +/ 45 and 378 +/- 83 for monoandrous butterflies compared to 403 +/- 40 and 611 +/- 46 for polyandrous butterflies. Bergström et al. (2002) also looked at polyandry and fecundity showing large differences in fecundity with the size of butterfly, but relatively weak effects of polyandry. In their experiment, polyandrous butterflies produced approximately 325 or 450 eggs on average (small vs. large), which fits well with the first-generation range found by Wiklund et al. (1993).

Implementation in the model

We assume that, under different constraints and conditions, different factors affect the size of the eggs produced and that adult size, adult age and temperature all play a part in a linear relationship for each. This can be expressed in the form of Equation 6, where e is the maximum possible egg mass, c, d and e are constants, a is the age of the female in degree days and T is the temperature.

megg=e-bmpupa-ca-dT Eq. 6

From Jones et al. (1982), that we have a range of egg sizes from 0.092 to 0.114 mg and from Wiklund et al. (1993), a similar lower bound, but higher upper bound of 0.155 (read from graph). If we assume these are close to the upper and lower possible bounds for egg size and that their recorded death limit of 400 day degrees is the limit of survival, then we can parameterise a linear model combining all three factors (Equation 6). However, in the study by Jones et al. (1982), they did not account for temperature. Therefore, we also need to adjust the maximum and minimum possible sizes to account for temperature variation.

In the study of Bergström et al. (2002), small females weighed 50.8 mg and large 76.9 mg. If we assume a linear relationship between these two points, this allows us to characterise the maximum number of eggs produced with changing weight (represented by x) as Total Fecundity = 4.7893x + 81.705 eggs. If we make the assumption that the results of the three experiments by Wiklund et al. (1993) were confounded by adult size, then we can use this relationship to determine the total fecundity of a female at birth.

Oviposition

Plant chemistry is a key factor in differential selection of egg-laying host sites (Chew and Renwick 1995). There is some plant overlap with P. rapae, but for most plants, P. napi strongly prefers plants avoided by P. rapae (Huang and Renwick 1993). The favoured host plants are natural crucifer species, with eggs laid on the leaves. Van Driesche et al. (2004) studied the butterfly in New England, finding that first generation butterflies laid their eggs in wooded habitats (95%), whereas second generation adults oviposited in meadows in full sun habitats. This indicates that host plants can be seasonally limiting and that prediction of the occurrence of suitable host plants is a critical part of the model.

Clarke (2022) list 59 species of plant used by Pieris napi larvae in Europe (Table 4). However, this list also includes less preferred hosts, including crop species. In Europe, the first generation of butterflies must use hosts with earlier flowering, such as Alliaria petiolata in woodlands and shady places, Raphanus raphanistrum and Sisymbrium officinale from disturbed land and hedgerows and Cardamine pratensis in damp open meadows. Later in the season, late-flowering crucifers must support the larvae of the second generation (e.g. Berteroa incana, Cakile maritima and Rorippa sylvestris). These later flowering crucifers are more typical of open grass or disturbed habitats rather than woodland. Thus, as described by Van Driesche et al. (2004) in New England, there is probably a shift in larval food-plant availability from woodlands and hedgerows for the first generation to open habitats, meadows and disturbed areas for the second-generation larvae in the later summer. However, as noted by Hardy and Dennis (2010), food-plants in typically non-habitat area can also be exploited.

Eggs are normally laid singly on host plants, but when host plants are scarce, multiple eggs can be laid and plants already having eggs on them may be selected for laying (Courtney 1988). In P. rapae, females with a high egg load will lay in clusters and will also lay in clusters in low plant-density habitats (Jones 1977). The increased egg-load per plant in low densities was attributed to the female finding the same host repeatedly. Thus, it appears that there is a tendency in Pieris to lay more eggs when egg loads are high, but P. napi will only lay in small clusters when host density is low.

Table 3.

Thermal performance curve parameters from von Schmalensee et al. (2024) for diapause termination rate and post-diapause developmental rates.

MODEL Parameter Estimate (Mode) 90% HPDI
Lower Upper
Diapause termination rate TPC (day−1) Tmin (males) −6.18 −8.07 −2.88
Topt (males) 0.683 −0.892 2.44
Tmax (males) 31.3 27.1 35.8
Ropt (males) 0.00746 0.00668 0.00816
Tmin (females) −4.11 −8.65 −1.19
Topt (females) 1.6 −0.851 2.78
Tmax (females) 30.9 27 35.8
Ropt (females) 0.00634 0.00582 0.00701
Post-diapause development rate TPC (day−1) Tmin 1.99 1.27 2.77
Topt 29.6 29.1 30
Tmax 36.9 35.5 39
Ropt 0.152 0.148 0.158
SD* 0.0852 0.0751 0.0959
Table 4.

The host plants recorded in literature for Pieris napi in Europe, collated by Clarke (2022).

Plant Name Authorities
Alliaria petiolata (M.Bieb.) Cavara & Grande
Arabidopsis arenosa (L.) Lawalrée
Arabidopsis thaliana (L.) Heynh.
Arabis alpina L.
Arabis ciliata Clairv.
Arabis hirsuta (L.) Scop.
Arabis sagittata (Bertol.) DC.
Arabis soyeri Reut. & A.Huet
Aurinia saxatilis (L.) Desv.
Barbarea intermedia Boreau
Barbarea verna (Mill.) Asch.
Barbarea vulgaris W.T.Aiton
Berteroa incana (L.) DC.
Biscutella laevigata L.
Brassica napus L.
Brassica nigra (L.) W.D.J.Koch
Brassica oleracea L.
Brassica rapa L.
Cakile maritima Scop.
Cardamine amara L.
Cardamine flexuosa With.
Cardamine heptaphylla (Vill.) O.E.Schulz
Cardamine hirsuta L.
Cardamine impatiens L.
Cardamine pentaphyllos (L.) Crantz
Cardamine pratensis L.
Cardamine trifolia L.
Diplotaxis tenuifolia (L.) DC.
Draba aizoides L.
Erucastrum gallicum (Willd.) O.E.Schulz
Erysimum cheiranthoides L.
Hesperis matronalis L.
Iberis saxatilis L.
Isatis tinctoria L.
Lepidium campestre (L.) W.T.Aiton
Lepidium coronopus (L.) Al-Shehbaz
Lepidium draba L.
Lepidium graminifolium L.
Lobularia maritima (L.) Desv.
Lunaria annua L.
Lunaria rediviva L.
Nasturtium officinale W.T.Aiton
Noccaea caerulescens (J.Presl & C.Presl) F.K.Mey.
Noccaea rotundifolia (L.) Moench
Pseudoturritis turrita (L.) Al-Shehbaz
Raphanus raphanistrum L.
Rorippa amphibia (L.) Besser
Rorippa palustris (L.) Besser
Rorippa sylvestris (L.) Besser
Sinapis alba L.
Sinapis arvensis L.
Sisymbrium irio L.
Sisymbrium loeselii L.
Sisymbrium officinale (L.) Scop.
Thlaspi arvense L.
Turritis brassica Leers
Turritis glabra L.
Reseda lutea L.
Tropaeolum majus L.

Implementation in the model

We will assume a habitat-specific density of food-plants, which will be seasonally variable. Woodland spaces and hedgerows will have a higher density early in the season, with woodlands being reduced sharply in the summer (Type A in Table 5). Open meadows will have a low-density early season which will increase towards summer (Type B in Table 5). Any disturbed habitats will increase larval food-plant density to peak in Summer (Type C in Table 5). All habitats will decline in food-plant density in later summer towards zero in October. This will provide three basic curves for larval plant density and all habitats will be assigned to one of these curves or to a zero curve (Type D excluded in Table 5).

The values for the curves will be developed under a future calibration step, but are expected to broadly be represented by the profiles suggested in Fig. 1. However, for a particular country or region, the habitat type will also need a scaling factor and variability factor. This means that, for a specific habitat instance in the model, the curve will be scaled by multiplying with a value (vp) drawn from a Poisson distribution around the mean expected maximum food plant density. This is a starting assumption creating stochastic variation, but, in future, the variability could be mechanistically generated if more data become available. These habitat factors will all be specific parameter inputs for any model run. Thus, the food plant density fp at month m will be the value of the specific curve type (A–C) c for habitat type h, multiplied by a variance factor v generated from a Poisson distribution specific to the habitat patch p (Equation 7).

fp=hc,mvp Eq. 7

The oviposition rate is assumed to be a declining function with age following a log-normal distribution such that the area under the curve is equivalent to the total number of eggs the female has when laid over her projected lifespan.

Table 5.

Habitat classification in terms of four curve types. Each column provides the list of ALMaSS habitat types assigned initially to each curve. The habitat types are defined as Types of Landscape Elements (Topping and Duan 2024b) using the names provided here preceded by tole. The tole-types not mentioned in this table are assumed to be in the category of zero to no larval food-plants.

Early Peak (A) Intermediate (B) Summer (C)
Hedges Marsh RoadsideVerge
WoodlandMargin Heath Railway
ForestAisle NaturalGrassWet FieldBoundary
Copse PermPastureTussockyWet PermPastureLowYield
PermPastureTussocky
PermanentSetaside
PermPasture
NaturalGrassDry
PitDisused
YoungForest
HedgeBank
BeetleBank
RoadsideSlope
HeritageSite
UnknownGrass
Wasteland
FarmYoungForest
OFarmYoungForest
GameCover
OPermPasture
OPermPastureLowYield
Figure 1. 

Generalised curves showing the proportion of maximum larval food-plant density assumed for each of the three types of non-zero habitat curves.

Dispersal

Both sexes fly by day (Carter 1984). Age influences butterfly flight, decreasing flight endurance with age. Male butterflies fly for a longer period than females and flight endurance increases with temperature in both sexes (Åhman and Karlsson 2009). There is some indication from the changing shape of the wings with season that flight capabilities might alter with the butterfly generation. Early forms may be specialised towards local searching, whereas later forms are better suited to dispersal (Kleckova et al. 2023).

Studies on P. rapae movement (Jones 1977; Jones et al. 1980), indicated that female butterflies fly about 700 m per day over the ground and ends the day 250–600 m from where they started. Each butterfly maintains a preferred direction throughout one day, but the direction changes unpredictably from day to day. If host plants are found in an area, then the rate of turning increased. Dennis and Shreeve (1997) classified butterflies into different groups depending on their dispersal abilities. In their study, P. napi and P. rapae only differed in the breadth of semi-natural habitats used (P. napi utilises more) and the fact that P. napi was designated as having intermediate movement, whereas P. rapae mobility was considered to be wide-ranging, out of three possible movement classifications. Wood and Pullin (2002) studied genetic isolation in P. napi in the UK and found no significant structuring and a lack of evidence linking genetic similarity and geographic proximity of populations. They suggest that dispersal ability is not the key factor determining distributions of P. napi at the spatial scale studied (the West Midlands). Therefore, it seemed that the species was restricted in its distribution by the availability of suitable habitat, rather than the ability of a species to disperse between available habitat patches.

Ohsaki (1980) studied P. napi movements in Japan and concluded that, due to the permanent nature of the habitats and the narrow range of host species, this butterfly could be characterised as relatively localised and sedentary (compared to P. rapae). Assuming movement observations can be used to parameterise European P. napi, then Ohsaki (1980) considered P. napi, for which there were few observations, to be rather similar to P. melete. P. melete had dispersal distances of approximately 500 m for the females and up to 750 m for males over a 2 to 4-week period. Both sexes did not migrate extensively, but rather drifted with time. The dispersal distances of P. rapae in the studies by Ohsaki were, however, much lower than those recorded by other studies (Jones 1977; Jones et al. 1980). However, Fric et al. (2006) found that P. napi differed morphologically between Spring and Summer generations, with better dispersal potential in later adults. Shkurikhin and Oslina (2015) found that the shape of the wings of P. napi in early and late generations differed. Later-generation adults had wing shapes that were more capable of longer and energy-consuming flights. Unfortunately, the actual dispersal distances are not recorded for either morphometric study.

Within a generation, age has a significant influence on a butterfly’s flight ability. Flight endurance in P. napi increases with temperature in both sexes, but decreases significantly with age (Åhman and Karlsson 2009). The change in flight endurance of female butterflies appears not to be easily explained by the breakdown of flight muscles with age (Åhman and Karlsson 2009).

Implementation in the model

We include the weather effects on flight by including an hourly assessment of the weather including two thresholds for lower temperature and rainfall. Those hours in the daylight hours of the day where rainfall is below the threshold and temperature is above the temperature threshold will be summed and the ratio of flight hours to non-flight hours (fh) will be used to scale the overall distance moved. The two threshold values will be parameters fitted from the calibration process.

It seems clear that there is a seasonal and age tendency for older and later butterflies to move differently. This may also relate to changing habitat distributions of the host plants, requiring greater dispersal later in the season. We assume the discrepancy between the dispersal distances from the studies in Japan and others is related to this difference and that dispersal in the first generation is naturally lower than in the second. Dispersal must also be related to the oviposition behaviour below, which could modify the actual distance moved. We can formulate the dispersal ability as dmax (distance moved) as a function of season (s), age (a) and available flight hours (h):

dmax=md+c×s×la×fh Eq. 8

where md is a minimum distance, c is a constant, s is a factor represented by two constants for early and late season and la is a fraction decreasing with increasing age (a), fh is the ratio of the flight hours to non-flight hours. dmax represents the movement allowed for foraging and oviposition behaviour in a day. It is, however, the displacement from the start location, not the distance moved in a day.

Two other constructs are needed to represent oviposition and foraging behaviour in the model. These are maps of forage resources and maps of larval food plants. We assume these to be available and updated on a daily basis at the scale of 1 m2.

We plan to represent the movement using the following rules:

  1. Determine whether the stomach contains nectar. If not, foraging for nectar will be triggered; otherwise, egg laying will be triggered if there are eggs to lay.
  2. If egg-laying, evaluate the location for egg-laying potential by retrieving a value for the host plant density.
  3. If host-plant density is enough for egg-laying, then lay an egg, otherwise determine whether the stomach contains nectar.
  4. If there is no nectar in the stomach or no eggs to lay, search for forage; otherwise, move to a new host area (or stay if egg-laying is successful).
  5. If there are enough eggs and movement steps left, search for a host plant or nectar source and repeat the behaviour.

The flow diagram for this behaviour is described in Fig. 2. This behaviour requires the specification of two functions, one for host-density acceptance and one to determine if the stomach is empty. The latter is presumed to be testing the stomach, which will be replenished by foraging and depleted by flying distance. For energy use, we assume a constant rate of sugar consumption per metre travelled and use the closest linear distance between points travelled to and from.

The function for determining egg laying at a given host density is given by a probability function between two thresholds such that no eggs are laid below Ll, a minimum host-plant threshold density (e.g. 0) and the probability of laying eggs then increases linearly to an upper threshold Lh. The values for the probabilities at Ll and Lh would be given as Llp and Lhp and, initially, would be calibrated to derive distributions of egg-laying patterns representative of the field situation. The host density value to drive this function will come from Equation 7.

Figure 2. 

Overview of behavioural decisions related to foraging and egg-laying movement. Each move uses energy which can be replenished by foraging. Each move also uses up the distance allocation for the day, which, when used, will result in stopping movement for the day.

Mortality

Parasitoids

Parasitoids are an important part of the ecology of P. napi. For instance, Benson et al. (2003) postulated that a Braconid parasite caused extinction of a bivoltine population of the butterfly due to poor survival of the second generation as a result Cortesia glomerata. Parasitoids are considered to be very efficient at controlling populations of Pieris butterflies. Mustata and Mustata (1999) identified 26 parasitoid species, belonging to the Hymenoptera (Braconidae, Ichneumonidae, Chalcididae, Pteromalidae) and Diptera (Tachinidae) controlling populations of Pieris spp. P. napi was recorded as suffering from 59.4% parasitisation rate. These parasitoids also have hyperparasitoids which limit their ability to control butterfly populations. Patriche and Andriescu (2005) looked at Pieris brassicae and P. rapae. in Moldavia, identifying nine primary parasitoids, eight of P. rapae and five for P. brassicae, with C. glomerata only found in P. brassicae. Of these primary parasitoids, there were two species of hyperparasitoids in the P. brassicae complex and five in the P. napae complex.

C. glomerata was studied by Brodeur and Vet (1995) who measured encapsulation rates of parasitoid eggs in Pieris spp., including P. napi, showing the rates increased with instar age. C. glomerata identify lower instar caterpillars and preferentially attacks these (Brodeur et al. 1991; Mattiacci and Dicke 1995). The parasitoid is also able to identify suitable patches of hosts at a distance using infochemicals (Geervliet et al. 1998).

Modelling relationships between parasitoids and their hosts is complicated by the spatial dynamics inherent in their ecology. To represent this, we need to know the likelihood of a parasitoid entering a patch with prey or staying or leaving to find another patch. This was modelled by Waage (1979) as a deterministic process based on the encounter rate with hosts. This model considers the time spent in a patch as a function of the density of hosts, the number of ovipositions in the patch and a declining reinforcement of the result of finding a host. It was updated to include a greater level of stochasticity by Pierre et al. (2012). However, applying this model to the agent-based simulation is complex because of the spatial structure in ALMaSS, which does not represent patches, but is more akin to a contour map of densities. The model also does not include the age distribution of the hosts, time to death of host, nor the population dynamics of the parasitoids. Therefore, we propose the development of a general parasitoid model to dynamically link parasitoid and host behaviour at the individual level in the model. This will also allow the specification of hyperparasitoids using the same model. This approach simplifies the mathematics considerably and is highly configurable.

Implementation in the model

The parasitoid model will be developed as an animal population within ALMaSS. The initial specification is assumed to follow the agent-based paradigm rather than invoking the subpopulation modelling approach (Duan and Topping 2024), although this decision will be revisited if computational speed becomes limiting.

The model is specified as a very simple parasitoid model with random searching for hosts. Superparasitism caused by different individuals is allowed and hyperparasitism is possible by creating a second parasitoid species targeting the first. Each parasitoid will need to be created with its unique class to manage the population. The general parasitoid class can be defined with the following attributes and behaviour:

Attributes:

HostSpecies – the target host species on which the eggs are laid.

HostSpeciesMaxAge – the age of the host at which eggs are no longer laid.

DispersalDistance – the maximum distance allowed for dispersal between foraging locations.

MaxNumberEggs – the number of eggs the parasitoid can lay (from a mean with distribution).

SearchEfficiency – the ability to find hosts in an area as a proportion of the hosts present.

SearchArea – the area searched for hosts during one day.

SearchingTime – the time needed to search, inversely related to host density.

HandlingTime – the period of time from finding a host until the next host is found without searching time.

DevelopmentalDegrees – the number of day degrees needed for completion of development.

DevelopmentalDegreesThreshold – the threshold above which day degrees count towards development.

DailyMortality – a daily chance of dying if in the adult phase and not hibernating.

HibernationSuccess – the probability of surviving hibernation and emerging the following year.

EmergenceTriggerDegrees – The number of day degrees about the EmergenceTempThreshold and after EmergenceDaylength before the adult will break hibernation.

EmergenceDaylength – the shortest day length (sunrise-sunset) that is possible post-hibernation.

EmergenceTempThreshold – the threshold above which day degrees are counted.

ReproductivePreparation – the number of days between emergence and the start of oviposition.

HibernationTrigger – the day length at which any newlyemerged adults enter hibernation.

Types of behaviour:

B_Hibernation – The adults will be in dormancy until the day length reaches EmergenceDaylength. They will then count day degrees above EmergenceTempThreshold until EmergenceTriggerDegrees is reached, at which point they will emerge from hibernation and start a pre-oviposition period. This method confers the flexibility to control the emergence by day length, temperature or both. On emergence, a mortality chance will be applied, based on the value of HibernationSuccess. This should also include the removal of the males from the population if any have been modelled through the juvenile stages.

B_PreOviposition – After emergence, the adult will wait ReproductivePreparation days before commencing reproduction. We assume mates are not limiting and can be found in this period. The value for MaxNumberEggs is assigned. This value will be drawn from a Gaussian distribution with a mean and variance specified (MEggs & MEggsVar).

B_Reproduction – the main activity of the adult. The behaviour starts with a mortality test, based on daily mortality probability DailyMortality; if the adult dies, it is removed immediately. The adult will move DispersalDistance and assess an area with radius of DispersalDistance/2 m for hosts. The density of hosts will then determine the number of eggs laid. This will be a function of host density H × SearchEfficiency. This value will be limited by the day length divided by HandlingTime + 1/H (SearchingTime). The result will be a lower rate of parasitisation at lower densities and a limit to the number of oviposition events at higher densities. We will assume no superparasitism as a result of multiple eggs of the same individual in a single host, but this will not prevent superparasitism from different individuals. Once all MaxNumberEggs are laid, if the adult survives this long, then it dies.

B_Development – once laid in a host, the development is modelled as a single stage to emergence from the host, based on DevelopmentalDegrees, which is the target number of degrees above DevelopmentalDegreesThreshold. Once reached, the adult emerges and the host is killed, as will be any younger parasitoids resulting from superparasitism. On emergence, the day length is assessed to determine if hibernation should start.

B_Hibernation – This is a threshold day length at which point adults are assumed to enter a diapause state. This test is only taken on emergence from the host. The assumption is that adults enter diapause at their current location and, in this state, are effectively suspended from the simulation until emergence.

Pesticide responses

Here, we follow the basic approach used in development of ApisRAM honey bee colony model (Duan et al. 2022). Pesticide exposure can occur through the consumption of contaminated nectar, contact and overspray. By considering these three exposure pathways, we provide an accurate representation of how model butterflies may encounter pesticides in contaminated environments. PB (t) is used to represent the accumulated pesticide burden for the individual until day t. PB (t) is updated by all three exposure paths and is used to calculate the stress caused by pesticide exposure. The new pesticide exposure on day t is represented by PN (t) which is calculated by Equation 9.

PN(t)=PI(t)+PC(t)+PO(t) Eq. 9

where PI (t) is the new contamination from intake, PC (t) is from the contact exposure and PO (t) is from the overspray exposure.

Intake exposure pathway

This pathway accounts for the pesticides that are ingested by larval and adult butterflies. Foraging butterflies may be exposed by collecting nectar from contaminated flowers, whilst larvae may eat contaminated leaves. The pesticide amount in the consumed resource is represented by PI (t) and is assumed to move directly to the butterfly's body (PN (t)).

Overspray and contact exposure pathways

The overspray exposure pathways are only relevant to any butterfly stage present at the moment of pesticide spraying. Contact exposure can occur when a mobile stage moves in a contaminated area. In both cases, a simple one-time absorption model is used here, i.e. we assume that the pesticide will be absorbed into the butterfly’s body once, which is controlled by a parameter for contact and overspray separately.

Contact exposure specifically pertains to foraging adults and occurs when a butterfly touches a contaminated plant surface during its foraging activities. The amount of pesticide transferred to the body through contact is controlled by two user-defined parameters as shown in Equation 10.

PC=scSPs Eq. 10

where PS is the pesticide amount per square metre on the plant surface and S is a constant representing the area of the butterfly that is in contact with the surface. Here, we might imagine that foot, abdomen and wing contact could occur. Sc is the absorption rate that determines the amount of pesticide PS that will be transferred to the body. This implementation does not take account of time of foraging, but both parameters can be adjusted to make this approach more or less conservative.

Overspray exclusively occurs when butterflies are active in a field where simultaneous pesticide spraying is taking place and/or where drift is happening in the adjacent area. Two parameters, starting time and ending time, for a pesticide spray event are used to control whether overspray can happen. Overspray can only occur if a butterfly is in the exact location and when the pesticide is being sprayed. When the overspray happens, the contaminated pesticide amount (PO) going to the butterfly’s body is calculated using Equation 11.

PO=soSa Eq. 11

where S represent a proportion of the body surface, as for Equation 10, a is the pesticide spraying application rate and so is the absorption rate for overspray.

Toxicology

To determine the effect of the pesticide body burden, the initial model only considers acute mortality. To do this requires three parameters. The first is a threshold for effect, Pt, above which there is a daily probability of mortality Pm. The third threshold D is the daily decay rate of the pesticide in the butterfly’s body. This approach will implement a first model for implementing pesticide effects and assumes that an ALMaSS or equivalent landscape model is available to generate the pesticide concentration in the nectar and vegetation, as well as overspray.

Other mortalities

There are a number of other mortalities considered in the model, the main one being a daily background mortality rate for all life stages. This mortality would be a fitting parameter specific to each stage (BMs for s = egg to adult). This mortality is a daily probability applied at the start of each daily time step.

Mortalities will also be associated with management events. We assume that all soil cultivation and harvest events result in 100% mortality for non-adult stages. Cutting of grass likewise would kill non-adult stages. Other managements could be considered on a case-by-case basis as needed when applying the model.

Finally, we assume mortality of non-diapause stages will occur with the onset of winter (here we assume 1 December) or with negative temperatures lower than a minimum temperature threshold (MinT), initially assumed to be -10 °C.

Discussion

Much of the Pieris model here is focused on development and, therefore, temperature relationships. We are particularly lucky that this species was the subject of detailed studies supporting the development of thermal performance curves. As such, we believe that this part of the model is quite robust. This is important because there is evidence from butterfly surveys in the Netherlands that the flight period starts earlier in recent years (van Sway, pers. comm.), even though the peak of activity appears not to have moved. The driver for this is likely the temperature-related plasticity. For example, Bauerfeind and Fischer (2014b) showed that, under realistic conditions of fluctuating temperatures, a moderately increased temperature of three degrees decreased larval and pupal development time by approximately 30%. This effect can work both ways. von Schmalensee et al. (2023) noted that populations of P. napi in Sweden show better overwintering survival than P. rapae. This manifests itself as a disproportionately larger first-generation population, which provides the necessary resilience for the population to be able to adapt to univoltinism. Thus, shorter and longer phonologies are possible depending on the temperature regime. We have attempted to build in this type of flexibility to the model by relying on temperate-related mechanisms as the major phenological driver. These observed phenology patterns under different weather conditions will also be extremely useful for pattern-oriented model fitting.

Two areas were particularly challenging when developing the formal model for P. napi. These were the interactions with parasitoids and movement. There are many studies of the parasitoids of P. brassicae and P. rapae. Still, the different types of behaviour of the three species, particularly in the different egg-laying behaviour, are likely to render the data from P. brassicae of limited use. For example, Brodeur et al. (1998) showed that the same parasitoid species reared on three different Pieris species had significantly different survival rates. Comparing Cotesia glomerata, a P. brassicase parasitoid, with C. rubecula, which is a parasitoid of P. rapae, Wiskerke and Vet (1994) found that female parasitoids used different spatial searching strategies related to the egg-laying pattern of the host. Therefore, we have opted to design a flexible and extensible model for parasitoid interactions that can be tuned, based on the empirical information for P. napi. This model may also be useful for other species with a generic parasitoid and hyper-parasitoid hierarchy.

The movement functions are difficult to design well, given that they must integrate a range of types of behaviour at a daily time step. The approach taken was to create a looped set of dispersal, foraging and egg-laying movement behaviour and to limit these, based on both energy and a maximum distance moveable. At the calibration stage, these types of behaviour will need to be carefully assessed since emergent patterns will come from the parameters implemented in the butterfly model in conjunction with the detail implemented in the underlying landscape model. The landscape model will provide the patterns of resources in space and time and, thus, the level of detail in the two model components needs to be comparable. This method does, however, not include longer-distance dispersal events. This could be included by adding a long low tail to the calculation of dmax from Equation 8, but given that the landscape scale used in ALMaSS is typically 10 × 10 km, this was not included in the current configuration.

As noted in framing the model, this formal model ignores males and includes a toxicological model designed to be extended in future iterations. The model, as it is described currently, would be detailed enough for the current regulatory approach of single product regulatory evaluation. However, the focus for future development is on creating a systems view for environmental risk assessment (Topping et al. 2024). Therefore, future versions should be adapted to include multiple pesticides, non-lethal effects and more detailed internal toxicological dynamics.

Acknowledgements

Thanks to Chris van Sway who provided input to an early draft of the paper. PollinERA receives funding from the European Union’s Horizon Europe Research and Innovation Programme under grant agreement No.101135005. Views and opinions expressed are those of the author(s) only and do not necessarily reflect those of the European Union (EU) or the European Research Executive Agency (REA). Neither the EU nor REA can be held responsible for them.

Additional information

Conflict of interest

The authors have declared that no competing interests exist.

Ethical statement

No ethical statement was reported.

Funding

PollinERA receives funding from the European Union’s Horizon Europe Research and Innovation Programme under grant agreement No.101135005.

Author contributions

Conceptualization: XD, CJT. Formal analysis: CJT, XD. Funding acquisition: CJT. Writing – original draft: CJT. Writing – review and editing: CJT, XD.

Author ORCIDs

Christopher John Topping  https://orcid.org/0000-0003-0874-7603

Xiaodong Duan  https://orcid.org/0000-0003-2345-4155

Data availability

All of the data that support the findings of this study are available in the main text.

References

  • Balachowsky AS, Mesnil LP (1936) Insect Pests of Cultivated Plants.
  • Bauerfeind SS, Fischer K (2014b) Simulating climate change: temperature extremes but not means diminish performance in a widespread butterfly. Population Ecology 56: 239–250. https://doi.org/10.1007/s10144-013-0409-y
  • Benson J, Van Driesche RG, Pasquale A, Elkinton J (2003) Introduced braconid parasitoids and range reduction of a native butterfly in New England. Biological Control 28: 197–213. https://doi.org/10.1016/S1049-9644(03)00058-6
  • Bergström J, Wiklund C, Kaitala A (2002) Natural variation in female mating frequency in a polyandrous butterfly: effects of size and age. Animal Behaviour 64: 49–54. https://doi.org/10.1006/anbe.2002.3032
  • Boggs CL, Ross CL (1993) The effect of adult food limitation on life-history traits in Speyeria mormonia (Lepidoptera, Nymphalidae). Ecology 74: 433–441. https://doi.org/10.2307/1939305
  • Brodeur J, Vet LEM (1995) Relationships between parasitoid host-range and host-defense – a comparative-study of egg encapsulation in 2 related parasitoid species. Physiological Entomology 20: 7–12. https://doi.org/10.1111/j.1365-3032.1995.tb00794.x
  • Brodeur J, Geervliet JBF, Fac L (1991) Host species affecting the performance of the larval parasitoids Cotesia glomerata and Cotesia rubecula (Hymenoptera, Braconidae) .1. preference for host developmental stage of Pieris (Lepidoptera, Pieridae). In: International Symp on Crop Protection. Ghent, Belgium, 543–545.
  • Brodeur J, Geervliet JBF, Vet LEM (1998) Effects of Pieris host species on life history parameters in a solitary specialist and gregarious generalist parasitoid (Cotesia species). Entomologia Experimentalis et Applicata 86: 145–152. https://doi.org/10.1046/j.1570-7458.1998.00275.x
  • Capinera J (2020) Handbook of vegetable pests. Academic press, 729 pp.
  • Carter DJ (1984) Pest Lepidoptera of Europe: With Special Reference to the British Isles. Springer Netherlands.
  • Courtney S (1988) Oviposition on peripheral hosts by dispersing Pieris napi (L.) (Pieridae). Journal of Research on the Lepidoptera 26(1–4): 58–63. https://doi.org/10.5962/p.266701
  • Daniel TL, Kingsolver JG, Meyhofer E (1989) Mechanical determinants of nectar-feeding energetics in butterflies – muscle mechanics, feeding geometry, and functional equivalence. Oecologia 79: 66–75. https://doi.org/10.1007/BF00378241
  • Dennis RLH, Shreeve TG (1997) Diversity of butterflies on British islands: ecological influences underlying the roles of area, isolation and the size of the faunal source. Biological Journal of the Linnean Society 60: 257–275. https://doi.org/10.1111/j.1095-8312.1997.tb01495.x
  • Dover JW (1989) The use of flowers by butterflies foraging in cereal field margins. Entomologist’s Gazette 40: 283–291.
  • Duan X, Topping CJ (2024) A general subpopulation model for the Animal Landscape and Man Simulation System (ALMaSS). Food and Ecological Systems Modelling Journal 5: e122467. https://doi.org/10.3897/fmj.5.122467
  • Ebert G, Rennwald E (1993) Die Schmetterlinge Baden-Württembergs / Band 1: Tagfalter I [Allgemeiner Teil: Systematik, Taxonomie und Nomenklatur, Faunistik und Ökologie, Gefährung und Schutz, Datenverarbeitung / Spezieller Teil: Papilionidae, Pieridae, Nymphalidae]. Vol. 1,Verlag Eugen Ulmer.
  • Emmet AM, Heath J (1990) The Moths and Butterflies of Great Britain and Ireland. Harley Books. [ISBN 0-946589-37-2]
  • Forsberg J (1987) Size discrimination among conspecific hostplants in 2 pierid butterflies – Pieris napi L and Pontia daplidice L. Oecologia 72: 52–57. https://doi.org/10.1007/BF00385044
  • Fric Z, Klimova M, Konvicka M (2006) Mechanical design indicates differences in mobility among butterfly generations. Evolutionary Ecology Research 8: 1511–1522.
  • Geervliet JBF, Ariëns S, Dicke M, Vet LEM (1998) Long-Distance Assessment of Patch Profitability through Volatile Infochemicals by the Parasitoids Cotesia glomerata and C. rubecula (Hymenoptera: Braconidae). Biological Control 11: 113–121. https://doi.org/10.1006/bcon.1997.0585
  • Greiser C, von Schmalensee L, Lindestad O, Gotthard K, Lehmann P (2022) Microclimatic variation affects developmental phenology, synchrony and voltinism in an insect population. Functional Ecology 36: 3036–3048. https://doi.org/10.1111/1365-2435.14195
  • Hardy P, Dennis R (2010) A butterfly exploiting the matrix: Pieris napi (Linnaeus, 1758) (Lepidoptera: Pieridae) ovipositing amongst mown grass in a city park. Entomologist’s Gazette 61: 155–158.
  • Jones R, Hart J, Bull G (1982) Temperature, Size and Egg Production in the Cabbage Butterfly, Pieris rapae L. Australian Journal of Zoology 30: 223–231. https://doi.org/10.1071/ZO9820223
  • Jones RE (1977) Movement patterns and egg distribution in cabbage butterflies. Journal of Animal Ecology 46: 195–212. https://doi.org/10.2307/3956
  • Jones RE, Gilbert N, Guppy M, Nealis V (1980) Long-Distance Movement of Pieris rapae. Journal of Animal Ecology 49: 629–642. https://doi.org/10.2307/4268
  • Karlsson B, Johansson A (2008) Seasonal polyphenism and developmental trade-offs between flight ability and egg laying in a pierid butterfly. Proceedings of the Royal Society B – Biological Sciences 275: 2131–2136. https://doi.org/10.1098/rspb.2008.0404
  • Kleckova I, Linke D, Rezende F, Rauscher L, Le Roy C, Matos-Maraví P (2023) Flight behaviour diverges more between seasonal forms than between species in Pieris butterflies. bioRxiv. https://doi.org/10.1101/2023.11.27.568806
  • Lees E, Archer D (1974) Ecology of Pieris napi (L.) (Lepidoptera, Pieridae) in Britain. Entomologist’s Gazette 25: 231–237.
  • Lehmann P, Pruisscher P, Posledovich D, Carlsson M, Käkelä R, Tang P, Nylin S, Wheat CW, Wiklund C, Gotthard K (2016) Energy and lipid metabolism during direct and diapause development in a pierid butterfly. Journal of Experimental Biology 219: 3049–3060. https://doi.org/10.1242/jeb.142687
  • Lehmann P, Van der Bijl W, Nylin S, Wheat CW, Gotthard K (2017) Timing of diapause termination in relation to variation in winter climate. Physiological Entomology 42: 232–238. https://doi.org/10.1111/phen.12188
  • Lehmann P, Pruisscher P, Koštál V, Moos M, Šimek P, Nylin S, Agren R, Väremo L, Wiklund C, Wheat CW, Gotthard K (2018) Metabolome dynamics of diapause in the butterfly Pieris napi: distinguishing maintenance, termination and post-diapause phases. Journal of Experimental Biology 221: jeb169508. https://doi.org/10.1242/jeb.169508
  • Lobry J, Rosso L, Flandrois J-P (1991) A FORTRAN subroutine for the determination of parameter confidence limits in non-linear models. Binary 3: 25.
  • Meyrick E (1927) A Revised Handbook of British Lepidoptera. Watkins and Doncaster, London, 914 pp.
  • Mustata G, Mustata M (1999) The parasitoid complex controlling Pieris populations in Moldavia-Romania. In: Entomologists Conference. Basel, Switzerland, 337–341.
  • Niitepold K, Perez A, Boggs CL (2014) Aging, life span, and energetics under adult dietary restriction in Lepidoptera. Physiological and Biochemical Zoology 87: 684–694. https://doi.org/10.1086/677570
  • Ohata M, Furumoto A, Ohsaki N (2011) Hatchling morphology and adaptation to host plants during juvenile stages in the butterfly Pieris napi. Ecological Research 26: 59–66. https://doi.org/10.1007/s11284-010-0758-3
  • Ohsaki N (1980) Comparative population studies of 3 Pieris butterflies, Pieris rapae, Pieris melete and Pieris napi, living in the same area.2. Utilization of patchy habitats by adults through migratory and non-migratory movements. Researches on Population Ecology 22: 163–183. https://doi.org/10.1007/BF02513543
  • Okamura Y, Tsuzuki N, Kuroda S, Sato A, Sawada Y, Hirai MY, Murakami M (2019) Interspecific differences in the larval performance of Pieris butterflies (Lepidoptera: Pieridae) are associated with differences in the glucosinolate profiles of host plants. Journal of Insect Science 19: 1–9. https://doi.org/10.1093/jisesa/iez035
  • Patriche G, Andriescu I (2005) The hyperparasitoid complex which limits the action of the primary parasitoids of the pieridae species (Insecta: Lepidoptera), defoliators in cabbage crops. Conference Proceedings. https://api.semanticscholarorg/CorpusID:54807481
  • Pierre J-S, Masson J-P, Wajnberg E (2012) Patch leaving rules: A stochastic version of a well-known deterministic motivational model. Journal of Theoretical Biology 313: 1–11. https://doi.org/10.1016/j.jtbi.2012.08.003
  • Posledovich D, Toftegaard T, Wiklund C, Ehrlén J, Gotthard K (2015) Latitudinal variation in diapause duration and post-winter development in two pierid butterflies in relation to phenological specialization. Oecologia 177: 181–190. https://doi.org/10.1007/s00442-014-3125-1
  • Ratkowsky DA, Reddy GVP (2017) Empirical model with excellent statistical properties for describing temperature-dependent developmental rates of insects and mites. Annals of the Entomological Society of America 110: 302–309. https://doi.org/10.1093/aesa/saw098
  • Ratte HT (1985) Temperature and Insect Development. In: Hoffmann KH (Ed.) Environmental Physiology and Biochemistry of Insects. Springer Berlin Heidelberg, Berlin, Heidelberg, 33–66. https://doi.org/10.1007/978-3-642-70020-0_2
  • Rosso L, Lobry JR, Flandrois JP (1993) An unexpected correlation between cardinal temperatures of microbial growth highlighted by a new model. Journal of Theoretical Biology 162: 447–463. https://doi.org/10.1006/jtbi.1993.1099
  • Shkurikhin AO, Oslina TS (2015) Seasonal phenotypic plasticity of the polyvoltine white butterfly Pieris napi L. (Lepidoptera: Pieridae) in the Southern Urals. Russian Journal of Ecology 46: 96–102. https://doi.org/10.1134/S1067413615010178
  • Shkurikhin AO, Oslina TS (2016) Seasonal variation of the forewing in polyvoltine whites Pieris rapae L. and P. napi L. (Lepidoptera: Pieridae) in the forest-steppe zone of the Southern Urals. Russian Journal of Ecology 47: 296–301. https://doi.org/10.1134/S1067413616030115
  • Stepanova LA (1962) An experiment in the ecological analysis of the conditions for the development of pests of cruciferous vegetable crops in nature. USSR Review of Entomology 41: 721–736.
  • Topping CJ (2022) The Animal Landscape and Man Simulation System (ALMaSS): a history, design, and philosophy. Research Ideas and Outcomes 8: e89919. https://doi.org/10.3897/rio.8.e89919
  • Topping CJ, Duan X (2024a) Managing large and complex population operations with agent-based models: The ALMaSS Population_Manager. Food and Ecological Systems Modelling Journal 5: e117593. https://doi.org/10.3897/fmj.5.117593
  • Topping CJ, Duan X (2024b) ALMaSS Landscape and Farming Simulation: software classes and methods. Food and Ecological Systems Modelling Journal 5: e121215. https://doi.org/10.3897/fmj.5.121215
  • Topping CJ, Hansen TS, Jensen TS, Jepsen JU, Nikolajsen F, Odderskaer P (2003) ALMaSS, an agent-based model for animals in temperate European landscapes. Ecological Modelling 167: 65–82. https://doi.org/10.1016/S0304-3800(03)00173-X
  • Topping CJ, Kondrup Marcussen L, Thomsen P, Chetcuti J (2022) The Formal Model article format: justifying modelling intent and a critical review of data foundations through publication. Food and Ecological Systems Modelling Journal 3: e91024. https://doi.org/10.3897/fmj.3.91024
  • Topping CJ, Bednarska AJ, Benfenati E, Chetcuti J, Delso NS, Duan X, Focks A, Laskowski R, Lombardo A, Marcussen LK, Metodiev T, Rubinigg M, Rundlöf M, Sgolastra F, Stoyanova C, Sušanj G, Williams JH, Ziółkowska EM (2024) PollinERA: Understanding pesticide-Pollinator interactions to support EU Environmental Risk Assessment and policy. Research Ideas and Outcomes 10: e127485. https://doi.org/10.3897/rio.10.e127485
  • Van Driesche RG, Nunn C, Pasquale A (2004) Life history pattern, host plants, and habitat as determinants of population survival of Pieris napi oleracea interacting with an introduced braconid parasitoid. Biological Control 29: 278–287. https://doi.org/10.1016/S1049-9644(03)00152-X
  • von Schmalensee L, Caillault P, Gunnarsdóttir KH, Gotthard K, Lehmann P (2023) Seasonal specialization drives divergent population dynamics in two closely related butterflies. Nature Communications 14: 3663. https://doi.org/10.1038/s41467-023-39359-8
  • von Schmalensee L, Hulda Gunnarsdóttir K, Näslund J, Gotthard K, Lehmann P (2021) Thermal performance under constant temperatures can accurately predict insect development times across naturally variable microclimates. Ecology Letters 24: 1633–1645. https://doi.org/10.1111/ele.13779
  • von Schmalensee L, Süess P, Roberts KT, Gotthard K, Lehmann P (2024) A quantitative model of temperature-dependent diapause progression. Proceedings of the National Academy of Sciences 121: e2407057121. https://doi.org/10.1073/pnas.2407057121
  • Waage JK (1979) Foraging for patchily-distributed hosts by the parasitoid, Nemeritis canescens. Journal of Animal Ecology 48: 353–371. https://doi.org/10.2307/4166
  • Wiklund C, Kaitala A, Lindfors V, Abenius J (1993) Polyandry and its effect on female reproduction in the green-veined white butterfly (Pieris napi L.). Behavioral Ecology and Sociobiology 33: 25–33. https://doi.org/10.1007/BF00164343
  • Wiklund C, Kaitala A, Wedell N (1998) Decoupling of reproductive rates and parental expenditure in a polyandrous butterfly. Behavioral Ecology 9: 20–25. https://doi.org/10.1093/beheco/9.1.20
  • Wiskerke JSC, Vet LEM (1994) Foraging for solitarily and gregariously feeding caterpillars – a comparison of 2 related parasitoid species (Hymenoptera, Braconidae). Journal of Insect Behavior 7: 585–603. https://doi.org/10.1007/BF01997434
  • Wood BC, Pullin AS (2002) Persistence of species in a fragmented urban landscape: the importance of dispersal ability and habitat availability for grassland butterflies. Biodiversity & Conservation 11: 1451–1468. https://doi.org/10.1023/A:1016223907962
  • Ziółkowska E, Bednarska AJ, Laskowski R, Topping CJ (2023) The Formal Model for the solitary bee Osmia bicornis L. agent-based model. Food and Ecological Systems Modelling Journal 4: e102102. https://doi.org/10.3897/fmj.4.102102.figure10
login to comment