Life history models and projection matrices

Matrix projection models (MPMs) are representations of the dynamics of apopulation across all life history stages deemed relevant, across themost relevant time interval (typically one year, and assumed to beconsistent within each analysis). They require a complete model of theorganism’s life history prior to construction, and this model mustexplicitly describe all life stages and all life history transitions tobe modeled. Each stage is mutually exclusive, meaning that an individualcan only be in a single stage at a given time. Each stage is representedin each projection matrix by a single column and a single row, and theintersection of each row and each column gives the probability or rateof a specific transition. Each transition takes exactly one full timestep, which needs to be defined consistently and should be equivalent tothe approximate amount of time between consecutive monitoring occasions.Matrix elements (\(a_{kj}\)) showeither the probability of transition for an individual in stage\(j\) at occasiont (along thecolumns), to stage\(k\) at occasiont+1 (along the rows), or the mean rate of production of newrecruits into stage\(k\) at occasiont+1 (along the rows) by individuals in reproductive stage\(j\) in occasiont (along thecolumns). Conceptually, each individual is in a particular stage in theinstance of monitoring or observation, and then either dies or completesa transition in the interval between that occasion’s observation and thenext occasion’s observation. Death is not an explicit life stage and sois not modeled as such, instead becoming a potential endpoint of eachtransition.

Stages are generally defined under the assumption that they are uniquemoments in an organism’s life. This uniqueness is assumed to extend tothe demography of the organism during its time within that stage. Insome cases, stages are defined as developmental stages, as happens withinsect instars. In other cases, stages are defined almost purely on thebasis of size, as occurs with perennial herbs that may exhibit adifferent number of stems, leaves, and flowers in each growing season.In still other cases, stages may be defined via other characteristics.For example, the vegetative dormancy stage of some perennial herbs isdefined by a lack of aboveground size, and so is characteristically notobservable(Shefferson et al. 2001). Instill other cases, age may also fall into a stage’s definition, as maythe maturity status, observability, etc.

The timing of monitoring relative to the reproductive season impacts thestructure of life history models. Life history models are typicallycategorized as eitherpre-breeding orpost-breeding. Here, “breeding” refers to theproduction of offspring, and so a pre-breeding model assumes thatmonitoring is conducted immediately before new recruits are born, whilea post-breeding model assumes that monitoring is conducted just afternew recruits are born (in both cases, breeding is assumed to beseasonal). In a pre-breeding model, fecundity equals the production ofnewborns multiplied by the survival probability of newborns to age 1,since fecundity must take place over a full time step and the timing ofthe census misses both the birth event itself and the time in which theorganism is a newborn(Kendall et al.2019). In a post-breeding model, fecundity equals the survival ofthe parent from the stage/age preceding reproduction to reproductionitself, multiplied by the number of newborns produced(Kendall et al. 2019).

In monitoring studies of plant populations, the typical life historymodel a pre-breeding model, in which fecundity is estimated as theproduction of seeds in a given year multiplied by their over-wintersurvival probability as seeds and their germination probability in thefollowing year. If the life history includes seed dormancy, then it ismodeled in the first instance by multiplying seed production by seedsurvival (which is given as the probability of maintaining seedviability from one year to the next), and then as seed survival forseeds already produced and dormant afterwards (seed survival may bemodeled to decrease with increasing age of dormant seed). Addedcomplexity can arise if there are multiple fecundity pathways, forexample when clonal reproduction is also possible, or if multiplepropagule stages exist, or if recruitment can be to seeds thatimmediately germinate as well as to seeds that are dormant for one ormore years. We urge users to be careful with this step, as properlydefining the life history model has important implications for allanalyses of population dynamics(Kendall et al.2019).

Figure 1.1 is a simple example of a stage-based life history model(technically alife cycle graph) and its associatedLefkovitch matrix for a terrestrial orchid species,Cypripediumcandidum(Shefferson et al. 2001).Here, we show each stage as a node, and each transition as auni-directional arrow (a). The rates and probabilities are shown asmathematical symbols (b), with\(S_{kj}\) denoting survival-transitionprobability from stage\(j\) inoccasiont to stage\(k\) inoccasiont+1, and\(F_{kj}\)denoting the fecundity of reproductive stage\(j\) into recruit stage\(k\) in this life history.

Figure 1.1. Simple life history model (a) and ahistorical MPM (b) forCypripedium candidum, a North American herbaceous plantspecies. Here,S is the dormant seed stage,J is the seedling stage,D is adultvegetative dormancy,V is adult vegetative sprouting,andF is flowering.

Ahistorical vs. historical matrix models

Figure 1.1b is an example of anahistorical MPM(ahMPM), which is a matrix projection model in which the future stage ofan individual is dependent only on its current stage, and not onprevious stages. In other words, individual history is not incorporatedinto ahMPMs. It may seem odd that individual history is not incorporatedinto the matrices given that demographic datasets are composed ofrecords of individual histories spanning several or even manyobservation events. However, construction of a typical ahMPM breaks upthese individual histories into pairs of consecutive stages across time,with each pair treated as independent of every other pair of consecutivestages. For example, if an individual is in stageA in occasion1, stageB in occasion 2, stageE in occasion 3, anddead in occasion 4, then its individual history is broken up intoA-B,B-E, andE-Dead, with each pair assumed to be independent. Theresulting projection matrix would be the same if these transitionsoriginated from different individuals or the same one - their order andrelationships do not have any impacts in ahMPM construction.

The matrix shown in figure 1.1 can be used to project populationdynamics by predicting the future numbers of individuals in each stage.For this example, we can project forward by multiplying the projectionmatrix by a column vector of the numbers of individuals in each stage,given in the same order as the order of stages corresponding to the rowsand columns of the matrix, in a particular time, as in

\[\begin{equation}\tag{1}\left[\begin{array}{rrrrr}n _{2,S} \\n _{2,J} \\n _{2,V} \\n _{2,D} \\n _{2,F}\end{array}\right] = \left[\begin{array}{rrrrr}S _{SS} & 0 & 0 & 0 & F _{SF} \\S _{JS} & S _{JJ} & 0 & 0 & F _{JF} \\0 & S _{VJ} & S _{VV} & S _{VD} & S _{VF} \\0 & S _{DJ} & S _{DV} & S _{DD} & S _{DF} \\0 & 0 & S _{FV} & S _{FD} & S _{FF}\end{array}\right] \left[\begin{array}{rrrrr}n _{1,S} \\n _{1,J} \\n _{1,V} \\n _{1,D} \\n _{1,F}\end{array}\right] \end{equation}\]

where\(n _{1,S}\) is the number ofindividuals in stage\(S\) in occasion1. The total population size in any given occasion is the sum of thenumbers of individuals in all stages in that occasion. The discretepopulation growth rate is calculated as the ratio of total populationsize in occasion 2 to the total population size in occasion 1, and theasymptotic growth rate\(\lambda\) canbe estimated as the dominant eigenvalue of this matrix.

The independence of consecutive stage-pairs reflects a centralassumption in ahMPM analysis: the stage of an individual in the nextoccasion is influenced only by its current stage. Conceptually, if anorganism’s stage in the next occasion is entirely determined by itscurrent stage, then its previous states do not influence thesetransitions. Thus, standard ahMPMs are two dimensional and reflect onlythe current and next immediate stage of individuals, given by thecolumns and rows, respectively. This is ultimately an extension of theiid assumption in statistics - that the states of individualsareindependent and originate fromidentically-distributed random variables.

MPMs have been estimated since the 1940s(Caswell2001). We have never done a meta-analysis of all of thesestudies, but nonetheless it is safe to say that studies consideringindividual history are rare. The typical MPM study uses ahMPMs, and soassumes independence of stage transitions across time even from the sameindividual. In fact, at the time of writing, we are aware of only fiveexamples of studies breaking these assumptions and using ahistorical approach, meaning that they incorporatedsome degree of individual history into matrix estimation and analysis(Ehrlén 2000; Shefferson, Warren II, and Pulliam2014; Shefferson, Mizuta, and Hutchings 2017; Shefferson et al. 2018;deVries and Caswell 2018).

Thehistorical MPM (hMPM) is an extension of the matrixprojection model that incorporates information on one previous occasioninto the determination of vital rates. Thus, the expectedsurvival-transition probability of an individual in stage\(j\) at occasiont to stage\(k\) at occasiont+1 depends notonly on its stage in occasiont but also on its stage inoccasiont-1. Population ecologists considering this problemanalytically might be inclined to add an extra dimension to the matrixto deal with this, thus creating a 3d array or cube. However, this ismathematically intractable, with many analyses becoming impossible(deVries and Caswell 2018). Instead, we utilizethe approach developed byEhrlén (2000),in which rows and columns represent life history stages paired inconsecutive occasions. Thus, columns now represent theFrom pair of stages (stage in occasionst andt-1), and rows now represent theTo pair ofstages (stage in occasionst andt+1), as in Figure1.2. This model is equivalent to the second-order model proposed bydeVries and Caswell (2018), although thelatter incorporates more history into transitions involving individualsjust born (we explain the difference further in this vignette).

Figure 1.2. Historical MPM forCypripedium candidum, a NorthAmerican herbaceous plant.

Here, we can refer to matrix elements as\(a_{kjl}\), where the element represents therate at which individuals transition to stage\(k\) in occasiont+1 after havingbeen stage\(j\) in occasiontand stage\(l\) in occasiont-1. If\(n_{kjl}\) is thenumber of individuals making this transition,\(n_{.jl}\) is the total number ofindividuals in stage\(j\) in occasiont and stage\(l\) in occasiont-1 regardless of stage or even status as alive or dead inoccasiont+1,\(m\) is thenumber of stages in the life history model, and\(d\) is the number of stages plus death inthe life history model (so\(d=m+1\)),then we have

\[\begin{equation}\tag{2}a _{kjl} = \frac{n _{kjl}}{n _{.jl}} = \frac{n _{kjl}}{\sum_{i=1}^{d} n_{ijl}}\end{equation}\]

These values can be used to determine the associated elements ofmatrices in ahistorical MPMs, as follows:

\[\begin{equation}\tag{3}a _{kj} = \frac{\sum_{l=1}^{m} n _{kjl}}{\sum_{l=1}^{m} \sum_{i=1}^{d} n_{ijl}}\end{equation}\]

Note that although one can use the actual numbers of individualsmaking historical transitions to estimate ahistorical MPMs, one cannotuse the matrix elements in a historical MPM to calculate the associatedahistorical MPM elements (and vice versa). Historical transitions mustbe weighted by the numbers of individuals corresponding to eachtransition from occasiont tot+1 that were in eachprevious stage in occasiont-1 to produce the correctahistorical matrix element. Historical matrix elements are equal only tothe ratios of numbers of individuals in these categories, rather than tothe numbers of individuals themselves, and contain no informationallowing the inference of the latter.

Historical MPMs can be projected forward in the same way that ahMPMscan. Here we see an example of a projection going forward one time step,using the hMPM shown in figure 1.2.

\[\begin{equation}\tag{4}\tiny\left[\begin{array}{rrrrrrrrrrrrrrrrr}n _{2,SS} \\n _{2,JS} \\n _{2,JJ} \\n _{2,VJ} \\n _{2,DJ} \\n _{2,FJ} \\n _{2,VV} \\n _{2,DV} \\n _{2,FV} \\n _{2,VD} \\n _{2,DD} \\n _{2,FD} \\n _{2,VF} \\n _{2,DF} \\n _{2,FF} \\n _{2,SF} \\n _{2,JF}\end{array}\right] = \left[\begin{array}{rrrrrrrrrrrrrrrrr}S _{SSS} & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0& 0 & 0 & 0 & 0 & 0 & 0 & F _{SSF} & 0\\S _{JSS} & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0& 0 & 0 & 0 & 0 & 0 & 0 & F _{JSF} & 0\\0 & S _{JJS} & S _{JJJ} & 0 & 0 & 0 & 0 & 0& 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & F_{JJF} \\0 & S _{VJS} & S _{VJJ} & 0 & 0 & 0 & 0 & 0& 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & F_{VJF} \\0 & S _{DJS} & S _{DJJ} & 0 & 0 & 0 & 0 & 0& 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & F_{DJF} \\0 & S _{FJS} & S _{FJJ} & 0 & 0 & 0 & 0 & 0& 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & F_{FJF} \\0 & 0 & 0 & S _{VVJ} & 0 & 0 & S _{VVV} & 0& 0 & S _{VVD} & 0 & 0 & S _{VVF} & 0 & 0& 0 & 0 \\0 & 0 & 0 & S _{DVJ} & 0 & 0 & S _{DVV} & 0& 0 & S _{DVD} & 0 & 0 & S _{DVF} & 0 & 0& 0 & 0 \\0 & 0 & 0 & S _{FVJ} & 0 & 0 & S _{FVV} & 0& 0 & S _{FVD} & 0 & 0 & S _{FVF} & 0 & 0& 0 & 0 \\0 & 0 & 0 & 0 & S _{VDJ} & 0 & 0 & S _{VDV}& 0 & 0 & S _{VDD} & 0 & 0 & S _{VDF} & 0& 0 & 0 \\0 & 0 & 0 & 0 & S _{DDJ} & 0 & 0 & S _{DDV}& 0 & 0 & S _{DDD} & 0 & 0 & S _{DDF} & 0& 0 & 0 \\0 & 0 & 0 & 0 & S _{FDJ} & 0 & 0 & S _{FDV}& 0 & 0 & S _{FDD} & 0 & 0 & S _{FDF} & 0& 0 & 0 \\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & S _{VFV}& 0 & 0 & S _{VFD} & 0 & 0 & S _{VFF} & 0& 0 \\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & S _{DFV}& 0 & 0 & S _{DFD} & 0 & 0 & S _{DFF} & 0& 0 \\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & S _{FFV}& 0 & 0 & S _{FFD} & 0 & 0 & S _{FFF} & 0& 0 \\0 & 0 & 0 & 0 & 0 & F _{SFJ} & 0 & 0 & F_{SFV} & 0 & 0 & F _{SFD} & 0 & 0 & F _{SFF}& 0 & 0 \\0 & 0 & 0 & 0 & 0 & F _{JFJ} & 0 & 0 & F_{JFV} & 0 & 0 & F _{JFD} & 0 & 0 & F _{JFF}& 0 & 0 \\\end{array}\right] \left[\begin{array}{rrrrrrrrrrrrrrrrr}n _{1,SS} \\n _{1,JS} \\n _{1,JJ} \\n _{1,VJ} \\n _{1,DJ} \\n _{1,FJ} \\n _{1,VV} \\n _{1,DV} \\n _{1,FV} \\n _{1,VD} \\n _{1,DD} \\n _{1,FD} \\n _{1,VF} \\n _{1,DF} \\n _{1,FF} \\n _{1,SF} \\n _{1,JF}\end{array}\right] \end{equation}\]

Here,\(n _{1,SS}\) is the number ofindividuals that were in stage\(S\) inboth occasions 0 and 1, and\(n_{2,DF}\) is the number of individuals that were in stage\(D\) in occasion 2 and stage\(F\) in occasion 1.

Historical MPMs involve much larger matrices, and so parameterize manymore matrix elements, than ahMPMs do. Figure 1.2 illustrates this issue,but the function-based matrix analysis in ourCypripediumcandidum vignette provides a more realistic example. In that casestudy, we use a life history with 54 life stages, of which the first isa dormant seed stage, the next four are immature stages, the sixth isadult vegetative dormancy, the next 24 stages are size-classifiednon-flowering adults, and the final 24 stages are size-classifiedflowering adults. A life history with 54 stages yields an ahistoricalmatrix with 54 columns and 54 rows, and so 2,916 elements. A fullhistorical matrix should have\(54 \times 54 =2,916\) columns and 2,916 rows, consisting of 8,503,056 elements,because the column and row stages represent pairs of stages.

In an ahistorical matrix, the only true zeros will be those elementscorresponding to biologically impossible transitions, or, in the rawmatrix case, those elements without any individuals taking therespective transition in a given time step. However, most elements in ahistorical MPM are structural zeros. This happens because transitionelements in hMPMs are only estimable if stage at occasiont isequal in the column and row stage pairs. For example, the transitionprobability between stage A in occasiont-1 and stage B inoccasiont (column stage), to stage C in occasiontand stage D in occasiont+1 (row stage) equals 0, because anorganism cannot be in both stages B and C in occasiont.Increasing the number of stages in a life history model causes thesestructural zeros to increase at a faster rate than the rate at which thenumber of truly estimable elements increases. In fact, if there are\(m\) stages in a life history,yielding\(m^2\) elements in anahistorical matrix, then although there will be\(m^4\) elements in the historical matrix,only\(m^3\) will bepotentially estimable while\((m-1)m^3\) are structural zeros (we saypotentially because some of the logically possible transitionsmay still be biologically impossible, or may have no data to yieldvalues other than 0). Thus, in theCypripedium candidumexample, all 2,916 elements of the ahMPM are potentially estimable(although some are biologically impossible), while only 157,464 of the8,503,056 elements of the associated hMPM are potentially estimable andthe rest are structural zeros.

Raw vs. function-based matrix models, including integral projectionmodels

Matrix projection models, whether historical or ahistorical, can also becategorized as eitherraw orfunction-based. The oldest and most common in theliterature is theraw MPM, in which transitionprobabilities in the matrix are estimated by counting all of theindividuals alive in a particular stage at occasiont, anddividing that number into the number of individuals from that set thatare still alive in each possible stage in occasiont+1(equations 2 and 3). For example, if 100 individuals are alive in stageD in year 2010, and 20 of these are alive in stage D in year 2011, asare a further 40 in stage V and a further 25 in stage F, then theassociated transition probabilities are 0.20, 0.40, and 0.25,respectively. The overall survival probability for stage D is the sum ofthese transitions, or 0.85.

Methods for estimating fecundity in raw matrices vary. In theCypripedium candidum vignette, we used counts of the number offruits produced by each individual in a given year, and multiplied bythe mean number of seeds per fruit and the mean germination probabilityin the next year. This study was conducted as a pre-breeding census -fecundity was estimated using a count of the seeds or other propagulesthat the plant produces (or a proxy thereof, such as the number offruits multiplied by the mean number of seeds per fruit) multiplied bythe probability of that propagule surviving to the next year in order tomake sure that fecundity reflects demographic processes across a fulltime step(Kendall et al. 2019). In othersystems in which offspring may be observed and tracked, counts of actualrecruits may be possible, such as in studies of nesting birds, and as inLathyrus vernus, an herbaceous perennial plant that is thesubject of some vignettes included in this package.

Infunction-based MPMs, matrix elements are populatedby kernels that link together functions representing key demographicprocesses governing each transition. Typically, demographic datasets areanalyzed via linear models to estimate demographic parameters such assurvival probability, the probability of becoming a certain sizeassuming survival, and fecundity. These linear models are developedusing nested subsets of the demographic data used in a study in order toestimate conditional probabilities reflecting all of the demographicprocesses that must occur over the span of a single time step. Matrixelements are then estimated as products of responses from these linearmodels set to particular inputs corresponding to stage at occasiont and occasiont+1, as well as any other parametersgoverning the construction of the matrix. Function-based MPMs are morerecent inventions than raw MPMs, and their strengths are making themincreasingly common in the literature.

Easterling, Ellner, and Dixon (2000)proposed a method of analyzing population dynamics called theintegral projection model (IPM), which modeled the lifehistory of an organism using a continuous size metric modeled on aGaussian distribution. At least theoretically, the aim was to produce amodel that was not contingent on a stage classification together with afunction-based kernel to model population dynamics. In practice,however, size was and is discretized into many fine-scale size classes.The result is that, in practice, IPMs are actually high resolutionfunction-based MPMs, breaking up a life history into many fine-scalesize classes using a continuous measure of size and then estimatingsurvival-transitions and fecundity according to these fine-scale sizeclasses(Ellner and Rees 2006). Packagelefko3 easily handles the creation and analysis of IPMs. Weinclude a vignette inlefko3 showing an example, and alsodevote a chapter of the free e-booklefko3: a gentleintroduction to the topic.

The example from Figures 1.1 and 1.2 may help to illustrate theprocedure of creating a function-based MPM. InCypripediumcandidum, the simple life history model shown in Figure 1 includesone adult stage that is not observable (D - vegetativedormancy), one adult stage that is observable but not reproductive(V - vegetative sprouting), and one adult stage that is bothobservable and reproductive (F - flowering). If\(\sigma _{.jl}\) is the probability ofsurvival from stage\(j\) in occasiont and stage\(l\) in occasiont-1 to any stage in occasiont+1,\(\rho _{.jl}\) is the probability ofsprouting in occasiont+1 after being stage\(j\) in occasiont and stage\(l\) in occasiont-1 and survivingto occasiont+1, and\(\gamma_{.jl}\) is the probability of flowering in occasiont+1after being stage\(j\) in occasiont and stage\(l\) in occasiont-1 and surviving to and sprouting in occasiont+1,then

\[\begin{equation}\tag{5}a _{DDD} = \sigma _{.DD} (1 - \rho _{.DD})\end{equation}\]

\[\begin{equation}\tag{6}a _{VDD} = \sigma _{.DD} \rho _{.DD} (1 - \gamma _{.DD})\end{equation}\]

\[\begin{equation}\tag{7}a _{FDD} = \sigma _{.DD} \rho _{.DD} \gamma _{.DD}\end{equation}\]

In both raw and function-based MPMs, most estimated elements aresurvival-transition probabilities. In ahMPMs, these give theprobability of an individual in stage\(j\) at occasiont surviving andtransitioning to stage\(k\) atoccasiont+1. Fewer elements are devoted to fecundity, althoughlife history models with more reproductive stages and more recruitstages will have more fecundity elements. In Figure 1.1, fecundity isshown in the top-right of the matrix, as the mean production of seedseither dormant or germinating in the next occasion. Ignoring thefecundity elements, we have asurvival-transitionmatrix (symbolized as eitherU orT), in which the column sums correspond to the expectedsurvival probabilities of individuals in each stage from occasiont to occasiont+1. Ignoring the survival-transitionterms, we have thefecundity matrix (symbolized asF), in which the column sums correspond to the expectedoverall time-specific fecundity of individuals in these respectivestages. The full MPM (symbolized asA) is the matrixsum of the survival-transition matrix and the fecundity matrix.

In principle, everything just noted about function-based ahMPMs alsoapplies to function-based hMPMs. However, there is a key difference interms of the distribution of elements within the matrix. The inclusionof stage or condition in occasiont-1 increases both the numberof survival-transition elements and the number of fecundity elements,and also has the effect of dispersing the locations of these elementsacross the matrix. Whereas estimated transitions in ahMPMs willgenerally occur in dense patches of estimated elements, in hMPMs, mostelements are structural zeros, with elements occurring eitherindividually surrounded by zeros, or in small patches of estimatedelements surrounded by zeros. This makes error-checking morechallenging, because the user needs to be aware of where elements aresupposed to occur in order to check their values properly.

Different parameterizations of the historical MPM

deVries and Caswell (2018) detaildifferent approaches to the development of hMPMs, differing in how stateor stage in occasiont-1 is incorporated into the matrix.Full prior stage dependence models deal with history byincorporating prior condition as the exact stage of an organism inoccasiont-1, yielding a matrix showing stage pairs in occasiont-1 andt along the columns of the matrix, and stagepairs in occasiont andt+1 along the rows of thematrix.Prior condition models deal with history bymaking the transition from staget to staget+1 afunction of stage at occasiont and condition in occasiont-1. Prior condition can be determined in the same way thatcurrent stage is determined (e.g. size classification), or a differentmeasure of condition can be used, such as growth (i.e. the change insize between occasionst-1 andt).

One key feature proposed bydeVries and Caswell(2018) is the addition of a new stage to account for the priorstatus of newborn individuals. This reflects a slightly differentinterpretation fromEhrlén (2000) of thehistorical transition. InEhrlén (2000),all matrix elements simply reflect the order of the events in a singletransition, including both fecundity and survival-transition events.However, indeVries and Caswell (2018),these matrix elements must also reflect the history of specificindividuals. In the latter case, the fact that a newborn in occasiont did not exist in occasiont-1 leads to the creationof a new prior stage to account for the time immediately prior to birth.This new stage only exists in the immediate prior occasion, and is onlyused for newborns in the first survival transition from birth. So, forexample, in a three stage MPM where the 1st stage is the newborn stageand only the 3rd stage is reproductive, we add a 4th stage to the priorportion of the row and column. Thus, we may start with the followinghMPM inEhrlén (2000) format:

\[\begin{equation}\tag{8}\tiny\left[\begin{array}{rrrrrrrrr}n _{2,AA} \\n _{2,BA} \\n _{2,CA} \\n _{2,AB} \\n _{2,BB} \\n _{2,CB} \\n _{2,AC} \\n _{2,BC} \\n _{2,CC}\end{array}\right] = \left[\begin{array}{rrrrrrrrr}S _{AAA} & 0 & 0 & S _{AAB} & 0 & 0 & S _{AAC}& 0 & 0 \\S _{BAA} & 0 & 0 & S _{BAB} & 0 & 0 & S _{BAC}& 0 & 0 \\S _{CAA} & 0 & 0 & S _{CAB} & 0 & 0 & S _{CAC}& 0 & 0 \\0 & S _{ABA} & 0 & 0 & S _{ABB} & 0 & 0 & S_{ABC} & 0 \\0 & S _{BBA} & 0 & 0 & S _{BBB} & 0 & 0 & S_{BBC} & 0 \\0 & S _{CBA} & 0 & 0 & S _{CBB} & 0 & 0 & S_{CBC} & 0 \\0 & 0 & S _{ACA} + F _{ACA} & 0 & 0 & S _{ACB} + F_{ACB} & 0 & 0 & S _{ACC} + F _{ACC} \\0 & 0 & S _{BCA} & 0 & 0 & S _{BCB} & 0 & 0& S _{BCC} \\0 & 0 & S _{CCA} & 0 & 0 & S _{CCB} & 0 & 0& S _{CCC} \\\end{array}\right] \left[\begin{array}{rrrrrrrrr}n _{1,AA} \\n _{1,BA} \\n _{1,CA} \\n _{1,AB} \\n _{1,BB} \\n _{1,CB} \\n _{1,AC} \\n _{1,BC} \\n _{1,CC}\end{array}\right] \end{equation}\]

This hMPM becomes as follows indeVries andCaswell (2018) format:

\[\begin{equation}\tag{9}\tiny\left[\begin{array}{rrrrrrrrrrrr}n _{2,AA} \\n _{2,BA} \\n _{2,CA} \\n _{2,AB} \\n _{2,BB} \\n _{2,CB} \\n _{2,AC} \\n _{2,BC} \\n _{2,CC} \\n _{2,AP} \\n _{2,BP} \\n _{2,CP}\end{array}\right] = \left[\begin{array}{rrrrrrrrrrrr}S _{AAA} & 0 & 0 & S _{AAB} & 0 & 0 & S _{AAC}& 0 & 0 & S _{AA,AP} & 0 & 0 \\S _{BAA} & 0 & 0 & S _{BAB} & 0 & 0 & S _{BAC}& 0 & 0 & S _{BA,AP} & 0 & 0 \\S _{CAA} & 0 & 0 & S _{CAB} & 0 & 0 & S _{CAC}& 0 & 0 & S _{CA,AP} & 0 & 0 \\0 & S _{ABA} & 0 & 0 & S _{ABB} & 0 & 0 & S_{ABC} & 0 & 0 & 0 & 0 \\0 & S _{BBA} & 0 & 0 & S _{BBB} & 0 & 0 & S_{BBC} & 0 & 0 & 0 & 0 \\0 & S _{CBA} & 0 & 0 & S _{CBB} & 0 & 0 & S_{CBC} & 0 & 0 & 0 & 0 \\0 & 0 & S _{ACA} & 0 & 0 & S _{ACB} & 0 & 0& S _{ACC} & 0 & 0 & 0 \\0 & 0 & S _{BCA} & 0 & 0 & S _{BCB} & 0 & 0& S _{BCC} & 0 & 0 & 0 \\0 & 0 & S _{CCA} & 0 & 0 & S _{CCB} & 0 & 0& S _{CCC} & 0 & 0 & 0 \\0 & 0 & F _{AP,CA} & 0 & 0 & F _{AP,CB} & 0& 0 & F _{AP,CC} & 0 & 0 & 0\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 &0 & 0 & 0 \\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 &0 & 0 & 0 \\\end{array}\right] \left[\begin{array}{rrrrrrrrrrrr}n _{1,AA} \\n _{1,BA} \\n _{1,CA} \\n _{1,AB} \\n _{1,BB} \\n _{1,CB} \\n _{1,AC} \\n _{1,BC} \\n _{1,CC} \\n _{1,AP} \\n _{1,BP} \\n _{1,CP}\end{array}\right] \end{equation}\]

Note that matrix 9 above can be reduced by two rows and two columns,since transitions from the newborn prior stage (P) can only be to thenewborn stage (A).

Packagelefko3 generally implements thefull priorstage dependence model approach, particularly in thedevelopment of raw matrices. However, in principle, function-basedmatrices developed withlefko3 can be seen as fallingwithin theprior condition model approach where priorcondition is determined in the same way as current stage. We alsoimplement both theEhrlén (2000) formatand thedeVries and Caswell (2018) format,and use the former as the default. We leave it to the user to decidewhich approach to use.

deVries and Caswell (2018) also suggestthat hMPMs can be organized differently to yield more intuitive models.In particular, matrices may be as organized as block matrices where thecomponent blocks are essentially matrices showing the transitions fromall stages in occasiont to stages in occasiont+1conditioned on individuals having been in the same stage in occasiont-1. We refer to these block component matrices asconditional matrices, since they show transitions fromoccasiont to occasiont+1 conditional on the sameprevious stage in occasiont-1. We have developed functioncond_hmpm() to take full hMPMs and decompose them intotheir associated conditional matrices.

Leslie (age-classified) MPMs

Some population biologists working withlefko3 may wish toproduce Leslie matrices, which classify an organism only by age. ALeslie matrix has the following general form

\[\begin{equation}\tag{10}\tiny\left[\begin{array}{rrrrrrr}n _{2, 0} \\n _{2, 1} \\n _{2, 2} \\n _{2, 3} \\\vdots \\n _{2, a-1} \\n _{2, a}\end{array}\right] = \left[\begin{array}{rrrrrrr}0 & F _1 & F _2 & F _3 & \cdots & F _{a-1} & F_a \\S _{1,0} & 0 & 0 & 0 & \cdots & 0 & 0 \\0 & S _{2,1} & 0 & 0 & \cdots & 0 & 0 \\0 & 0 & S _{3,2} & 0 & \cdots & 0 & 0 \\\vdots & & & & & & \vdots \\0 & 0 & 0 & 0 & \cdots & 0 & S _{5+B,4C} \\0 & 0 & 0 & 0 & \cdots & S _{a,a-1} & S_{>a,a} \\\end{array}\right] \left[\begin{array}{rrrrrrr}n _{1, 0} \\n _{1, 1} \\n _{1, 2} \\n _{1, 3} \\\vdots \\n _{1, a-1} \\n _{1, a}\end{array}\right] \end{equation}\]

Here,\(n _{2, 0}\) refers to thenumber of newborn (age = 0) individuals in time 2,\(a\) is the final age in the model,\(S _{1, 0}\) is the survival probability ofa newborn individual (age = 0) to the next age,\(F _1\) is the fecundity of an individual inage 1, and\(S _{>a, a}\) is thesurvival of individuals in the final modeled age to the next age (andlater). Note that in Leslie models, survival transitions are always justbelow the diagonal or in the final lower-right element of the diagonal,and fecundity transitions are always in the top row. Note also thatage-classified models are inherently ahistorical, and so no historicalmodels are possible. Packagelefko3 estimates both raw andfunction-based forms of these.

Age-by-stage MPMs

Packagelefko3 can also handle the estimation and analysisof age-by-stage MPMs. These matrices incorporate both stages (whethersize-classified or not), as well as the ages in which they occur. Thestandard approach uses a block matrix to repeat relevant stages by eachage. For example, in a three stage life history with a single newbornstage and two adult stages, only one of which is reproductive (as inequations 1, 4, 8, and 9), and in which age impacts survival andfecundity only for the first 4 years, beyond which survival andfecundity transitions remain the same, we have the following ahistoricalage-by-stage matrix.

\[\begin{equation}\tag{11}\tiny\left[\begin{array}{rrrrrrr}n _{2,1A} \\n _{2,2B} \\n _{2,2C} \\n _{2,3B} \\n _{2,3C} \\n _{2,4B} \\n _{2,4C}\end{array}\right] = \left[\begin{array}{rrrrrrr}0 & 0 & F _{1A,2C} & 0 & F _{1A,3C} & 0 & F_{1A,4+C} \\S _{2B,1A} & 0 & 0 & 0 & 0 & 0 & 0 \\S _{2C,1A} & 0 & 0 & 0 & 0 & 0 & 0 \\0 & S _{3B,2B} & S _{3B,2C} & 0 & 0 & 0 & 0 \\0 & S _{3C,2B} & S _{3C,2C} & 0 & 0 & 0 & 0 \\0 & 0 & 0 & S _{4B,3B} & S _{4B,3C} & S _{5+B,4B}& S _{5+B,4C} \\0 & 0 & 0 & S _{4C,3B} & S _{4C,3C} & S _{5+C,4B}& S _{5+C,4C} \\\end{array}\right] \left[\begin{array}{rrrrrrr}n _{1,1A} \\n _{1,2B} \\n _{1,2C} \\n _{1,3B} \\n _{1,3C} \\n _{1,4B} \\n _{1,4C}\end{array}\right] \end{equation}\]

The above matrix is ahistorical. Packagelefko3 does notcurrently support historical age-by-stage matrices, although they arecertainly possible.

Packagelefko3 also supports the construction of age-hybridMPMs, which are simple age-by-stage MPMs in which a singlenon-age-classified stage is added to what is otherwise a Leslie MPM.Situations in which such an approach is useful include an age-classifiedplant model with a dormant seed stage, among others.

Step 1: Life history model development

The basic workflow to analysis with packagelefko3 startswith the development of a life history model that encapsulates all ofthe appropriate life stages relevant to population dynamics. We do notpresent a complete discussion on this step, instead focusing on a sketchof the process focused on the life cycle graph and noting key issuesrelated to this package. We encourage users to explore the literature onthis subject, but particularly note chapters 3 and 4 inCaswell (2001) for a good treatment useful forbeginners. Stages need to be defined carefully, and should stronglyaccount for variability in vital rates. Statistical approaches such asJenks natural breaks optimization can help determine ideal“natural” breaks in the distribution of size data to act as bordersbetween size-classified stages in raw MPMs(Jenks1967). In contrast, stages within function-based MPMs are morestraightforward to define, with regular breaks across the size spectrum.We also noteKendall et al. (2019) as agood reference detailing common problems in MPM construction and how toavoid them through proper development and operationalization of lifehistory models.Beissinger and Westphal(1998),Wardle (1998), andSalguero-Gómez and Casper (2010) provide gooddiscussions of the proper application of life history models tounderstand population dynamics through the MPM approach.

Here we utilize alife cycle graph approach to build alife history model. A life cycle graph shows the life history of anorganism as a series of nodes and arrows. The nodes represent uniquestages that the individual may go through during development, andindividuals are noted as occurring in these stages during monitoringsessions conducted over a roughly or strictly regular frequency, such asat the start of the growing season every year, or over specific dateseach year. Arrows represent transitions between stages acrossconsecutive occasions. These may be interpreted either as probabilitiesof survival when reflecting the individual’s stage across time, or asrates when reflecting the production of offspring. A good life cyclegraph, and in general a good life history model, must explicitly includerepresentations of all stages possible for the organism to enter, andall biologically plausible transitions between these transitions. Figure1a shows an example of a life cycle graph.

The construction of a good life history model is not as easy as it mightappear, particularly because it can be influenced by factors seeminglyexternal to the life history of the organism. Some of the most importantconsiderations outside of the species’ biology include the size of thedataset, and the decision of whether to create a raw MPM or afunction-based MPM. All else being equal, a larger dataset and afunction-based MPM will allow more stages to be used in the lifehistory. Assuming a dataset of fixed size, raw MPMs will include moreand more zeros as the number of stages in the life history modelincreases. This happens because the chance of having any individualsactually move through a specific transition decreases as the number ofpossible transitions increases. These functional zeros that becomeincreasingly common with increasing stages can cause problems inanalysis, because they suggest that some transitions are impossiblewhen, in fact, individuals moving through them were missing just bychance. Even a transition with extremely high probability may be missingin some years just by chance. One impact of this would be the reductionof the mean transition within an element-wise mean matrix to anunrealistic level, causing the predicted population growth rate to beartificially low.

Function-based MPMs are better able to handle large numbers of matrixelements because the risk of overparameterization is reduced byparsimonious model selection when vital rate models are determined.However, the size of the dataset influences the statistical power ofthese vital rate models, with smaller datasets yielding larger errorassociated with matrix elements and hence a loss of variability withpotentially important factors. Loss of statistical power in vital ratemodels can lead to the loss of important process variance, and can makethe population appear more static than it really is. Further, developingfunction-based MPMs for life history models that include stages rarelyor never actually observed can yield problematic inferences, such asmatrices suggesting that some stages have survival probabilities greaterthan 1.0. Chapter 6 inMorris and Doak(2002) provides a good description of techniques that may be usedto define stages and to determine the exact number to use given thecontext imposed by a given dataset of a given size.Ellner and Rees (2006) andDoak et al. (2021) provide some discussion ofthis issue within the context of integral projection model development.

Once a life history model is developed, we need to characterize the mainlife stages in a way that is relevant to the dataset. We do this withthesf_create() function, which creates astageframe object. This term is a combination of thetermslife history stage anddata frame, and theobject itself is a data frame that describes all of the life historystages. This object is created in a way that shows how stage relates tosize, reproductive status, observation status, propagule status,maturity status, entry status (i.e. a stage’s ability to serve as thefirst stage of an individual’s life), and presence in the dataset, andin which the description of each life stage is unique. A good stageframedescribes stages in ways that matrix-estimating functions utilize tocompute elements accurately, and also allows the user to stipulate extradescriptive information on each stage as text. If the user wishes tocreate an IPM or a function-based MPM with many stages, then thisfunction also allows stage definition to be automated.

Step 2a: Dataset organization

Once the stageframe has been created, the demographic data needs to beorganized into a format thatlefko3 can use. The properformat is one that we termhistorically-formatted verticalformat, orhfv format for short. We havedeveloped two functions for this purpose, depending on the originalformat of the data. If the data is arrangedhorizontally,meaning that individual life histories are recorded in a spreadsheet inwhich each unique individual’s data is organized within a single row,where columns correspond to descriptor variables and to condition atdifferent times, then theverticalize3() function canstandardize this dataset into the proper format. If the dataset is invertical, ahistorical format, in which an individual’scondition across time is recorded across rows in a spreadsheet, with asingle row corresponding to either a single observation or to a pair ofconsecutive observations, then thehistoricalize3()function can standardize the dataset properly. These functions both canalso add some new variables that might be of use. For example, bothfunctions can estimate radial density if Cartesian coordinates for eachdatum and a radius are provided. Both functions are powered by binariesbuilt in C++ using single control loops to propagate each data frameproperly and efficiently, and so run very quickly (typically under 1s ona 2019 MacBook Pro). The case studies included inlefko3illustrate the usage of these functions on real data. The functionsummary_hfv() can be used to provide quick details aboutthe resulting standardized datasets.

Step 2b: Provide supplemental information for matrix estimation

MPMs are often estimated only partially from available demographicdatasets. Some transitions are parameterized using information gatheredfrom other studies, whether through direct input in the matrix orthrough the development of kernels contingent on external information.Other transitions might also be estimated via proxy transitionselsewhere in the matrix. At other times, users may wish to manipulatespecific transitions either through addition or multiplication byspecific amounts. Inlefko3, this information can beprovided through thesupplement table, which can bedeveloped using thesupplemental() function. This functionallows users to create a data frame detailing:

1. specific transitions to set as the main fecundity transitions;
2. specific transitions to set as constants, with values entered by theuser;
3. specific transitions to alter additively, with values entered by theuser;
4. specific transitions to estimate via proxy transitions;
5. specific multipliers for fecundity in cases in which fecundityestimates resulting from linear modeling must be modified tocharacterize the full transition, or for transitions estimated viaproxy; and
6. specific multipliers for survival in cases where survival isincorporated by proxy and must be set at a reduced or elevated levelrelative to the original, directly estimated transition.

Examples might illustrate where this approach is useful. If we lack dataon subterranean juvenile stages in a plant species, but we haveestimates of survival for those stages from other studies, then we mightuse those estimates as constants in the MPM. We can also set sometransitions explicitly to 0 if we believe that there is a biologicalneed for such an override. Further, if we lack demographic data on thedevelopment of germinated seeds to the seedling or earliest adult stage,but have reason to believe that the survival probabilities should besimilar to survival within an observed stage of seedlings or smalladults, then we can use the latter survival-transitions as proxies andeven decrease them via an offset or multiplier. Finally, if fecundity isa function of seed production, survival to the next year, andgermination probability, then germination probability might be estimatedvia a separate field germination study. We might wish to incorporatethis germination probability both as a constant transition, and as amultiplier on estimated fecundity. Supplement tables provide the meansto include all of this information.

Step 3: Assess whether historical MPM is justified, and developmodels of demographic parameters

Historical MPMs are not necessarily justified in all cases. We urge thereader to use the principle of parsimony in deciding this. For example,users may analyze whether lagged effects of previous stage appear toinfluence vital rates. Packagelefko3 includes one methodto do this based on linear modeling. A number of other methods exist,mostly focused on what are known asmemory models(Brownie et al. 1993; Pradel, Wintrebert, and Gimenez2003; Cole et al. 2014). Here, we focus only on the mainlefko3 method, and leave it to the reader to explore otheroptions.

Inlefko3, functionmodelsearch() allows usersto explore whether history influences component vital rates, and todevelop linear models of these vital rates. The results can be used todecide not only whether a historical MPM is justified, but also todevelop function-based hMPMs and ahMPMs themselves. Packagelefko3 can estimate linear models to estimate up tofourteen different vital rates:

1. Survival probability - This is the probabilityof surviving from occasiont to occasiont+1, giventhat the individual is in stage\(j\)in occasiont (and, if historical, in stage\(l\) in occasiont-1). Inlefko3, this parameter may be modeled as a function of upto three size metrics, reproductive status, patch, year, age, and anumber of individual or environmental covariates in occasionstandt-1, and individual identity. This parameter is required inall function-based matrices.

2. Observation probability - This is an optionalparameter denoting the probability of observation in occasiont+1 of an individual in stage\(k\) given survival from occasiontto occasiont+1. This parameter is only used when at least onestage is not observable. For example, some plants are capable ofvegetative dormancy, in which case they are alive but do not necessarilysprout in all years in which they are in this stage. In these cases, theprobability of sprouting may be estimated as the observationprobability, and its complement is the probability of becoming dormantassuming survival. Note that this probability does not refer to observereffort, and so should only be used to differentiate completelyunobservable stages where the observation status refers to an importantbiological phenomenon, such as when individuals may be alive but have asize of 0. Inlefko3, this parameter may be modeled as afunction of up to three size metrics, reproductive status, patch, year,age, and a number of individual or environmental covariates in occasionst andt-1, and individual identity.

3. Primary size transition probability - This isthe probability of becoming size\(k\)in occasiont+1 assuming survival from occasiont tooccasiont+1 and observation in that time. Inlefko3, this parameter may be modeled as a function of upto three size metrics, reproductive status, patch, year, age, and anumber of individual or environmental covariates in occasionstandt-1, and individual identity.

4. Secondary size transition probability - This isthe probability of becoming size\(k\)in occasiont+1 assuming survival from occasiont tooccasiont+1 and observation in that time, within a second sizemetric used for classification in addition to the primary metric. Inlefko3, this parameter may be modeled as a function of upto three size metrics, reproductive status, patch, year, age, and anumber of individual or environmental covariates in occasionstandt-1, and individual identity.

5. Tertiary size transition probability - This isthe probability of becoming size\(k\)in occasiont+1 assuming survival from occasiont tooccasiont+1 and observation in that time, within a third sizemetric used for classification in addition to the primary and secondarymetrics. Inlefko3, this parameter may be modeled as afunction of up to three size metrics, reproductive status, patch, year,age, and a number of individual or environmental covariates in occasionst andt-1, and individual identity.

6. Reproduction probability - This is an optionalparameter denoting the probability of reproducing in occasiont+1 given survival from occasiont to occasiont+1, and observation in that time. Note that this should beused only if the researcher wishes to separate breeding fromnon-breeding mature stages. If all adult stages are potentiallyreproductive and no separation of reproducing from non-reproducingadults is desired, then this parameter should not be estimated. Inlefko3, this parameter may be modeled as a function of upto three size metrics, reproductive status, patch, year, age, and anumber of individual or environmental covariates in occasionstandt-1, and individual identity.

7. Fecundity rate - Under the default setting, thisis the rate of successful production of offspring in occasiontby individuals alive, observable, and reproductive in that time, and, ifdesired and sufficient information is provided in the dataset, thesurvival of those offspring into occasiont+1 in whateverjuvenile class is possible. Thus, the fecundity rate of seed-producingplants might be split into seedlings, which are plants that germinatedwithin a year of seed production, and dormant seeds. Alternatively, itmay be given only as produced fruits or seeds, with the survival andgermination of seeds provided elsewhere in the MPM development process,such as within a supplement table. An additional setting allowsfecundity rate to be estimated using data provided for occasiont+1 instead of occasiont. Inlefko3,this parameter may be modeled as a function of up to three size metrics,reproductive status, patch, year, age, and a number of individual orenvironmental covariates in occasionst andt-1, andindividual identity.

8. Juvenile survival probability - This is anoptional parameter that is used to model the probability of survivingfrom juvenile stage\(j\) in occasiont to occasiont+1. It is used only when the userwishes to model juvenile vital rates separately from adults, assumingdifferent relationships with explanatory variables than those used toclassify adults. Inlefko3, this parameter may be modeledas a function of up to three size metrics, patch, year, and a number ofindividual or environmental covariates in occasionst andt-1, and individual identity.

9. Juvenile observation probability - This is anoptional parameter denoting the probability of observation in occasiont+1 of an individual in mature stage\(k\) given survival from a juvenile stage inoccasiont to occasiont+1. It is used only when theuser wishes to model juvenile vital rates separately from adults,assuming different relationships with explanatory variables than thoseused to classify adults. Inlefko3, this parameter may bemodeled as a function of up to three size metrics, patch, year, and anumber of individual or environmental covariates in occasionstandt-1, and individual identity, and all other caveats notedin (2) above apply.

10. Juvenile primary size transition probability -This is an optional parameter denoting the probability of becomingmature size\(k\) in occasiont+1 assuming survival from juvenile stage\(j\) in occasiont to occasiont+1 and observation in that time. It is used only when the userwishes to model juvenile vital rates separately from adults, assumingdifferent relationships with explanatory variables than those used toclassify adults. Inlefko3, this parameter may be modeledas a function of up to three size metrics, patch, year, and a number ofindividual or environmental covariates in occasionst andt-1, and individual identity.

11. Juvenile secondary size transition probability -This is an optional parameter denoting the probability of becomingmature size\(k\) in occasiont+1 assuming survival from juvenile stage\(j\) in occasiont to occasiont+1 and observation in that time, in a secondary size metric inaddition to the primary size metric. It is used only when the userwishes to model juvenile vital rates separately from adults, assumingdifferent relationships with explanatory variables than those used toclassify adults. Inlefko3, this parameter may be modeledas a function of up to three size metrics, patch, year, and a number ofindividual or environmental covariates in occasionst andt-1, and individual identity.

12. Juvenile tertiary size transition probability -This is an optional parameter denoting the probability of becomingmature size\(k\) in occasiont+1 assuming survival from juvenile stage\(j\) in occasiont to occasiont+1 and observation in that time, in a tertiary size metric inaddition to the primary and secondary size metrics. It is used only whenthe user wishes to model juvenile vital rates separately from adults,assuming different relationships with explanatory variables than thoseused to classify adults. Inlefko3, this parameter may bemodeled as a function of up to three size metrics, patch, year, and anumber of individual or environmental covariates in occasionstandt-1, and individual identity.

13. Juvenile reproduction probability - This is anoptional parameter denoting the probability of reproducing in maturestage\(k\) in occasiont+1given survival from juvenile stagej in occasiont tooccasiont+1, and observation in that time. It is used onlywhen the user wishes to model juvenile vital rates separately fromadults, assuming different relationships with explanatory variables thanthose used to classify adults. Inlefko3, this parametermay be modeled as a function of up to three size metrics, patch, year,and a number of individual or environmental covariates in occasionst andt-1, and individual identity, and all othercaveats in (4) apply.

14. Juvenile maturity probability - This is anoptional parameter denoting the probability of becoming mature inoccasiont+1 given survival from juvenile stagej inoccasiont to occasiont+1. This is a differentparameter than (13) above and is particularly used to weight transitionprobabilities to juvenile vs. mature stages properly in matrices withdefined juvenile stages. It is used only when the user wishes to modeljuvenile vital rates separately from adults, assuming differentrelationships with explanatory variables than those used to classifyadults. Inlefko3, this parameter may be modeled as afunction of up to three size metrics, patch, year, and a number ofindividual or environmental covariates in occasionst andt-1, and individual identity, and all other caveats in (1)apply.

Of these fourteen vital rates, most users will estimate at leastparameters (1) survival probability, (3) size transition probability,and (7) fecundity. These three are the default set for functionmodelsearch(). Parameters (2) observation probability and(6) reproduction probability may be used when some stages are includedthat are completely unobservable (and so do not have any size), or thatare mature but non-reproductive, respectively. Parameters (8) through(14) should only be added if the dataset contains juvenile individualstransitioning to maturity, and these juveniles live essentially as asingle stage before transitioning to maturity, or before transitioningto a stage that is size-classified in the same manner as adult stagesare.

Functionmodelsearch() handles the entire modeling processusing an exhaustive model building process combined with informationtheoretical model selection using model AICc and the number of estimatedparameters in each model. It handles the development of global models,exhaustive model building, and the selection of best-fit models. Usersshould provide this function with information about the following:

1. Individual history - Are the matrices to bebuilt historical or ahistorical? If the former, then the state of theindividual in occasiont-1 will be included inmodeling.

2. Modeling approach - Should the models beestimated as generalized linear models (GLMs) or mixed linear models?While most function-based matrix models are estimated as the former, thelatter approach can account for repeated observations of the sameindividual by including individual identity as a random factor. Mixedmodels also allow time and patch to be treated as random variables,which allows for broader and more theoretically soundinference.

3. Suite of factors - Should both size andreproductive status be tested as causal factors? Should two-wayinteractions be included? Alternatively, should models be estimated asconstants?

4. Suite of vital rates - Which adult demographicparameters should be estimated? The defaults are (1) survival, (3)primary size, and (7) fecundity. Should (2) observation status, (6)reproductive status, or (4) secondary size or (5) tertiary size also bemodeled?

5. Juvenile vital rate estimation - Should juvenileparameters (8) through (14) also be modeled separately from adultparameters?

6. Best-fit criterion - If a model with fewerparameters exists within 2.0 AICc units of the model with the lowestAICc, then should this model be used as the best-fit model (thedefault), or should the model with the lowest AICc always bechosen?

7. Size distribution - Should size be modeled as acontinuous variable under a Gaussian or Gamma distribution, or as acount variable under either a Poisson or negative binomial distribution?If a count variable, then should the distribution be zero-inflated toaccount for excess zeros, zero-truncated to account for a lack of zeros,or left unaltered?

8. Fecundity distribution - Should fecundity bemodeled as a continuous variable under a Gaussian or Gamma distribution,or as a count variable under either a Poisson or negative binomialdistribution? If a count variable, then should the distribution bezero-inflated to account for excess zeros, zero-truncated to account fora lack of zeros, or left unaltered?

9. Continuous distribution estimation method - Ifusing a Gaussian or Gamma distribution for size, then should theprobability of transition be estimated via the cumulative densityfunction (CDF) method (the default), or the midpoint method? SeeDoak et al. (2021) for a discussion of therelative merits of both.

10. Timing of fecundity - Functionmodelsearch() assumes that linear models of fecundity use ametric counted or measured in occasiont as the response. Thisapplies well with most herbaceous plant datasets, where flowers or seedsproduced in one year might be the fecundity response measured. However,users not wishing to follow this default behavior can use thefectime option to stipulate a fecundity metric measured inoccasiont+1.

11. Age - Is an age or age x stage MPM the maingoal? Currently only age MPMs and ahistorical age x stage MPMs areoffered inlefko3.

12. Individual and annual covariates - Shouldindividual or annual covariates be tested as causal factors on vitalrates? If so, then should they be included as fixed, quantitative termsor as random, categorical terms?

13. Stage groups - If stage groups are noted in thestageframe, then should they also be used as fixed, independentcategorical predictors in modeling?

14. Missing values for occasions and patches -Should occasions or patches with coefficients that drop to 0 in modelingbe estimated as random draws from the corresponding distributions ofoccasion or patch coefficients?

15. Censoring - Should data points marked asquestionable be used or eliminated?

16. Variable names - The names of all relevantvariables in the dataset need to be specified. Note that the defaultbehavior assumes variable names output via thehistoricalize3() orverticalize3() functions,which produce standardized historically-formatted vertical (hfvformat) datasets.

Once all decisions are made and associated input terms are provided,modelsearch() goes to work. The result is alefkoMod object, which is a list in which the first 14elements are the best-fit models developed for each vital rate. Theseare followed by an equivalent number of elements showing the full modeltables developed and tested, followed by a table detailing the names anddefinitions of terms used in modeling, an element detailing the best-fitcriterion used, and ending on quality control data showing the number ofindividuals and the number of unique transitions used in the estimationof each model, as well as the overall accuracy of the models. Dependingon user choices, linear modeling is handled through thelm() andglm() functions in thestats package, theglm.nb() function in theMASS package, thelmer() andglmer() functions inlme4(Bates et al. 2015), theglmmTMB()function inglmmTMB(Brooks et al.2017), theglm.nb() function in packageMASS, thevglm() function in packageVGAM(Yee and Wild 1996; Yee2015), or thezeroinfl() function in packagepscl(Zeileis, Kleiber, and Jackman2008). Exhaustive model building proceeds through thedredge() function in packageMuMIn(Bartoń 2014). Model selection is handledthrough assessment of AICc and the number of parameters estimated permodel (see 6.Best-fit criterion above).

Ifmodelsearch() is set for historical analysis(historical = TRUE, the default), then the decision ofwhether to develop a historical MPM can be made on the basis of whetherany best-fit vital rate model includes size or reproductive status inoccasiont-1. If at least one vital rate does, then ahistorical MPM is justified. Otherwise, it is not. Regardless, theoutput can be used to create a function-based MPM in the next step.

Advanced users oflefko3 may wish to create their ownmodels without the package’s automated model selection function. In thatcase,lefko3’s matrix functions can accommodate singlemodels developed using the functions and packages mentioned before.

Step 4: Estimate matrices

MPM creation can be accomplished with nine different functions:rlefko2() for raw ahistorical stage-classified MPMs,rlefko3() for raw historical stage-classified MPMs,arlefko2() for raw ahistorical age-by-stage MPMs,rleslie() for raw age-classified MPMs,flefko2() for function-based ahistorical MPMs and IPMs,flefko3() for function-based historical MPMs and IPMs,aflefko2() for function-based ahistorical age-by-stageMPMs,fleslie() for function-based age-classified MPMs, andmpm_create() as a general, all-purpose MPM and IPM creationfunction. These functions incorporate binary kernels developed to handlethe estimation of matrix elements quickly and efficiently. A single runofflefko3(), for example, should be able to yield allannual matrices for all patches for theCypripedium candidumdataset provided withlefko3 in under a half-minute on mostmachines (14s or so on RPS’ 2019 MacBook Pro with 2.3 GHz 8-Core IntelCore i9). Parallel computing should not be necessary even with theslowest of current machines, provided that the machine is current enoughto handle at least R 3.6.3.

The output for each of these functions is alefkoMatobject, which is an S3 object (i.e. a list with pre-defined elements).These lists include elementsA,U, andF, which are themselves lists of complete projectionmatrices, survival-transition matrices, and fecundity matrices,respectively. For example, code such asmatobject$A[[1]]would access the first complete projection matrix in alefkoMat object namedmatobject. They arefollowed by elements referred to asahstages,hstages, andagestages, which provide dataframes describing the stages (technically showing the stageframe, thoughedited by the function for consistency), the historical stage-pairs, andthe age-stage pairs shown in the order in which they occur within thematrices. Thelabels element is a data frame giving adescription of each matrix in the order in which it occurs within theA,U, andF elements, includingthe population, patch, and monitoring time designations. The finalelements are quality control outputs that vary in content depending onwhether the output matrices are raw or function-based, but nonethelessalways include at least the numbers of individuals andindividual-transitions used for estimation.

Users wishing to check the output MPMs for errors or simply understandmore about them may utilize functionssummary() andcond_hmpm(). Functionsummary() summarizeseachlefkoMat object used as input, and also checks thesurvival of stages, stage-pairs, or age-stages by assessing the columnsums for eachU matrix and providing warnings if anysurvival probabilities fall below 0 or are greater than 1. Functioncond_hmpm() decomposes hMPMs into MPMs showing transitionsfrom occasiont to occasiont+1, but conditioned onthe same stage at occasiont-1. Thus, a single hMPM with fivestages and 25 stage pairs would be broken down into five conditionalMPMs, one for each stage at occasiont-1, with five columns andfive rows denoting stage at occasionst andt+1,respectively. These can be examined individually, and the survival ofstages from occasiont to occasiont+1 can be assessedas a function of stage in occasiont-1 by taking the associatedcolumn sums of these conditional matrices.

Step 4a. Importing matrices and IPMs

Some population ecologists will wish to utilize the power andfunctionality oflefko3 on already existing matrices.Although many analytical functions work with simple lists of matrices,users may import sets of matrices aslefkoMat objects inorder to utilize the most powerful analyses. Thecreate_lM() function allows users to take a list of squarematrices of equal dimensions, and import them into alefkoMat object. Users wishing to do this will also need todevelop a stageframe, and to input the order of populations, patches,and years associated with the order of matrices in the input list.Functioncreate_lM() runs several quality checks on theinput data, and assuming that there are no obvious inconsistencies,yields the finallefkoMat object for further use inanalysis.

Users may also import IPMs. If the characteristics of the kernelspropagating an IPM are properly described in a publication, and thekernel is composed of generalized linear models or generalized linearmixed models, then functionvrm_import() allows thisinformation to be input intolefko3 and used to recreatethe IPMs in discretized form. Please seelefko3: a gentleintroduction for more information, or see the function example.

Step 4b. Edit matrices

Once the MPM is created, users may wish to make modifications to some orall of the matrices. While thesupplemental() function wasused to fix certain transition values and transition multipliers acrossall matrices within the MPM creation function, theedit_lM() function can be used to make similar sorts ofchanges and updates to individual matrices or groups of matrices. Themethodology is similar to that used in functionsupplemental() (seeStep 2b: Provide supplementalinformation for matrix estimation), but with a few extraoptions thrown in for greater ability to customize matrices themselves.

Step 5. MPM analyses

Packagelefko3 includes a number of functions to aidanalysis once matrices are created. These may be of greater utility insome circumstances than established functions such aseigen(), because our functions are made to handle evenextremely large, sparse matrices. Currently, we include functions toestimate element-wise arithmetic mean matrices, discrete and stochasticpopulation growth rate, stable and long-run stage structure,reproductive value, deterministic and stochastic sensitivity andelasticity, and deterministic and stochastic life table responseexperiments (LTRE and sLTRE). We also include a basic projectionfunction that can work under deterministic, stochastic, cyclical, anddensity dependent assumptions. Bootstrapping is possible, with functionbootstrap3() providing bootstrapped replicates of the coredemographic dataset, all MPM generation functions yielding bootstrappedMPMs based on these bootstrapped data, and most MPM analysis functionscapable of handling the resulting bootstrapped MPMs. Plans are in theworks for further analysis functions in the future.