Experimental conditions

The biotic conditions, particularly those related to the experimental organisms in both controlled and open-field environments, are critical in shaping host–parasitoid interactions. A fundamental consideration is the selection of realistic host densities, as functional response analysis depends entirely on host availability. For instance, egg parasitoids often encounter host patches with egg numbers that exceed their daily parasitism capacity. Telenomus busseolae, a specialist pro-ovigenic egg parasitoid, provides an example of this challenge. As a specialist, it relies heavily on its target host, Sesamia spp. (Lepidoptera: Noctuidae). It efficiently detects and parasitizes host eggs, with minimal self-superparasitism due to host-marking pheromones. Moreover, due to its high reproductive potential, it is capable of laying more than 50% of its total eggs (approximately 80 under laboratory conditions) within the first two days of its lifespan (Agboka et al., Reference Agboka, Schulthess, Chabi-Olaye, Labo, Gounou and Smith2002; Olaye et al., Reference Olaye, Schulthess, Shanower and Bosque-Pérez1997). Given that Sesamia spp. egg batches can reach up to 200–300 eggs during peak oviposition periods, far exceeding the early-life parasitism capacity of T. busseolae (Dimotsiou et al., Reference Dimotsiou, Andreadis and Savopoulou-Soultani2014), experiments should be tailored to reflect these conditions occurring in field conditions (Jamshidnia et al., Reference Jamshidnia, Kharazi-Pakdel, Allahyari and Soleymannejadian2010). A similar principle applies to generalist egg parasitoids, such as Trichogramma spp., which parasitize the eggs of a wide range of lepidopteran host species across varying host densities in natural conditions. Studies indicate that using ecologically relevant host densities during the parasitoids’ peak parasitism periods improves the accuracy of functional response characterisation (de Oliveira and Reigada, Reference de Oliveira and Reigada2023; Tonğa, Reference Tonğa2024; Tonğa et al., Reference Tonğa, Erkek, Ali, Fathipour and Özder2024). These studies typically resulted in either a Type II or Type III functional response. Focusing solely on host densities that do not surpass daily parasitism capacity, under the assumption that handling time is negligible and does not constrain host attack constants may lead to incomplete interpretation of functional response models (Fellowes et al., Reference Fellowes, van Alphen, Shameer, Hardy, Wajnberg and Jervis2023). For example, the parasitism rates of Trichogramma species were better explained by a Type I functional response as host egg densities submitted to them did not exceed their daily parasitism potential (Kalyebi et al., Reference Kalyebi, Overholt, Schulthess, Mueke, Hassan and Sithanantham2005). The authors argued that while a Type I functional response may not fully represent the parasitoid’s theoretical potential for controlling pests with high reproductive rates, it is still valuable for evaluating the control potential of these parasitoids against pests with smaller egg batches (Kalyebi et al., Reference Kalyebi, Overholt, Schulthess, Mueke, Hassan and Sithanantham2005). This suggests that both lower and higher host densities can reveal important aspects of host–parasitoid dynamics, depending on the targeted pest species (Fellowes et al., Reference Fellowes, van Alphen, Shameer, Hardy, Wajnberg and Jervis2023). For generalist egg parasitoids, variation in densities of different host species may provide insights into their host–parasitoid interactions. In contrast, for specialist parasitoids, functional response models derived from experiments with both lower and higher host egg densities can perform better. Therefore, inclusion of widely spanning host densities helps interpret more realistic Type II and Type III models. However, a clear distinguish between Type II and Type III models require sufficient number of replications in lower host densities since underestimation of the effect of lower host densities may mislead decision-making between Type II and III models (Uszko et al., Reference Uszko, Diehl and Wickman2020).

The same considerations extend to Braconidae, parasitoids of aphids, which encounter colonies with varying host densities under natural conditions. In such cases, experimental constraints on host density may not be necessary (Byeon et al., Reference Byeon, Tuda, Kim and Choi2011; Talebi et al., Reference Talebi, Kazemi, Rezaei, Mirhosseini and Moharramipour2022). However, for parasitoids attacking larval, pupal, or adult stages of lepidopteran pests, host availability differs significantly. These life stages are often solitary, in small groups, concealed, or mobile. Consequently, since synovigenic parasitoids may not reach their full ovigeny potential during foraging, ovigeny alone may not determine functional response outcomes (Dannon et al., Reference Dannon, Tamò, van Huis and Dicke2010). Instead, the rate of egg maturation and, in some cases, adult host feeding can also influence the results, indicating that host availability is the key component of functional response evaluations and should be carefully considered. Accordingly, host aggregation patterns also influence experimental designs. For example, Dipteran larvae often exhibit aggregative behaviour and tend to pupate in close proximity, increasing the probability of parasitoids encountering a broad range of host densities (Heaton et al., Reference Heaton, Moffatt and Simmons2018; Khan et al., Reference Khan, Khuhro, Awais, Memon and Asif2020).

The majority of functional response experiments use a 24-hour duration for practical reasons. This timeframe is often sufficient for comparative assessments of idiobiont and some koinobiont parasitoids, but it does not fully capture field conditions, particularly for koinobionts with longer lifespans. Thus, a critical factor is the time required for the development of an adequate number of parasitoid eggs in response to host density. Traditional functional response experiments often use short time frames, which can introduce bias, especially for synovigenic parasitoids which require time to mature their eggs before they are fully prepared for effective foraging (Aguirre et al., Reference Aguirre, Logarzo, Triapitsyn, Diaz-Soltero, Hight and Bruzzone2024; Griffen, Reference Griffen2021). In contrast, proovigenic parasitoids, which emerge with their full complement of eggs, are immediately capable of foraging, making them suitable for short-term experiments (24 h). However, many parasitoids exhibit intermediate reproductive strategies, being neither fully synovigenic nor proovigenic (Jervis et al., Reference Jervis, Heimpel, Ferns, Harvey and Kidd2001). In such cases, as well as in synovigenic species, additional factors, such as non-foraging time budgets, must be incorporated into functional response models to ensure accurate evaluations (Aguirre et al., Reference Aguirre, Logarzo, Triapitsyn, Diaz-Soltero, Hight and Bruzzone2024). Extending experiments beyond 24 h can help egg maturation processes, but only if the risks of host depletion and superparasitism are carefully managed (Tonğa, Reference Tonğa2024; Tonğa et al., Reference Tonğa, Erkek, Ali, Fathipour and Özder2024).

Another challenge in functional response experiments is how to address host depletion and superparasitism, both of which can skew estimates of parasitism dynamics. In egg parasitoids, host quality often deteriorates over time, while larval growth can alter the suitability of hosts for larval parasitoids. One approach has been to replace hosts to maintain constant densities, but this can introduce artificial conditions. Modern analytical tools allow declining host densities to be modelled directly (Bolker, Reference Bolker2008; Rosenbaum and Rall, Reference Rosenbaum and Rall2018), reducing the need for replacement to avoid complications during host replacement. For generalist parasitoids, counting host depletion without replacement is recommended to account for parasitism dynamics accurately (Hassell et al., Reference Hassell, Lawton and Beddington1977), whereas specialist parasitoids often rely on host-marking pheromones to avoid superparasitism, making them ideal model organisms for functional response studies particularly compared to those that do not account for host replacement (Agboka et al., Reference Agboka, Schulthess, Chabi-Olaye, Labo, Gounou and Smith2002; Bruce et al., Reference Bruce, Schulthess, Makatiani and Tonnang2021). This is primarily because the lower frequency of re-encountering hosts during a parasitoid’s search can be more accurately represented by models that incorporate host replacement dynamics (Juliano, Reference Juliano, Scheiner and Gurevitch2001). Whether host replacement produces more realistic evaluations of parasitoid behaviour remains an open question that should be considered carefully.

The examples discussed above serve as guidance for considering host density, reproductive strategy, and species of parasitoids in functional response studies. However, the complexity of host–parasitoid interactions means that experimental designs will continue to evolve, particularly with advances in data science and modelling approaches. Further refinement in experimental methodologies will be essential to improve comparative analyses and predictive modelling.

Model selection

The first step in functional response analysis is model selection, which requires statistical validation by investigators (Fellowes et al., Reference Fellowes, van Alphen, Shameer, Hardy, Wajnberg and Jervis2023). Traditionally, researchers characterised functional response types based solely on curve shapes. Statistical significance tests, such as logistic regression with p-values, remain widely used and provide a consistent framework for classification. At the same time, modern approaches increasingly emphasise model-based inference, including more probabilistic approaches fitting Hill exponents directly, which are advantageous when significance testing falls short in application (DeLong, Reference DeLong2021; DeLong et al., Reference DeLong, Coblentz and Uiterwaal2025). A strategic starting point is to begin model selection fitting polynomial logistic regression which estimates key parameters from the proportion of parasitized hosts as a function of host density: intercept (/constant, P0), linear (P1), hyperbolic (/quadratic, P2), and sigmoid (/cubic, P3) (Eq. 1 in table 1) (Juliano, Reference Juliano, Scheiner and Gurevitch2001).

Table 1. The equations along with their definitions and characteristic functions mentioned in current study

where N_a is the number of parasitized host eggs, N ₀ is the initial host egg density and N_a/N ₀ is the proportion of parasitized host eggs, P ₀–P ₃ are the estimates corresponding to intercept, linear, quadratic and cubic parameters, respectively, T is the experimental duration, a is the attack rate, b is the attack coefficient, h is the handling time.

The assessment procedure in this manuscript is performed in R statistical software (R Core Team, 2025). The estimates of polynomial logistic regression as a function of parasitism proportion (binary scale from 0 to 1) at different densities are obtained using a generalised linear model (GLM) with a binomial distribution and logit link function using MASS package (Hardy and Smith, Reference Hardy and Smith2023; Ripley et al., Reference Ripley, Venables, Bates, Hornik, Gebhardt and Firth2023; Venables and Ripley, Reference Venables and Ripley2002). The significant P1 and P2 parameters allow characterisation of functional response type (Juliano, Reference Juliano, Scheiner and Gurevitch2001).

One criterion to separate Type II and III functional responses by analysing the proportion of parasitized host is to test for significant positive or negative linear coefficients in the expression fit by the method of binomial GLM to data on proportion parasitism versus N_O (initial host density). The significant negative linear coefficient (P1) revealed Type II (the proportional parasitism decreases as the host egg density increases) functional response whereas a significant positive linear (P1) coefficient accompanied by a negative quadratic (P2) coefficient revealed Type III (initially increasing proportional parasitism is followed by a decreasing tendency as the host egg density increases) functional response (fig. 1) (Juliano, Reference Juliano, Scheiner and Gurevitch2001).

Figure 1. Typical characterisation of functional response curves for Type II (red curve) and Type III (blue) depicting differential relationship in host–parasitoid models as a function of proportional host parasitism.

Several examples of previously published datasets are used to guide researchers in the evaluation of the functional response. The example codes and model summaries included in this paper are derived from the statistical evaluation of various parasitoid strains and species, as documented in our earlier publications (Tonğa, Reference Tonğa2024; Tonğa et al., Reference Tonğa, Erkek, Ali, Fathipour and Özder2024).

The first example performs a clear functional response model with a Type II curve as the linear coefficient is significantly negative (Box 1, P < 0.001).

Box 1: Statistics of polynomial logistic regression established as a generalised linear model with a binomial distribution clearly demonstrating a Type II functional response.

In another example (Box 2), the model generates a negative linear parameter which could easily be referred to as Type II functional response as well while quadratic and cubic parameters are insignificant (Box 2, P > 0.05). In this case, a researcher could explore opportunities for reduced models. When the cubic parameter was discarded from the model and the analysis was re-performed, the model approvingly provided a negative linear parameter (Type II) (Box 2, reduced model, P < 0.001). The consistency in AIC value also confirms the performance of the model.

Box 2: Statistics of polynomial logistic regression established as a generalised linear model with a binomial distribution clearly demonstrating a Type II functional response after model reduction. The outputs for ‘model’ and ‘reduced_model’ are presented in this box.

In another example, a dataset resulting in complete set of significant coefficients including negative quadratic parameter accompanied by a positive linear parameter, reveals a Type III functional response model (Box 3, P < 0.001).

Box 3: Statistics of polynomial logistic regression established as a generalised linear model with a binomial distribution clearly demonstrating a Type III functional response.

The non-significant P values of coefficients in full polynomial logistic regression models can introduce bias, indicating that reduced models could be valuable alternatives while evaluations require careful consideration. Although models with lesser numbers of coefficients (restricted to linear and quadratic terms) are considered promising, experts clearly advise a more focused investigation in data modelling (Pritchard et al., Reference Pritchard, Paterson, Bovy and Barrios-O’Neill2017). As a promising approach, the frair_test() function from the package frair generates the sign and significance of first-order coefficient (linear) and second-order coefficient (quadratic) terms in the polynomial logistic regression equation as typical automatically reduced models (Pritchard et al., Reference Pritchard, Paterson, Bovy and Barrios-O’Neill2017). This approach has been promising to provide a basic statistical test to differentiate whether a dataset fits to a Type II or a Type III functional response. The authors of the package frair suggest that this approach offers a phenomenological approach to response shape characterisation, but mechanistic models remain necessary for deeper ecological insights (Barrios‐O’Neill et al., Reference Barrios‐O’Neill, Dick, Emmerson, Ricciardi and MacIsaac2015; Pritchard et al., Reference Pritchard, Paterson, Bovy and Barrios-O’Neill2017). In an alternative approach, a generalised form of functional response involving a scaling exponent (q) parameter can help distinguish between Type II and III models, where q = 0 represents a Type II and q > 0 refers to a Type III functional response (Kalinkat et al., Reference Kalinkat, Rall, Uiterwaal and Uszko2023; Pritchard et al., Reference Pritchard, Paterson, Bovy and Barrios-O’Neill2017; Rosenbaum and Rall, Reference Rosenbaum and Rall2018). Therefore, we recommend that researchers explore detailed methodologies to avoid over-reliance on direct reduced models. Such reliance can influence the direction of model selection, potentially introducing constraints in functional response data analysis. Reduced models may impact the model selection process and introduce bias in subsequent steps, particularly in the exploration of relevant parameters. In cases where full vs. reduced models present ambiguities, Akaike Information Criterion (AIC) comparisons can aid in selecting the most appropriate functional response model. A thorough understanding of this approach is essential for accurate functional response estimation.

Parameter estimation

Estimating the functional response curve provides valuable insights into parasitoid–host interactions. However, once the type of functional response is determined and similar curves are observed across different parasitoid groups, further investigation is necessary to elucidate the specific differences and similarities between these groups. A Type II functional response, characterised by a decelerating parasitism rate, can destabilise host–parasitoid population dynamics by inducing inverse density-dependent host mortality (Hassell, Reference Hassell1978). In contrast, a Type III functional response, which involves density-dependent host mortality, may contribute to stabilising these dynamics (Fellowes et al., Reference Fellowes, van Alphen, Shameer, Hardy, Wajnberg and Jervis2023; Murdoch and Oaten, Reference Murdoch and Oaten1975).

The second step in functional response analysis involves data modelling and parameter estimation which can be carried out using various fitting algorithms. For example, non-linear least square has been widely used (Trexler et al., Reference Trexler, McCulloch and Travis1988). Bayesian approaches are quite advanced allowing full posterior distributions of parameters and model selection (Papanikolaou et al., Reference Papanikolaou, Kypraios, Moffat, Fantinou, Perdikis and Drovandi2021). More frequently, maximum likelihood estimation (MLE) is applied to better quantify uncertainty and exploit likelihood‐based inference. This method optimises the relationship between the number of parasitized hosts and the initial host density at the start of the experiment (table 2; Bolker, Reference Bolker2008; Bolker et al., Reference Bolker, Team and Giné-Vázquez2023; Juliano and Williams, Reference Juliano and Williams1987). At this stage, this optimisation on the basis of arbitrary probability distributions, allows parameter estimation in both Type II and Type III functional responses, depending on the selected model. This process is implemented using the mle2 function of the bbmle package (table 2; Bolker et al., Reference Bolker, Team and Giné-Vázquez2023). For a Type II functional response, where the proportion of parasitized hosts declines monotonically with increasing host density, the attack rate (a) and handling time (h) parameters were estimated using the Roger’s random predation model (Eq. 2 in tables 1 and 2). This model accounts for host depletion, assuming no host replacement during the experimental period (Rogers, Reference Rogers1972).

Table 2. R codes employed to derivate the respective parameters of equations listed in Table 1 as well as their definitions and functions

Hymenopteran parasitoids predominately exhibit Type II and Type III functional responses, requiring appropriate model selection. For Type II functional response, Holling’s disc equation is the classical model that estimates the attack rate (a) and handling time (h) under the assumption that prey (or hosts) are continuously replaced (Holling, Reference Holling1959). However, in experiments where hosts are not replaced, as is typical in parasitoid studies where each host can only be parasitized once or the host matures, host depletion occurs over time. In such cases, the Rogers random predator equation is more appropriate, as it explicitly accounts for the declining number of available hosts during the exposure period (table 2, Eq. 2). The solutions to these equations were derived using the Lambert W function below (table 2, Eq. 3; Bolker, Reference Bolker2008; Rosenbaum and Rall, Reference Rosenbaum and Rall2018). This approach allows for an explicit solution of the functional response equation without requiring iterative numerical methods, facilitating faster and more accurate parameter estimation in empirical studies.

For Type III functional responses, characterised by a sigmoidal pattern, where parasitism initially increases at low host densities before reaching an asymptote at higher densities, the attack coefficient is treated as a function of initial host density. To ensure consistency of theoretical approach with Type II functional response models, the Hassell Type III model was employed as an alternative to classical Type III model that describes that host depletes over time and provides an analytical correct solution for estimated parameters (Eq. 4 in table 1; Hassell et al., Reference Hassell, Lawton and Beddington1977; Rosenbaum and Rall, Reference Rosenbaum and Rall2018).

Which model fits best?

Selecting the best-fitting functional response model can be challenging, as multiple equations may be applicable. Visual inspection of data points, as a subjective method, unsystematically fitted curves, limited host availability, and high variability in observations, in particular at lower host densities, can all complicate the decision-making process (Kalinkat et al., Reference Kalinkat, Rall, Uiterwaal and Uszko2023). For example, in case of high variance, a Type II functional response could be more advantageous than generalised models with three parameters while increasing the number of observations can improve statistical power and model support to establish a Type III model (Marshal and Boutin, Reference Marshal and Boutin1999; Novak and Stouffer, Reference Novak and Stouffer2021). Moreover, although a model may exhibit a high R², this metric alone may not fully reflect model’s ability to capture the parasitoid’s responses at low and high host densities. These effects at different density levels require a deeper understanding of host-handling mechanisms (Okuyama, Reference Okuyama2012, Reference Okuyama2024). To address these challenges, additional tools may be used. For example, a recent study has used the Deviation Information Criterion (DIC) index for model selection, which can help identify well-fitting models by balancing explanatory power and parameter complexity (Aguirre et al., Reference Aguirre, Logarzo, Triapitsyn, Diaz-Soltero, Hight and Bruzzone2024). Yet, there are further approaches to find the best fitting equations. The Akaike Information Criterion (AIC) is particularly useful for comparing models with similar functional response types (Tonğa, Reference Tonğa2024). Both AIC and DIC, belonging to a family of information-theoretic indices, evaluate the trade-off between model fit and complexity, helping decide which of several candidate equations (e.g., different parameterisations of Type II or Type III models) provides the best balance between explanatory power and parsimony. A combined approach using above-mentioned tools to choose the best-fitting model seems convenient, although its effectiveness depends on the evaluator’s statistical expertise.

Parameter comparison

Parameter estimation in functional response analysis produces two key parameters. First, it should provide an attack rate or attack coefficient depending on the type of functional response which has recently been explained as the space (area or volume) containing host that is effectively cleared of the host by the attacker and named as space clearance rate for predators (DeLong, Reference DeLong2021). In Type II functional responses for parasitoids, attack rate (a = b) represents the number of hosts successfully parasitized as a function of host density. In Type III functional responses, the attack coefficient (b) is density-dependent and derived as a/N₀ = b with addition of a Hill exponent (DeLong, Reference DeLong2021). Since these parameters represent fundamentally different processes, direct comparisons between Type II and Type III parasitoids are not valid (Fellowes et al., Reference Fellowes, van Alphen, Shameer, Hardy, Wajnberg and Jervis2023). However, recent advances in the field demonstrates that the derivation of attack rate (space clearance rate therein) for an alternative Type III model, having the same biological meaning for Type II models, is possible with integration of a variety of functions that produce an asymptotic a (DeLong, Reference DeLong2021). Instead, comparisons are generally restricted to the parameters within the same functional response type. For example, the attack rate should be compared between parasitoids exhibiting Type II responses while the attack coefficient should be compared between parasitoids exhibiting Type III responses. Confidence intervals can be used to assess differences – non-overlapping CI’s indicate significant differences in attack parameters, but overlapping 95% CI’s may not necessarily refer to lack of statistical significance (Cumming, Reference Cumming2009). Handling time can be compared across different functional response types, as it remains consistent across models (Box 4).

Box 4: Maximum likelihood estimation model output for attack rate and handling time, demonstrating the use of confidence intervals in distinguishing parasitoid groups. The non-overlapping confidence intervals significantly separates two parasitoid groups in handling time, while overlapping confidence intervals may not significantly differ two parasitoid groups in attack rate.