Potential-Growth Indicators Revisited: More Merits of Indication

: The notion of potential-growth indicator came to being in the field of matrix population models long ago, almost simultaneously with the pioneering Leslie model for age-structured population dynamics, albeit the term has been given and the theory developed only recent years. The indicator represents an explicit function, R ( L ), of matrix L elements and indicates the position of the spectral radius of L relative to 1 on the real axis, thus signifying the population growth, or decline, or stabilization. Some indicators turned out useful in theoretical layouts and practical applications prior to calculating the spectral radius itself. The most senior (1994) and popular indicator, R 0 ( L ), is known as the net reproductive rate, and we consider two more ones, R 1 ( L ) and R RT ( A ), developed later on. All the three are different in what concerns their simplicity and the level of generality, and we illustrate them with a case study of Calamagrostis epigeios , a long-rhizome perennial weed actively colonizing open spaces in the temperate zone. While the R 0 ( L ) and R 1 ( L ) fail respectively because of complexity and insufficient generality, the R RT ( L ) does succeed, justifying the merit of indication.


Introduction
The concept of Potential-Growth Indicator (PGI) was developed in the theory of matrix population models (MPMs) for the dynamics of discrete-structured biological populations [1,2]. This kind of model represents the basic tool in the mathematical demography of plant and animal populations structured with regard to a certain classification trait such as the age, size, or developmental stage of individuals in a local population of a given species [1,2]. Mathematically, the MPM is a system of difference equations, x(t +1) = L(t)x(t), t = 0,1, 2, ..., (1) where the vector of population structure, x(t)  , belongs to the positive orthant of the n-dimensional vector space and L(t) is a nonnegative nn matrix called the population projection matrix (PPM) [1,2]. Each component of x(t) is the (absolute or relative) number of individuals in the corresponding status-specific group at time moment t, while the elements of L, called vital rates [1], carry information about the rates of demographic processes in the population.
The pattern of (i.e., the allocation of zero/nonzero elements in) matrix L corresponds to the associated directed graph [3], which is called the life cycle graph [1] (LCG) as it represents graphically the biological knowledge of life histories involved into the model and the way the population structure is observed in the field or laboratory (Figure 1 gives an example). The LCG may be strongly connected [3], signifying certain integrity of the individuals' life history and providing for the PPM being irreducible gt, terminally generative plants, the stages being distinguishable in the field. Solid arrows indicate transitions occurring for one year (no transition, in particular); dashed arrows correspond to the annual recruitment [7].
In theoretical layouts and practical applications, it turned out useful to consider the PPM as the sum, of its parts responsible respectively for the transitions (T) between individual statuses and population recruitment (F) [8,1]. In particular, we have for the LCG shown in Figure 1. When calibrated reliably from data, matrix L(t) gives rise to a rich repertoire of qualitative properties and quantitative indices characterizing the population under study at the place where and the time when the data were mined [1,2]. In particular, its dominant eigenvalue 1(L) > 0, a positive eigenvalue coincident with (L) (the spectral radius of L) and existing by the classical Perron-Frobenius theorem [4][5][6], gives a quantitative measure of how the local population is adapted to its environment [7], thus serving as an efficient tool of comparative demography [1] and enabling a forecast of population viability [9]. This ability ensues from the dynamics of trajectory x(t) as t tends to + when (a primitive [5]) L(t) = L does not change with time: where x* is a positive eigenvector corresponding to 1 [6]. Therefore, we have so that the location of 1 relative to 1 is of principal importance.
Of biological importance is one more characteristic, namely, the net reproductive rate Preprints (www.preprints.org) | NOT PEER-REVIEWED | Posted: 1 June 2021 doi:10.20944/preprints202106.0029.v1 which turned out meeting all of the following conditions [8,1,10]: where  means 'if and only if'. Matrix L belongs a set of PPMs, i.e., it is representable as the sum (2) where T is substochastic (none of the column sums exceeds 1) but not pure stochastic (each of the column sums equals 1). Remark 1. Since expression (6) suggests the convergence of the power series in T, it requires a stronger restriction on L than just T being substochastic. The convergence requires which mathematically is equivalent to (T) < 1 and biologically means the lack of immortal individuals in the population. . Equivalence (7) was later called indicator property [11] (as it indicates the location of 1 relative to 1). Any function R(L) that meets conditions (7), or equivalently was called a potential-growth indicator (PGI), and an indicator R 1 (L) was proposed much simpler than R0(L) (6) [11]: where I denotes the identity matrix.
As an immediate consequence from the theory of rank-1 corrections of nonnegative matrices [12], the indicator property of R 1 (L) was proved for the class of PPMs in which the fertility matrix F has rank 1. A particular case is cited in (3), meaning a single stage where the population recruitment appears (Figure 1), or F may have only one nonzero column meaning a single reproductive stage in the life cycle ( Figure 2 and matrix (3) in [13]).
These two particular cases are however not exhaustive among PPMs, and since R 1 (L) has revealed certain "merit of indication" [2] motivated also by the current modeling practice [14], the challenge is to expand R 1 (L) to a wider class of PPMs or to propose another PGI simple enough in calculation. Aimed at responding this challenge, we present here a construction of a more general PGI by means of a general algorithm developed recently by V.N. Razzhevaikin and E.E. Tyrtyshnikov [15,16] for nonnegative matrices. We expose this algorithm in Section 2. Section 3 reports on a case study of Calamagrostis epigeios, a perennial herbaceous graminoid, a traditional object to study in the Russian botanical school [17], with the PPMs beyond the class of rank-1 F. "More Merits of Indication" by PGIs are thereafter discussed in Section 4.

Razzhevaikin-Tyrtyshnikov's Algorithm to Construct a PGI
Note first that R1(L) (9) cannot be expanded for an L with rank(F) = 2, which is illustrated in Section 3.3. Therefore, the task to develop a more general PGI is really urgent, and a more general framework for that [16] is presented below. The above notations are standard, and we introduce some original ones below. Definition 1. We call renumbering in a square matrix any permutation of its rows and the same columns.
Corollary 1. Renumbering in a matrix does not affect its spectrum, hence the spectral radius, too. Matrix A being irreducible is then equivalent to the lack of a renumbering that might transform A to a block-triangular form with two or more nonempty diagonal blocks.

Definition 2.
A submatrix B of matrix A is called a principal submatrix if it is formed of the elements that are symmetric about the principal diagonal. There is a renumbering in A that places B into the lower right position. When the renumbering is needless, the principal submatrix is called the lower principal submatrix.
Definition 3. Let f1 and f2 be two functions defined on the same domain  and taking values in the spaces of real matrices sized n1n1 and n2n2, respectively. We call f1 and f2 to be R-equivalent if the equality holds true for any   .
Remark 2. Domain  determines certain sets within the spaces defined by the elements of f1 and f2. In particular, if those elements depend (explicitly or implicitly) on some parameters, then  can be regarded as a set in the space of the parameters.
Remark 3. Definition 3 reduces to the PGI definition (7) when  is the domain where R(L) does possess the PGI property and n2 = 1.
The following definition generalizes the PGI notion to a general set . When  is a cone of nonnegative matrices Definition 4 reduces to the PGI notion (7) by the classical Perron-Frobenius theorem for nonnegative matrices [5].
Theorem 1. Let the identity function, fk: The proof is given in [16], as well as those of the other theorems providing for the construction of a stability indicator for an arbitrary nonnegative matrix.
If the theorem conditions are met at each successive step of the dimension-reducing chain: then function f1 = a 11 1 can serve as a stability indicator for matrix A ≥ 0. In what follows, we refer to this indicator as RRT(A).
The following obstacles may occur to perform the next step in the chain (12) when applied beyond the terms of Theorem 1.
2) The next a kk k = 1 and a ii k ≠ 1 for some i < k. In this case, a proper renumbering in matrix Ak relocates a ii k to the lower right angle and, if a ii k < 1, we can continue the chain, but if a ii k > 1, we come to point 1) above.
3) All a ii k = 1, i ≤ k. Two subcases are possible here.

3.1)
There is an irreducible principal submatrix B of Ak having dimension greater than one. In this subcase, (Ak) ≥ (B) > 1, so that (An) > 1 by virtue of Theorem 1.

3.2)
There exists a renumbering in Ak that brings it to a triangular form with units on the principal diagonal. In this case, (Ak) = 1, so that (An) = 1 by virtue of Theorem 1.
When the terms of Theorem 1 are met at every step of chain (12), the resulting indicator has the form of a multi-storey fraction of the matrix elements. This fraction can be simplified by means of the following downgrade procedure: we replace the fraction gained at the current step in the form of +1. This procedure eventually enables getting an indicator in the form of a polynomial in the matrix elements.
It turns out (see [16]) that this polynomial is equal, in the above case, to 1 -p(1; A) = R 1 (A) . The aforementioned alternatives allow a similar consideration applied to the lower principal ll submatrices Bl of matrix A used to calculate the elements a kk k with k = n + l -1 to be compared with 1 at the l-th step in the chain (12). In the robust case, i.e., for a kk k > 1 and a jj j < 1 for j < k, the inequality 1 -p(1; Bl) > 1 holds true, while the inequalities a kk k < 1 holding for all k ≤ n imply that RRT(A) < 1. Here, we have where B runs the chain of n lower principal submatrices Bj, j ≤ n, of matrix A (beginning with the matrix itself). Thus, expression (13) is appropriate as a stability indicator in case if the strict inequalities hold. Since renumbering in a matrix does not affect its spectrum, one can also use a chain of upper principal submatrices in (13).
Problematic cases arising for a fixed numbering of rows and columns in the construction of a stability indicator of matrix A due to the equality a kk k = 1 holding true at the next step force us to use the same expression instead of expression (13), albeit with B that runs through all possible principal submatrices of matrix A (there are 2 n -1 such Bs, see [16] for greater detail).

Material and Method
Reed grasses of the Calamagrostis genus are perennial herbaceous polycarpic graminoids growing in a wide range of ecological conditions, mainly in forests and meadows, a traditional research object in the Russian school of population botany. Their ontogenies are well studied including that of Calamagrostis epigeios [17,19], where the ontogenetic stages of individual plants are detectable from the morphology of the above-ground parts. Moreover, the chronological age (in years) of the plant can also be determined from the annual increments in its tillering zone [19].
As a pioneer species of plant succession, C. epigeios invades open plant communities, such as spruce clear-cut areas, due to extensive population growth by means of vegetative propagation. However, the mosaic of its colonies (systems of clonal plants) has not yet formed a continuous cover by the initial years of colonization, so that the colonies have distinguishable boundaries [21]. Annual observations (late August, when the herbs complete their development) show that different plants spend different numbers of years at the same stage of ontogeny, i.e., there are individuals of different ages among the plants at the same stage. The diversity of individual pathways of development among the plants within a single colony can be represented in the form of LCGs shown in Figure 2. These graphs are constructed on a two-dimensional integer-valued "lattice" of all possible states that an individual plant may have in terms of age-stage. Botanists call this diversity polyvariant ontogeny and consider it as a major mechanism of plant adaptation to the environment [22,17,19,23]. The numerical values of ontogenetic transition rates assigned to the solid arrows in Figure 2 were obtained from the data of "identified individuals" type [1, p. 134], while the reproduction rates are shown on the dashed arcs of the LCG as symbolic parameters to be calibrated rather as a tribute to the tradition of reproductive uncertainty, i.e., inability to detect the status of parents for each recruited plant, hence inability to calibrate the status-specific reproduction rates in a unque way [24,19]. In fact, the 2017 project was aimed to get rid of the reproductive uncertainty inherent in the analysis of aboveground vegetation [ibidem], and for this purpose the whole reed colony (two of them in fact) was digged up, with preserving the entire system of rhizome parent-offspring links. The fragments of the linkage system related to each parent plant (Figure 2 in [21]) were analyzed and summerized over the colony, thus resulting in a total picture of paret-offspring links in a colony. Figure 3 represents the diagrams of offspring survival for all the age-statge-specific parental groups in each C. epigeios colony [21]. The integer-valued parameters a, b, ..., n, o are the numerators of the reproduction rates that are shown on the arcs of the LCG (Figure 2), and the denominators are the (absolute) sizes of the corresponding parental groups in 2014, i.e., one year before the excavation ( Table 2 in 20]). Thus, the vector-matrix equation (1) takes on the form of

Calibration of the C. epigeios PPMs
which is satisfied for each of 11 (for Colony I) or 14 (for Colony II) components of vector x(t) with matrix L presented in Table 1. Unlike equation (1), which describes the "projection" of the population

II I
structure one step into the future, equation (14) describes the transition to the current state of the colony from its last year's state, retrospectively restored based on the excavation data. But this circumstance can no way affect the properties of the model, and in particular, the dominant eigenvalue, 1(Lrec), of the reproductive-core submatrix of matrix L (outlined by dashed lines in Table   1, with deleting unnecessary rows and columns) still serves as a quantitative measure of the colony's adaptability to environmental conditions at a given time step. For comparison purposes, this measure is modified by a "scaling factor" S0/S1 accounting for the expansion of colony area in one year [21]. Calibrated without any reproductive uncertainty, the C. epigeios PPMs give rise to population characteristics of certain biological interest [21]. In particular, the values of 1 much greater than 1 ( Table 1) witness the high speed of population growth in the young colonies.
As we can see from Table 1, the fertiity matrices have rank 2 and 3 for Colonies I and II respectively, and we therefore cannot apply the PGI R1(L). Expression (13) for the RRT(L) produce explicit formulae for Colonies I and II (Appendix A, (A6), (A7)), which result in the numeric values shown in Table 1 just to illustrate the PGI property.

Calibration from the above-ground data alone
The excavation of colonies is a much more laborious method of field study than the above-ground observations of marked individuals on permanent sample plots bearing reproductive uncertainty [7,20]. Though revealing valuable information on the clonal rhizome system, the excavation destroys the field object, thus ceasing its further study. Therefore, to develop reliable methods to cope with the reproductive uncertainty [25,7] is still an urgent challenge. The question is how the PGI (if found in a relevant form) could help respond this challenge.
It follows from Figure 1 that the reproductive-core submatrices for Colonies I and II are respectively Although the ontogenetic-transition matrices T1 and T2 were calibrated from the excavation data [21], they could also be calibrated, in a unique way, from the above-ground data of the "identified individuals" type [1,20]. These data are however unable to determine the status-specific reproduction rates (the elements of F), leaving them as uncertain parameters, which have only to satisfy certain linear constraints. These constraints express the observation that population recruitment at a status-specific group consists of the contributions by parental plants at several generative statuses, i.e., { + + + + = 163, + + ℎ = 22 (17) for Colony I and { + + + + + = 207, + ℎ + + + = 54, for Colony II.
It is known that during the vegetative propagation of long-rhizome grasses, the young rhizomes of the first-year virginal plants develop most actively [26] and this activity decreases over the years.
The following hierarchy of inequalities expresses this expert knowledge for Colony I: where C(…) denotes the total number of alive daughter plants produced by the parents of this status-specific group. A similar hierarchy for Colony II looks as follows: Preprints (www.preprints.org) | NOT PEER-REVIEWED | Posted: 1 June 2021 doi:10.20944/preprints202106.0029.v1 . (20) In terms of our integer-valued parameters, these hierarchies take on the following forms: for Colony I, and for Colony II. Now, the basic equation (14) is satisfied with the 8-xparameter family of PPMs constrained by (17), (21) for Colony I, and with the 13-parameter family constrained by (18), (22) for The integer-valued feature of these families lead us to Diophantine systems, which may have a lot of solutions in our case, following the number of ways to split a given sum into a given number of summands (diminished by the additional constraints) multiplied with similar numbers for the second given sum and for the third. But the principal point is that this number is finite (no matter how great), which means that all these solutions can be found by enumerating all the feasible options. Every solution thus found, i.e., a set of nonnegative integer numbers a, b, ..., n, o, satisfying the system of equations and inequalities, forms its own reproductive-core submatrix Lrec (as shown in Table 1) with the corresponding value of 1(Lrec), and we find those ones among (a finite number of) such matrices that give 1min and 1max, the minimum and maximum value of 1, respectively. Thus, we obtain a certain range [1min, 1max] for uncertain 1(Lrec), or the fitness bounds (lower and upper) for the local population of a given clonal species with polyvariant ontogeny ( Table 2).
To illustrate the relationships of 1s with various potential-growth indicators, Table 2 shows also the corresponding PGI values. In particular, we see that 1 < 1(Lrec) < R0(Lrec) in both Colonies, thus illustrating Theorem 3.3 of Li and Schneider [10]. The negative values of R1 demonstrate the PGI inability of R 1 (L) when L has a rank-2 or rank-3 fertility matrix F. Two extreme columns to the right show that 1 and RRT attain their extreme values at quite different parameter points.

Discussion
When models look complicated, like those shown in Figures 2, 3, the additional attributes like PGIs are of logical use, rather than attempting to simplify the model, for instance, by aggregating several variables. However, the complexities shown in Figures 2, 3 reflect the reality of empirical data, while simplifications may distort the real picture. For example, although the LCG shown in Figure 1 is fairly simple, it bears what we call "reproductive uncertainty" [27, p. 43], i.e., the inability to calibrate the status-specific reproduction rates, a and b, from field data because only the total contribution to the population recruitment by the two generative groups, g and gt, is observed in the field. As a result, we could only obtain a certain range of uncertain 1 values rather than a particular number [28]. Among various, poorly efficient, ways to get rid of reproductive uncertainty, there was an attempt to aggregate the two generative groups into one, thus eliminating the very basis of uncertainty and providing for the PPM calibration in a unique way [29]. A logical was to obtain the 1 of the aggregated model within the range of the original one, and this happened indeed for several observation years. But there were also a pair of years (2013, 2014, Table 3, ibidem) that gave the aggregated 1 outside the range. Moreover, it was greater than 1, whereas the range was located entirely to the left of 1 (Table 2, ibidem), i.e. the aggregation resulted in a principally wrong forecast of population viability.
integer-valued calibration (Section 3.3) turns out resource-consuming, too. Here, we expected another merit of indication for saving computer time when solving the Diophantine systems (17), (21) and (18), (22) by means of enumeration. Computer algorithm to find all the solutions for Colony I (8 unknowns, Appendix A) poses no problem, but that for Colony II (13 unknowns, Supplementary Materials) does. The corresponding 13 nested cycles contain 57×14×29×9×4×24×50×25×9×4×24×6×6 = 7.773910 14 repetitions of the inmost check, and a crude estimate suggested 14 to 24 days of continuous operation with a 2.3 GHz Intel Core i9 processor. Since the Colony was growing fast, there was an idea to preliminary check for RRT(L2rec) > 1 and exit otherwise in order to accelerate the inmost cycle. The idea has however failed once we have found that RRT(L2rec) > 1, hence 1(L2rec) > 1, already at the lower bounds of all the parameter ranges (see (A7) in Appendix A). Since 1 is a monotone function of matrix elements [5], it follows that 1(L2rec) > 1, hence RRT(L2rec) > 1 over the total ranges of parameters, and this makes the check useless in each particular case.
However, the idea was fruitful, in this regard, when the task was to investigate the effects that uncertain parameters of the seed bank exert on the model outcome and to find the bounds of uncertainty in a case study of Androsace albana, another alpine short-lived perennial [14]. The R 1 (L) turned out linear w.r.t. the uncertain parameters, and this linearity caused a great 'merit of indication' while checking whether  1 (L) be less or greater than 1.
Similarly, the calculation of RRT(A) has turned out fruitful in the evolutionary optimality problems [33] where (L) and RRT(L) were considered as selection criteria and an additional condition, (L) = RRT(L) = 1 held true. The condition guaranteed the optimal values being reached at the same sets of parameters -unlike the extremal points in our current problem ( Table 2) -thus reducing the optimal search to that for a polynomial function.
Also, the use of simple algorithms was highly efficient to solve stability problems for highly sparse nonnegative matrices [34].

Conclusions
The history of PGIs in matrix population models, which began with a simplest expression for the Leslie matrix three quarters a century ago [Les45], has now come to its logical end with the explicit formulas for R0(L), R1(L), and RRT(L), where L is a PPM representable as L = T + F (2). The theme has been developing along two directions, the first direction tending to simplicity, the second one to generality.
While R0(L) (6), the net reproductive rate, confines only to the convergence of a series in the definition (6), hence (T) < 1, thus leading (till 2020) in the generality, it has been loosing in simplicity as the spectral radius represents a much more difficult computational task than the determinant in R1(L) (9) or a polynomial of matrix elements in RRT(L) (cf. (A4-5) and (A6)). While R1(L) provides for the maximal simplicity, its PGI property holds true only when rank(F) ≤ 1 [12], whereas practical applications may suggest rank(F) ≥ 2 (cf. (15)(16)); [20]). The simplest PGI thus has the weakest generality.
The dimension-reducing algorithm by Razzhevaikin and Tyrtyshnikov [15,16] has reached the ceiling in generality as it can be applied to any square nonnegative matrix. Although, it fails in the monotone PGI property (when inequality R(L1) > R(L2) being equivalent to 1(L1) > 1(L2)), so that the extremal points over a compact set are different for R(L) and 1(L) (see e.g., Table 2), the indicators can still be useful for solving constraint optimization problems [7].
To conclude, the current theory of PGIs offers three explicit functions -R0(L), R1(L), and RRT(L)of matrix elements, which complement the mathematical toolbox of matrix population models. At certain levels of generality and simplicity, these PGIs are able to cover a vast variety of particular situations. Particular problem formulation, be it a calibration task under parameter uncertainty, or stability study in a certain class of PPMs, or anything else with nonnegative matrices, should prompt a pertinent PGI to get 'more merits of indication' in solving the problem. Acknowledgments: Programming and calculations were implemented in MATLAB_R2020b.

Conflicts of Interest:
The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, or in the decision to publish the results.
Either of expressions (A4), (A5) is obviously much more complicated than (A1) at least because of square-root terms. It is also more complicated than any component of the Razzhevaikin-Tyrtyshnikov PGI expression obtained by means of the dimension-reducing algorithm (13) We have eliminated a pair of original rows as some current ones evidently majorized them.
The expression for RRT(L2rec) as a polynomial function of 13 parameters is similar to (A6) but cumbersome enough (Supplementary Materials) to be not reproduced in a publication. However, its value at the lower bounds of the parameter ranges is hence 1(L2rec) > 1 due to the PGI property, and it is so everywhere within the ranges due the the monotone property of 1 [5].