Predator size affects the intensity of mutual interference in a predatory mirid

Abstract Interference competition occurs when access to an available resource is negatively affected by interactions with other individuals, where mutual interference involves individuals of the same species. The interactive phenomena among individuals may be size‐dependent, since body size is a major factor that may alter prey consumption rates and ultimately the dynamics and structure of food webs. A study was initiated in order to evaluate the effect of mutual interference in the prey‐specific attack rates and handling times of same size class predators, incorporating variation in consumer size. For this purpose, laboratory functional response experiments were conducted using same age predators, that is, newly hatched (first instar) or mature (fifth instar) nymphs of the polyphagous mirid predator Macrolophus pygmaeus preying on Ephestia kuehniella (Lepidoptera: Pyralidae) eggs. The experiments involved four predator density treatments, that is, one, two, three, or four predators of same age, that is, either first‐ or fifth‐instar nymphs, which were exposed to several prey densities. The Crowley–Martin model, which allows for interference competition between foraging predators, was used to fit the data. The results showed that mutual interference between predator's nymphs may occur that affect their foraging efficiency. The values of the attack rate coefficient were dependent on the predator density and for the first‐instar nymphs were significantly lower at the highest predator density than the lower predator densities, whereas for the fifth‐instar nymphs in all density treatments were significantly lower to that of the individual foragers' ones. These results indicate that mutual interference is more intense for larger predators and is more obvious at low prey densities where the competition level is higher. The wider use of predator‐dependent functional response models will help toward a mechanistic understanding of intraspecific interactions and its consequences on the stability and structure of food webs.


| INTRODUC TI ON
Competition for resources is common in food webs (Kratina et al., 2009;Zhang et al., 2015), driving species distribution and evolution (Abrams, 2000;Craine & Dybzinski, 2013). Competitive interactions between individuals are assumed to be of two discrete types (Begon et al., 1996): exploitation competition, when species compete for the same limited resource, and interference competition in which one individual actively interferes with another individual's access to a resource. Both types of competition have negative effects mainly on the fitness of the weaker competitors (Lang & Bendow, 2013;McCormick & Weaver, 2012). Furthermore, mutual interference involves direct interactions between individuals of the same species that adversely may affect their searching efficiency, that is, the rate that a predator search for its prey, and ultimately its feeding rate (Beddington, 1975;DeLong & Vasseur, 2011;Hassell, 1978).
Recently, Papanikolaou et al. (2016) showed that mutual interference may alter the feeding rate of aphidophagous coccinellid species, depending on predator density exposed in a patch.
Understanding the relationship between per capita predator consumption in increasing prey density, that is, the functional response (Solomon, 1949), is crucial for describing and estimating the feeding interactions between a consumer and a resource (see, e.g., Abrams, 1989;Englund et al., 2011;Jeschke et al., 2002;Sentis, Hemptinne, & Brodeur, 2012 and essentially the cornerstone for most quantitative ecological studies (Moffat et al., 2020;Zhang et al., 2018). Holling (1959a) proposed various types of functional response; in type I and II, the prey consumption is assumed to increase linearly with prey density or increase asymptotically, respectively, reaching a plateau, while in a type III functional response, the prey consumption is supposed of a sigmoid form as prey density increases.
The well-known disk equation, which describes type II functional responses, assumes that a single predator spends its time on two kinds of activities: prey searching and prey handling (Holling, 1959b); prey searching is depicted on the attack rate coefficient of the disk equation, representing the per capita prey eaten by a predator at low prey densities, that is, the densities that predator is not satiated, and determining the initial slope of the functional response curve; prey handling refers to the handling time coefficient, that is, the time a predator spends pursuing, subduing and eating a prey item, indicating the per capita maximum feeding rate at a given time. However, the disk equation is considered as a purely prey-dependent functional response model because it did not incorporate parameters to assess behavioral interactions between foragers that may affect per capita functional responses due to competition. Nevertheless, intraspecific interactions act on predator search efficiency, and thus, predator dependence may facilitate predation or may reduce the potential of the interacting predators to locate and consume a prey.
For example, Kratina et al., (2009) reported high encounter rates between individuals of the benthic flatworm Stenostomum virginanum when fed on Paramericum aurelia which reduced the available foraging time for the predator.
Extending Holling's assumptions, more realistic functional response modeling approaches have been developed to predict predation rates when the predator and the prey density vary. These models are considered as predator-dependent (e.g., Beddington, 1975;Crowley & Martin, 1989;DeAngelis et al., 1975).
Although the viability of these models has been shown (Skalski & Gilliam, 2001), their application remains underutilized in the ecological literature, and thus, our knowledge on how density dependence affects the linkages between predators and their prey resources remains little elaborated (e.g., Kratina et al., 2009;Papanikolaou et al., 2016).
Apart of predator density, the interactive phenomena among predators may be size-dependent, since body size is a major factor that may alter prey consumption rates and ultimately the dynamics and structure of food webs (Christensen, 1996;Fuiman, 1989;Gravel et al., 2013;Mittlebach, 1981;Tsai et al., 2016;Yodzis & Innes, 1991). Similarly, body size may alter predator foraging capacity and may lead to changes in competitive ability and predator-prey dynamics (Hin & de Roos, 2018;Yodzis & Innes, 1991). For example, Thorp et al. (2018) reported that predator body size alters the type of functional response of the African clawed frog, Xenopus laevis (Anura: Pipidae), to Culex pipiens (Diptera: Culicidae); small predators exhibited a type II response, while medium and large predators type III responses. Therefore, evaluation of density-dependent effects from conspecific individuals on size-dependent predation may offer valuable information in understanding the factors involved in predator-prey interactions.
Omnivorous heteropteran predators are well known as high value naturally occurring natural enemies in several agroecosystems (i.e., Schaefer & Panizzi, 2000). Their important contribution to integrated pest management and biological control strategies is generally recognized, due to their predation efficiency and distribution in several habitats (Coll & Guershon, 2002;Fantinou et al., 2009;Molla et al., 2011;Perdikis et al., 2011;Pérez-Hedo & Urbaneja, 2015). Μacrolophus pygmaeus, (Rambur) (Hemiptera Miridae) a generalist predator of whiteflies, aphids, mites, and several lepidopteran species including Tuta absoluta (Meyrick) (Lepidoptera Gelechiidae), is native to the Mediterranean area and is commonly used in pest management (Biondi et al., 2018;Chailleux et al., 2017;Perdikis & Lykouressis, 2000). Fantinou et al. (2008), Fantinou et al. (2009) showed that Μ. pygmaeus exhibited a type II functional response on each of the nymphal instars of M. persicae. Moreover, the presence of floral resources reduced the number of aphids consumed and therefore the functional response plateau (Maselou et al., 2014).  Lampropoulos et al. (2013) reported a negative, antagonistic effect between two 5th-instar nymphs of M. pygmaeus at intermediate but not at high levels of Trialeurodes vaporariorum (Westwood) (Hemiptera: Aleyrodidae) nymph's availability. However, foraging effort and hunting behavior of the predator were not strongly affected by the presence of another conspecific at intermediate densities of Myzus persicae nymphs (Sulzer) (Hemiptera: Aphididae) (Maselou et al., 2014). Michaelides et al. (2018) showed that interactive effects on consumption of eggs of T. absoluta of two 5th-instar nymphs of M. pygmaeus were more intense at high prey densities. Adult body weight of predatory mirids was significantly reduced when nymphs raised in groups of 2, 4, 8, and 16 individuals, in the presence of prey, rather than reared in isolation (Arvaniti et al., 2019). Therefore, the above evidence supports that interference may exist between nymphs of M. pygmaeus.
In this study, we searched whether mutual interference interacts with predator body size. For this purpose, laboratory functional response experiments were conducted to test the effect of various predator-prey density combinations on feeding rate of similar sized individuals of either the 1st-or the 5th-instar nymphs of M. pygmaeus. To achieve this, the effect of predator density, prey density, and predator body size on conspecific interference was explored by utilizing specialized functional response models to estimate parameters that show the trends of interference and its intensity among the variable traits tested.

| Functional response data
The experimental setup consisted of Petri dishes (9 cm diameter, 1.5 cm height) with a mesh-covered hole in the lid (3 cm diameter) to reduce the accumulation of humidity. A tomato leaflet (6.45 ± 0.06 cm in length and 3.34 ± 0.13 cm in width) was placed, abaxial surface up, on a layer of water-moistened cotton wool on the bottom of each Petri dish. In all of the experiments, 1st-or 5th-instar nymphs of Μ. pygmaeus were used, <24 hr of age. Therefore, we expect interference interactions to be symmetric, that is, it affects all individuals the same. The 1st-instar nymphs were collected from stems with eggs and 5th-instar nymphs were obtained from first-instar nymphs that were transferred from wood-framed rearing cages to cages with potted tomato plants with food at 25°C, 65 ± 5% RH, and 16-hr light per day, and left to develop until the 5th instar. Then, 1st-or 5th-instar nymphs were introduced to caged tomato plants and were deprived of prey for 24 hr prior the beginning of the experiments to exclude the influence of variable hunger levels. Single and multiple predator treatments were conducted under controlled conditions of 25°C, 65 ± 5% RH and 16-hr light per day. One, two, three, and four predators from each size class were introduced into a dish with the tomato leaflet and the prey. As prey, eggs of E. kuehniella were used (EPHEScontrol ® , Agrobio S.L.,). In each of 1st-instar predator treatment, the prey densities were instars; however in the case of 5th instars, an additional prey density was included. In all cases, the prey/predator ratio was kept unchanged among the treatments. After 12 hr, the predators were removed from the dishes and the number of eggs consumed was recorded. There were performed 10-16 replicates at each prey density.

| Statistical analyses
As the dependent variable is binomial (each egg prey was eaten or not at the end of the experiments), a logistic regression analysis on the proportion of egg prey eaten as a function of the initial egg density was conducted to determine the shape of the functional response (Trexler et al., 1988). A polynomial function (Juliano, 2001) was fitted using the method of maximum likelihood: , where N e denotes the total number of eggs eaten, N 0 the initial egg number available, and P 0 , P 1 , P 2 , and P 3 are the intercept, linear, quadratic, and cubic coefficients. If the linear coefficient is significantly negative, the proportion of egg prey eaten declines monotonically, suggesting a type II functional response; a significantly positive linear coefficient followed by a significantly negative quadratic coefficient suggesting a type III functional response (Juliano, 2001).
In all treatments, the linear coefficient was significantly negative (Tables 1, 2). Therefore, as the data indicated a type II functional response, the Crowley-Martin model (Crowley & Martin, 1989) was used to fit the data: where N denotes the prey density; P the predator density; a the predator's attack rate, that is, the per capita prey mortality at low prey densities; T h the handling time which reflects the time, a predator spends on pursuing, subduing, eating, and digesting its prey; b the rate of encounter of a single predator with other predators; and t w the time wasted on an encounter. Although there are several models that allow for interference competition, we choose the Crowley-Martin model as in our recent study (Papanikolaou et al., 2016) it fitted several functional response data slightly better compared to the Beddington-DeAngelis (Beddington, 1975;DeAngelis et al., 1975), Holling type II (Holling, 1959b) and Hassell-Varley (Hassell & Varley, 1969) models.
The parameters b and t w are structurally nonidentifiable since they always appear as a product. Therefore, we grouped them into one parameter c which is a positive constant describing the magnitude of interference among predators (Kratina et al., 2009;Papanikolaou et al., 2016;Skalski & Gilliam, 2001). Hence, the Crowley-Martin model become: In case of a single predator (P = 1), the Crowley-Martin model is equal to the disk equation. Fitting was performed using the Bayesian approach suggested by Papanikolaou et al. (2016). In this task, a prior distribution is assigned to the parameters to express our prior beliefs about them before seeing the data. This information is then updated in the light of experimental data using Bayes theorem by multiplying it with the likelihood leading to the posterior distribution which contains all the information about the parameters. We used vague priors for most of the model parameters; typically, exponential distributions with very high variance (e.g., 1,000,000) to reflect our ignorance about them and allowed the posterior distribution of the parameters to be mostly informed by the data. However, we used the posterior distribution of the attack rate parameter, obtained from fitting the Holling model and by approximating it with a Gamma distribution (using the method of moments) as an informative prior, in order to overcome issues of practical nonidentifiability.
Our inference procedure was based upon Markov chain Monte Carlo methodology. In particular, we employed a random-walk Metropolis algorithm to draw samples from the posterior distribution of the parameters. The variance of the proposal distribution was tuned in order to achieve an acceptance rate of 25% (Roberts et al., 1997). All statistical analyses were conducted in R (R Core Team, 2019).

| RE SULTS
Macrolophus pygmaeus first-instar nymphs exhibited type II functional response, as the linear coefficient of the polynomial function was significant negative in all predator combinations (Table 1)    In case of fifth-instar nymphs, the data indicated anew a type II functional response for all the experiments (Table 3)   three, and four predators were used, the mean handling times were 0.542, 0.507, 0.605, and 0.534 hr, respectively.
The magnitude of interference (parameter c) among predators did not include zero values in 95% credible intervals in all cases.
This further indicates that predation of Macrolophus pygmaeus follows the predator-dependent assumptions of the Crowley-Martin model.

| D ISCUSS I ON
Our results showed density dependence of the functional response of M. pygmaeus. Density increase of the predator demonstrated a decrease in the estimated attack rate, and therefore the per capita searching predator efficiency, which indicates a limitation of its predation ability at low prey densities. An aspect that has received considerable attention in functional response studies is the limit to the number of prey attacked at high prey densities. It has been considered that handling time is the limiting factor for parasitoids and satiation for predators (Getz & Mills, 1996;Jeschke et al., 2002).
Digestion and handling prey are discrete biological processes, as digestion influences the predators hunger level and its probability of searching for prey (Jeschke et al., 2002;Papanikolaou et al., 2014).
It is likely that M. pygmaeus is a digestion-limited predator, that is, digestion limits its predation efficiency. Thus, at low prey densities M. pygmaeus is not satiated, so that digestion breaks do not exist, and predators constantly search for prey. Therefore, we anticipate that at these densities, frequent encounters among predators distract individuals from prey searching activity. In contrast, at high prey densities the predator searching time is relatively low, and consequently the time loss due to encounters is much reduced. In addition, as individuals become satiated at these high prey densities, the time lost during digestion breaks may partially or fully accommodate the cost of interference; if digestion breaks proceed during interference interactions, the time cost of interference may be negligible (van Gils & Piersma, 2004).
Our data showed that M. pygmaeus first-instar nymphs were less susceptible to interference competition than 5th instars as shown from the comparable values of attack rate to that of individual foragers. In the case of fifth-instar nymphs, interference was evident since the values of attack rate at all foraging densities were significantly lower than that for single predators, and furthermore, significantly reduced with the increase of predator density. Thereafter, these results indicate a much stronger mutual interference among foragers of the fifth than the first nymphal instar. These results support the hypothesis that increase of predator size intensifies the phenomena of mutual interference among conspecific predators.
Larger insect nymphs tend to show higher mobility (Hagstrum & Subramanyam, 2010), increasing the probability for encounter with competitors, so that they may be disoriented by the foraging activity.
In the outcome of mutual interference, besides predator density and size, arena/prey patch size may be a key component determining the intensity of interactions (Papanikolaou et al., 2016).
We expect that the encounter rate among conspecifics is reversely associated with the patch size (Dostalkova et al., 2002;Kindlmann & Dixon, 2001). Therefore, a predator density threshold may exist above which conspecific interactions become significant for a given patch size. This hypothesis is supported by comparing our results in experiments conducted in same size arenas for both instar nymphs.
Therefore, a predator dependence effect was recorded only at the highest density of 1st-instar nymphs (i.e., 4 individuals), but in contrast, this effect emerged already when two 5th-instar nymphs were enclosed together. In a previous study, the feeding rate of a predatory coccinellid was modified at low prey densities but only above a critical predator density (Papanikolaou et al., 2016).
A critical aspect of the evaluation of biological control agents is their response to increasing prey density. Most studies are focused on individual's functional response, which do not account for competition between individuals (e.g., Cabral et al., 2009;Jalali et al., 2010;Lee & Kang, 2004;Papanikolaou et al., 2011). In fact, the application of predator-dependent functional response models may offer valuable outcomes, however is underutilized so far in the ecological literature, partly due to their complexity. Fitting such models to experiments can be challenging due to issues of practical identifiability, which arise when the experimental data do not provide enough information about the parameters of interest (Papanikolaou et al., 2016). The Crowley-Martin model shows statistically significant improvement over the purely prey-dependent Holling's disk equation, as the 95% credible intervals of the parameter that account for interference competition did not include zero. This is in accordance with Skalski and Gilliam (2001), as they suggest the use of predator-dependent functional response models over the Holling's disk equation. Their wider use in predator-prey interactions will help toward a better understanding of intraspecific interactions and its consequences on the stability and structure of food webs. In this task, our study shows that mutual interference is dependent on predator body size, with further impact on the attack rate. We expect this size dependence of the searching efficiency of foraging predators to affect predator-prey dynamics (e.g., Aljetlawi et al., 2004;De Roos et al., 2003;Hassell & May, 1973). The differences in competitive ability between different sized predators and the importance of the attack rate on stabilizing predator-prey dynamics have been highlighted by Persson et al. (1998), as it may result in different types of dynamics (stable, chaotic, etc.). Therefore, the quantitative assessment of mutual interference in predator-prey systems could lead to a further understanding of their dynamics.
The results of the current work deliver useful insights to be considered in the mass rearing of predators such as M. pygmaeus and their use in biological control. Firstly, prey densities should be kept at high levels particularly for larger nymphs rearing; otherwise, competition is increased reducing prey consumption rates with possible negative effects in the development and reproduction of the predator. Although much smaller in size than 5th-instar nymphs, first-instar nymphs' predation rates are negatively affected by conspecific interactions, too.
Apart from prey risk reduction, intraspecific competition can lead predator individuals to prey switching by using other food resources for which less competition exists for their use (Svanbäck & Bolnick, 2007). Considering that optimal foraging theory supposes that predator may add alternative prey to its diet as a result of prey depletion (Stephens & Krebs, 1986), we expect that mutual interference can drive predatory mirids to alternative diets. According to our results, this should occur at low prey densities, and its intensity should increase with predator density, as a result of the competitive interactions among individuals. In such situations, phenotypic variation is likely to force individuals of a cohort to alternative food resources (Svanbäck & Bolnick, 2007). Searching these effects in case of omnivorous predators, in addition to alternative prey, alternative plant food resources must also be considered, as they may alter intraspecific interactions (Maselou et al., 2015).
In conclusion, our study highlights the importance of mutual interference on the predator functional response of a predatory omnivorous predator. It also provides evidence that predator size regulates mutual interference. These findings indicate that studying intraspecific interactions needs to utilize functional response modeling and take realistic levels of predator dependence into account. Future studies should investigate these effects under more complex environments, for example, field resembling conditions, or study the intensity of these effects in the presence of alternative food resources.

ACK N OWLED G M ENTS
NEP acknowledges financial support from the Foundation for Education and European Culture (2016-2019).

CO N FLI C T O F I NTE R E S T
The authors declare no conflicts of interest.