Erickson and Nadon: Stepwise stochastic simulation for distributions of missing life history parameter values 89 
variability is interspecific genetic differences that are 
a result of natural selection. This source of variability 
was modeled explicitly by using parameter estimates 
for individual species as data points in our regression 
models. A certain amount of intraspecific geographical 
variability in both genetics and habitat conditions also 
contributes to variability in estimates of life history 
parameters. Consequently, there may be some concern 
that we pooled data from across the world and there- 
fore ignored geographical differences within species. 
Unfortunately, given that published data for most spe- 
cies came from only one location, it was impossible for us 
to directly control for geographical heterogeneity in life 
history parameters. Geographical differences in values 
of life history parameters for species simply add to the 
overall variability of the stepwise parameter estimates 
and should not cause any specific bias. 
Another source of variability that was directly mod- 
eled are regression model coefficients. Coefficient 
variability was included by randomly sampling the mul- 
tivariate normal distribution, as defined by the variance— 
covariance matrix. Similarly, uncertainty in the starting 
Lynax Parameter was directly included by sampling from 
a distribution or by bootstrapping raw L,,,,, observations. 
Furthermore, the lack of access to raw age-at-length and 
length-at-maturity data, coupled with the fact that most 
published reports from life history studies do not include 
a parameter variance—covariance matrix, prevented us 
from directly including uncertainty from fitted values 
(i.e., L.., K, and L,,,,) in our models. This source of uncer- 
tainty, although not modeled explicitly, likely adds to the 
global residual error. 
Finally, a few minor sources of uncertainty were not 
directly included because they are unlikely to have a major 
effect on output variability. Specifically, variability in Ay is 
hard to measure given the sensitivity of this parameter 
to juvenile age-at-length data points. In our approach, we 
simply fixed this parameter to the overall mean of —0.6. 
Further, we did not include uncertainty in the relation- 
ships of fork or standard length to total length that we 
used to convert all parameter data to total lengths. Not 
including uncertainty from length conversions is unlikely 
to be a major source of error, and the studies that publish 
these conversion factors rarely provide their associated 
variability. Finally, we did not specifically include vari- 
ability in the M estimates derived from A,,,, (Equation 8). 
The relationship between M and longevity depends on the 
assumption of S. We fixed S at 0.05 to convert our raw 
observations of A,,,, to M in order to model the M~K rela- 
tionship. This assumption can be relaxed in the R package 
StepwiseLH by specifying an S value and generating M 
estimates under this new value. Uncertainty in the values 
of the S parameter could be included by researchers as a 
separate step by using Monte Carlo simulations on this 
parameter. 
Despite the various sources of uncertainty described in 
this section, it is important to note that use of the stepwise 
approach still results in reasonably accurate parameter 
estimates, as demonstrated in our tests. 
Parameter estimation 
In addition to producing results from the 4 models that 
relate key life history traits, we calculated the following 
ratios: Lyjax to La, Lmat tO Liamax> Lmat to L.., and M to K, 
all Beverton—Holt invariants. Results from compari- 
son of these ratios from our study with those previously 
published indicate that the invariants differ among tax- 
onomic groups in our study (groupers, grunts, wrasses, 
and lamniform and carcharhiniform sharks). Our findings 
are consistent with previous investigations of life history 
invariants that found similar variance in relationships of 
life history traits between taxonomic groups of fish (Char- 
nov and Berrigan, 1991; Nadon and Ault, 2016; Thorson 
et al., 2017). 
The ratio of M to K is of particular interest because 
the relationship between M and K determines the over- 
all shape of a growth trajectory (Hordyk et al., 2015). In 
our study, we used the criteria from Prince et al. (2015) 
to define a ratio of M to K that was >1.0 as indeterminate 
growth and a ratio of M to K that was <1.0 as determinate 
growth. Species with determinate growth are limited to a 
Maximum size after which growth mostly ceases (Lincoln 
et al., 1982), and species with indeterminate growth are 
not limited to a maximum size but typically see growth 
slow with size (Sebens, 1987). Sharks were the only group 
with a ratio of M to K that was >1.0 (1.53), indicating inde- 
terminate growth. This observation for sharks is consistent 
with the findings of indeterminate growth in the literature 
(Pardo et al., 2013; Heupel et al., 2014). All other families 
in the study described here had determinate growth pat- 
terns with ratios of M to K <1.0. 
On average, species with indeterminate growth will 
reach only a fraction of L,, at A,,,,, With an increasingly 
smaller proportion at higher ratios (Hordyk et al., 2015). 
Consequently, the ratio of M to K has a strong effect on 
the divergence between the parameters L,, and Lamy. AS 
this ratio increases, the curvilinearity of the growth trajec- 
tory decreases and L,, can increase drastically and become 
more of a fitting parameter than one related to Lg,,,,. An 
increase in the ratio of M to K ultimately reduces the ratio 
of Linat to L.., while the ratio of L,,,¢ to Lamax Stays similar. 
The described interactions between life history ratios can 
also lead to L,, being much higher than L,,,,, an outcome 
that may seem surprising but is to be expected for spe- 
cies with high ratios of M to K. This effect of the stepwise 
approach on the ratio of M to K was observed in our study: 
the groups with the highest ratios of M to K (groupers and 
sharks) had the greatest differences between the ratios of 
Lat to L., and of Lys tO Lamax- Subsequently, as in the 
study of Nadon and Ault (2016), we used the L,,4~Lamax 
relationship as the preferred model for obtaining L,,,t- 
Our decision was further validated by r* values being 
higher when L,,,,x, Was used as the predictor than when 
L,. was used. 
To varying degrees, the relationships estimated with 
stepwise simulation in this study were similar to those 
established by Nadon and Ault (2016) with positive 
linear relationships between parameters in both the 
