Stage-Dependent Structured Discrete-Time Models for Mosquito Population Evolution with Survivability: Solution Properties, Equilibrium Points, Oscillations, and Population Feedback Controls

: This paper relied on the investigation of the properties of the stage-structured model of coupled larvae and adult mosquito populations’ evolution when parameterized, in general, by time-varying (or stage-dependent) sequences. In particular, the investigated properties were the non-negativity of the solution under non-negative initial conditions, the boundedness of the sequence solutions under any ﬁnite non-negative initial conditions, the equilibrium points, and the convergence conditions to them in the event that the parameterizing sequences converge to ﬁnite limits. Some further properties that were investigated relied on deriving the oscillation conditions of the solutions under certain conditions of the parameterizations. The use of feedback controls to decrease the foreseen numbers of alive mosquitoes in future evolution stages is also proposed. The proposed control actions are exerted on the birth rate and / or the maximum progression rate sequences. Some illustrative examples are also given.


Introduction
It is well-known from the related background literature that mosquitoes undergo complete metamorphosis following a life cycle of four stages, namely, egg, pupae, larva, and adult. See, for instance, [1] for an easily comprehensible description of those stages and the duration period of each of them. Mosquitoes typically stay in an aquatic environment during the first three life stages and in an aerial environment during their adult stage. Mathematical models are useful to describe the evolution of the mosquito population. The introduction of sterile mosquitoes in the environment in order to reduce the fertility of the whole population and control their populations under acceptable levels of tolerance has been proposed in the background literature. See, for instance, [1][2][3], and some of the references therein. Also, the study of releases of sterile mosquitoes has been investigated in [4], with a Beverton-Holt type model for survivability, as well as in [5]. It is well-known that Beverton-Holt type models are very popular to study the evolution of species that reproduce by eggs, [4,[6][7][8][9][10][11][12] such as 1 + η 1n x n ; y n+1 = γ n x n 1 + η 2n x n ; ∀n ∈ Z 0+ (1) where Z 0+ = Z + ∪ {0} and α n and γ n are the birth rate and the maximum progression rate sequences, respectively, given by α n = α n (x n , y n ) = f n k 1n ; γ n = γ n (x n , y n ) = g n k 2n (2) and define the survival probabilities of larvae and adults as s 1n = s 1 (x n , y n ) and s 2n = s 2 (x n , y n ), given by: s 1n = s 1 (x n , y n ) = k 1n 1 + η 1n = α n f n (1 + η 1n ) ; s 2n = s 2 (x n , y n ) = k 2n 1 + η 2n = γ n g n (1 + η 2n ) whose respective maxima are given by k 1n = k 1n (x n , y n ) ; k 2n = k 2n (x n , y n ) (4) and f n = f (x n , y n ) ; g n = g(x n , y n ) (5) are the number of off-springs produced per adult and the progression rate of larvae or the adult emergency rate, and η 1n = η 1n (x n , y n ) ; η 2n = η 2n (x n , y n ) (6) are density-dependent factors. Note that k 1n and k 2n are, respectively, merged into α n and γ n . Note from (1) that total extinction of the larvae in finite time occurs in two steps ahead from the extinction of the adult if the adult extinguish, since x n = 0 ⇒ y n+1 = 0 ⇒ x n+2 = 0 . However, x n+1 = α n y n 0 if y n 0. However, the extinction of the adult mosquitoes happens one step ahead of the extinction of the larvae, as expected, since y n+1 = 0 if x n = 0.

Remark 1.
Note from (2) that α n can be reduced by external controls on either f n or k 1n , or on both, by reducing, respectively, the off-springs or the survival probability of the larvae, for instance, by spreading aquatic insecticide in lagoons and unsuitable static waters, by drying unnecessary or unsuitable artificial or natural water reservoirs, or also by introducing in the environment sterile mosquitoes. In the same way, γ n can be reduced by external controls on either g n or k 2n , or on both, by reducing, respectively, the survival probability of the adult mosquitoes, for instance, by spreading aerial insecticide in urban drains and bush water lands or by drying unsuitable artificial or natural water depots.
It is assumed through the paper that the sequences {α n } ∞ n=0 , γ n ∞ n=0 , η 1n ∞ n=0 , and η 2n ∞ n=0 are positive and bounded. Note that it is also possible to combine and rearrange the two above equations under the form of Beverton-Holt evolution population structures of the larvae and adult growing rules, while defining "ad hoc" carrying capacity-like and intrinsic growth rate-like parameterizations. The addendum "like" used reflects the feature that the carrying capacity and intrinsic growth rate are artificially introduced here, while having the purpose of describing the evolution of the species as it were related to a single stage, that is, a self-evolution description of either the larvae or the adult population stages. So, with such an equivalent description, the larvae and adult evolution sequences are self-generated by embedding one population stage in the other one. In this way, the mutual couplings are deleted at the expense of joining two evolution steps in each evolution equation. The information of the other species' evolution stage in each case is contained in a carrying capacity-like parameter, which depends on two previous stages, as expected. The larvae sequence self-evolution rule is as follows: x n+1 = µK n x n K n + (µ − 1)x n = α n γ n−1 x n−1 (1 + η 1n x n )(1 + η 2,n−1 x n−1 ) ; ∀n ∈ Z + .
The following result stands for the non-negativity of the carrying capacity-like parameter of the larvae sequence, so as to endow it with a physical and biological sense: (8) is a non-negative real sequence if and only if x n x n−1 ≥ α n γ n−1 (1 + η 1n x n )(1 + η 2,n−1 x n−1 )µ (11) for all n ∈ Z + . If the intrinsic growth rate-like sequence is redefined in (10) as K n = K n (µ n , x n , x n−1 ) via an intrinsic growth rate-like sequence µ n ∞ n=0 by replacing µ → µ n , then the result still holds if µ n ≥ max 1, α n γ n−1 x n−1 (1 + η 1n x n )(1 + η 2,n−1 x n−1 )x n (12) for all n ∈ Z + .
The proof of Proposition 1 is given in Appendix A.
The subsequent result stands for the boundedness, non-negativity, and eventual extinction of the adult mosquitoes and the conditions to get them: Proposition 2. The following properties hold: (i) Assume that µ ≥ 1 and x 0 ≥ 0 are finite. Then, {x n } ∞ n=0 is bounded with x n+1 ≤ x n ; ∀n ∈ Z + if and only if the sequence { K n } ∞ n=0 is such that {x n − K n } ∞ n=0 is non-negative, which is equivalent to α n γ n−1 ≥ (1 + η 1n x n )(1 + η 2,n−1 x n−1 )x n /x n−1 ; ∀n ∈ Z + . A necessary condition for that is that The above results keep valid if µ is stage-dependent with µ n ≥ 1 for all n ∈ Z + .
The sequence {x n } ∞ n=0 is strictly decreasing if and only if {x n − K n } ∞ n=0 is positive, that is, if and only if α n γ n−1 < (1 + η 1n x n )(1 + η 2,n−1 x n−1 )x n /x n−1 ; ∀n ∈ Z + . A necessary condition for that is α n γ n−1 < (1 + η 1n x n )(1 + η 2,n−1 x n−1 ); ∀n ∈ Z + . (ii) Assume that µ > 1 and x 0 > 0 are finite. Then, The proof of Proposition 2 is given in Appendix A. The idea used in the proof of Proposition 2 (ii) consisted of using the population's inverses to build the evolution of the population's inverse sequence was proposed by Stevic [6] when proving the so-called Cushing-Henson conjecture with a new short proof. Basically, note that although the Beverton-Holt equation is nonlinear, its inverse is linear, which facilitates the calculation of the population evolution. The idea was later on exploited in [7][8][9] to introduce control theory tools for such an equation. Basically, the carrying capacity inverse can sometimes be modified by convenience in a control context by using feedback to improve the expected economic results in fisheries, aquaculture exploitations, etc. by acting, for instance, on the water temperature or food delivery in a convenient way. Further results which relate the carrying capacity sequence related to the larvae population evolution as well as solution boundedness properties are discussed in the next result.
, for which sufficient conditions have been given in Proposition 2, and that µ ≥ 2. Then, { K n } ∞ n=0 is a nonnegative real sequence and { x n } ∞ n=0 is bounded if its initial value is non-negative and finite. If µ > 2 then { K n } ∞ n=0 is a positive sequence.
The proof of Proposition 3 is given in Appendix A. Note that Proposition 3 basically establishes that the larvae numbers at each evaluation are non-less than the carrying capacity-like parameter, as it can be expected from the parallel results in usual models based on of the Beverton-Holt equation. See, for instance, [6][7][8][9][10][11] and some references therein.

Stability and Convergence of the Solution
The next result proves the obvious expected fact due to the nature of the problem, that the sequences of populations of larvae and adults jointly extinguish at the same rate in the event that extinction occurs. First note from (1) that which is being used as an auxiliary equation in the proof of the subsequent result. Such a result establishes that the populations of larvae and adults are bounded if they converge asymptotically to zero (extinction). In the event that extinction occurs, it has to jointly happen for both populations of larvae and adult, as expected by heuristic considerations.
n=0 → 0 , while they are also bounded as a result.
The proof of Proposition 4 is given in Appendix A.
The subsequent result establishes that, under nonnegative finite initial conditions, the sequences generated by the two equations of (1) for larvae and adult, respectively, are bounded without invoking related properties on the alternative description of the carrying capacity-like sequence, as it was done in Proposition 2.

Proposition 5.
Assume that x 0 ≥ 0 and y 0 ≥ 0 are finite. Then, {x n } ∞ n=0 and y n ∞ n=0 are non-negative and bounded. In particular, x n+1 ≤ α n γ n−1 η 2,n−1 and y n+1 ≤ γ n η 2n ; ∀n ∈ Z 0+ . The proof of Proposition 5 is given in Appendix A. The subsequent technical result is supported by some calculations of Appendix A and it is of interest to obtain further results later on. In particular, the results emphasize both self-evolution descriptions of the larvae sequence progress and its evolution via couplings to the adult stage evolution. Proposition 6. The following result stands to rewrite both equation (1) with equivalent expressions depending only on one of the larvae or the adult populations.
The proof of Proposition 6 is given in Appendix A. The subsequent result discusses the existence and stability and the nonexistence of a positive equilibrium point. See also [1] and [12] for similar equilibrium studies under invariant parameterizations in this kind of model and [13,14] for studies of convergence of classes of more general discrete sequences to their equilibrium points.
n=0 → η 2 with those limits being positive. Then, the following properties hold: There exists a unique globally stable non-negative equilibrium point (x e , y e ) of (1), if and only if αγ ≥ 1 implying extinction, that is (x e , y e ) = (0 , 0), if and only if αγ = 1. The nonzero equilibrium point is If the equilibrium point is zero then it is not locally asymptotically stable. The positive equilibrium point is locally asymptotically stable if and only if 1 A sufficient condition for Property (iii) to hold is 1 and a necessary condition for such a condition to hold is The proof of Proposition 7 is given in Appendix A. Proposition 7 is now revisited under the removal of the assumption Then, the populations of larvae and adult mosquitoes converge asymptotically to zero at an exponential rate. The result folds if condition 2 is replaced with the stronger one: lim sup n→∞ α n − α n−1 γ n−2 x n−2 (1 + η 1n x n )(1 + η 2,n−1 x n−1 ) γ n−1 x n−1 (1 + η 1,n−1 x n−1 )(1 + η 2,n−2 x n−2 ) ≤ 0.
The proof of Proposition 8 is given in Appendix A. Since the sequences {x n } ∞ n=0 and y n ∞ n=0 are non-negative and bounded, then there exist non-negative real constants X m , X M (≥ X m ), Y m , and Y M (≥ Y m ), such that X m ≤ x n ≤ X M and Y m ≤ y n ≤ Y M ; ∀n ∈ Z 0+ . Those inequalities exhibit boundedness and stability of the populations of larvae and adults, and are also compatible with their potential convergence to the equilibrium points and with eventual presence of bounded oscillations discussed later on. Note also that the constants X m , X M , Y m and are related as follows: Again note from (1) that If g n = x n+1 /x n , one gets: A refinement of the above bound can be obtained by writing x n = g n−1 x n−1 leading to The combined use of the two above upper bounds leads to g n ≤ g nM = α n X m min Y M , γ n−1 X M 1 + η 2,n−1 X m min 1 1 + η 1n X m , 1 1 + η 1n g n−1 X m .
Also, one can get lower bounds as follows: Note also that we can use the expression x n = n−1 i=k g i x k for any k ∈ Z 0+ . This yields the sets of upper bounds and lower bounds: In the following result, one considers potential constant sequences as the trivial case of nonstrict oscillatory sequences in the most general case. For other interesting investigations of oscillatory solutions in discrete sequences, one can refer to [10,12]  (i) Assume that for each given n ∈ Z 0+ , there exist finite k 1 = k 1 (n) ∈ Z 0+ and k 2 = k 2 (k 1 (n)) = k 2 (n)(> k 1 ) ∈ Z + , such that either Then, the mosquito larvae evolution solution {x n } ∞ n=0 is oscillatory.
The solution {x n } ∞ n=0 is strictly oscillatory if in each of the above double conditions at least one of the inequalities is strict.
The oscillatory solution {x n } ∞ n=0 is periodic of, in general, a time-varying (or stage-dependent) oscillation period, if and only if there exists a strictly increasing sequence {k n } ∞ n=0 with k n+1 − k n ≤ k s < +∞, such that (ii) A sufficient condition for the solution {x n } ∞ n=0 to be oscillatory is that for any fixed non-negative integer number or any set of non-negative integer numbers (depending on the available computational possibilities) k(i) ≤ i − 2 for i = n, n + 1, . . . . , n + k 2 and each n ∈ Z 0+ : (iii) Different alternative sufficiency-type conditions for the solution {x n } ∞ n=0 to be oscillatory are obtained with any of the following replacements in the sufficiency-type conditions of Property (ii): g iM (k(i), i) → g 0iM , g 1iM , g 2iM or by the alternative use of any of the combination of minima appearing in their definitions, g im (k(i), i) → g 0im , g 1im , g 2im or by the alternative use of any of the combination of maxima appearing in their definitions.
It turns out that an oscillatory solution of the larvae solution sequence does not necessarily imply an oscillatory solution of the adult sequence.

Worked Examples
Example 1. Oscillatory larvae solution of period two. From (1), one can write equivalently: x n+1 α n ; y n+1 = γ n x n 1 + η 2n x n so that y n+1 y n = α n γ n x n g n so that if x n+1 = g n x n and y n+1 = f n y n then: (a) If g n ≥ 1 then Note that an oscillation in the larvae solution sequence does not necessarily imply an oscillation in the adult solution sequence of the same period.
Assume for instance that g 2n > 1 and g 2n+1 < 1 for any n ∈ Z 0+ . It is possible that So, in this case, the larvae solution has an oscillatory solution of period two, while the adult solution has no oscillatory solution of such a period.

Example 2.
Limit cycle of two levels of the larvae solution for model time-varying parameterization. Assume that for some n 0 ∈ Z 0+ and any n(≥ n 0 + 1) ∈ Z + , one has a solution x n+1 = x n−1 = a x n−2 = x n = b for the larvae evolution then one has from Equation (A2) of Appendix A that and then Therefore, the birth rate subsequence generates a cycle of a period of two stages in the larvae solution sequence on [n 0 + 1 , ∞), given by the If a = b then this subsequence of the solution is, in particular, constant. The above periodic oscillatory sequence for the larvae evolution does not imply that the solution for the adults' evolution necessarily also has a periodic solution subsequence. To this end, additional constraints are needed. Now, assume that for some n 01 ∈ Z 0+ and any n(≥ n 01 + 1) ∈ Z + , one has a solution y n+1 = y n−1 = c y n−2 = y n = d for the adult mosquito evolution. Then, one can conclude from (A.3) of Appendix A that x n+1 = α n y n 1 + η 1n x n ; y n+1 = γ n x n 1 + η 2n x n ; and then γ n = γ n−1 +(η 1,n−1 −η 2,n−1 (1+η 2n α n−1 c))d+η 2n γ n−1 α n−1 c . Therefore, the progression rate subsequence generates a cycle of a period of two stages in the adult solution sequence on [n 01 + 1 , ∞) given by the subsequence y n Note that both cycles are achievable in an independent fashion, each through its respective sequence of gains. If the birth and the progression rates are generated by (17) and (18), then both solutions are periodic of period two stages on the interval max(n 0 , n 01 ) + 1. On the other hand, if the birth rate (17) does not hold for n ≥ n 0 , but instead one has: then, the solution of the larvae evolution {x n } ∞ n 0 +1 converges to a limit cycle. A similar conclusion applies for the solution sequence of the adult mosquitoes if (18) does not hold, while it is replaced with a limit condition for the progression rate.
Example 3. Limit cycle of two levels of the larvae solution for model constant parameterization. If a stage-independent parameterization is considered in Example 2 with α n = α, γ n = γ, η 1n = η 1 , η 2n = η 2 ; ∀n ∈ Z 0+ , then one can conclude from the two equations of (1) in Example 2 that in order to guarantee the non-negativity and boundedness of the oscillation if it exists. Furthermore, direct calculations with the various equations in (19) conclude that Note that: (a) At the equilibrium point, a = b so that η 1 and η 2 can be either identical or distinct. (b) If a cycle of a two-stage period being distinct of an equilibrium point occurs, then a b so that η = η 1 = η 2 .
In the above second case, Equation (19) becomes: Thus, which implies in addition that for the existence of non-negative cycle values, a has to be small enough to satisfy the constraint a < αγ−1 η . Also, As a result, any two-stage period cycle, for both non-negative populations of larvae and adult mosquitoes, that is not coincident with the equilibrium point, has to satisfy η 1 = η 2 , αγ > 1, a < αγ−1 η , and for a given a, b(a), c(a), and d(a) are given by (22)- (24). The extinction case is obtained as a particular case, implying that the cycle coincides with the zero equilibrium point under a nonstrict upper bounding constraint of a leading to 0 = a ≤ αγ−1 η , implying αγ = 1 and b = c = d = 0 from (21) and (22).

Monitored Control of the Mosquito Populations
It has been verified that the intensity of transmission can be determined through the vectorial capacity, equivalent to the basic reproduction ratio of a disease. It describes the total number of potentially infectious bites that would occur from all the mosquitoes in an area biting a single infective human along a single day. The so-called McDonald model is very sensitive to interventions focused on adult mosquitoes. In fact, such interventions cause a reduction in both the probability of vector survival and the ratio of vector to humans. The following items have to be properly identified and fixed for the monitored control intervention: (a) The physical space of intervention; (b) The intensity and periodicity of the intervention, taking into account general details and information, such as seasonality influencing temperature, density, and larval food amounts availability; (c) The appropriate period of the day of intervention, for instance, at nights, for the case of Anopheles, where mosquitoes are more active; (d) The use of large allowed amounts of insecticide compatible with the ecosystem preservation and limited health influence damage.
Details on potential control strategies and protocols to follow can be found in [24,25] and some references therein. The proposed controls are performed by manipulation of the progression rates which depend, in general, on temperature, larval diet and density effects, and, of course, on the type of mosquito [26,27].
The use of controls, subject to feedback information, for the monitored reductions of larvae and adult mosquito populations is developed in this section. Its practical implementation has to take into account the constraints generated by the considerations of the above items a) to c) by the insecticide company's planning.
Assume that for some real sequences α n Then, the interpretation of the necessary and sufficient condition of Proposition 2 guarantees that Since 0 ≤ αγ ≤ 1 and 0 ≤ α n γ n ≤ 1, then the above inequality is trivially true, since any nonpositive real number is less than or equal to some non-negative real number.
The density factors η 1n and η 2n could decrease as the populations decrease by the action of the controls, leading to an increase of the survival probabilities. However, this feature is neglected in the model.

Remark 2.
Note from the subsequent equilibrium point equations obtained in Proposition 7: that the equilibrium point under the controls has a smaller value related to the control-free situation if α is replaced with α c < α and γ is replaced with γ c < γ, being smaller limits of the corresponding gain sequences, assumed to be convergent, through the use of insecticide dropping or the introduction or sterile mosquitoes, under the kept constraint αγ ≥ 1, compared to the control-free case.
Now, Equation (7) becomes modified as follows: Mathematics 2019, 7, 1181 13 of 29 The decreasing of the larvae population at the (n + 1) stage by the control action related to the uncontrolled case is . From the inverses of Equation (1), and provided that x n 0 , y n 0, one obtains the following evolutions of the inverse populations of larvae and adults: Note that if x 0 0 and y 0 0, then x n 0 , y n 0; ∀n ∈ Z 0+ and the above equations are well-posed. The use of the inverses sometimes facilitates the calculation for several stages ahead due to the nonlinear forms of (1), which makes their inverses have a linear evolution, [4,6,9,11]. Note that the above equations can be written more compactly as follows: ; ∀n ∈ Z 0+ so that along two evolution stages, one has: (27) It turns out that for n even, the matrix of the dynamics of the above inverse system is diagonal and for n odd it is antidiagonal. In the following, the solution for n even is addressed a follows. Define the state sequence vector of the inverse system as x I n = x −1 n , y −1 n T , parameterized by: Then, By using recursive computation, one gets: By replacing the subscript n → 2n in (27), one gets: Assume that for some K 2n = K 11 2n K 12 2n K 21 2n K 22 2n , the following constraint holds for some α 2n+1 , γ 2n+1 ; ∀n ∈ Z 0+ or equivalently, such that the following identity holds for some prefixed suitable matrix sequence A * 2n ∞ n=0 : subject to the constraint lim in f n→∞ ∞ n=0 ρ * 2n > 0, which avoids the uncontrolled, unbounded growing of the mosquito population, which is equivalent to the asymptotic convergence to zero of the larvae and adult populations in an infinite number of evolution stages.
Note that A * 2n ∞ n=0 defines the suited evolution of the dynamics of the inverse system for the selection of the maximum progression gains α 2n+i , γ 2n+i ; i ∈ Z + . However, the sequence A * 2n ∞ n=0 is not, in practice, of full design freedom. In particular, it is needed to prevent the environment against untolerated damage due to abusive amounts of insecticide drops. It also has to accommodate for the temperature (a key factor in the abundance of mosquitoes in a certain environment), larval diet supply, and density effects. Therefore, the intensity of the intervention has to take into account the seasonal issues and the control of environment damage under admissible bounds. Note also that if lim n→∞ ∞ n=0 ρ * 2n = 0, then the state sequence of the population's inverse converges to zero so that the populations of larvae and adult diverge to infinity, which should be avoided if possible. Now, define the subsequent subsets of Z 0+ : N 1 = n ∈ Z 0+ : ρ 2n > 1 ; N 2 = n ∈ Z 0+ : ρ 2n < 1 ; N 3 = n ∈ Z 0+ : ρ 2n = 1 and N 3 ∅ if N 1 = ∅ and 2 ∅. Thus, the condition lim in f N→∞ N n=0 ρ * 2n > 0 might be characterized alternatively as follows which is proven in Appendix A: 2n > 0 if and only if there is a real constant ε > 0 such that: A worked example about population controls follows: Example 4. One and two stage ahead monitoring controls. The manipulation of the birth rate and the maximum progression rate sequences as external controls to govern the population growth can be programmed stage-by-stage in the life cycles of mosquitoes or several stages ahead. Assume that: n ≥ 0; ∀n ∈ Z 0+ such that α 0 n and γ 0 n are the control-free nominal values, which can be constant α = α 0 and γ = γ 0 or not. Note that such a kind of maximum saturating constraints α (2) n , γ (2) n are needed to avoid unsuitable, intolerable, excessive environment damage (for instance, to water lands, water reservoirs, nearby bush, agricultural camps, or plantations, wildlife etc.), for instance, due to distribution or overdoses of insecticide. The lower constraints are due to the lack of possibility of physically extinguishing the whole population of mosquitoes, except in extremely small isolated environments.
(a) One-stage-ahead minimization population monitoring: One gets from (1) that (26) holds. So, given the current values of x 2n and y 2n , one concludes that their minimum values at the next stage (i.e., the respective maxima of their inverses) for α 2n ∈ α 2n are the subsequent ones: 2n ] 2n α 0 2n x 2n 1+η 2,2n x 2n (b) Two-stage ahead minimization population monitoring: One concludes from (30) and (31) that, given x 2n , then x 2n+1 ) and then, for already given x 2n and x 2n+1 , one has However, note that 2n ] x 2n+1 As a result, it turns out that the minimization of x 2n+2 given the larvae and adult populations at the stage 2n is obtained via a maximization of the adult population at the stage 2n + 1, which is not coincident with the one-stage-ahead minimization of such a population. One can conclude that the one-stage-ahead minimization strategy does not lead to a sustained minimization through future stages.

Simulation Examples
In this section, three numerical simulation examples illustrating some of the theoretical developments presented in the previous Sections are considered. The examples discuss the existence and location of equilibrium points, the existence of oscillatory solutions to the proposed model, and the use of controls to reduce the populations of larvae and adult mosquitoes.

Non-Negativeness, Boundedness, and Location of Equilibrium Points
Consider the stage-structured model (1) parameterized by the following sequences: with positive initial conditions given by x 0 = 0.2 and y 0 = 0.1. These sequences are generation-dependent and their order of magnitude is similar to that considered in [1]. A periodic behavior has been selected for α n and γ n in order to mathematically describe a seasonal behavior of the birth and maximum progression rates, which is typical in the reproduction cycles of many animals [27]. Thus, these values were selected as a periodic variation around a positive value in order to account for this biological fact. According to Proposition 5, the sequences {x n } ∞ n=0 and y n ∞ n=0 remain non-negative and bounded with an upper-bound given by x n+1 ≤ α n γ n−1 η 2,n−1 and y n+1 ≤ γ n η 2n , respectively. This fact is illustrated in Figures 1 and 2, where the evolution of {x n } ∞ n=0 and y n ∞ n=0 is displayed, respectively. In these figures it can be observed that the elements of both sequences are positive at all times while remaining bounded. In addition, the upper bounds for both sequences are also depicted in Figures 1  and 2       Furthermore, Proposition 7 discusses the location of the equilibrium points along with the conditions under which these exist. In order to numerically check the results claimed in Proposition 7, we consider the following parameterization for the model sequences: Now, all the sequences converge to the positive constant values: Notice that under these values, the condition holds. Consequently, according to Proposition 7 there exists a positive equilibrium point with coordinates: which is locally asymptotically stable. The following Figure 3 depicts the evolution of both populations with the above considered parameterization.
The final values for both populations obtained from simulation are given by x end = 3.3327 and y end = 0.889, which are in great accordance with the theoretical values presented in Proposition 7 and calculated above. Furthermore, it can be observed that the populations are positive and remain bounded as stated in Proposition 5. Moreover, Figures 4 and 5 show the evolution of both populations for different initial conditions given by: Case 1: x 0 = 3 and y 0 = 1; Case 2: x 0 = 4 and y 0 = 3; and Case 3: x 0 = 1 and y 0 = 0.2. As it is shown in Figures 4 and 5, all the trajectories converge to the same equilibrium point, which is locally asymptotically stable, as it is claimed in Proposition 7(iii).        Finally, we now change the value of γ n ∞ n=0 to γ n = 0.2 + 0.1e −0.3n cos(2π0.1n) so that the condition αγ = 1 holds. The initial conditions are x 0 = 3 and y 0 = 1. The Figure 6 displays the evolution of both populations in this case.  In Figure 6 it is shown that both populations converge to the extinction equilibrium point corresponding to (x e , y e ) = (0 , 0). If we now complete a simulation for a large number of generations and make a zoom on the last part of it, we obtain Figure 7. In this Figure it is shown that the trajectories do not asymptotically converge to the extinction point, but they remain close to it since this point is not locally asymptotically stable as, concluded in Proposition 7(ii). Therefore, this Section illustrates some of the theoretical results contained in Section 2 regarding the non-negativity, boundedness, and existence of equilibrium points of the stage-structured model (1).

Oscillatory Solutions
In this Section, a numerical example concerning the oscillatory character of the solution of the stage-structured model (1) is presented. Thus, consider the sequences given by:
In this Section, a numerical example concerning the oscillatory character of the solution of the stage-structured model (1) is presented. Thus, consider the sequences given by:    (1) with this parameterization will be oscillatory. Figure 9 displays the solution of system (1), corroborating the result concluded from Proposition 9. From Figure 8 we can observe that there exist intervals where the sequence takes values smaller than unity, while on other intervals the values are larger than unity. Thus, it is easy to deduce that we can find two constants k 1 , k 2 such that we are in the position to apply Proposition 9 and conclude that the solution to the stage-structured model (1) with this parameterization will be oscillatory. Figure 9 displays the solution of system (1), corroborating the result concluded from Proposition 9.

Use of Controls to Reduce the Larvae and Adult Mosquito Populations
In this last subsection, a simulation example concerning the use of controls to reduce the larvae and adult mosquito populations is considered. Thus, consider the sequences given by:

Use of Controls to Reduce the Larvae and Adult Mosquito Populations
In this last subsection, a simulation example concerning the use of controls to reduce the larvae and adult mosquito populations is considered. Thus, consider the sequences given by: 1n sin(2π0.1n), γ 0 n = 0.4 + 0.1e −0.3n cos(2π0.1n), η 1n = 0.15 1 + e −0.1n , η 2n = 0.15 1 + e −0.15n where α 0 n and γ 0 n are referred to as the nominal control-free values. The evolution of larvae and adult mosquito populations under this parameterization are displayed in Figure 3. As it is observed in Figure 3, both populations converge to finite positive values, which are the equilibrium points. The objective of using controls is to reduce the values of these populations. To this end we follow Example 4 guidelines. Thus, we take α n ∈ α (1) n , α (2) n and γ n ∈ γ (1)     Figure 10 compares the control-free evolution of populations with the controlled case.  As it is observed in Figure 10, the control of the birth and maximum progression rate sequences allows for reducing the number of larvae and adult mosquitos, since the equilibrium point attained by the model is given by lower values than with the control-free situation. Furthermore, assume that the natural dynamics of the system for nominal parameters given by α 0 n → α 0 and γ 0 n → γ 0 make the system converge to the equilibrium point given by so that the achieved values are not biologically acceptable. We would like to reduce them to the values x e < x 0 e and y e < y 0 e . To this end, we would calculate the necessary values of the birth and progression rates α, γ from the above equations, so as to get the desired new equilibrium point (x e , y e ) as: x e y e (1+η 2 x e ) (2η 1 η 2 x e +(η 1 +η 2 )) 2 −(η 1 +η 2 ) 2 4η 1 η 2 The last step is to calculate the amount of insecticide that must be sprayed in the atmosphere to achieve these birth and progression rates. The relationship between the amount of insecticide and the attained values for the rates depends on the spraying method and other environmental characteristics. Once the amount of insecticide is calculated, the effect of overdose and its environmental impact can be specifically considered. In this way, the formulation developed in Section 3 allows for designing feasible values for these sequences. Overall, the presented numerical examples corroborate some of the theoretical results discussed in the previous Sections 2 and 3.

Conclusions
In this paper, the mathematical formulation of a stage-structured model of larvae and adult mosquito populations has been presented. The model is described by time-varying parameters in order to account for their potentially seasonal dependence. This situation is accepted to be of practical importance due to the typical seasonal behavior of reproduction cycles in many animal populations. The work also investigated the conditions under which the model remains non-negative, possesses attainable equilibrium points, and is bounded and/or oscillatory. These conditions are derived in terms of the sequences parameterizing the system or in the terms of its asymptotic limits, when they exist. Consequently, the paper provides useful conditions to determine when equilibrium points exist and where they are located. The value of equilibrium points inform of the steady-state larvae and mosquito populations. When these values are not acceptable for environmental or human reasons, the paper proposes the way to change the birth and progression rates in order to perform control on them. This information is necessary to calculate the amount of insecticide that should be sprayed. It has also to be taken into account the seasonal temperature, the diet stocks for larvae, and the local density effects in order to accommodate the saturated values of the maximum progression rates to be used in the control feedback implementation. follows. From the Beverton-Holt-like equation (7) under the carrying-like capacity (8) and the given intrinsic growth rate-like parameter, one concludes that: ; ∀n(≥ 2) ∈ Z + with x −1 = 0, x 0 > 0 and x 1 = α 0 x 0 1+η 10 x 0 . The above equations delete the first stages n = 0, 1 in order to avoid spurious calculations in the recursive equations, since x −1 = 0 and K −1 0 = +∞. Then, since µ > 1, one has the claimed result of Property (ii).

Proof of Proposition 3.
From Proposition 1, (1+η 1n x n )(1+η 2,n−1 x n−1 )xn . Note that by recalculating α n γ n−1 x n−1 from (8) as dependent on K n and from Proposition 2 that Thus, the conditions of Proposition 1 for In addition, from Proposition 2, {x n } ∞ n=0 is bounded. If the inequality is strict then the carrying capacity is positive for any n ∈ Z + if x 0 is finite and non-negative.
Since αγ > 1 for the nonzero equilibrium point x e = x e1 to be non-negative, one gets that the equilibrium point is positive and locally stable if and only if: and, again, since for the positivity of the equilibrium point it is again needed that αγ > 1, a sufficient condition from the above constraints to guarantee the local asymptotic stability of the equilibrium together with its positivity is: so that a necessary condition for the above sufficient condition to hold is: 1 > √ αγ − √ η 1 η 2 √ αγ−1 η 1 , since αγ > 1. Properties (ii) to (iv) have been proven.
It turns out that the asymptotic extinction of the larvae also implies the extinction of the adult mosquitoes.
The strict oscillatory case and the periodic one follow as appropriate in this particular case. Property (i) has been fully proven. Properties [(ii)-(iii)] are trivial modifications to guarantee the existence of an oscillatory solution by using the various obtained upper bounds and lower bounds of the sequences g n ∞ n=0 .