On a Coupled Time-Varying Beverton–Holt Model with Two Habitats Subject to Harvesting, Repopulation, and Mixed Migratory Flows of Populations

,


Introduction
Discrete-time Beverton-Holt equation-based models, based on the logistic equation, are a common tool for modelling the evolution of species which reproduce by eggs such as insects, birds, and fshes [1][2][3][4].Te most elementary Beverton-Holt equation is parameterized by the carrying capacity of the environment and the species intrinsic growth rate.Te above parameters can be generalized to sequences in order to consider diferent conditions on the population evolution in diferent seasonal periods.Additional parameters such as the independent consumption, associated with eventual disturbances, and harvesting, associated with fshing or hunting, are also included in generalized versions of the equation.Te elementary invariant Beverton-Holt equation possesses two equilibrium points, namely, the extinction one and the carrying capacity level.Under typical evolution conditions of intrinsic growth rates exceeding unity, the extinction equilibrium point is unstable while the nonextinction one is globally asymptotically stable [4,5], implying also the population persistence.If the Beverton-Holt equation is subject to periodic carrying capacities and intrinsic growth rates of the same period then the averaged periodic sequence of population numbers lies below the average of the carrying capacities.Te above, socalled, Cushing-Henson conjecture has been rigorously proved to be true by Stevic [2].Some of their extensions for the q-diference Beverton-Holt equation have been described and proved in [3].On the other hand, control theory issues for Beverton-Holt have been investigated and discussed in [4][5][6].In particular, control actions have been proposed consisting of an online design of the environment carrying capacity for monitoring its suitable sequence of values necessary to match a prescribed suitable population evolution sequence.Tis method is claimed to be potentially applicable in semiopen environments such as some fsheries.It can be pointed out that the intrinsic growth rate and the environment carrying capacity can be related to each other, in real situations, as discussed in [7].On the other hand, it has also been proposed an impulsive extended competition Beverton-Holt model among species from the stability point of view [8].Beverton-Holt models as well as other mathematical related models like, for instance, Ricker's-type models are often used in Maritime Biology for population evolution studies and population management.See, for instance [9].In [6][7][8][10][11][12][13], harvesting actions have been investigated in the framework of extended Beverton-Holt models.Typically, harvesting is related to legal regulated fshing or hunting but illegal poaching might also be included [6,8].Periodicity of the solutions of the Beverton-Holt equation due to seasonality has been focused on in [14,15].On the other hand, a global dynamics analysis of extended models [16], as well as the appearance of bifurcations [17,18] or the presence of resonances [19], has also been investigated.More recently, an extended Beverton-Holt model defned on isolated time scales has been discussed in [20] while an extended Beverton-Holt equation that considers discrete delays has been described in [21].Also, more sophisticated Beverton-Holt-based models are often utilized for monitoring fshing stocks and migration fows in order to elucidate the recommended top of captures [22,23].Te so-called independent consumption, that is, the unforeseen changes of populations due to perturbations can be characterized also under stochastic formulations.See, for instance [24,25] and some of the references therein.Recently, a dynamic version of the discrete Ricker and Beverton-Holt model has been proposed and discussed in [26].Te population evolves according to a density-dependent population contribution and limited replenishment of the stock.Te persistence of the population and the stability properties has been also investigated.Also, it is known that the exchanges of populations are useful for recovery actions of the habitats [27], what is of interest for its regeneration and its eventual exploitation.
On the other hand, evolutionary frameworks for the description and mathematical characterization of the species evolution have received important attention from long time ago.Te pioneering work later on leading to those frameworks in more recent times was undoubtedly Darwin´s book "On the origin of the species" [28], which describes the fundamental evolution principle based on three axioms of natural selection related to variation, inheritance and competition.See also, for instance [29,30].Variation refers to diferent traits of phenotype of individuals belonging to the same population.Roughly speaking, inheritance is the mixture of both parents traits transmitted to the new-borns.Competition primarily refers to the fact that the best transmitted genetic inheritance to the environment make that individuals endowed with them can better survive related to the others.See, for instance [28][29][30][31], and some related references there in.Competition can also refer to the case of various species sharing a certain habitat which are competing for food, refuge for defence of predators, etc. and to the case of predator-prey competition for survival.It can be pointed out that, in the background literature, it has been also usual to develop Ricker's-type models for the characterization of population evolution in the frameworks of interspecifc and intraspecifc competition related, respectively, to competition between individuals or groups of the same or diferent species.See, for instance [32,33], although the related literature is very rich.
Also, Beverton-Holt-type evolution models are also useful for those purposes.In the above context, a discrete evolutionary Beverton-Holt model has been proposed and studied in [29], supported by the evolutionary game framework by considering the strong Allee efect related to predation saturation.Also, Neimar-Sacker bifurcation with chaotic behaviour has been detected in the model of this research, while its asymptotic stability is studied from fxed point theory.Also, a new approach for a specifc class of discrete time evolution models has been proposed in [30] for a class of new single-and multispecies evolutionary competition models by considering the model as a nonautonomous diference equation.
In this paper, two time-varying Beverton-Holt models are proposed and studied which consider that individuals of the same species are living and reproducing in two diferent habitats which are not mutually isolated since they admit populations interchanges.Tose habitats can have diferent characteristics of intrinsic growth rates and environment carrying capacities due, for instance, to be subject to disposal, diferent characteristics of temperature, humidity, available refuge, etc. Mutual fows of population in-between both habitats are considered as well as the existence of harvesting actions, basically, legal hunting or fshing, and/or illegal poaching, and also eventual repopulation actions, those ones being formulated as negative harvesting.Te migrations between habitats as well as the harvesting actions are considered to be proportional to the numbers of individuals while the whole balances of changes of population they generate with respect to the standard population variations, being governed by the intrinsic growth rates and the carrying capacities, are considered as a generalized harvesting contribution.It can be pointed out that exchanges of populations could be useful, in particular, in the context, for instance, of habitat recovery actions, which are receiving certain interest in part due to the efects of climate changes.See, for instance [27], and some of the references therein.Te frst proposed model considers that the generalized harvesting takes place on juvenile individuals so that they do not contribute to the standard population growth at the reproductive stage.In typical real cases, this problem view could describe, for instance, repopulations with either fngerlings or juvenile incorporations to the habitats or by removal of some of them for later exploitations in fsheries.In this sense, the generalized harvesting is viewed as an "a priori" action on the nonadult population which can coexist with the reproductive stage performed by adults.Te second model considers that the generalized harvesting takes place on the adult populations after the reproduction cycle they have performed has ended.It can be viewed as an "a posteriori" generalized harvesting action to the reproductive stage which is able to modify the stocks of population.
Tis paper presents a contribution to the scientifc literature by introducing two interconnected time-varying Beverton-Holt models with the inclusion of generalized harvesting actions alongside classical reproduction stages.Trough a rigorous analysis, we establish the well-posedness of both models, thoroughly explore the stability of equilibrium points, and delve into the oscillation properties of the solutions.Beyond its theoretical signifcance, this research holds practical promise for the management of production fsheries operating under dynamic and diverse conditions, ofering valuable insights that can inform more efective decision-making in this critical domain.
Te paper is organized as follows.Section 2 establishes both models based on two time-varying Beverton-Holt coupled equations being, respectively, eventually subject to generalized "a priori" and "a posteriori" generalized harvesting actions, that is, in parallel to the standard reproduction stage or after such a stage has ended.Section 2 also gives sufciency-type nonnegativity conditions for the solution of both models as well as the characterization of the equilibrium points of both models in stationary regime.Te local stability and the oscillation conditions of Model 2 are also investigated.On the other hand, the characterization and the local stability properties of the equilibrium points and some oscillation properties of Model 1 are investigated in Section 3. Te equilibrium points are shown to be identical in both models although their stability conditions difer.Te boundedness of the solutions, under reasonable conditions on the time-varying sequences of parameters which parameterize the models, is also proved in Sections 2 and 3. Tis fact concludes that, under such conditions, global stability is also achieved.Some numerical simulations are given and discussed in Section 4, and, fnally, conclusions end the paper.
Z + and Z − denote, respectively, the sets of positive and negative integer numbers.R 0+ and R 0− denote, respectively, the sets of nonnegative and nonpositive real numbers.R + and R − denote, respectively, the sets of positive and negative real numbers.
M ≻ 0 and M ≻≻ 0 stand, respectively, for a matrix M with at least one positive entry, respectively all positive entries.

Problem Statement and Main Results
Te parameterizations of Beverton-Holt equations can be given by sequences of intrinsic growth rates and environment carrying capacities in a more general context that periodic sequences associated, for instance, with seasonality issues [26].In the sequel, two time-varying discrete Beverton-Holt models are proposed and studied which consider that numbers of members of the same species are living and reproducing in two diferent habitats or environments.Tose habitats can have diferent characteristics of intrinsic growth rates and environment carrying capacities due, for instance, to diferent characteristics of temperature, humidity, available refuge, etc. Mutual fows of populations inbetween both mentioned habitats are considered as well as the existence of harvesting actions, basically, legal hunting or fshing or illegal poaching, and also eventual repopulation actions being formulated as negative harvesting.It includes also eventual independent consumption efects on the populations.Te migrations between habitats as well as the harvesting actions are considered proportional to the numbers of individuals and the whole balances of changes of population they generate with respect to the standard population variations are considered as a generalized harvesting.Te frst model considers that the generalized harvesting takes place on juvenile individuals so that they do not contribute to the standard population growth at the evolution stage where the standard reproduction takes place.In real cases, this could describe, for instance, repopulations with fngerlings or juvenile incorporations to the habitats or removal of some of them for later exploitations in fsheries.Te second model considers that the generalized harvesting takes place on adults after the reproduction cycle has ended.

Model 1 Subject to "
A Priori" Generalized Harvesting on Juvenile Individuals.Te population of the two groups of the same species evolves according to the following time-varying coupled Beverton-Holt equations: Discrete Dynamics in Nature and Society ∀k ∈ Z 0+ under initial conditions x i (0) ≥ 0 for i � 1, 2. Te populations in the habitats i � 1, 2 at the k-th sampling instant are, respectively, x 1 (k) and x 2 (k).Te sequences are the intrinsic growth rates and the environment carrying capacities, respectively.Te terms h 1,2 (k) involve eventual couplings between both populations and describe the generalized harvesting terms which can be positive, negative, or null.Te coefcients 1 ≥ λ 1,2 (k) ≥ 0 describe the fractions of the respective populations which leave each of the habitats, respectively, 1, 2 to migrate to the other habitat, respectively, 2, 1. Te coefcients δ 1,2 (k) ≥ 0 represent the fractions of immigrated individuals to each habitat from the other habitat, some of them could not reach their fnal new habitat so that the constraints δ 1 (k) ≤ λ 2 (k) and δ 2 (k) ≤ λ 1 (k) hold.Te terms ω i (k)x i (k); i � 1, 2, with positive coefcients, describe the standard harvesting (i.e., fshing/hunting) including the legally regulated one and the eventual illegal poaching.If a coefcient ω i (k) < 0 for some i � 1, 2, and since the generalized harvesting are allocated with minus signs in (1) and (2), that means that the eventual repopulation dominates the eventual standard harvesting action of fshing or hunting.
Note that it is considered, for the sake of simplicity, that the eventual reproduction takes only place over of the existing numbers of individuals in the habitat.In this way, if there is no population in one of the habitats, the repopulation strategy is not activated and, if it is activated, then the numbers of new individuals do not exceed the already existing population.Tus, . Also, note from (1)-( 4) that the populations are nonnegative for any nonnegative initial conditions if the generalized harvesting sequences are nonpositive.Tis can happen as related to repopulation.However, if the generalized harvesting sequences have positive members, those should be upperbounded to keep the nonnegativity of the populations.Te following assumption is made in the sequel related to this concern.
∞ k�0 ⊆ R + are bounded sequences, and which can be expressed equivalently as Note that Assumption 1 is equivalent to. equivalently, Te following result addresses the nonnegativity of the populations under nonnegative initial conditions.It is used that, for all k ∈ Z 0+ and i � 1, 2, it holds that Proposition 2 (Nonnegativity of the solution of Model 1).

Both populations of Model 1 are nonnegative for all samples under nonnegative initial conditions if for any
4 Discrete Dynamics in Nature and Society under the necessary condition which always holds, irrespective of Proof.From (3) and ( 4) in (5): equivalently, which yield directly the proof.Te total generalized harvesting in both habitats is given by k) so that at least one of the ω i (k) is negative implying a dominating repopulation action over alternative population interchanges or fshing/ hunting in the corresponding habitat.Te inverse populations evolve according to where Discrete Dynamics in Nature and Society Note that, for i � 1, 2, h i (k) � 0 ⟺ h I i (k) � 0 as expected.Note also that, as expected as well, if h i (k)≠0 then sgn(h I i (k)) � sgnh i (k) so that a positive (respectively, negative) generalized harvesting leads to a negative (respectively, positive) contribution on the corresponding inverse population dynamics for any while the respective contributions of h i (k + 1) and h I i (k + 1) to x i (k + 1) and x − 1 i (k + 1) are defned with opposed signs in the corresponding formulas.
A further assumption for the Beverton-Holt equation is the following one which is inspired in a usual parallel one for the case of one habitat for the species.For that single habitat case, an intrinsic growth rate exceeding the unity threshold guarantees the species nonextinction.

□ Assumption 3 (Intrinsic growth rate exceeding the unity threshold). Assumption 1 holds with μ
In summary, the standard nonnegativity conditions for Model 1 are 5) to be fulflled by the generalized harvesting sequence.
Note that by weakening the above conditions, one gets that if One also concludes that the (k + 1)− -th population is zero, implying extinction, if K i (k) � 0. Note that in the proposed model the carrying capacities and intrinsic growth rates can be, in general, distinct for both habitats.
Remark 4. Note that the generalized harvesting contributions to the (k + 1)-th sample from the k-the sample in Model 1 assume implicitly that the involved individuals in hunting or fshing of migration or immigration do not contribute to reproduction at the k-th sample.Typical reasons are that they are moving due to migration, they are being eliminated from the population due to hunting or fshing and/or they are introduced being juvenile individuals for stock repopulation but they are not still valid for the species reproduction.Note also that the model relies on the fact that the sampling instants are not instantaneous while they have a certain potentially nonnull duration ε k , that is, the sampling instants are t k � kT + ε k ; k ∈ Z 0 + while the intersample period fulfls T ≫ max k∈Z 0 + ε k .

Model 2 Subject to "A Posteriori" Generalized
Harvesting on Adult Individuals.Compared to Model 1, it is now assumed that the generalized harvesting occurs at the (k + 1)-th sample.Tat is fshing/hunting and/or migrations which can alter the nominal recruitment (i.e., that evolving for null harvesting) occur after the species reproduction period.Te harvesting function now, depends on the (k + 1)-th sample in order to accurately depict the harvesting occurring on adult individuals after the reproduction stage.Tus, if the Beverton-Holt equation represents the efect of the reproduction at the k-the sample on the at the k-th sample, and since harvesting occurs thereafter, this sequence must be defned at (k + 1)-th sample since it does not afect to the reproduction stage at the k-th sample.Contrarily in Model 1, harvesting has been assumed on juvenile individuals so that it would afect to the reproduction stage at the next (k + 1)-th.Terefore, harvesting has been scheduled at the k-th sample for Model 1. Tus, the proposed Model 2 becomes: ∀k ∈ Z 0 + under initial conditions x i (0) ≥ 0 for i � 1, 2. Te above equations can be rearranged as follows: where Remark 5.It is assumed in Model 1 that the harvesting actions are performed on juvenile elements which migrate of their habitat or either they are introduced, for instance, via repopulation, or removed, for instance, under controlled recovery of juvenile elements for later breeding in fsheries or for selling for consumption as it happens with fshing actions on young eels.So, these numbers do not contribute actively to reproduction through the parameterization of the intrinsic growth rate.In Model 2, generalized harvesting, such as fshing/hunting, repopulation, or migration actions, take place on adults after the reproduction period took place.
For presentation convenience, it is frst discussed Model 2 along the next section.

Nonnegativity of the Solution of Model 2. Assumption 1
or its stronger version Assumption 3 given on Model 1 are also kept for Model 2 while in addition, one assumes for Model 2 that: Assumption 6 implies that negative harvesting, for instance, repopulation (of adults) after the reproduction period at some sampling instants requires for its modulus to be less than 1 Te upper-bounding constraint (5) on the generalized harvesting for nonnegativity of the solution now takes the form.
Te nonnegativity of the solution of Model 2 for fnite initial nonnegative conditions is characterized in the subsequent result.It is also discussed the conditions from previous samples to the current to be nonnegative provided that the solution at the current sample is nonnegative.

Proposition 7 (Nonnegativity of the solution of Model 2). Te following properties hold:
(i) Assume that, for each k ∈ Z 0+ , any of the two following sets of constraints holds: (a) Migrations between both habitats with harvesting actions consisting of repopulations:
Proof.First, note that the coefcient matrix M(k + 1) of the left-hand-side of ( 21) is nonsingular if and only if Ten, x 1 (k+1), where where then Note that for each k and ω 1,2 (k + 1) > 0 is positive harvesting (i.e., fshing or hunting) while ω 1,2 (k + 1) < 0 indicates repopulation.If δ i (k) < λ j (k); i, j(≠i) � 1, 2 then a part of the migrated individuals from each of the habi1tat does not reach successfully the other habitat.Ten, (a) Conditions (25) imply that, for any k ∈ Z 0 + , x 1,2 (k) > 0 imply that x 1,2 (k + 1) > 0 with repopulations (ω 1,2 (k + 1) < 0) and nonnull migration/immigration.In this case, adjM(k + 1)≻≻0, i.e., the adjoint matrix has strictly positive entries) and det (26) imply that, for any k ∈ Z 0 + , x 1,2 (k) > 0 imply that x 1,2 (k + 1) > 0. In this case, adjM(k + 1)≻0, i.e., the adjoint matrix has nonnegative entries with the diagonal ones being positive) and det Property (i) has been proved.Property (ii) follows since one has that N(k) is nonsingular and Tis implies that for any k ∈ Z 0+ and i � 1, 2, Note that the constraint of Property (ii) of Proposition 7 implies that Proof.Proceed by contradiction arguments.Note from ( 22) that { }, then one concludes in both cases that 0 { } is not feasible so that if one of the populations diverges the other one diverges too so that their solutions are jointly unbounded sequences.Now, it follows from (32) that the matrix sequence converges to zero as k ⟶ ∞,from (32), while it is less than (1 − ρ) for any prefxed real constant ρ ∈ (0, 1) and all nonnegative integer number k ≥ k 0 for some fnite k 0 � k 0 (ρ).Terefore, if (x 1 (0), x 2 (0)) T is fnite then, from (29), For any given fnite initial conditions.Hence, one reaches a contradiction to Next, the eventual equilibrium points of Model 2 are characterized in the case when the parameterization sequences converge to constant values.It will be seen that the discussion of the feasibility of the nonextinction equilibrium points, in terms of their positivity, is not a trivial task since, in the general case of existence of migration and harvesting, they are the nonunique solutions of a cubic algebraic equation.By introducing some reasonable constraints on the limits of the parameterizing sequences, the existence of at least a feasible nonextinction equilibrium point is proved, and this fact is then corroborated in the numerical simulations described in Section 4.
Assume that all the parametrical sequences in Model 2 are constant and denote the eventual components of the equilibrium point as x 1 and x 2 which are the equilibrium points of each of the two habitats.Note the following features: (a) From ( 17)-( 20), (x 1 , x 2 ) � (0,0) is the joint extinction point of both populations which implies also that null generalized harvestings h 1 � h 2 � 0 as well as null population interchanges due to mutual migrations as well as null fshing/hunting actions.(b) From ( 17)-( 20), (x Discrete Dynamics in Nature and Society Tus, (b1) if δ 1 ≠ 0 then the extinction of x 1 implies that of x 2 and, if δ 2 ≠ 0, then the extinction of x 2 implies that of x 1 (b2) if δ 1 δ 2 ≠ 0 (that is, if there are joint mutual migrations from each habitat to the other habitat) then x 2 � 0 if x 1 � 0 (b3) if δ 1 � 0 then the extinction of x 1 does not necessary imply that of x 2 , but the extinction of x 2 can happen, and, if δ 2 � 0, then the extinction of x 2 does not necessarily imply that of x 1 although such an extinction can happen (c) Te combination of ( 17)-( 20) for a stationary solution (x 1 , x 2 ) yields the following constraints: Now, if δ 1 δ 2 ≠ 0 then x 2 � 0 if x 1 � 0 and, equivalently x ≠ 0 if x 2 ≠ 0. Ten, if both x 1 and x 2 are nonnull then.
with m � x 2 /x 1 so that x > 0 implies that A similar discussion can be applied for x 2 > 0. Te case of numerator and denominator of the second expression in (40) being negative is not feasible since it leads to the contradiction if m � x 2 /x 1 > 0 If both numerator and denominator of that expression are positive then x 2 > 0 under the constraint: Note that if δ i � λ i + ω i then x i � K i for i � 1, 2. Tis means that the equilibrium points in both habitats equalize the respective environment carrying capacities if the mutual interchanges of population coming from the other habitat equalize.Equation (38) may be rewritten as Since m > 0 and x 1 x 2 > 0 for the nonextinction equilibria then.
(44) Discrete Dynamics in Nature and Society By zeroing the frst member of the above relations, one obtain By zeroing the second member of the above relations, one obtains to the necessary condition: Te given derivations are compacted as follows: Proposition 9. Te jointly nonextinction stationary values for both species satisfy the constraint: Another related result for the steady states follows: Proposition 0 (Constraints for nonzero equilibrium points).Te following properties hold: (i) Te jointly nonextinction stationary values for both species satisfy the constraint: (ii) If δ 1 � 0, then, subject to λ 1 + ω 1 < μ 1 − 1, and If δ 2 � 0, then, subject to λ 2 + ω 2 < μ 2 − 1, and (iii) Te components of the nonzero, and, in general nonunique, nonextinction equilibrium point satisfy the constraints: Proof.Te equilibrium conditions might be expressed as A nonzero equilibrium vector is subject to the condition that the above coefcient matrix be rank defective, that is, δ 1 δ 2 � ∆.Property (i) has been proved.
To prove Property (ii), note that, if δ 1 � 0 then 1 and then it has to be solved in x 1 for given x 2 the following equation: whose positive root is the nonzero equilibrium point of the second habitat if and proceeding in a similar way as above, we get the equilibrium point component x 1 if δ 2 � 0 as the positive zero of the second order equation in x 1 for the given equilibrium component x 2 : which completes the proof of Property (ii).
To prove Property (iii), note that the extinction point satisfes (50) and also nonzero equilibrium points have to satisfy (50) so that for some nonnegative scalars α 1 and α 2 with α 1 α 2 � δ 1 δ 2 , the subsequent constraints have to hold in order for the null space of the coefcient matrix not to be identically zero: which lead to with. 1 2 � σ for some real σ > 0. However, note that if α 1 ≠ δ 1 or α 2 ≠ δ 2 while α 1 α 2 � δ 1 δ 2 ≠ 0 then the solution sequence would not be unique so that (61) has to be fulflled for x 1,2 calculated from (61) with α i � α i (x i ) � σ i δ i ; i � 1, 2 which leads to (55).Some more general results follow below without considering the constraint that δ 1 δ 2 � 0 of Proposition 10. □ Proposition .Let x 1 x 2 ≠ 0 and m � x 2 /x 1 .Ten, for any nonzero real ρ, one has: Proof.Note by equalizing to zero the two frst amounts of (38) that the following constraints holds for nonzero x 1 and x 2 : and Discrete Dynamics in Nature and Society Tus, in order for ( 64)-( 65) to have the same solutions in x 1 , the coefcients of x 1 and the independent term have to be mutually proportional with the same nonzero proportionality constant λ satisfying.

□
Proposition 2 (Nonextinction consensus).Te equilibrium point components are identical and nonzero if subject to max(0,1 and it holds Proof.On gets from (38) that if x 1 � x 2 ≠ 0 then which is equivalent to which is compacted as (68) subject to max(0, 1

and (69).
Te local stability of the equilibrium points is now discussed via the convergence properties of the Jacobian matrix at such points.

□ Proposition 3 (Local stability of the equilibrium points). Te equilibrium point
equivalently, if the intrinsic growth rates are lower-bounded, related to the remaining stationary parameters, accordingly to: 14 Discrete Dynamics in Nature and Society In particular, the extinction equilibrium point (0,0) is locally asymptotically stable if the intrinsic growth rates are upper-bounded, related to the remaining stationary parameters, accordingly to min max 1≤i,j(≠i)≤2 Proof.Consider ( 21)- (23), which lead to (32), with � M, for a constant parameterization at the asymptotic limits of the parameterizing sequences Te limit evolution equation ( 22) for the limits of the parametrizing sequences can be expressed as Te Jacobian matrix at the equilibrium point (x 1 , x 2 ) becomes Discrete Dynamics in Nature and Society which leads to after replacing the values of the equilibrium points from (55).Since the spectral radius of a square matrix is less than or equal to any of its norms then, by taking the l ∞ and the l 1 norms in (77), the spectral radius of the Jacobian satisfes equivalently, if (72) or its equivalent (73) hold.Te conditions (74)-(75) are got as particular cases of local asymptotic stability if x 1 � x 2 � 0 since in that case the Jacobian matrix becomes to be so that 16 Discrete Dynamics in Nature and Society guaranteed by (74).

□
Remark 14 (Oscillation condition).Te oscillation condition for a period of N samples is easy to characterize from (32).Note that if there exists a set of N values of the solution , then such values defne an oscillation of the solution where and I 2 is the second-order identity matrix.

Further Stability and Oscillation Properties of Model 1
Te equilibrium points of Model 1 are the same as those of Model 2 since they satisfy also the same stationary constraints Ten, both generalized harvesting sequences take identical steady-state values at the equilibrium points.However, their evolution dynamics are distinct and the local asymptotic stability conditions of the equilibrium points difer as well.Terefore, there are some variations on the sufciency-type conditions related to guaranteeing the basic properties.Note that, for Model 1, equation (32) becomes modifed as follows: x 1 (k+1) Proposition 5 (Local stability of the equilibrium points).Te equilibrium point (x 1 , x 2 ) of Model 1 is locally asymptotically stable if and only if the spectral radius of the Jacobian matrix at (x 1 , x 2 ) defned by Discrete Dynamics in Nature and Society is less than one while the evolution sequence of populations is nonnegative, that is, if In particular, the extinction equilibrium point (0,0) is locally asymptotically stable if Proof.Te equilibrium point (x 1 , x 2 ) is locally asymptotically stable if the spectral radius of the Jacobian matrix is less than one while the evolution sequence of populations is nonnegative.By simple inspection, both conditions hold if (84) holds or if (85) holds under similar arguments as in the proof of Proposition 13.Te extinction equilibrium point (0,0) is locally asymptotically stable if
Remark 22 (Global stability).Note that Proposition 8, related to Model 2, and Proposition 18, related to Model 1, conclude that, under reasonable conditions of the timevarying sequences which parameterize the models, the global stability is also achievable since the solution sequences are globally bounded for all samples.

Numerical Simulations
Tis section contains three examples that illustrate some of the main theoretical results proved and discussed in the previous sections: , for i � 1, 2, are, respectively, generated by means of the following diference equations: related to the habitat 2. Te following initial conditions are considered for obtaining the evolution of the sequences generated by (91).Namely, the resulting sequences are for i � 1, 2 and ∀k ∈ Z 0+ .In this way, for instance, the sequence μ 1 (k)   ∞ k�0 starts with the value μ 1 (0) � 1.5 and it converges to the value ρ μ1 /(1 − ε μ1 ) � 4 as k tends to infnity.Moreover, the sequence μ 1 (k)   ∞ k�0 is monotonically nonincreasing since μ 1 (0) < ρ μ1 /(1 − ε μ1 ).Te same result can be said about all the sequences defned by (91) from the choice of the values for the parameters in such equations and the initial conditions for μ One can see in Figure 2 that q i (k) ≤ p i (k), for i � 1, 2, ∀k ∈ Z 0+ which guarantees the nonnegativity of the solution of Model 1 as Proposition 2 establishes provided that μ i (k) > 1∀k ∈ Z 0+ .In summary, Figure 1 shows that x i (k) ≥ 0, for i � 1, 2, ∀k ∈ Z 0+ , accordingly to the result of Proposition 2 since the conditions q i (k) ≤ p i (k) and μ i (k) > 1, for i � 1, 2, and ∀k ∈ Z 0+ , are satisfed.

Example 2. Tis example illustrates the results of Proposition 7. Te sequences μ
, for i � 1, 2, are, respectively, generated by means of the following diference equations: with the following values for the parameters: related to the habitat 1 and the parameters: related to the habitat 2. Te following initial conditions are considered for obtaining the evolution of the sequences generated by (96). Figure 3 shows the evolution of the species populations x i (k)   ∞ 0 , for i � 1, 2, if the population are initially x 1 (0) � 200 and x 2 (0) � 250 while Figure 4 displays the evolution of s One can see in Figure 4 that: Which guarantees the nonnegativity of the solution of Model 2 as Proposition 7(a) establishes.As a consequence, the fact that x i (k) ≥ 0, for i � 1, 2, ∀k ∈ Z 0+ is shown in Figure 3 accordingly to the result of Proposition 7.

Example 3.
Tis example illustrates the results of Proposition 10.All the parametrical sequences associated to both habitats are considered constant, namely:  17)-( 20) by taking into account that x i (k + 1) � x i (k) � x i for i � 1, 2 at the equilibrium point.In this way, one obtains that From (101), with i � 1 and j � 2, one obtains by direct calculations that From (101), with i � 2 and j � 1, and introducing (102), one obtains that the population of the habitat 1 at the potential equilibrium point is one of the roots of the following cubic equation: Discrete Dynamics in Nature and Society 23 , Only the real and strictly positive roots r i , for i � 1, 2, 3, of (103), if they exist, have signifcance in a population model.Ten, the potential equilibrium point really exits if r i , for i � 1, 2, 3, is real and positive and the population for the habitat 2 obtained from (102) by substituting x 1 � r i , with r i being real and positive, in (102) is also real and positive.Such a fact is fulflled if: By considering the values for the parameters in (100) the roots of the polynomial (103) are the following: (106) Only the value r 1 satisfes the condition (105) since (K 1 (μ 1 − 1 − λ 1 − ω 1 ))/((1 + λ 1 + ω 1 )(μ 1 − 1)) � 1571.4.Ten, there is only an equilibrium point where the population in each habitat is, respectively, given by x 1 � 8164.7 and x 2 � 7390.4.Such a result can be seen in Figure 5 where the evolution of the populations of the habitats 1 and 2 are displayed in Figure 5 when the initial condition is x 1 (0) � 400 and x 2 (0) � 100.Tis result is in accordance with the result in Proposition 10(iii) related with the convergence of the model 2 to a nonzero equilibrium point.

Conclusions
Tis paper has proposed and discussed two discrete coupled time-varying Beverton-Holt which consider that individuals of the same species are living and reproducing split into two nonisolated diferent habitats, in fact, admitting mutual populations interchanges, standard harvesting (fshing, hunting, illegal poaching, etc.), and negative harvesting being associated to repopulation actions.Both habitats can have diferent parameterizing characteristics of intrinsic growth rates and environment carrying capacities due, for instance, to their diferent characteristics of temperature, humidity, available refuge, etc. Te migrations between habitats, as well as the various harvesting actions, are considered to be proportional to the available numbers of individuals while the whole balances of changes of population they generate with respect to the standard population variations, being governed by the intrinsic growth rates and the carrying capacities, are jointly considered as a generalized harvesting contribution which can have a positive, null, or negative global balance at each sample.It can be pointed out that exchanges of populations could be useful for preservation or exploitation, in particular, in the context, for instance, of habitat recovery actions, which are receiving certain interest in the last years.Te frst proposed model relies on a generalized harvesting which takes place on juvenile individuals so that they do not contribute to the standard population growth at the reproductive stage.In typical real cases, this approach might describe, for instance, repopulations with either fngerlings or juvenile incorporations to the habitats or by removal of some of them for later exploitations, for instance, in fsheries.In this sense, the generalized harvesting is viewed as an "a priori" action on the nonadult population which can coexist with the reproductive stage performed by adults.Te second model considers that the generalized harvesting takes place on the adult populations after each reproduction cycle they have performed has ended.It can be viewed as an "a posteriori" action to the reproductive stage which is able to modify the stocks of population.Te equilibrium points of the stationary solution in the presence and absence of harvesting action have been formally characterized as well as their local asymptotic stability properties in the case of intrinsic growth rate exceeding unity and eventual execution of harvesting actions.Te equilibrium points are identical for both models although their respective stability conditions can be distinct.It has been found that three nonextinction equilibrium points can exist of which, at least one, is feasible under reasonable constraints of the immigration levels related to the remaining model parameterizing parameters.Some numerical examples have also been discussed.Te restriction to two habitats of the proposed models has been done for exposition clarity but the formal extension of the results to N > 2 habitats is direct.

Figure 1
shows the evolution of the species populations x i (k)   ∞ k�0 , for i � 1, 2, if the population are initially x 1 (0) � 200 and x 2 (0) � 250 while Figure 2 displays the evolution of q
Boundedness of the Solution of Model 2. It is now discussed the boundedness of the solution of Model 2 under some weak reasonable constraints: ≥ ε > 0 for i � 1, 2. Assume that Assumption 1, with eventually its stronger version Assumptions 3 and 6 hold.Assume, furthermore, that both habitats are mutually coupled in the sense that δ 1,2 (k) > 0 for infnitely many samples.Ten, the sequences x 1 (k)