A Stochastic Differential Equation Model for Predator-Avoidance Fish Schooling

This paper presents a system of stochastic differential equations (SDEs) as mathematical model to describe the spatial-temporal dynamics of predator-prey system in an artificial aquatic environment with schooling behavior imposed upon the associated prey. The proposed model follows the particle-like approach where interactions among the associated units are manifested through combination of attractive and repulsive forces analogous to the ones occurred in molecular physics. Two hunting tactics of the predator are proposed and integrated into the general model, namely the center-attacking and the nearest-attacking strategy. Emphasis is placed upon demonstrating the capacity of the proposed model in: (i) discovering the predator-avoidance patterns of the schooling prey, and (ii) showing the benefit of constituting large prey school in better escaping the predator's attack. Based on numerical simulations upon the proposed model, four predator-avoidance patterns of the schooling prey are discovered, namely Split and Reunion, Split and Separate into Two Groups, Scattered, and Maintain Formation and Distance. The proposed model also successfully demonstrates the benefit of constituting large group of schooling prey in mitigating predation risk. Such findings are in agreement with real-life observations of the natural aquatic ecosystem, hence confirming the validity and exactitude of the proposed model.


Introduction
Swarm dynamics, one of the most commonly observed phenomena in the real world, have attracted the interest of many researchers from diverse fields including biology, mathematics, and computer engineering [1,2,3,4,5,6].As one manifestation of collective behavior in nature, swarm dynamics offer many advantages to the associated individuals: higher foraging success [7,8,9], improve protection upon predators [8,10,11,12,13,14], increase mating chance [8,10], and more efficient energy consumption due to hydrodynamic or aerodynamic advantages [8,9,10].Based on observational and empirical investigations of the interaction between congregating animals and their neighbors, researchers have constructed artificial life [15], for example, mathematical models of schooling prey in the presence of a predator.
As far as mathematical modeling of predator-prey interaction with schooling behavior upon the associated prey is concerned, two different approaches of modeling are routinely considered.The first approach examines the predator-prey system from the standpoint of tracking population densities of the associated entities.In this framework, mathematical models are typically constructed by integrating functional responses into the general description (see, for example, [16,17,18,19,20,21,22,23]).The second approach of modeling examines the predator-prey system as selfpropelled units in an artificial ecosystem.Here, the associated entities are treated as agents (or particles) that move around and interact with each other perpetually in the prescribed domain.Hereinafter, the discussion concerning mathematical modeling of the predator-prey system follows the latter approach.
Let us review some existing studies.Oboshi et al. [24] proposed a mathematical model to describe predator-prey interaction in an aquatic ecosystem.Therein, they designed an artificial habitat where schooling prey fish coexist with solitary predator.
They used a local rule that each prey fish in the school chooses one way of action among four possibilities according to the closest mate to construct a model of difference equations for anti-predator schooling behavior.Their model is able to return single predator-avoidance maneuver that fits real observation at the Port of Nagoya Public Aquarium.Such an evasive maneuver is also observed in [12,25,26,27,28].
Later on, the same authors [29] extended their assessment by incorporating three different conditions for prey-predator encounters in the artificial habitat, namely predator faces center, side, and vicinity of the schooling prey fish.By assimilating such configurations into their model, they managed to obtain two additional anti-predator patterns.
Nishimura [30] investigated mechanisms for a predator to select an individual prey from a group.Therein, the author introduced a priority function as pertinent parameter determining the selectivity of predator upon a particular prey in the group.
The study is later extended [31] by including three different priority functions, each corresponds to the selection of nearest, peripheral, and split individuals among the schooling prey.It is shown that, from the predator's viewpoint, the selection of peripheral individual is the best of the three available options: such a strategy yields the highest average number of successfully captured prey.
Zheng et al. [32] considered the evasive maneuver of schooling prey fish as complex collective traits resulting from the combination of schooling, cooperative, and selfish behavior.Three mathematical models of prey fish behavior are constructed, each corresponds to a fundamental trait.A combination of such models yields one antipredator pattern which is similar to the one obtained by Oboshi et al. [24].Lee et al. [33] implemented molecular dynamics (MD) to model the dynamical behavior of schooling prey in response to predator's attack.They retrieved one predator-avoidance pattern of the schooling prey, namely the crescent-shaped maneuver.In a separate study [34], the transitional regimes associated with such a maneuver are elaborated further.
Demšar and Lebar Bajec [35] studied three hunting tactics of the predator (attack center, nearest, and isolated prey) in association with two distinct behavior of the prey (schooling and solitary) using a mathematical model based on fuzzy logic.They found that compared to individualistic behavior, schooling is the optimal defense mechanism to avoid predation.Later on, Demšar et al. [14] extended the assessment by including composite hunting tactics of the predator.Two composite hunting tactics were proposed: (i) the predator may choose one of three simple hunting strategies (attack nearest, center, and peripheral prey) in successive attacks based on probability, and (ii) the predator was set to initially disperse the schooling prey and then pursue an isolated individual.They concluded that predator's confusion plays an important role upon the evolution of the composite tactics.Moreover, they conveyed that when confusion is taken into account, the latter composite tactic is the better predation mode that provides favorable outcome for the predator.
Despite the vibrant research activity in the field of mathematical modeling of predator-prey system, up to our knowledge, there is still no mathematical model that provides not only one but all real predator-avoidance patterns mentioned above.Moreover, there is also no mathematical model which shows the benefit of making a large school of fish in surviving predator's attack.In this paper, we aim to construct a mathematical model based on stochastic differential equations (SDEs) for predatoravoidance fish schooling which can show not only simulated patterns adapting to the real patterns but also the benefit of making a large school of fish in surviving predator's attack.The SDE model is chosen as it is well-suited for capturing complex movements and interaction rules of individuals in diverse environment.Furthermore, an SDE model is amenable for undertaking numerical simulations in computer platform.
In [36,37], we introduced a system of stochastic differential equations to describe the process of fish schooling in a free, unbounded environment based on a set of behavioral rules of fish schooling proposed by Camazine et al. [38].The rules are stated as follows: (a) The school has no leaders and each fish follows the same behavioral rules.
(b) Each fish uses some form of a weighted average of the position and orientation of its nearest neighbors to decide where to move.
(c) There is a degree of uncertainty in the individual's behavior that reflects both the imperfect information-gathering ability of a fish, and the imperfect execution of the fish's actions.
However, these rules do not consider obstacle-avoiding mechanisms and the foraging process of a school of fish.In accordance with this, in [6,39], we proposed behavioral rules for both phenomena, respectively.Furthermore, using the model in [36], we constructed two mathematical models for obstacle-avoiding mechanisms and the foraging behavior of fish schools.As a result, four obstacle-avoiding patterns are discovered.In addition, the probability of foraging success is observed to escalate up to an optimal value of school size, then gradually decreases as the school size gets larger.
In order to construct the mathematical model for predator-avoidance fish schooling in the present paper, we integrate two hunting strategies of a predator into our model in [36].The first strategy is that the predator always moves to the center of the prey fish school and catches the neighboring fish while moving.The second one is that it chases the "nearest" prey fish.For the latter strategy, the predator chases all prey fish in the group with imposed weighted coefficients that are inversely proportional to the distance from it to the corresponding individuals.Therefore, in this study, the "nearest" attacking strategy means that the predator pursues all the prey fish in the school while giving utmost attention to its nearest prey.
The paper is organized as follows.In Section 2, we introduce our SDE model for predator-avoidance fish schooling and explain how we construct mathematical expressions for the two hunting strategies and integrate them into the general model.
We also explain "being eaten" conditions for the prey fish.Section 3 provides four predator-avoidance patterns obtained from the model which fit real-life antipredatory maneuvers.We call them Split and Reunion, Split and Separate into Two Groups, Scattered, and Maintain Formation and Distance.In addition, pseudocodes to obtain these patterns are also given.Section 4 demonstrates the capability of our model in describing the benefit of making a large school of fish for better surviving the predator's attack-an advantage that is particularly well-recognized in real-life situation-by means of probability of survival of prey fish and the average time of prey being alive when the school size increases.Finally, Section 5 poses the concluding remarks of the present study.

Model Description
In this section, we newly introduce an SDE model to describe the dynamics of schooling prey under the attack of a predator in an obstacle-free domain R d (d = 2 or 3).In constructing the model, we adopt the particle-like approach; namely, each of the associated units in the predator-prey system is treated as if they were particles moving in space R d where interactions among them are manifested in the forms of fundamental forces analogous to the ones occurred in molecular physics [5,40].The proposed model is formulated as follows: (2.1) Here, N depicts the population of the schooling prey; x i (t) and v i (t) (i = 1, 2, . . ., N ) respectively denote the position and velocity of the i-th prey at time t; y(t) and v(t) correspondingly represent the position and velocity of the predator at time t, while • specifies the Euclidean norm of a vector.
As can be seen in the system (2.1), the proposed model consists of two parts: (i) a group of expressions that describe the dynamics of the schooling prey (the first two equations) and (ii) a group of expressions that depict the behavior of the predator (the last two equations).Let us discuss each of these groups of the model in more detail.

The Schooling Prey
The dynamics of the schooling prey are described by the first two equations in (2.1).
These expressions are extended versions of the ones developed in our previous study [36].The first equation is an SDE that tracks the position x i (t) of the i-th prey fish.It contains a stochastic term σ i dw i (t), considered as a white noise that reflects uncertainty in an individual's behavior within the school due to some kind of defect in the process of information-gathering among the individuals and imperfect execution of an individual's action.
Subsequently, the second equation in (2.1) solves the velocity v i (t) of the i-th prey fish.Therein, α and β are positive coefficients; parameters 1 < p < q < ∞ denote positive exponents; and r > 0 is the critical distance between two individuals in the school.Together with the Euclidean distance x i − x j , the overall term − r p x i −x j p − r q x i −x j q (x i − x j ) determines interaction between two prey individuals i and j in the form of attraction or repulsion.If the distance exceeds r, namely x i −x j > r, then they would move towards each other due to the induced attractive force between them.Contrarily, if x i − x j r, then the interaction would be the repulsive force so as to avoid collision.
To account for the effect of the predator's presence upon the behavior of the schooling prey, the term H(x i , y) is introduced within the expression of v i (t).Overall, such a remark manifests repulsive interaction between each individual prey and the predator.In this study, we propose the following mathematical formula for Here, R 1 > r is a constant that reflects the threshold distance where an individual prey begin to display the anti-predation maneuver to avoid the approaching predator; δ and θ 1 are positive coefficients.The strength of repulsion between the schooling prey and the predator depends upon the distance between them, x i −y .The closer the distance, the stronger the induced repulsive force between them.

The Predator
Let us now turn to discuss the expressions that represent the behavior of the predator in our model.Analogous to the prey counterpart, the position y(t) of the predator is described by the third equation in (2.1).Therein, the term σdw t represents uncertainty in predator's action due to the imperfect decision-making during hunting activity.
The velocity of predator v(t) in (2.1) is formulated based on generalized Newton's second law of motion.Therein, the resultant force term F (x i , v i , y, v) can be regarded as an attractive driving force that draws the predator towards the schooling prey.
Apart from its conventional definition in Newton's second law of motion, in our model, the force term F (x i , v i , y, v) is specifically designed to mimic the hunting strategy adopted by the predator.In this regard, we propose two hunting tactics for the predator, namely (i) the center-attacking and (ii) the nearest-attacking strategy.
Each of these tactics occupied its own mathematical formula for F (x i , v i , y, v).

Hunting Strategy I: Predator Attacks the Center of the Schooling Prey
In an aquatic ecosystem, this hunting strategy is mainly exemplified by swordfish (Xiphias gladius) [14,41].Observation [41] revealed that a similar hunting tactic is also adopted by other predatory fish species such as tuna (Thunnus spp.), barracuda (Sphyraena spp.), marlin (Makaira spp.), and sailfish (Istiophorus platypterus).
The associated mathematical expression of F (x i , v i , y, v) for this hunting strategy is devised as: Here, R 2 > r > 0 reflects the threshold distance where the predator begins to display predation maneuver to capture the nearby prey; θ 2 , γ 1 , and γ 2 are positive coefficients; the center position x c and velocity v c of the schooling prey are defined respectively as the average positions and velocities of the individuals constituting the school: The strength of attractive force depends upon the distance between the predator and the central positional point of the schooling prey, y − x c .The closer the distance between them, the stronger the triggered attractive force that draws the predator towards the center of the schooling prey.

Hunting Strategy II: Predator Attacks the Nearest Prey
The second hunting strategy is devised in such a way that allows the predator to pay more attention to the nearest (visible) prey.Put in other words, the predator focuses to hunt the most vulnerable target from the schooling prey.In the marine ecosystem, this particular predation tactic is mainly exhibited by the camouflage-fish species such as stonefish (Synanceia verrucosa) [42], scorpionfish (Iracundus signifier ) [43], frogfish (Antennarius commerson) [44], and the Atlantic stargazer fish (Uranoscopus scaber ) [45].
For this particular hunting strategy, the corresponding mathematical expression for F (x i , v i , y, v) is given as: In (2.5), all the pertinent parameters occupy similar definitions as the ones given for (2.3).The striking discrepancy with the former hunting strategy, however, lies in the fact that in the current hunting tactic the predator's behavior is affected by the movements of all individual prey constituting the school.The tendency of the predator to put special attention upon the nearest prey is depicted by the expression 1 y−x j θ 2 .Such a remark serves as the weighted coefficient for each of the associated prey.It implies that the closest prey would trigger the strongest attraction to the predator.

"Being Eaten" Condition for the Prey
Having discussed the proposed model, this section proceeds to elaborate the condition when one or more prey is captured by the predator during its attack.We use the term "being eaten" to describe such a situation.Mathematically, the condition of the i-th prey "being eaten" is defined as: where 0.01 m 1 specifies the multiplier constant.Put it in words, a particular prey is declared to be eaten by the predator when its Euclidean distance with the approaching predator is less than the critical distance between each individual prey in a school r multiplied by a positive constant m.From this definition, it is therefore possible for the predator to capture multiple prey simultaneously.Such a situation can be specified in the model by tuning the multiplier constant m.

Predator-Avoiding Patterns
As outlined in the Introduction, one of our aims in this study is to implement the proposed predator-prey model to discover the predator-avoiding patterns exhibited by the schooling prey.To fulfill such an objective, we carry out numerical simulations upon our predator-prey model.The associated computer codes are constructed in the form of programming functions thereby allowing us to better organize the assessment.
In the beginning, numerical simulations are carried out to establish schooling condition for the associated prey.Each of the prey is initially set to occupy a random position in the domain R d .The first task is therefore to bring all these scattered prey together in a manner that follows certain behavioral rules such that they manage to form a school.In the light of such an objective, the predator is (deliberately) excluded from the system.For the purpose of generating schooling prey, we adopt the model from our previous studies [36,37]: This system is similar in form to the first two equations in (2.1) except for the last term in the velocity expression.Therein, the friction force −kv i supersedes the term depicting contribution from the predator.Parameter k denotes a constant that manifests the magnitude of such a force.All the other quantities in (3.1) occupy similar definitions as the ones given for the first two equations in (2.1).Scheme 1 outlines the pseudocode for undertaking simulation to generate the schooling prey.
Therein, σ w denotes the magnitude of the white noise; ∆t is the time step; t max is the maximum simulation time; and v max is the maximum allowable velocity of individual prey within the school.
Subsequently, the outcomes of the schooling prey simulation in terms of prey's position x i and velocity v i are passed on as initial condition to the predator-prey function.By doing so, the predator-prey simulation can be carried out under the condition of the associated prey already formed a school.Scheme 2 summarizes the sequential steps taken for executing the predator-prey simulation to acquire the anti-predator maneuvers of the schooling prey.
In the computer program, we do not implement the conventional loop technique to carry out the successive calculations.Instead, we use the vectorized adapted parallel computation to accelerate the calculation process.In our experience during computer Transform x i such that the center point lies at coordinate (0, 0) (for d = 2) or (0, 0, 0) (for d = 3); return x i , v i 13: end function Scheme 1. Pseudocode to generate the schooling prey.
programming, such a practice improves the calculation speed significantly than the traditional loop technique thereby allowing efficient execution of the associated codes.
Furthermore, for the purpose of taking full advantage of the available computational resources, we arrange the computer codes to be amenable for execution both in the central processing unit (CPU) or in the more advanced graphical processing unit (GPU) platform.To search out for the evasive maneuver traits exhibited by the schooling prey to avoid predation, we run several predator-prey simulations in two-dimensional space R 2 under the condition of moderate prey's school size N = 30; the simulations are carried out for both hunting strategies.In this regard, we prescribe the simulations such that the predator only attacks the schooling prey once.Consequently, the total simulation time for all runs is fixed at t max = 3, 000 while the simulation time step is set to be ∆t = 5 × 10 −3 .Furthermore, appropriate tuning upon key model parameters is undertaken as well.During the assessment, once a particular avoiding pattern is discovered, we then repeat the simulation 100 times under the same parameter settings.Such an approach is carried out to examine the capacity of our model in returning a particular avoidance pattern consistently, considering the stochastic nature of our model.Although minor variation in the detailed structure of a particular evasive pattern is observed, the general characteristic of such a pattern is consistently obtained from all 100 runs.Such a condition indicates that our model occupies excellent robustness in returning a particular predator-avoidance pattern.
Based on the results of numerical simulations, we identify four evasive patterns of the schooling prey as defensive responses to the approaching predator.We label them as: (i) Split and Reunion, (ii) Split and Separate into Two Groups, (iii) Scattered, and (iv) Maintain Formation and Distance.Figure 1 shows the corresponding predator-avoidance profiles resulted from our model.Table 1 summarizes the adopted model parameter values to obtain the concomitant patterns, whereas Table 2 outlines the associated predator's hunting strategy to retrieve such patterns.
Pattern I (Split and Reunion) is the first evasive maneuver trait of the schooling prey obtained from our model (see Figure 1a).As can be seen in the figure, the schooling prey is initialized as a single composite unit with a circular shape in two dimensional domain, consistent with the prescribed initial condition of our numer- Table 2.The adopted hunting strategy in the predator-prey simulation for the attainment of predator-avoiding patterns.
Predator-Avoiding Pattern The Associated Hunting Strategy Adopted in the Simulation ical simulations.An observation study [46] identified such a circular structure as the defensive, look-around schooling mode of the associated prey.As the predator approaches, the schooling prey reacts accordingly by exhibiting a polarized state: all individual members of the school are oriented in a parallel arrangement, facing similar direction, and move together as an egalitarian unit (the leftmost profile of Figure 1a).As the predator arrived in the vicinity of the schooling prey, they readily displayed an anti-predator maneuver by expanding their formation at the right angles to the direction of predator's attack.Such a maneuver created some kind of vacuole in the middle-portion of the school (the second profile from the left of Figure 1a).Following this trait, the schooling prey momentarily broke the unitary formation into two smaller groups on the left and right side of the advancing predator (the third profile from the left of Figure 1a) before immediately rejoined behind the predator as a single school unit (the rightmost profile of Figure 1a).It appears that this particular pattern emerged due to the strong interactions among individual prey within the school invoked in the model (α = 15 and p = 4).In natural ecosystem, this particular evasive pattern closely resembles the fountain effect [11,41,47,48] (also called the F -maneuver).Such a term was introduced by Potts [47] to describe a particular anti-predation behavior exhibited by the schooling prey under mediumintensity of predation threat: the school splits into two smaller groups and rejoin behind the attacking predator [47,48].
Pattern II (Split and Separate into Two Groups) generally exhibits similar avoidance characteristics as Pattern I in regard to schooling prey's responses during early and mid-time predator attack (see the leftmost until the third profile from the left of Figure 1b).In this particular evasive maneuver, however, after the school splits into two smaller groups, the latter formation is maintained until the end of simulation: the schooling prey does not rejoin to form a unitary structure (the rightmost profile of Figure 1b).It is apparent that this particular anti-predator pattern emanated due to the weak interactions among the individual prey invoked in the simulation (α = 1).Based on the observed similarity in behavioral traits with Pattern I, it is therefore conceivable to treat Pattern II as an adaptive form of the conventional F -maneuver pattern due to some kind of limitations imposed by the natural habitat of the living organisms (for example, the presence of coral reef in the vicinity that limits the evasive maneuver of the schooling prey).
In Pattern III (Scattered), the schooling prey broke the unitary formation and displayed a disordered trait.In this case, the predator's attack managed to cause total disruption of the school, rendering the individual prey to move radially outward in such an abrupt manner (see Figure 1c).In regard to our model, it seems that this particular evasive maneuver appeared due to the large value of the fright factor (δ = 5) imposed in the simulation.A large value of δ represents the panic condition of the prey.Such a situation is in agreement with observational studies [11,47]: the schooling prey exhibits a tendency to break school formation when the school's members become panic due to the predator's attack.Potts [47] associated such a trait with a strong fright response of the schooling prey due to the high-intensity of predation threat.
The final evasive maneuver obtained from our model is Pattern IV (Maintain Formation and Distance).Here, the schooling prey demonstrates vigilant behavior: the school's members kept attention to the approaching predator and react accordingly by maintaining a safe distance from the predator.The schooling prey is readily polarized and moved in harmony with each other to avoid predation (see Figure 1d).
The school formation is therefore preserved until the end of the simulation period.In the framework of our model, it seems likely that this particular pattern was induced by the large value of γ 2 employed in the simulation (γ 2 = 10): such a condition allows the predator to move in alignment with the schooling prey.This particular pattern shares similar traits with those of low-intensity predation threat defined by Potts [47].
Let us now discuss the characteristics of the discovered predator-avoidance patterns in more detail.Figure 2 displays the characteristics of the retrieved antipredation patterns in terms of schooling prey's diameter (d schooling prey ) and standard deviation of prey's velocity (s v i ) as a function of simulation time, respectively.Both parameters are defined accordingly as: As can be seen in Figure 2a, each of the discovered predator-avoidance patterns exhibits different characteristics in terms of diameter of the schooling prey.For Pattern I, as the school's formation breaks into two smaller groups, the value of  the fact that in the corresponding pattern, the prey breaks the school and move radially outward in an abrupt manner.The maintenance of single school formation in Pattern IV is clearly demonstrated by a constant value of d schooling prey , after a small decrease at the initial stage of simulation.The characteristics of s v i for each corresponding evasive patterns reflect the tendency of the schooling prey (in the associated pattern) to move synchronously in the presence of predation threat.As we can see in Figure 2b, Pattern I and II share relatively similar traits of s v i .Here, both patterns exhibit a sharp increase in s v i in the middle portion of the simulation period, indicating an increase variability in the individual prey's velocity in association with the breaking of a unitary school formation into two smaller groups.After the prey splits formation and moves away from the attacking predator, they become less agitated and correspondingly reestablish the uniform characteristics at the end of the simulation period (shown by a decreasing trend of s v i followed by reestablishment of a new "stable" condition of s v i ).For Pattern III, the panic condition of the prey is clearly demonstrated in the sharp increase of the corresponding s v i at the early stage of the simulation.This is further supported by the fact that such a high s v i value prevails until the end of the simulation.A relatively constant value of s v i is observed for Pattern IV, indicating the maintenance of the synchronous individual movements within the school.
We would like to recall the unique feature that occurs in Pattern IV.Therein, it is observed that the school's diameter d schooling prey decreases during the initial stage of the simulation.Such a condition reflects important characteristic of the schooling prey under imminent predation threat, namely the polarization mode of the schooling prey.Partridge [11] defined such a trait as a parallel arrangement of the school members: when the school is under imminent predation threat, the individual members of the school move closer to one another and establish synchronous movements with their neighbors [11].In this regard, our model occupies excellent capability in demonstrating the occurrence of the polarization mode.To demonstrate the appearance of such a trait, we carry out numerical simulations of Pattern IV under different school size (N = 15, 20, 30, 50, and 80).Next, we would like to point out the capacity of each of the discovered predatoravoidance patterns in providing protection to the corresponding prey.To do so, we calculate the number of eaten prey at the end of the simulation period (t max = 3, 000) for every observed patterns.For each evasive pattern, we carry out the numerical simulations 1, 000 times. Figure 4 outlines the results of such a procedure.Therein, the largest number of eaten prey is exhibited by Pattern III (on average, up to 11 prey is captured by the predator).From Figure 4, it can therefore be inferred that Pattern III offers the least protection mode for the prey to avoid predation.It seems that the prey becomes more vulnerable to the predator when they are agitated and broke the unitary school formation.Such a finding is consistent with observational studies of Pavlov and Kasumyan [41] and Shaw [48].On the other hand, the smallest number of eaten prey is demonstrated by Pattern IV (the average number of captured prey = 0).Apparently, when the prey manages to keep the unitary school formation and maintains safe distance from the predator, the latter entity can never successfully catch an individual prey until the end of the prescribed simulation period.It is therefore conceivable to deem Pattern IV as the most effective evasive mode to avoid predation.Furthermore, our model suggests Pattern II as the second-best evasive mode of the schooling prey to mitigate predation risk (on average, only 2 prey is captured by the predator).Pattern I, however, does not come up as the favorable mode for the schooling prey as the average number of eaten prey is observed to be four times as large as Pattern II.
As a final annotation of this section, let us point out a few relevant examples in real ecological environment related to the discovered evasive patterns from our model.
By doing so, we sought to demonstrate the excellent conformity of the simulated predator-avoidance patterns with the real-life evasive patterns that occurred in the natural ecosystems.In this regard, the examples are taken mainly from predatorprey interaction in aquatic environment-including those of piscivorous birds-that include schooling behavior in the associated prey.Evasive maneuver in the form of Pattern I is well-exhibited by a school of grouper fish in avoiding the attack of the red bass predatory fish (Lutianus bohar ) [47].A similar trait is also displayed by a school of sand-eels (Ammodytes sp.) in avoiding the incursion of mackerel [49].On the other hand, the defensive trait of Pattern II is observed to be adopted by a school of hardyheads (Pranesus capricornensis) in escaping the attack of avian predators, such as tern (Thalasseus bengalensis) [28] and Australasian gannet (Morus serrator ) [26].Other observational studies [25,50] reported similar behavior exhibited by a school of herring (Clupea harengus) in dodging the attack of a killer whale (Orcinus orca).Anti-predation maneuver of Pattern III is reported to be observed in a school of Atlantic silversides (Menidia menidia) in escaping the predation of a predatory black sea bass fish (Centropristes striata) [41].Finally, the evasive maneuver of Pattern IV is found to be performed by a school of snappers (Lutianus monostigma) in escaping the attack of a barracuda or a black-tipped shark [11,47,48].

Benefit of Prey's Schooling in Surviving Predator's Attack
In this section, we investigate the contribution of schooling prey's population size N upon the chance of survival (not being eaten) of each individual prey in the school as well as their corresponding living time.
First, let us define a mathematical expression for the probability of individual prey to be eaten by the predator: Here, P eaten (t) and P eaten (t − 1) denote the probability of individual prey being eaten at current and previous time steps, respectively while N survived (t) and N survived (t − 1) are the number of survived prey at current and previous time steps, respectively.
P eaten (t) is therefore calculated recursively with initial condition set as P eaten (0) = 0.
The primary interest here is examining the characteristics of P eaten at the allotted final simulation time t max , namely P eaten (t max ).sufficiently large number of repetitions.Here, the repetitions (sampling trials) are carried out 3, 000 times so as to obtain adequate data to capture inherent variability induced by the stochastic nature of our model.Clearly, such an operation involves a large number of simulations.We therefore implement a parallel computational process.Such a configuration allows the simulations to be executed up to ten times Let us first discuss the characteristics of P eaten (t max ) as a function of school size N .As can be seen in Figures 5(a) and 6(a), in general, the values of P eaten (t max ) for They found that for all predator-prey systems, increasing the number of prey from 1 to 6 and then to 20 reduced the probability of success of the predator in capturing the prey.The functionality of prey group size in reducing the risk of predation is also reported in [11,48,53].
As a final annotation, we would like to relate our findings with the established theoretical framework of predator-prey interaction.In this regard, mitigation of predation risk due to increasing group size in the schooling prey species has been acknowledged as the dilution effect [9,12].Pitcher [12] argued that such an effect should be considered alongside the number of predator's attack upon the associated schooling prey, resulting in the so-called attack-abatement effect.In our simulation study, we realize that the advantage of schooling ceases as the prey's population size N reaches the third segment (120 N 200): as N increases above 120, there is no additional reduction in the probability of individual prey being eaten by the predator.
Such a condition may be attributed to the fact that we limit the number of predator's attack in our simulation: only single attack of the predator is allowed.Therefore, in the light of theoretical framework, our model occupies excellent capability in demonstrating the dilution effect.Consideration of the attack-abatement effect is, however, not within the scope of this study.Such a point of view is left for future work.
population in the school increases, the probability of an individual prey to be eaten by the predator decreases.This result is consistent with observations in the real-life aquatic ecosystem.
(4) Our model is consistent with the established theoretical framework of predatorprey system: the model possesses excellent capability in demonstrating the occurrence of the polarization mode of the schooling prey and the dilution effect.

Figure 1 .
Figure 1.Predator-avoidance patterns in two-dimensional space obtained from the proposed predator-prey model: (a) Split and Reunion, (b) Split and Separate into Two Groups, (c) Scattered, and (d) Maintain Formation and Distance.In the figure, the blue dots represent the schooling prey while the red dot denotes the solitary predator; the length and direction of the arrows attached to each of the dots designate the magnitude of velocity and direction of movements of each corresponding entity at the particular time, respectively.
Pattern I: Split and Reunion II Pattern II: Split and Separate into Two Groups I Pattern III: Scattered II Pattern IV: Maintain Formation and Distance I schooling prey increases.As the two smaller groups get closer, it decreases accordingly until it returns to the original value, indicating that the prey rejoins into a unitary school formation.Distinct characteristic is apparent for Pattern II where the value of d schooling prey increases gradually until the end of the simulation period, indicating that the two smaller groups of prey do not rejoin into a single formation.For Pattern III, a sharp increase of d schooling prey is observed.This condition is consistent with

Figure 2 .
Figure 2. Characteristics of the discovered predator-avoidance patterns in terms of: (a) Diameter of the schooling prey (d schooling prey ) and (b) Standard deviation of prey's velocity (s vi ) as a function of simulation time, respectively.

Figure 3
shows the corresponding characteristics of d schooling prey and s v i , respectively.As shown in Figure3a, the first attribute of the polarization mode (namely, the individuals in the school move closer to one another) is clearly displayed by the decrease of d schooling prey during 0 t 500 for all the corresponding cases.Such a condition reflects that as soon as the schooling prey realizes that they are under imminent predation threat, they immediately become more vigilant and positioned themselves closer to each other.The second attribute of the polarization mode (namely, the individuals in the school align themselves synchronously with their neighbors) is shown in Figure3b.Therein, after a sharp

Figure 3 .
Figure 3. Polarization mode displayed by the schooling prey (represented by the evasive maneuver of Pattern IV (Maintain Formation and Distance)), outlining the characteristics of: (a) Diameter of the schooling prey (d schooling prey ) and (b) Standard deviation of prey's velocity (s vi ) as a function of simulation time, respectively.

Figure 4 .
Figure 4.The total number of eaten prey for each corresponding predator-avoidance patterns at the final simulation time t max = 3, 000.

Figure 5 .
Figure 5. Numerical results of center-attacking strategy: (a) Profiles of P eaten (t max ) as a function of school size N , (b) Profiles of talive as a function of school size N .

Figure 6 .
Figure 6.Numerical results of nearest-attacking strategy: (a) Profiles of P eaten (t max ) as a function of school size N , (b) Profiles of talive as a function of school size N .