{'fraction_non_alphanumeric': 0.0287730593041971, 'fraction_numerical': 0.008617060789897239, 'mean_word_length': 4.822664359861592, 'pattern_counts': {'":': 0, '<': 8, '<?xml version=': 0, '>': 8, 'https://': 0, 'lorem ipsum': 0, 'www.': 1, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Recognizing and prevention of probable regime shift in density regulated and Allee type stochastic harvesting model with application to herring conservation

An ecological system with multiple stable equilibria is prone to undergo catastrophic change or regime shift from one steady-state to another. It should be noted that, if one of the steady states is an extinction state, the catastrophic change may lead to extinction. A suitable manual measure may control the prevention of catastrophic changes of different species from one equilibrium to another. We consider two stochastic models with linear and nonlinear harvesting terms. We inspect either density regulation or Allee type density regulated models [Saha et al., Ecological Modelling, 2013], which have substantial applications in the herring fish population's viability study. Both the deterministic models we consider here contain bi-stability under certain restrictions, and in that case, one of the stable states is the extinction state. We assume that the dynamical system under consideration is closed, i.e., immigration and emigration are absent. The demographic noise is introduced in the system by substituting an ordinary differential equation with a stochastic differential equation model, where the birth and death rates of the deterministic process are used to obtain the instantaneous mean and variance in the stochastic differential equation. Our study reveals that, the catastrophic changes can be avoided manually by a suitable choice of handling time that will eventually help to prevent the sudden extinction of the harvested population. The entire study is illustrated through the herring population size data obtained from the Global Population Dynamics Database (GPDD) and simulation experiment.


Introduction
The exploitation of biological renewable resources and harvesting of these resources are common practices in fisheries, forestry, etc. A suitable way is essential to carry out the adequate requirement of the human need for a long time. The interaction between human activities and the ecosystem has a significant impact on nature at local to global scales [14]. Indiscriminate harvesting of biological resources may hamper the sustainability of the different rare and endangered species. Hence, harvesting should be carried out scientifically to obtain the maximum benefit without hampering the ecological balance. Improper harvesting may also cause a catastrophic change in the population under consideration in the ecosystem. Many theoretical models suggest that the ecosystem may switch instantly from a stable state to an alternative stable state which is termed as regime shift [17,29,30]. It may happen that if a regime shift occurs, it will be impossible to restore the stock to previous levels [26]. Moreover, regime shift is an important feature of the physical environment of ecosystems and has the capability to make an effect on stock productivity. The shift refers to low frequency, high magnitude fluctuation in the marine ecosystems presuming changes in community composition, population abundances, and trophic structure. These changes happen in the abundance of both exploited and unexploited populations. An alternative stable state for regime shift is an important ecological phenomenon in the environment. When an alternative stable state occurs in the ecology, multiple stable equilibria exist in its underlying deterministic skeleton, and many exciting phenomena can occur [25]. A suitable example of regime shift phenomena in ecology is lake eutrophication, where water is polluted by high concentrated nutrients such as nitrogen and phosphorus [11,36].
Generally, the Allee effect manifests a population bi-stability in a deterministic setup [21], whereas the theta-logistic model does not exhibit the same stability criteria. Although, in the case of the density regulated growth process, one can not rule out the possibility of the existence of multiple stable states in the presence of harvesting other than linear harvesting. Most management measures are directed at the individual stock of single species and do not take into account species interaction, such as prey-predator relationship [23]. This motivates us to consider the important single-species models with density regulated parameters in the Allee type phenomena. A study on individual population dynamics by Sibly et al. [27] concludes that most of the density and per capita growth rate (pgr) relationship of different species is concave in nature. This relationship can be well defined by the density-dependent theta logistic model. In many cases, the rapid depletion of population size may reduce the species' fitness in a concave pattern with additional exposure to the Allee effect. Hence the density regulated Allee model is appropriate to explain such growth phenomena. Sau et al. extended the Allee model with density regulation around carrying capacity under harvesting. The author introduced this model in studying the extinction status of one of the economically beneficial fishes, viz. Atlantic herring [28]. We believe that, exploring the harvesting issue in the theta-logistic and Allee model must be an exciting research area to be explored as it can cover the growth profile of a wide range of species.
It is worthy to mention that, if the growth process of species follows the Allee mechanism, it enhances the morality of species [12]. Incorporating the harvesting issue in this model can be interpreted as additional mortality. In addition to this, when the harvesting is incorporated into the system, it may act in favour of the regime shift from a higher stable state to a lower stable state. Harvesting is associated directly with the reproduction process of the species. If the reproduction success is more, the chance of effective harvesting is naturally increased. In many research articles on harvesting, authors introduced stochasticity through white noise instead of the birth-death process [15,22]. Harvesting is an important issue for the closed system as immigration and emigration do not occur directly in this case. Thus, extinction is a purely absorbing state under the birth-death set up, but this is not the case for white noise. Note that, in the second case, the population can recover from an extinction state which is not possible for a closed system. This justifies our argument that the association of demographic stochasticity through the white noise is not pertinent to the harvesting model. This inconsistency can be avoided if the demographic stochasticity is incorporated through the birth-death process in place of any other commonly used noises. Thus, in this work, we frame a suitable birth-death process using the deterministic model to incorporate the stochasticity in the system [6,33,34].
In this article, we consider two stochastic models viz. theta-logistic and Allee type density regulated models with linear and nonlinear harvesting terms. Both the linear and nonlinear harvesting strategies are considered, although the emphasis is given on nonlinear harvest rate as it is more realistic compared to the linear one [23]. The nonlinear harvest rate involves an additional parameter known as handling time. We can search for an optimum choice of density regulated parameter and the handling time so that we can restrict the species from undergoing regime shift from a higher stable state to a lower one (may be an extinction state). This leads the species under consideration to stay in a stable and safe mode instead of the risk of extinction. Similar stochastic models were also analyzed in Sau et al. [28]. The extract of the paper was mainly associated with the species' probability of extinction and the expected time to extinction. We believe that these two measures are insufficient to properly explain species' extinction threat since these two measures cannot capture the ensuing catastrophic change or possibility of regime shift, which may lead to species extinction. This may happen for a system having multiple stable equilibria. For example, suppose the estimated probability of extinction is found to be small or the expected time to extinction is large. It may happen if the species maintain its stable state at a moderately high equilibrium state. But this does not give us any clue about the possible regime shift to a lower stable equilibrium (may be an extinction state) if the system is perturbed by the noise of sufficiently large strength. This entire part was missing in the work of Sau et al. [28].
We demonstrated the theoretical findings with the help of time-series data of the herring fish population, obtained from the Global Population Dynamic Database (henceforth GPDD, http://www.sw.ic.ac.uk/cpb/cpb/dpdd.h with GPDD ID-1741 and GPDD ID-1772. We use the early warning signal toolbox in R to generate generic early warning signals of the time series data. The value of the control parameter, handling time, is estimated under the suggestion of the proper harvesting policy so that the possibility of extinction can be minimized. We derive the stationary distribution using the Kolmogorov forward equation. The infinitesimal mean and variance are obtained from the expressions of birth rate and death rate. In a stochastic environment, the sustainability of the resource population can be viewed with the help of its statistical distribution or by an effective potential function. The number of local minima of a potential function gives some significant input in the nature of stability and the possibility of the regime shift from one stable state to another stable state.
The paper is organized as follows. The section 2 contains some ideas of both the density regulated and density regulated Allee type model considering the linear and nonlinear harvesting. Then we consider our models in both deterministic and stochastic frameworks. In section 3, we provide the sketch of the stationary distribution for both linear and nonlinear harvesting. We examine sustainability through the potential function, regime shift phenomena, etc. In the next section 4, we simulate our result through numerical technique, and in section 5, we interpret our result through real-world data. Finally, we end with a conclusion in section 6.

Model formulation 2.1 The deterministic model
The Allee effect is a biological phenomenon identified by a correlation between population density and the mean individual fitness of a population. An Allee effect may be pervasive when the population size becomes low, and a critical threshold value exists underneath the per-capita growth rate becomes negative. The population may be vulnerable to extinction due to this critical depensation level [12]. For example, many aquatic animals such as fishes move from one region to another in a school, as schooling is generally considered a defense policy against predator attacks. Again the indiscriminate harvesting in dense school is a major reason for the decline of the population into danger [16]. In such cases, the strength of density dependence can have potential implications for the preservation of animal populations [10]. Note that, apart from the Allee threshold, the strength of the density regulation parameter θ has an important role. To capture Allee phenomena together with density regulation around carrying capacity, we consider the equation according to Saha et al. [10] dx where x(t) is the biomass at time t, r is called the intrinsic growth rate, K is the carrying capacity, θ is called the density regulation around carrying capacity. We call this model as Allee-Saha model (ASM) [28]. However, in reality, most of the species follow θ−logistic process [27], and hence we also consider the θlogistic growth process to encompass a large number of species. In the θ-logistic model, the per-capita growth rate (pgr) decreases monotonically with the increase in population size. The theta logistic model is widely used to clarify the density dependence in the real population as well as the time-series data. The general theta-logistic model can be represented as [9] dx dt (2. 2) The parameters are described earlier. The pgr growth profile of the species is convex or concave depending on the value of θ. The relationship is convex or concave according as the value of θ is greater or less than 1 respectively [27]. When θ is equal to 1, the equation (2.2) become classical logistic model. The sustainability of different rare and endangered species may greatly depend on harvesting strategy [28]. We consider two types of harvesting strategies viz. proportional or linear harvest rate and nonlinear harvest rate. In general, the linear harvesting strategy is determined by the catch-per-unit effort hypothesis [13] and written as H(x, E) = qEx, where q is the catchability coefficient and E is the harvesting effort. Note that, this harvesting strategy has some defect in the real-life situation, which may be described as follows [19]: 1. It considers the random search for every resource population.
2. It assumes that every individual to the resource population is equally likely to be captured.
3. H is unbounded with respect to E for a fixed x.

4.
H is unbounded with respect to x for a fixed E.
Handling time is one of the important parameters to avoid the above four unrealistic features. Hence, the functional form of the harvesting policy is written as H(x, E) = qEx aE + Lx [19], in which a is the degree of competition among the boats, fishermen and other technology used in fishing [4]. L is the product of capture rate and handling time [3]. If the capture rate is assumed to be constant, L is dependent only on handling time, and hence, L can be treated as a representative of handling time. This improved version of harvesting exhibit a saturation effect with respect to both resource size and effort. It is clearly observed that, h → q a x as E → ∞ when x is fixed and h → q L E as x → ∞ when L is fixed [19]. Henceforth, we call this modified harvest rate as non-proportional or nonlinear harvest rate. Introducing the harvesting phenomena in the models (2.1) and (2.2), we have where H(x) is the harvesting component. In case of our model, we consider both the linear and nonlinear harvesting.

Stochastic model
In many cases, the system of Ito stochastic differential equations are used to study a random dynamical system. A commonly used technique to develop SDE is by studying possible changes in the system components during a small interval of time. In this case, the transitions in the system are studied over a small time interval, and a differential equation is formed by considering the time interval approaches zero [5,6]. As an example, the birth-death process may be applicable to introduce stochasticity in the system. A stochastic process is addressed as a birth-death process where jumps from a particular state (number of individuals, lineages, cells, etc.) are only allowed to neighboring states. If the number of individual or similar quantities is increased by one, it represents the birth, while a decline to the left-hand part constitutes death. These characteristics of the birth-death process help us clarify the mathematical analysis, but we apply this process in different real-world mathematical modeling. This type of model permits us to deal with any questions regarding the transition or state probabilities of the process, mean, variance, stationary distribution, the meantime of existence, probabilities of extinction, etc. The outcomes procure with these models can be compared with real-world data.
The presence of the Allee effect invites the possibility of regime shift, which may lead to the extinction of species. Harvesting is associated directly with the reproduction process of the species. Demographic stochasticity is an inherent property of the species which is directly related to the reproduction success. In comparison, environmental stochasticity is a sudden event that affects species reproduction by certain chance factors. In the present investigation, we have focused on the effect of demographic stochasticity in the system. Environmental stochasticity, where random fluctuation in the environment is effective to the entire population, is not considered here.
Introducing demographic stochasticity through white noise has certain limitations. This type of stochastic perturbation biologically represents disturbances, which occurs independently of current population densities but depends on the immigration and emigration of the species [2]. Harvesting is such a system where immigration and emigration do not occur directly in the model. So the incorporation of demographic stochasticity through white noise is not appropriate. This complication can be avoided if the demographic stochasticity is incorporated through the birth-death process instead of any other commonly used demographic noise [6,33,34]. In our study, both the time and state space are continuous variables. Let p(x, t) indicates the probability density function (p.d.f)for the continuous random variable X(t). Then, the mathematical notation of p.d.f is written as It is a well known established formula [6] that, the probability density function for the general birth and death process satisfies the Kolmogorov forward differential equation are the birth and death rates.
The sample paths, X(t) of a stochastic process satisfy the succeeding Ito stochastic integral equation [6,18,37] Here, the first integral part of the (2.6) is the deterministic part, i.e., Riemann integrable, and the second integral part is a stochastic Ito integral part. Here W is known as Wiener process and defined as △W (t) = W (t+△t)−W (t) and it satisfies normal distribution, N (0, △t). For simplicity, this equation can be expressed as the SDE [6] In a population model, the birth and death function should follow the following properties. For ASM, we assume that there exist real numbers M , N and A * such that, In addition to the above four properties b(x) and d( For θ-logistic model, we assume that there exist two real numbers M and N with M < N such that [6] 1 In addition to the above four properties, b(x) and d(x) are C 2 ([0, N ]). The feasibility of these regularity conditions is well discussed in Appendix A.
There are multiple ways to split up the single growth rate into the difference of a birth and death rate for the θ-logistic and ASM. Among many considerable forms, we have chosen the most ecological meaningful choice for our model considered here. Basically, in the θ-logistic model, the per capita growth rate (pgr) is high at low population density, and it is maximum at zero population size. The pgr diminishes with the increase in population size as more populations are competing for the same limited resources. Hence, the crowding effect, which is commonly considered to be the reflection of competition is present in the species which follows the θ-logistic growth process. On the other hand, a species the growth process of which follows the Allee mechanism generally adopt the cooperation at low density to avoid becoming extinction [1]. Thus the competition or, in other words, the crowding effect can be avoided in the Allee model.
After careful consideration of these facts and the above four properties, we obtain the expressions for the birth and death rates corresponding to (2.3) and (2.4) according to [38]. We obtain the following expressions for birth and death rates respectively, Using two sets of birth and death rates, the ASM and the θ-logistic model with harvesting term are governed by the following stochastic differential equations For the stochastic differential equations 2.10 and 2.11, X = 0 is an exit boundary point [24]. It implies that X(t) converges in finite time to the absorbing state 0. Thus the only stationary distribution is the absorbing state at 0 with respect to the Dirac delta measure.

Existance of bistability under nonlinear harvesting:
The existence of the bistability of ASM is already discussed in [28]. Here 0 is always a stable equilibrium point, and the nontrivial solutions are stable or unstable depends on certain conditions. But for the θ-logistic model, the trivial solution is not always stable, and it also depends on some conditions. We have The linearization around K will yield the following equations The equilibrium points of the above system are 0, x 1 and x 2 where x 1 , x 2 are given by If (q − ra) > 0 and L − aE K > 0 then x 2 exists. But when (q − ra) > 0, 0 is stable [28]. This shows that when (q − ra) > 0, the two stable equilibrium 0 and x 1 is separated by the unstable equilibrium point x 2 . On the other hand if (q − ra) < 0, x 2 does not exist and 0 becomes unstable.

Stationary distribution and possibility of regime shift
Stationary distribution is one of the important measures to visualize the model's behavior in a long time. It is the probability distribution that remains unchanged as time progress. The potential function is an alternative measurement of the stationary distribution, which is equally capable of determining the possibility of the regime shift of a system. The increase in the asymmetry in the potential function is the indication of catastrophic change or regime shift [20,31]. Moreover, the local minima of a potential function indicate the location of the stable equilibrium points of the deterministic counterpart of the model. Consider {X(t), t 0} be a continuous stochastic variable with transition density p(t, x, y) for t 0, which satisfies the Kolmogorov forward equation

Stationary distribution: θ-logistic model
Where µ(x) is the drift coefficient or infinitesimal mean and σ 2 (x) is infinitesimal variance. Hence, there exist stationary density f (.) (may or may not exist), which satisfy the following condition f (x) = f (y)p(t, x, y) dy for all t 0.
In the case of θ-logistic harvesting model (2.4), where the parameters are represented as earlier, the stochastic differential equation corresponding to the deterministic model (2.4) is given by (2.11). We find the stationary distribution for the stationary state by the Kolmogorov forward equation and equating ∂p ∂t = 0.
We define the probability density function at stationary state as P s (x). In the case of the model (2.5), Where N is the normalizing constant and dx , is the effective potential function [20]. In particular for model (2.11) where, N (E) is the normalizing constant and defined as

Case 2: Nonlinear harvest rate
In the previous section, we give the expression of stationary distribution when the harvest is linear. Traditionally, the harvesting policy is generalized by the catch-per-unit effort hypothesis [13]. So the handling time is sensible for the bounded functional response. In this section, we consider the nonlinear harvesting Similar by previous technique, the stationary probability density function is given by, Here N is the normalizing constant and U (x) is called the effective potential function, which is defined as [20] U Infinitesimal variance σ 2 (x) is defined in the previous section and σ ′ (x) = dσ dx . In this case, the θ-logistic model with nonlinear harvest rate takes the form, From the equation (3.5), we can write the birth and death terms as: The Kolmogorov forward equation in this case can be written as, So from (3.4) we obtain the following expression, The exact mathematical expression of U (x) can not be established due to the presence of nonlinearity in the population biomass. To resolve the problem, we linearize the term 1 − ( x K ) θ . The approximation of 1 − ( x K ) θ can be made using binomial expansion as follows: Substituting the value in the above expression of U (x), we get, , C = L(r + 2), D = qE + aE(r + 2). Note that, this form of U(x) is valid for 0 < x < 2K.

Linear and nonlinear harvest rate
In this section, we will derive the stationary distribution under nonlinear harvest rate and the corresponding stationary distribution in case of linear harvest rate can be easily obtained as a special case of the former distribution with L=0. Here, we will give some outline or sketch of the stationary distribution for ASM. The derivation of the analytic expression will be difficult in this case, so we mainly emphasize the numerical technique to analyze the potential function. For the ASM with harvesting, the forms of birth and death terms are . The harvest term may be linear or non linear according to our choice. Considering nonlinear harvest rate, the mean and variance can be written as Substituting these values in the equation (3.4), we will obtain the potential function. This expression is not solvable analytically. So we will analyze this properly through numerical technique, which is discussed in section 4. The potential function and hence the stationary distribution in case of proportional harvesting can be achieved by assigning L = 0.

Stochastic perturbation on density regulated Allee type model
In nature, the dynamical systems are dependent on external interference. Many inconveniences such as pest outbreaks, changes in weather conditions, rainfall, fire, etc., may be harmful to the populations in multiple ways, such as exterminate a portion of it. If the system has a single steady-state, the system will settle back to the same state after such a perturbation. Whereas, if the system has multiple stable states, such as the Allee effect, a sufficiently large disturbance of the ecological state may put the ecosystem into the basin of attraction of another stable state. The plausibility does not depend only on the perturbation but also on the size of the basin of attraction. In terms of stability landscapes, if the valley is small, a small perturbation may be sufficient to displace the ball far enough to push it over the hill, resulting in a shift to the alternative stable state. For example, in fisheries, increased harvesting on piscivores can result in a shift from high piscivores low-planktivore to the state low-piscivores high planktivore regime [35]. Several fisheries are suspected of having suffered this type of transition [32,8].
In general, when the growth process of the species follows ASM in the absence of harvesting, the species has two stable equilibrium points, 0 and K, separated by an unstable equilibrium point at Allee threshold A.
Here we introduce the nonlinear harvesting strategy, which is obtained incorporating handling time during the harvesting process.
Since the exact expression for non zero equilibrium points is intractable, we rely on the bifurcation diagram in figure 1. When the handling time is small, the extinction state is the only stable equilibrium point. However, when we increase the handling time, there exists another stable equilibrium point, which is near the carrying capacity i.e, two stable equilibrium is separated by an unstable equilibrium. This observation establishes the fact that, if handling time (L) lies below a certain threshold value corresponding to the bifurcation point (as indicated in figure 1), extinction is inevitable. This bifurcation figure is also discussed in the work of Sau et al. [28]. In the long run, the distribution of population size can be represented by the stationary distribution of the population size. However, the computation of stationary distribution for ASM is very complex, as we observe in section 3. To avoid this complexity, we compute the potential function of the population size, which possesses just an opposite behavioral characteristic of the stationary probability density function. The potential function attains minima where the stationary distribution has mode and vice-versa. In this case, we observe that if L increases, the minimum value of the potential function at the lower equilibrium point remains the same, although it decreases at the higher equilibrium point. It indicates that the population is more likely to stay around a higher equilibrium point if L increases ( figure 2(a)). In addition, an increase in the minimum value at a higher equilibrium point indicates a higher chance of regime shift from upper stable equilibrium to lower stable equilibrium, which is the extinction state. Furthermore, it is observed that if L is increased, the basin of attraction for the upper stable equilibrium point increases ( figure 3). This shows that a proper handling time (sufficiently large) can prevent possible catastrophic change and thus restrict population extinction. This argument is validated by numerical simulation in data analysis.
A similar result can be observed for θ also. If θ increases the intra-species competition decreases (since (x/K) θ decreases for x < K) and naturally, the population is more likely to stay at higher equilibrium density (see figure 2 (b)).

Density regulated harvesting model under stochastic environment
In the previous section, we demonstrate the effect of different relevant parameters on the possible regime shift for ASM, which possesses bi-stability. In this section, we consider the θ-logistic harvesting model, which may possess multiple stable equilibrium points under a nonlinear harvesting strategy. We also analyze the effect of different relevant parameters on the stationary density, which is obtained in the range of population size [1, ∞). Figure 5(a) shows that, for linear harvest rate, the resource population becomes more likely to stay around the carrying capacity as θ increases. This phenomenon occurs because, increase in θ diminishes the intra-species competition, which in turn increases the species birth rate. It should be noted that, when θ is large, the species such as mammals generally likes to stay around carrying capacity and are more stable around it [27]. When a linear harvesting strategy is applied, for small values of θ, the probability of staying near a low population size is high, which is alarming for the population. This problem can be overcome if nonlinear harvesting is applied. We observe that, if a suitable amount of handling time is applied, the probability of staying near population size 1 is negligible, although density regulation parameter θ is small ( figure 5(b)). Furthermore, it can be observed ( figure 5 (c)) that if the handling time (L) is increased keeping all other parameters fixed, the population size has a tendency to stay near carrying capacity with high probability. Thus, we may enhance the possibility of survivability of species whose density regulation around carrying capacity is small by adopting nonlinear harvesting instead of linear. Density regulation around carrying capacity (θ) is an inherent property of the species and for most of the species θ < 1 [27]. So we may conclude that nonlinear harvesting is beneficial for the sustainability of species.
A dynamical system faces a catastrophic shift when it has multiple stable equilibria separated by unstable equilibria. The theta-logistic model with nonlinear harvesting may possess multiple stable equilibria under certain choices of the parameters. To avoid mathematical complexity, we consider the linearized model as discussed in the section 2.3. In this case, we observe that if (q − ra) > 0, there exist another two nonzero equilibrium points x 1 , x 2 and x 1 is stable. This shows that under the constraint (q − ra) > 0, the two stable equilibrium 0 and x 1 are separated by an unstable equilibrium x 2 . The exact expression for the potential function is also obtained for this linearized model. The potential function gradually loses its symmetry when the system enters into the unstable region from the stable region. The onset of asymmetry in the potential function indicates a high chance of catastrophic change in the system [20]. This may lead the system to become extinct. The change in equilibrium density of the resource population with respect to L and θ is depicted in figure 6, which is also given in [28]. The equilibrium density marked with red color represents the stable equilibrium, and those marked in blue color are considered as unstable points. As the choice of the parameters is shifted from the stable to the unstable region, the potential function loses its symmetry, which is evident in figure 7. In addition, it can be seen in figure 8 that an increase in L leads to an increase in the symmetry of the potential function, which establishes the fact that increased handling time reduces the chance of catastrophic regime shift. As the parameter theta is an inherent feature of the concerned species, we do not have any control over it. We can adjust handling time to avoid catastrophic change when a system contains bistable equilibrium points, and the species maintains its stable existence at the upper stable equilibrium point. We observe from figure 4 that if L is increased, the basin of attraction of upper stable equilibrium also increases. This means that the population can withstand a larger depletion in population than that when L is small.
The linear (proportional) harvesting strategy is a special case of nonlinear harvesting when L is chosen to be 0. Hence, in every simulation activity, L may be chosen to be 0 to incorporate the proportional harvesting policy. It is observed that, when the parameters L, θ are chosen from the stable region (marked in red in figure 6), the population size is almost normally distributed with a mean near carrying capacity ( figure 9(a)). If the chosen parameters are away from the stable region, the said distribution becomes asymmetric, and the clustering of population near 0 increases (figure 9).

Data analysis
The extract of the study in the previous sections altogether reveals that, when a species remain in the upper stable equilibrium, we can prevent the species from going to the zero stable equilibrium by controlling the handling time in the nonlinear harvesting. This is one of the important mechanisms to control species sustainability. We validate our numerical results through population time-series data of herring fish population in two specific regions available in the GPDD. We consider the population time-series data of 16 years (1972-1987) with GPDD ID-1741 and 41 years  for GPDD ID-1772. For the data set GPDD ID-1741, the experiment was performed in Baltic sea areas covering 28 and 29 S area, and for the data set GPDD ID-1772, it was done in Prince Robert district, British Columbia. Sau et al. [28] fitted these population size data (in Metric Tonnes) to the θ-logistic and ASM considering both the linear and nonlinear harvest rates. The fitting through the RGR model using a grid search technique was performed for the species GPDD ID-1741, and ASM with nonlinear harvest rate is found to be the bestfitted model. For the species with GPDD ID-1772, they used the size modeling approach, where the raw population time-series data can be directly applied for estimation. They found that the θ-logistic model with nonlinear harvesting is the best-fitted curve in this model.
In this article, we generate sample paths using the Gillespie algorithm. We have simulated the sample paths for different handling times; the remaining parameter values are obtained from Sau et al. [28]. The simulated time series for the herring fish population corresponding to GPDD ID 1741 shows that, a catastrophic regime shift may occur near the time point 1500 if the harvesting strategy is linear. The possible regime shift can be delayed effectively if nonlinear harvesting is applied. In figure 10(b) it is observed that a small increase in the handling time is sufficient to delay this shift up to 2000 time points.
The sample path for GPDD id 1772 for different values of L is depicted in figure 11. For this species, the θ-logistic harvesting model with a nonlinear harvest rate is found to be the best-fitted model. As in the previous case, here also we observe that the sufficient increase in the handling time can delay the catastrophic change and thus can prevent population extinction (see figure 11 (b).) This observation conveys to the management authority that non-linear harvesting with suitable handling time can be permissible to prevent this species from extinction.

Generic early warning
Different early warning signals can be useful to detect impending regime shift directly from the time series data. In the previous section, for the GPDD ID 1741 and ID 1772, we observe the possibility of a catastrophic shift in the fish population density from an adequately high stable population density to an extinction state in the long run when handling time is low ( figure 10, 11). To find the generic early warning signals, we restrict our observation on the time series data for the population size of last 500 time points for both the species since the early warning signals may fail to predict the regime shift for a long-term time series [7] . In both cases, the standard deviation, autocorrelation increase where as the return rate decreases prior to the catastrophic shift ( figure 12, 13). This analysis supports the incidence of occurrence of regime-shift. Thus the generic early warning signals are found to be robust in this case.

Conclusion
One of the major challenging issues in ecology is maintaining an adequate abundance of different essential biological resources to prevent its extinction. Marine fisheries are salient examples of such an ecological system as the fishes are being depleted rapidly due to indiscriminate harvesting [10]. In this context, an appropriate harvesting strategy is important to maintain the sustainability of the species. Besides this, the marine system and similar ecological systems are highly exposed to catastrophic changes, particularly when they have multiple or alternate stable states. The system with multiple stable states may be prone to extinction if one of the stable states is the extinction state.
In this article, we consider a stochastic harvesting model where the growth process of the species follows theta-logistic and Allee mechanism [28]. Stochasticity is introduced through the concept of the birth-death process. The extinction state is one of the stable states for both models. Our findings in the present work reveal that the possibility of regime shift or catastrophic changes can be restricted by the suitable choice of handling time. Note that, if the handling time (L) is increased, the chance of maintaining the abundance of resource population at a higher stable equilibrium point is more, and this can minimize the chance of species extinction.
The density regulation parameter (θ) is one of the important structural properties of the species, and it is fixed for a particular location. In comparison, the handling time is controlled through the experimenter. The species with low theta value is more threatened, as our study suggests. Hence, in this case, it is more significant to choose proper handling time. To be more precise, in this work, we have identified the parameter region which allows bi-stability. So, for a given data set, if the estimated values of the parameters fall in this region, the biological resource manager should raise the flag to detect the possibility of a regime shift. This detection becomes urgent if one of the stable states is the extinction state. The resource management in particular the fishery management for marine or lake systems, should detect the early warning signal by applying a generic early warning toolbox or identifying the basin of attraction, or observing the asymmetry of potential function. If anyone or all of the measures indicates the possibility of regime shift, our present work can help to take proper harvesting policy by choosing suitable handling time so that the regime shift can be averted.
We establish our findings through the population time-series data of herring fish population at the Baltic sea and Prince Robert district, British Columbia. The numerical simulation suggests that there may exist a regime shift from the non-zero state to the extinction state for both species in the near future, if handling time is absent or very low. We observe that this regime shift can be avoided if handling time is taken sufficiently large. Let us consider β = −m (say),m > 0. Therefore (x − β) > 0, ⇒ (x − α) < 0, ⇒ x < α. We can assume Remark: We can apply the same process to find M and N for proportional harvesting also.