BASIC WORKFLOW
Our dataset is organized in horizontal format, with rows correspondingto unique individuals and columns corresponding to individual conditionin specific monitoring occasions (which we also refer to asyears). The original spreadsheet file used to keep the datasethas a repeating pattern to these columns, with each year having asimilarly arranged group of variables.
This dataset includes information on 1,119 individuals arrangedhorizontally, by row. There are 38 variables, by column. The first twocolumns identify each individual. This is followed by four sets of ninecolumns, each namedVolumeXX,lnVolXX,FCODEXX,FlowXX,IntactseedXX,Dead19XX,DormantXX,Missing19XX,andSeedlingXX, whereXX corresponds to theyear of observation and with years organized consecutively. Thus,columns 3-11 refer to year 1988, columns 12-20 refer to year 1989, etc.This strictly repeating pattern allows us to manipulate the originaldataset quickly and efficiently vialefko3. There are fouryears of data, from 1988 to 1991. Ideally, we should also have arrangedthe columns in the same order for each year, with years in consecutiveorder with no extra columns between them (this order is not required,provided that we input all variable names in correct order whenstandardizing the dataset).
Step 1. Life history model development
We will now create astageframe, which is a data framethat describes all stages in the life history of the organism in a wayusable by the functions in this package. It links the dataset to thelife cycle graph used to model the organism’s life history. Itshould include complete descriptions of all stages that occur in thedataset, with each stage defined uniquely. Since this object can be usedfor automated classification, all sizes, reproductive states, and othercharacteristics defining each stage in the dataset need to be accountedfor explicitly. The final description of each stage occurring in thedataset must not completely overlap with any other stage also found inthe dataset, although partial overlap is allowed and expected.
Here, we create a stageframe namedlathframe based on theclassification used inEhrlén (2000). Webuild this by creating vectors of the characteristics describing eachstage, with each element always in the same order within the vector,using thesf_create() function. The most important settingshave to do with the size bins for the stages. Two methods are currentlyallowed: 1) a vector of representative sizes in the option calledsizes along with the associated minimum and maximum sizesfor each stage in the options calledsizemin andsizemax, or 2) a vector of stage size bin midpoints in thesizes option along with the half-width of each size bin inthebinhalfwidth option. Note that the representativesizes, whether midpoint or otherwise, are not actually used in anycalculations - the size bin minima and maxima are the most importantdata used in calculations. Here, we will use approach 2), focused on theuse of size bin half-widths.
If size values are not to be binned, then narrow half bin-widths can beused (the default is 0.5, and all entries must be positive). Forexample, in this dataset, vegetatively dormant individuals necessarilyhave a size of zero, and so we can set thehalfbinwidth forthis stage to 0.5 provided that the resulting size bin does not overlapwith any other size bin matching the other characteristics of vegetativedormancy. We can also set this manually with thesizeminandsizemax options, if we prefer.
Here we use a single size variable, but up to three size variables maybe used. Vectorstagenames must include unique names only.Vectorsrepstatus,obsstatus,matstatus,immstatus,propstatus,andindataset are binomial vectors referring to status as areproductive stage, status as an observed stage, status as a maturestage, status as an immature stage, status as a propagule stage, andstatus as a stage occurring within the user-supplied dataset,respectively. The combination of these characteristics must becompletely unique for each stage. The final vector, calledcomments, holds strings of stage descriptions.
sizevector<-c(0,100,13,127,3730,3800,0)stagevector<-c("Sd","Sdl","Tm","Sm","La","Flo","Dorm")repvector<-c(0,0,0,0,0,1,0)obsvector<-c(0,1,1,1,1,1,0)matvector<-c(0,0,1,1,1,1,1)immvector<-c(1,1,0,0,0,0,0)propvector<-c(1,0,0,0,0,0,0)indataset<-c(0,1,1,1,1,1,1)binvec<-c(0,100,11,103,3500,3800,0.5)comments<-c("Dormant seed","Seedling","Tiny vegetative","Small vegetative","Large vegetative","Flowering","Vegetatively dormant")lathframe<-sf_create(sizes = sizevector,stagenames = stagevector,repstatus = repvector,obsstatus = obsvector,propstatus = propvector,immstatus = immvector,matstatus = matvector,indataset = indataset,binhalfwidth = binvec,comments = comments)#lathframeCare should be taken in assigning sizes to stages, particularly whenstages occur withsize = 0. In most cases, a size of zerowill mean that the individual is alive but not observable, such as inthe case of vegetative dormancy. However, a size of zero may have othermeanings. For example, if the size metric used is a logarithm of themeasured size, then observable sizes of zero and lower may be possible.These situations may impact matrix construction and analysis,particularly when dealing with function-based MPMs, such as IPMs.
Step 2a. Data standardization
Next, we will standardize the dataset intovertical format.Vertically formatted datasets are structured such that each rowcorresponds to the state of a single individual in two (if ahistorical)or three (if historical) consecutive monitoring occasions. Functionverticalize3() will create ahistorically formattedvertical dataset, orhfv dataset. We can use one of twoapproaches for this task, both using this function. In the most generalcase, we may input all column names from the original dataset in eachinput option in the same order, corresponding to monitoring occasion.This allows any format of dataset to be used. An easier approach may beused in datasets with strictly repeating sets of columns. In the lattercase, most of the inputs should constitute the names of the firstvariables coding for particular states. For example,Seedling1988 is the first variable coding for status as ajuvenile,Volume88 is the first variable coding for themain size metric, andIntactseed88 is the first variablecoding for fecundity. Functionverticalize3() can determinethe locations of variables for later monitoring occasions assuming arepeating pattern of such variables organized as blocks for each year(noted asblocksize, which here equals 9 because there arenine such variables per year in the same order). There are fourmonitoring occasions (noted asnoyears). We also have arepeated censor variable, the first of which isMissing1988, and we note here that we wish to censor thedata and to keep data points withNA values in the censorterm. Thepatchidcol andindividcol terms arevariables denoting which patch/subpopulation and individual each row ofdata belongs to, respectively. The settingNAas0 = TRUEtells R that missing values in size and fecundity should be interpretedas zeros, which allows us to infer incidents of vegetative dormancy ascases where size is equal to zero. Finally,stageassignties the dataset to the correct stageframe. Note that prior tostandardizing this dataset, we need to create a new individual identityvariable,indiv_id, composed of the identity of the patchthat the individual is in as well as its patch-specific ID, becauseotherwise individual identity will be shared by plants in differentpatches.
lathyrus$indiv_id<-paste(lathyrus$SUBPLOT, lathyrus$GENET)lathvert<-verticalize3(lathyrus,noyears =4,firstyear =1988,patchidcol ="SUBPLOT",individcol ="indiv_id",blocksize =9,juvcol ="Seedling1988",sizeacol ="Volume88",repstracol ="FCODE88",fecacol ="Intactseed88",deadacol ="Dead1988",nonobsacol ="Dormant1988",stageassign = lathframe,stagesize ="sizea",censorcol ="Missing1988",censorkeep =NA,censorRepeat =TRUE,censor =TRUE,NAas0 =TRUE)summary_hfv(lathvert)#>#> This hfv dataset contains 2527 rows, 54 variables, 1 population,#> 6 patches, 1053 individuals, and 3 time steps.The resulting reorganization has dramatically changed the dimensions ofthe dataset, which started with 1119 rows and 38 variables, and now has2527 rows and 54 variables. Theverticalize3() functionincludes error-checking measures designed to find instances whereindividual characteristics do not match those assigned to stages in theassociated stageframe. In those instances, warning messages will bedisplayed and the instances will be markedNoMatch instage1,stage2, orstage3. Thesummary_hfv() function allows us to quickly summarize themain characteristics of our resultinghfv dataset, and we cansee a longer summary by addingfull = TRUE.
If we removefull = FALSE from thesummary_hfv() call, then the extended portion of oursummary will reveal that theverticalize3() function hasautomatically subset the data to only those instances in which theindividual is alive in occasiont (please look at the variablealive2). Knowing this can help us interpret othervariables. For example, the mean value foralive3 suggestsvery high survival to occasiont+1 (92.2%). Further, theminimum values for all size variables are zero, suggesting thatunobservable stages occur within the dataset (these are instances ofvegetative dormancy).
Theverticalize3() function has made stage assignments forall individuals at each time. This can be seen in the subset summary inthestage2index column, which shows all individuals thatare alive and have a size of 0 in occasiont to be in theseventh stage. Typelathframe and enter at the prompt tocheck that the seventh stage in the stageframe is really vegetativedormancy.
A data subset summary can teach us how reproduction is handled here. Thereproductive status of flowering adults is set as reproductive inlathframe, and a subset summary of just those individualsflowering in occasiont shows all of those individuals set toreproductive (see the distribution of values forrepstatus2in the output below). However, fecundity ranges from 0 to 66 (seefeca2, which codes for fecundity in occasiont),meaning that some flowering adults do not actually produce anyoffspring. In fact, only 44.2% of flowering plants produced seed inoccasiont. This happened because our reproductive statusvariable,FCODE88, notes whether these individuals floweredbut not whether they produced seed. Since this plant must be pollinatedby an insect vector, some flowers should yield no seed. This issue doesnot cause problems in the creation of raw matrices, but it might causedifficulties in the creation of function-based matrices under someconditions. It helps to consider whether the definitions used for stagesare appropriate, and so whether reproductive status must necessarily beassociated withsuccessful reproduction or merely the attempt.Here, we associate it with the latter, but in other vignettes we willreconsider this assumption.
summary(lathvert[which(lathvert$stage2=="Flo"),c(25:30)])#> sizea2 size2added repstra2 repstr2added feca2 fec2added#> Min. : 98.4 Min. : 98.4 Min. :1 Min. :1 Min. : 0.000 Min. : 0.000#> 1st Qu.: 732.5 1st Qu.: 732.5 1st Qu.:1 1st Qu.:1 1st Qu.: 0.000 1st Qu.: 0.000#> Median :1141.8 Median :1141.8 Median :1 Median :1 Median : 0.000 Median : 0.000#> Mean :1388.5 Mean :1388.5 Mean :1 Mean :1 Mean : 4.793 Mean : 4.793#> 3rd Qu.:1758.0 3rd Qu.:1758.0 3rd Qu.:1 3rd Qu.:1 3rd Qu.: 6.000 3rd Qu.: 6.000#> Max. :7032.0 Max. :7032.0 Max. :1 Max. :1 Max. :66.000 Max. :66.000summary(lathvert[which(lathvert$stage2=="Flo"),c(31:36)])#> censor2 juvgiven2 obsstatus2 repstatus2 fecstatus2 matstatus2#> Min. :0 Min. :0 Min. :1 Min. :1 Min. :0.0000 Min. :1#> 1st Qu.:0 1st Qu.:0 1st Qu.:1 1st Qu.:1 1st Qu.:0.0000 1st Qu.:1#> Median :0 Median :0 Median :1 Median :1 Median :0.0000 Median :1#> Mean :0 Mean :0 Mean :1 Mean :1 Mean :0.4424 Mean :1#> 3rd Qu.:0 3rd Qu.:0 3rd Qu.:1 3rd Qu.:1 3rd Qu.:1.0000 3rd Qu.:1#> Max. :0 Max. :0 Max. :1 Max. :1 Max. :1.0000 Max. :1Step 2b. Provide supplemental information for matrix estimation
Now we will createsupplement tables, which provideexternal data for matrix estimation not included in the main demographicdataset. Specifically, we will provide the seed dormancy probability andgermination rate, which are given as transitions from the dormant seedstage to another year of seed dormancy or to the germinated seedlingstage, respectively. We assume that the germination rate is the sameregardless of whether the seed was produced in the previous year or hasbeen in the seedbank for longer. We will incorporate both terms asconstants for specific transitions within our matrices, and as constantmultipliers for fecundity, since fecundity will be estimated as theproduct of seed produced and either the seed germination rate or theseed dormancy/survival rate. The fecundity multipliers will also serveto tell R which transitions are the fecundity transitions.
lathsupp2<-supplemental(stage3 =c("Sd","Sdl","Sd","Sdl"),stage2 =c("Sd","Sd","rep","rep"),givenrate =c(0.345,0.054,NA,NA),multiplier =c(NA,NA,0.345,0.054),type =c("S","S","R","R"),stageframe = lathframe,historical =FALSE)#> Warning: NA values in argument multiplier will be treated as 1 values.#lathsupp2The supplement table above will only work with ahistorical andage-by-stage MPMs, while the next supplement table will work only forhistorical MPMs. The primary difference is the incorporation of stage intimet-1. Going back a further step to look at timet-1 shows us that there are two ways that a dormant seed intimet can be a dormant seed in timet+1 - it couldhave been a dormant seed in timet-1, or it could have beenproduced by a flowering adult in timet-1. Likewise, a seed inoccasiont could have become a seedling in occasiont+1 after being produced by a flowering adult in occasiont-1, or remaining in the seed bank in occasiont-1.So, we will enter four given transitions in the historical case, ratherthan two transitions as in the ahistorical case. We will also designatethat transitions from seed in occasiont-1 and seedling inoccasiont to mature stages in occasiont+1 shouldoriginate from theNotAlive stage in occasiont-1,because the seed stage does not occur in the dataset and so functionrlefko3() cannot estimate these transitions withoutassuming an appropriate proxy.
lathsupp3<-supplemental(stage3 =c("Sd","Sd","Sdl","Sdl","Sd","Sdl","mat"),stage2 =c("Sd","Sd","Sd","Sd","rep","rep","Sdl"),stage1 =c("Sd","rep","Sd","rep","npr","npr","Sd"),eststage3 =c(NA,NA,NA,NA,NA,NA,"mat"),eststage2 =c(NA,NA,NA,NA,NA,NA,"Sdl"),eststage1 =c(NA,NA,NA,NA,NA,NA,"NotAlive"),givenrate =c(0.345,0.345,0.054,0.054,NA,NA,NA),multiplier =c(NA,NA,NA,NA,0.345,0.054,NA),type =c("S","S","S","S","R","R","S"),type_t12 =c("S","F","S","F","S","S","S"),stageframe = lathframe,historical =TRUE)#> Warning: NA values in argument multiplier will be treated as 1 values.#lathsupp3These two supplement tables show us that we have survival-transitionprobabilities (type = 1, whereas fecundity rates would begiven astype = 2) and fecundity multipliers(type = 3), that the second and fourth transitions involvereproduction from occasiont-1 to occasiont while theothers involve survival (type_t12 = 1 for survival betweenoccasionst-1 andt, andtype_t12 = 2 forfecundity), that the given transitions originate from the dormant seedstage (Sd) in occasiont (and seeds orreproductive stages in occasiont-1 in the historical case),and the specific values to be used in overwriting and to multiplyfecundity values by:0.345 and0.054. If wewished, we could have used the values of transitions to be estimatedwithin this matrix as proxies for these values, in which case we wouldenter the stages corresponding to the correct transitions in theeststageX columns, and thegivenrate columnwould be blank. This is what we did with one set of transitions - namelytransitions from the seedling class to mature stages. Finally, we havealso included fecundity multipliers for newly produced seed in thebottom two rows.
Note that thesupplemental() function allowed us to isolatespecific transitions to alter, and to use shorthand to identify largegroups of transitions (e.g. usingmat,immat,rep,nrep,prop,npr,obs,nobs,groupX, orall to signify all mature stages,immature stages, reproductive stages, non-reproductive stages, propagulestages, non-propagule stages, observable stages, non-observable stages,stages within groupX, or simply all stages, respectively).Functionsupplemental() also allowed us to tell R whichtransitions to treat as primary reproductive transitions, which wasaccomplished by identifying these transitions with reproductivemultipliers.
Step 3. Tests of history
We have chosen to build both ahistorical and historical MPMs in thisvignette. However, in a typical analysis, it is most parsimonious totest whether history influences the demography of the populationsignificantly first, using historical MPMs only if the test supportsdoing so. A number of methods exist to conduct these tests, and werecommendBrownie et al. (1993),Pradel, Wintrebert, and Gimenez (2003),Pradel (2005), andColeet al. (2014) for good discussions and tools to help with this.
Inlefko3, functionmodelsearch() can be usedto test for the effects of history directly from the historical verticaldataset. This function estimates best-fit linear models of the key vitalrates used to propagate elements in function-based MPMs. There are up to14 different vital rates possible to test, of which seven are adultvital rates and seven are juvenile vital rates. The standard vital ratesthat we may wish to test are survival (marked assurv invitalrates), primary size (size), andfecundity (fec), which are the default vital rates assumedby the function. Here, we also test observation status(obs), which can serve as a proxy for sprouting probabilityin cases where plants do not necessarily sprout, and reproductive status(repst), which assesses the probability of reproduction incases where reproduction is not certain. This dataset also includesjuveniles whose vital rates we wish to estimate. We designate this bysettingjuvestimate to the correct juvenile stage in thedataset. Because size varies in juveniles, we also setjuvsize = TRUE.
To test history with functionmodelsearch: 1) use thehistorical vertical dataset as input, 2) sethistorical = TRUE, 3) input the relevant vital rates toestimate, 4) set the suite of independent factors to test (size andreproductive status in occasionst andt-1 and allinteractions, or some subset thereof), 5) set the name of the juvenilestage (if juvenile vital rates are to be estimated and such stages occurin the dataset), 6) set the proper distributions to use for size andfecundity, and 7) note which variables code for individual identity(used to treat identity as a random factor in mixed linear models),patch identity (if there are multiple patches and vital rates should beestimated with patch as a factor), and observation occasion. We also setquiet = "partial" to limit the amount of text output whilethe function runs.
histtest<-modelsearch(lathvert,historical =TRUE,suite ="main",vitalrates =c("surv","obs","size","repst","fec"),juvestimate ="Sdl",sizedist ="gaussian",fecdist ="gaussian",indiv ="individ",year ="year2",juvsize =TRUE,quiet ="partial")#>#> Developing global model of survival probability...#>#> Global model of survival probability developed. Proceeding with model dredge...#>#> Developing global model of observation probability...#>#> Global model of observation probability developed. Proceeding with model dredge...#>#> Developing global model of primary size...#>#> Global model of primary size developed. Proceeding with model dredge...#>#> Developing global model of reproduction probability...#>#> Global model of reproduction probability developed. Proceeding with model dredge...#>#> Developing global model of fecundity...#>#> Global model of fecundity developed. Proceeding with model dredge...#>#> Developing global model of juvenile survival probability...#>#> Global model of juvenile survival probability developed. Proceeding with model dredge...#> Warning: Juvenile maturity status in time t+1 appears to be constant (1). Setting to constant.#>#> Developing global model of juvenile observation probability...#>#> Global model of juvenile observation probability developed. Proceeding with model dredge...#>#> Developing global model of juvenile primary size...#>#> Global model estimation failed. Dropping individual identity term.#>#> Redeveloping global model of juvenile primary size...#>#> Global model of juvenile primary size developed. Proceeding with model dredge...#> Warning: Juvenile reproductive status in time t+1 appears to be constant, and so will be set to constant.#>#> Finished selecting best-fit models.#summary(histtest)The summary generated is quite long, and likely resulted in a series ofwarnings from the model-building functions utilized. For our purposes,we need only look at elements corresponding to the best-fit models forour tested vital rates, which are the elements markedsurvival_model,observation_model,size_model,repstatus_model,fecundity_model,juv_survival_model,juv_observation _model,juv_size_model,juv_reproduction_model, andjuv_maturity_model. The line beginningFormula: in each of these sections shows the best-fitmodel, in standard R formula notation (i.e. \(y = ax + b\) is given asy ~ x). The independent terms tested include variablescoding for size in occasionst andt-1, given assizea2 andsizea1, and reproductive status inoccasionst andt-1, given asrepstatus2andrepstatus1, respectively. Several vital ratesincorporatesizea1 orrepstatus1 in thebest-fit models (notably adult survival, primary size, reproductionprobability, and fecundity), and so individual history does have asignificant impact on the demography of this population. The qualitycontrol information in elementqc is also of interest,showing that some of our models are quite strong, with accuracy ofroughly 1.0 for adult survival and 0.903 for observation, while othersare less accurate, most notably primary size and fecundity.
Step 4. MPM estimation
Now let’s create some raw Lefkovitch MPMs based onEhrlén (2000). We have seen that history shouldbe included in these analyses, which justifies creating only historicalmatrices. However, to introduce these functions in greater depth anddetail, we will first create ahistorical MPMs. The functions we will useto build these MPMs includerlefko2() for raw ahistoricalMPMs, andrlefko3() for raw historical MPMs.
Ehrlén (2000) shows a mean matrix coveringyears 1989 and 1990 as occasiont. We will utilize the entiredataset instead, covering 1988 to 1991. Note that we will not creatematrices for subpopulations in this case (to include them, add theoptionspatch = "all", patchcol = "patchid" to the inputbelow).
ehrlen2<-rlefko2(data = lathvert,stageframe = lathframe,year ="all",stages =c("stage3","stage2"),supplement = lathsupp2,yearcol ="year2",indivcol ="individ")#ehrlen2The output from this analysis is alefkoMat object, whichis list with the following elements:
A: a list of full population projection matrices, inorder of population, patch, and year (order given inlabels)
U: a list of matrices where non-zero entries arelimited to survival-transition elements, in the same order asA
F: a list of matrices where non-zero entries arelimited to fecundity elements, in the same order asA
hstages: a data frame showing the order of pairedstages (given if matrices are historical; otherwiseNA)
agestages: a data frame showing the order ofage-stage pairs (given if matrices are age-by-stage ahistorical;otherwiseNA)
ahstages: the stageframe used in analysis, withstages potentially reordered and edited as they occur in the matrix
labels: a table showing the order of matrices bypopulation, patch, and year
matrixqc: a short vector used insummary statements to describe the overall quality of eachmatrix (used insummary() calls)
dataqc: a short vector used insummarystatements to describe key sampling aspects of the dataset (used insummary() calls)
Input options for functionrlefko2() includeyear = "all", which can be changed toyear = c(1989, 1990) to focus just on years 1989 and 1990,as inEhrlén (2000), oryear = 1989 to focus exclusively on the transition from1989 to 1990 (the year entered is the year in occasiont).Matrix-estimating functions inlefko3 have a defaultbehavior of creating matrices for each year in the dataset except thefinal year, rather than lumping all years together to produce a singlematrix. However, patches or subpopulations will only be separated if apatch ID variable is provided as input.
We can understandlefkoMat objects in greater detailthrough thesummary() function.
We learn that three fullA matrices and theirU andF decompositions were estimated, andthat they are ahistorical. This is expected given that there are fourconsecutive years of data, yielding three time steps, and an ahistoricalmatrix requires two consecutive years to estimate transitions. Thefollowing line notes the dimensions of those matrices. The third,fourth, fifth, and sixth lines of the summary show how many survivaltransition and fecundity elements were actually estimated, both overalland per matrix, the number of populations, patches, and time stepscovered by the MPM, and the number of individuals and transitions thematrices are based on (these last numbers can be used to understandmatrix quality). The next section shows summaries of the column sumsfrom the survival-transition (U) matrices in thislefkoMat object. These column sums correspond to thesurvival probabilities of the different life stages, and so thesummaries must show numbers ranging from 0.0 to 1.0. If any matrices hadelements below 0.0, and if any survival matrix columns summed to anumber greater than 1.0, then we would see warnings to that effect inthe summary output. Finally, the very last line of the summary notesthat the first matrix appears to have a broken life cycle, meaning thatsome of the transitions expected based on our assumed life cycle areequal to 0.0.
The matrix creation functions in this package sort the stages providedin the stageframe according to a standardized rubric, and so the orderof stages in theahstages element of alefkoMat object may differ from the original input.Particularly, they sort propagules first, followed by immature stages,followed by non-reproductive adult stages, followed by reproductiveadult stages. Let’s see the order in our matrices. We see that the keydifference from the original input is that the order of the floweringstage and the vegetative dormancy stage have been flipped.
Now we’ll estimate historical matrices. Because of the size of thesematrices, we will only show thelefkoMat summary. We willonly create Ehrlén-format matrices - users wishing to createdeVries-format hMPMs should addformat = "deVries" to theinput options (the resulting matrices will be bigger, but contain thesame number of estimated elements).
ehrlen3<-rlefko3(data = lathvert,stageframe = lathframe,year ="all",stages =c("stage3","stage2","stage1"),supplement = lathsupp3,yearcol ="year2",indivcol ="individ")#summary(ehrlen3)The summary output shows several differences from the ahistorical case.First, there is one less matrix estimated in the historical case than inthe ahistorical case, because raw historical matrices require threeconsecutive occasions to estimate each transition instead of two.Second, these matrices are larger than ahistorical matrices, with thenumbers of rows and columns generally equaling the number of ahistoricalrows and columns squared (although deVries-format increases thedimensions through the addition of a prior stage for newborns, and inboth Ehrlén and deVries formats the numbers of rows and columns may bereduced under some circumstances as shown in the next block). Finally, amuch greater proportion of each matrix is composed of zeros in thehistorical case than in the ahistorical case, although there arecertainly more non-zero elements as well. This sparseness results fromhistorical matrices being composed primarily of structural zeros. As aresult, the historical matrices in this example have non-zero entries in(79 + 7) / 2401 = 3.5% of matrix elements, while the equivalentahistorical matrices have non-zero entries in (24.67 + 2) / 49 = 54.4%of matrix elements.
We can see the impact of structural zeros by eliminating some of them inthe process of matrix estimation. This can be done by settingreduce = TRUE, which tellsrlefko3() toeliminate stage pairs in which both column and row are zero vectors.Here, we now have matrices with 19 fewer rows and columns, and (79 + 7)/ 900 = 9.6% of elements as potentially non-zero.
ehrlen3red<-rlefko3(data = lathvert,stageframe = lathframe,year ="all",stages =c("stage3","stage2","stage1"),supplement = lathsupp3,yearcol ="year2",indivcol ="individ",reduce =TRUE)#summary(ehrlen3red)Next we will create the element-wise mean ahistorical matrix. Functionlmean() creates alefkoMat object retainingall descriptive information from the originallefkoMatobject. The full output not in any order includes the composite meanmatrix (shown in elementA), the mean survival-transitionmatrix (U), the mean fecundity matrix (F), adata frame outlining the definitions and order of historical pairedstages (hstages, shown asNA in this casebecause the matrices are ahistorical), a data frame outlining theage-by-stage combinations used (agestages, shown asNA because these are not age-by-stage matrices), a dataframe outlining the actual stages as outlined in thestageframe object used to create these matrices(ahstages), a data frame outlining the definitions andorder of the matrices (labels), and two quality controlvectors used in output for thesummary() function(matrixqc anddataqc).
Now we will estimate the historical element-wise mean matrix. We willshow only the top-left corner of the rather large matrix (a sectioncomprised of the first twenty rows and eight columns of the 49 x 49matrix).
The prevalence of zeros in this matrix is normal because most elementsare structural zeros and so cannot equal anything else. This matrix isalso a raw matrix, meaning that transitions that do not actually occurin the dataset cannot equal anything other than zero. Now let’s take alook at thehstages object associated with this meanmatrix.
There are 49 pairs of life history stages, corresponding to the rows andcolumns of the historical matrices. The pairs are interpreted so thatmatrix columns represent stage pairs in occasionst-1 andt, and rows represent stage pairs in occasionst andt+1. For an element in the matrix to contain a number otherthan zero, it must represent the same stage at occasiont inboth the column stage pairs and the row stage pairs. For example, Theelement [1, 1] represents the transition probability from dormant seedat occasionst-1 andt (column pair), to dormant seedat occasionst andt+1 (row pair) - the occasiont stages match, and so this element is possible. However,element [1, 2] represents the transition probability from seedling inoccasiont-1 and very small adult in occasiont(column pair), to dormant seed in occasiont and in occasiont+1 (row pair). Clearly [1, 2] is a structural zero because itis impossible for an individual to be concurrently a dormant seed and avery small adult.
Error-checking is more difficult with historical matrices because theyare typically one or two orders of magnitude bigger than theirahistorical counterparts, but the same basic strategy can be used hereas with ahistorical matrices. In these cases we can usesummary() function to assess the key quality controlcharacteristics of the mean hMPM, such as the distribution of survivalprobability estimates for historical stages.
Let’s try another approach, looking at some conditional historicalmatrices. Conditional matrices are matrices of ahistorical dimensionthat show transitions from stage at occasiont to occasiont+1, conditional on all individuals having been in the samestage in occasiont-1. They are calculated from historicalMPMs, and the output below shows all conditional matrices developed fromthe firstA matrix ofehrlen3mean. The firstmatrix, for example, shows all transitions involving individuals thathad been in the dormant seed stageSd in occasiont-1, while the last matrix shows transitions involvingindividuals that had been vegetatively dormant in occasiont-1.
Quick scans will show many transitions missing, because each stage hasonly certain stages that can transition from it, and to which it cantransition. Further transitions are missing because the MPM is raw, andso some transitions are not parameterized because no individuals madethem.
One last error-checking technique before we analyze our MPMs: matrixvisualization plots. These plots provide a relatively easy way tounderstand the “spatial spread” of values throughout a matrix. Packagelefko3 includes functionimage3(), whichprovides an easy way to make these images. Let’s start off by looking atthe ahistorical mean matrix.
The resulting image shows non-zero elements as red spaces, and zeroelements as white spaces. Rows and columns are numbered, and we can seethat this matrix is reasonably dense.
Now let’s take a look at the mean historical matrix. The historical meanmatrix has many more rows and columns, and has more of both zero andnon-zero elements. However, it has become a sparse matrix, and it turnsout that increasing the numbers of life stages will increase thesparsity of historical projection matrices.
Step 5. MPM analysis
Packagelefko3 includes functions to conduct some analysesof population dynamics. We will start by estimating the asymptoticpopulation growth rate (\(\lambda\))and the stochastic population growth rate (\(a= \text{log} \lambda _{S}\)) from the ahistorical MPMs, includingboth the annual MPM and the mean. For the stochastic case, we will setthe seed forR’s random number generation to make ouroutput reproducible. Note that each\(\lambda\) estimate includes a data framedescribing the matrices in order (given as thelabelsobject within the output list). Here is the set of ahistorical annual\(\lambda\) estimates, followed by\(\lambda\) for the mean matrix, andthe stochastic population growth rate (\(a =\text{log} \lambda _{S}\)).
# Deterministic ahistoricallambda3(ehrlen2)#> pop patch year2 lambda#> 1 1 1 1988 0.8952585#> 2 1 1 1989 0.9235493#> 3 1 1 1990 1.0096490# Deterministic mean ahistoricallambda3(ehrlen2mean)#> pop patch lambda#> 1 1 1 0.9574162# Stochastic ahistoricalset.seed(42)slambda3(ehrlen2)#> pop patch a var sd se#> 1 1 1 -0.04490197 0.03154986 0.1776228 0.001776228We can also estimate standard errors using bootstrapped data. Packagelefko3 makes this relatively easy, using thebootstrap3() function with the previously used functions.
lathboot<-bootstrap3(lathvert,reps =100)ehrlen2_boot<-rlefko2(data = lathboot,stageframe = lathframe,year ="all",stages =c("stage3","stage2"),supplement = lathsupp2,yearcol ="year2",indivcol ="individ")ehrlen2_boot_means<-lmean(ehrlen2_boot)e2bm_lambdas<-lambda3(ehrlen2_boot_means)lambda_se<-sd(e2bm_lambdas$lambdas)/sqrt(100)lambda_se#> [1] 0.0007775337We will now analyze the historical MPM. There are several differences inthe output in addition to the lower growth rate estimates. First,because there are four years of data, there are three ahistoricaltransitions possible for estimation: year 1 to 2, year 2 to 3, and year3 to 4. However, in the historical case, only two are possible: fromyears 1 and 2 to 3 (technically, from years 1 [t-1] and 2[t] to years 2 [t] and 3 [t+1]), and fromyears 2 and 3 to 4. Second, historical matrices cover more of theindividual heterogeneity in a population by splitting ahistoricaltransitions by stage in occasiont-1. This heterogeneity mayreflect many sources, for example the impacts of trade-offs operatingacross years(Shefferson and Roach 2010).One particularly common trade-off is the cost of growth: an individualthat grows a great deal in one time step due to great environmentalconditions in that year might pay a large cost of survival, growth, orreproduction in the next if those environmental conditions deteriorate(Shefferson, Warren II, and Pulliam 2014;Shefferson et al. 2018).
# Deterministic historicallambda3(ehrlen3)#> pop patch year2 lambda#> 1 1 1 1989 0.8863920#> 2 1 1 1990 0.9855435# Deterministic mean historicallambda3(ehrlen3mean)#> pop patch lambda#> 1 1 1 0.9182809# Stochastic historicalset.seed(42)slambda3(ehrlen3)#> pop patch a var sd se#> 1 1 1 -0.08872396 0.0177556 0.1332501 0.001332501We can also examine the stable stage distributions, as follows for theahistorical case.
ehrlen2mss<-stablestage3(ehrlen2mean)ehrlen2mss#> matrix stage_id stage ss_prop#> 1 1 1 Sd 0.29261073#> 2 1 2 Sdl 0.04579994#> 3 1 3 Tm 0.22875637#> 4 1 4 Sm 0.18625613#> 5 1 5 La 0.07716891#> 6 1 7 Dorm 0.05588715#> 7 1 6 Flo 0.11352077The data frame output shows us the stages themselves(stage, and associated number instage_id),which matrix they refer to (matrix), and the stable stagedistribution (ss_prop). Interpreting these values, we findthat the mean matrix suggests that, if we project the population forwardindefinitely assuming the population dynamics are static and representedby this matrix, we will find that approximately 29% of individualsshould be dormant seeds (suggesting a large seedbank). A further 23% and19% should be very small and small adults, respectively, and 11% shouldbe flowering adults. Almost 6% of the population should eventually becomposed of vegetatively dormant adults.
We can estimate the stable stage distribution for the historical case,as well. Because the historical output for thestablestage3() function is a list with two data frames,let’s take a look at each of these data frames in turn. The first willbe the stage-pair output.
This data frame is structured in historical format, and so shows thestable stage distribution of stage pairs. We may wish to see which stagepair dominates, in which case we might look at the row with the maximumss_prop value.
ehrlen3mss$hist[which(ehrlen3mss$hist$ss_prop==max(ehrlen3mss$hist$ss_prop)),]#> matrix stage_id_2 stage_id_1 stage_2 stage_1 ss_prop#> 17 1 3 3 Tm Tm 0.2404459About 24% of the population is expected to be composed of tiny adultsmaintaining themselves as tiny adults.
The longer format of the historical stable stage output makes it harderto read. However, historical values can also be combined by stage atoccasiont (stage_2) to estimate thehistorically-corrected stable stage distribution shownin theahist element, which allows comparison to a stablestage distribution estimated from a purely ahistorical MPM. Notice belowthat thess_prop column shows values that are differentfrom the ahistorical MPM, suggesting the influence of individualhistory.
ehrlen3mss$ahist#> matrix stage_id stage ss_prop#> 1 1 1 Sd 0.26275496#> 2 1 2 Sdl 0.04112686#> 3 1 3 Tm 0.29476860#> 4 1 4 Sm 0.19794066#> 5 1 5 La 0.06553412#> 6 1 7 Dorm 0.04645369#> 7 1 6 Flo 0.09142110To see the impact of history on the stable stage distribution, let’splot the ahistorical and historically-corrected stable stagedistributions together. We will also include the stochastic long-runstage distribution in our output, which is estimated with the samefunction but using thestochastic = TRUE option. This willallow us to see the impact of random temporal variation.
ehrlen2mss_s<-stablestage3(ehrlen2,stochastic =TRUE,seed =42)ehrlen3mss_s<-stablestage3(ehrlen3,stochastic =TRUE,seed =42)ss_put_together<-cbind.data.frame(ehrlen2mss$ss_prop, ehrlen3mss$ahist$ss_prop, ehrlen2mss_s$ss_prop, ehrlen3mss_s$ahist$ss_prop)names(ss_put_together)<-c("d_ahist","d_hist","s_ahist","s_hist")rownames(ss_put_together)<- ehrlen2mss$stagebarplot(t(ss_put_together),beside=T,ylab ="Proportion",xlab ="Stage",ylim =c(0,0.35),col =c("black","orangered","grey","darkred"),bty ="n")legend("topright",c("det ahistorical","det historical","sto ahistorical","sto historical"),pch =22,col ="black",pt.bg =c("black","orangered","grey","darkred"),bty ="n")Let’s take the deterministic portion first. Accounting for individualhistory increased the prevalence of tiny and small adults, but decreasedthe prevalence of dormant seeds, and dormant, large, and floweringadults. Now when we also take in the impact of temporal stochasticity,we can see differences in the proportion of dormant seeds, seedlings,and all adult stages, with the greatest differences in dormant seeds andtiny adults.
Let’s take a look at the reproductive values now, in similar order tothe stable stage distribution case. Initially, we will create all setsof reproductive value objects, and then we will plot them. The structureof these objects is the same as that of the stable stage structureoutputs. Because the four vectors are all standardized such that thefirst non-zero reproductive value is set to 1.0, they are on differentscales, and so we will make them comparable for plotting purposes bystandardizing them against their vector sums.
ehrlen2mrv<-repvalue3(ehrlen2mean)ehrlen3mrv<-repvalue3(ehrlen3mean)ehrlen2mrv_s<-repvalue3(ehrlen2,stochastic =TRUE,seed =42)ehrlen3mrv_s<-repvalue3(ehrlen3red,stochastic =TRUE,seed =42)rv_put_together<-cbind.data.frame((ehrlen2mrv$rep_value/sum(ehrlen2mrv$rep_value)), (ehrlen3mrv$ahist$rep_value/sum(ehrlen3mrv$ahist$rep_value)), (ehrlen2mrv_s$rep_value/sum(ehrlen2mrv_s$rep_value)), (ehrlen3mrv_s$ahist$rep_value/sum(ehrlen3mrv_s$ahist$rep_value)))names(rv_put_together)<-c("det ahist","det hist","sto ahist","sto hist")rownames(rv_put_together)<- ehrlen2mrv$stagebarplot(t(rv_put_together),beside=T,ylab ="Relative reproductive value",ylim =c(0,0.4),xlab ="Stage",col =c("black","orangered","grey","darkred"),bty ="n")legend("topleft",c("det ahistorical","det historical","sto ahistorical","sto historical"),pch =22,col ="black",pt.bg =c("black","orangered","grey","darkred"),bty ="n")Both deterministic and stochastic analyses show that flowering adultshave the greatest reproductive value in both ahistorical and historicalanalysis, while dormant seeds have the least. Historical MPMs suggestlower contributions of seedlings and vegetative dormancy, but largercontributions of flowering adults.
Now we will look at the deterministic sensitivities of\(\lambda\) to the ahistorical mean matrixelements.
ehrlen2sens<-sensitivity3(ehrlen2mean)#> Running deterministic analysis...print(ehrlen2sens$ah_sensmats[[1]],digits =3)#> [,1] [,2] [,3] [,4] [,5] [,6] [,7]#> [1,] 0.0182 0.00284 0.0142 0.0116 0.00479 0.00347 0.00705#> [2,] 0.2061 0.03225 0.1611 0.1312 0.05434 0.03936 0.07994#> [3,] 0.2565 0.04016 0.2006 0.1633 0.06766 0.04900 0.09953#> [4,] 0.4110 0.06432 0.3213 0.2616 0.10838 0.07849 0.15943#> [5,] 0.5639 0.08826 0.4408 0.3589 0.14871 0.10770 0.21877#> [6,] 0.3067 0.04801 0.2398 0.1952 0.08089 0.05858 0.11899#> [7,] 0.7221 0.11302 0.5645 0.4596 0.19043 0.13791 0.28014# The highest sensitivity value is:max(ehrlen2sens$ah_sensmats[[1]])#> [1] 0.7220817# This occurs in element:which(ehrlen2sens$ah_sensmats[[1]]==max(ehrlen2sens$ah_sensmats[[1]]))#> [1] 7# The highest sensitivity value among biologically plausible elements is:max(ehrlen2sens$ah_sensmats[[1]][which(ehrlen2mean$A[[1]]>0)])#> [1] 0.4596282# This occurs in element:which(ehrlen2sens$ah_sensmats[[1]]==max(ehrlen2sens$ah_sensmats[[1]][which(ehrlen2mean$A[[1]]>0)]))#> [1] 28The highest sensitivity value is associated with a biologicallyimpossible transition - dormant seeds (stage/column 1) cannot transitionto flowering (stage/row 7) (element 7). Among biologically plausibleelements, the highest sensitivity is associated with element 28, whichis the transition from small adult (stage/column 4) to flowering(stage/row 7).
We will now look at the sensitivity of\(\lambda\) to elements in the historicalmean MPM.
The first element produced in this analysis ish_sensmats,which is a list composed of sensitivity matrices of the historicalmatrices, in order (we have only one in our mean matrixlefkoMat object). These matrices are the same dimensions asthe historical matrices used as input, and so can be quite huge. This isfollowed byah_sensmats, which is a list composed ofhistorically-corrected sensitivity matrices of corresponding ahistoricalmatrix elements (calculated using the historically-corrected stablestage distribution and reproductive value vector produced inehrlen3mss andehrlen3mrv, respectively).These values are different from the sensitivities estimated from theahistorical matrices themselves, but are formatted similarly. Next,h_stages andah_stages give the order ofpaired stages and life history stages used in the historical andhistorically-corrected sensitivity matrices, respectively. Finally, wehave the originalA,U, andFmatrices used as input.
Our historical matrices are large and full of zeros. So, we will lookfor the highest sensitivity associated with a biologically plausibleelement (i.e. non-zero matrix elements). Then, we will assess thehighest biologically plausible sensitivity in the historically-correctedsensitivity matrices, to compare against the ahistorical sensitivityanalysis.
# The highest sensitivity value among biologically plausible elements:max(ehrlen3sens$h_sensmats[[1]][which(ehrlen3mean$A[[1]]>0)])#> [1] 0.3148525# This value is associated with element:which(ehrlen3sens$h_sensmats[[1]]==max(ehrlen3sens$h_sensmats[[1]][which(ehrlen3mean$A[[1]]>0)]))#> [1] 802# Highest historically-corrected sensitivity value among plausible elements is:max(ehrlen3sens$ah_sensmats[[1]][which(ehrlen2mean$A[[1]]>0)])#> [1] 0.8591236# This occurs in element:which(ehrlen3sens$ah_sensmats[[1]]==max(ehrlen3sens$ah_sensmats[[1]][which(ehrlen2mean$A[[1]]>0)]))#> [1] 20The maximum biologically plausible sensitivity value in the historicalmatrix is element 802, which is associated with column 17 (tiny adult inoccasionst-1 andt) and row 18 (tiny adult inoccasiont to small adult in occasiont+1). Thistransition is from tiny adult in occasionst-1 andtto small adult int+1. The historically-corrected sensitivityanalysis finds that\(\lambda\) is mostsensitive to element 20, which is the transition from tiny adult todormant. This is a little different from the ahistorical MPM, whichsuggested element 28 (small adult to flowering adult).
We can perform the same analyses as above with stochastic sensitivities.In that circumstance, simply run thesensitivity3()function with thestochastic = TRUE argument set, as below.
ehrlen2sens_s<-sensitivity3(ehrlen2,stochastic =TRUE)#> Running stochastic analysis...ehrlen3sens_s<-sensitivity3(ehrlen3,stochastic =TRUE)#> Running stochastic analysis...A complementary approach to sensitivity analysis is elasticity analysis.Elasticities are easier to interpret because projection matrix elementsequal to 0 produce elasticity values also equal to 0, thus eliminatingbiologically impossible transitions from consideration. Elasticityvalues are also scaled to sum to 1.0, making elasticities of survivaltransitions easier to compare to those of fecundity. However, they arealso interpreted differently, because while sensitivity analysis showsthe impact of a tiny butabsolute change to a matrix element on\(\lambda\), elasticity analysis showsthe impact of a tiny butproportional change to a matrixelement on\(\lambda\). In fact, bothsensitivities and elasticities are essentially local slopes, and so arenot unit free. It is therefore not unusual for sensitivity andelasticity analysis to yield different inferences.
Let’s look at the elasticity of\(\lambda\) to matrix elements in theahistorical mean matrix.
ehrlen2elas<-elasticity3(ehrlen2mean)#> Running deterministic analysis...print(ehrlen2elas$ah_elasmats,digits =3)#> [[1]]#> [,1] [,2] [,3] [,4] [,5] [,6] [,7]#> [1,] 0.00655 0.00000 0.000000 0.0000 0.00000 0.00000 0.0116#> [2,] 0.01162 0.00000 0.000000 0.0000 0.00000 0.00000 0.0206#> [3,] 0.00000 0.02784 0.151678 0.0164 0.00052 0.00415 0.0000#> [4,] 0.00000 0.00127 0.041102 0.1586 0.02210 0.02184 0.0167#> [5,] 0.00000 0.00000 0.000604 0.0290 0.04851 0.01744 0.0532#> [6,] 0.00000 0.00315 0.007180 0.0223 0.00896 0.00628 0.0107#> [7,] 0.00000 0.00000 0.000000 0.0353 0.06862 0.00887 0.1673# The maximum elasticity value:max(ehrlen2elas$ah_elasmats[[1]])#> [1] 0.167339# This value is associated with element:which(ehrlen2elas$ah_elasmats[[1]]==max(ehrlen2elas$ah_elasmats[[1]]))#> [1] 49Elasticity analysis exhibits strong differences from sensitivityanalysis. In particular, we find that\(\lambda\) is most strongly elastic inresponse to changes in element 49, which is the stasis transition forflowering adults (stage 7). We can sum the columns of the elasticitymatrix to see which stages\(\lambda\)is most and least elastic in response to, as below.
print(colSums(ehrlen2elas$ah_elasmats[[1]]),digits =3)#> [1] 0.0182 0.0323 0.2006 0.2616 0.1487 0.0586 0.2801Here we see that\(\lambda\) is moststrongly elastic in response to changes in transitions associated withflowering adults, followed by transitions involving small adults.Dormant seeds and seedlings have the smallest impact on\(\lambda\), and the impacts of fecundity(shown in the top-right corner of the elasticity matrix) appear quitesmall.
Now on to elasticity analysis of the historical MPMs. Once again, wewill not output the matrices. Typeehrlen3elas at theprompt to see these matrices.
ehrlen3elas<-elasticity3(ehrlen3mean)#> Running deterministic analysis...# The highest deterministic elasticity value:max(ehrlen3elas$h_elasmats[[1]])#> [1] 0.1718096# This value is associated with element:which(ehrlen3elas$h_elasmats[[1]]==max(ehrlen3elas$h_elasmats[[1]]))#> [1] 801The highest elasticity appears to be associated with the element 801,which is at row 17, column 17. This corresponds to the stasis transitionof tiny adults (tiny int-1 to tiny int to tiny int+1).
Elasticities are often treated as additive, making the calculation ofhistorically-corrected elasticity matrices easy. These are stored in theah_elasmats element ofelasticity3() outputoriginating from a historical MPM. Eyeballing this matrix, it appearsthat the historically-corrected elasticity matrix supports stasis in thetiny adult stage as the transition that\(\lambda\) is most elastic to, with stasisin vegetative dormancy a close second.
print(ehrlen3elas$ah_elasmats,digits =3)#> [[1]]#> [,1] [,2] [,3] [,4] [,5] [,6] [,7]#> [1,] 0.00699 0.000000 0.00000 0.0000 0.00e+00 0.01161 0.00000#> [2,] 0.01161 0.000000 0.00000 0.0000 0.00e+00 0.00000 0.00000#> [3,] 0.00000 0.010736 0.19877 0.0344 8.26e-05 0.00000 0.00219#> [4,] 0.00000 0.000292 0.04500 0.1615 2.56e-02 0.02781 0.01907#> [5,] 0.00000 0.000000 0.00000 0.0306 4.54e-02 0.05418 0.00920#> [6,] 0.00000 0.000000 0.00000 0.0366 6.09e-02 0.16395 0.00497#> [7,] 0.00000 0.000583 0.00238 0.0162 7.42e-03 0.00889 0.00315Next, we will create a barplot of the elasticities of life historystages from ahistorical vs. historically-corrected analyses. We willalso incorporate stochastic elasticity analysis here to assess theimportance of temporal environmental stochasticity on population growth.
ehrlen2elas_s<-elasticity3(ehrlen2,stochastic =TRUE)#> Running stochastic analysis...ehrlen3elas_s<-elasticity3(ehrlen3,stochastic =TRUE)#> Running stochastic analysis...elas_put_together<-cbind.data.frame(colSums(ehrlen2elas$ah_elasmats[[1]]),colSums(ehrlen3elas$ah_elasmats[[1]]),colSums(ehrlen2elas_s$ah_elasmats[[1]]),colSums(ehrlen3elas_s$ah_elasmats[[1]]))names(elas_put_together)<-c("det ahist","det hist","sto ahist","sto hist")rownames(elas_put_together)<- ehrlen2elas$ahstages$stagebarplot(t(elas_put_together),beside=T,ylab ="Elasticity",xlab ="Stage",col =c("black","orangered","grey","darkred"),bty ="n")legend("topleft",c("det ahistorical","det historical","sto ahistorical","sto historical"),pch =22,col ="black",pt.bg =c("black","orangered","grey","darkred"),bty ="n")Historical analyses generally find that population growth rate is lesselastic in response to seedlings, and flowering adults, and more elasticto tiny and dormant adults than in ahistorical analyses. Stochastic anddeterministic population growth rates are barely elastic in response todormant seeds and seedlings. Note the dramatic impact of environmentalstochasticity and individual history combined on the impact of tinyadults.
Let’s now look at the elasticities of different kinds of transitions. Wewill use thesummary() function, which outputs data framessummarizing elasticity sums by the kind of transition. First, we willcompare ahistorical against historically-corrected transitions.
ehrlen2elas_sums<-summary(ehrlen2elas)ehrlen3elas_sums<-summary(ehrlen3elas)ehrlen2elas_s_sums<-summary(ehrlen2elas_s)ehrlen3elas_s_sums<-summary(ehrlen3elas_s)elas_sums_together<-cbind.data.frame(ehrlen2elas_sums$ahist[,2], ehrlen3elas_sums$ahist[,2], ehrlen2elas_s_sums$ahist[,2], ehrlen3elas_s_sums$ahist[,2])names(elas_sums_together)<-c("det ahist","det hist","sto ahist","sto hist")rownames(elas_sums_together)<- ehrlen2elas_sums$ahist$categorybarplot(t(elas_sums_together),beside=T,ylab ="Elasticity",xlab ="Transition",col =c("black","orangered","grey","darkred"),bty ="n")legend("topright",c("det ahistorical","det historical","sto ahistorical","sto historical"),pch =22,col ="black",pt.bg =c("black","orangered","grey","darkred"),bty ="n")We see similar patterns across types of elasticity analysis.Particularly, population growth rate is most elastic in response tochanges in stasis transitions, and least elastic to changes infecundity. Environmental stochasticity appears to exacerbate thesedifferences.
Packagelefko3 also includes functions to conduct manyother analyses, including deterministic and stochastic life tableresponse experiments, and general projection including quasi-extinctionanalysis and density dependent analysis. Users wishing to conduct theseanalyses should see our free e-manual calledlefko3: agentle introduction, as well as through the resources found onthe{r}evolutionarydemography website.