As outlined in my earlier posts, particularly blog post 4 on Growth Rate Issues, choosing and fitting a good model for growth is not straightforward. I mentioned in a very early post that I had difficulty producing reliable growth parameter estimates for abalone in the 1980’s because the estimates varied according to the data subsets I used (varying length classes and truncated or expanded times at liberty).
Knight (1968), Francis (1988a,1988b), Sainsbury (1990) and Haddon et al (2008), among many others, have addressed the difficulties. In particular, Wang & Thomas (1995) posited that there were no currently available models which accommodated the individual variability apparent in most datasets and that most existing models should be used with caution. Haddon (2021) has very clearly worded explanations of the issues.
I am hoping that this time around I will make better progress with turban snail data and the currently available software than I did with abalone in the 1980’s.
Two key problems with the von Bertalanffy model are the high correlation between the two key parameters, K and Linfinity (see Pilling et al, 2002 and Figures 5 to 7 below) and the ambiguity surrounding asymptotic growth. Most forms of the von Bertalanffy model predict linearly decreasing growth increments with increasing size at tagging. An alternative Gompertz model predicts increasing increments for juveniles then decreasing increments for adults. The inverse logistic model suits linear juvenile growth which then decreases in adulthood but may also flatten out.
My simple take on this, not having the mathematical smarts of the authors mentioned above, is that growth rates of juveniles, and of larger (older) snails near a supposed asymptotic length, are important for model selection and interpretation. Accordingly, my approach will be to use multiple models, and with each one, apply multiple filters so that I understand what is in the data, and how this might relate to the biology of these snails, which is more my focus than the maths.
For the asymptotic growth issue, in an earlier post I mentioned the likelihood that while a snail lives, it grows (albeit slowly). My largest Lunella shell, at 125mm, is much larger than the estimated Linfinity values from data modelling (85-95mm). Similarly, I have several Turbo militaris shells much larger than my early Linfinity estimates of around 110-120mm. Haddon et al (2008) also raised the possibility of indeterminate growth in which the dynamics would permit growth to continue unabated (possibly very slowly) until each animal died.
A third problem is individual variability, and I assume there are two research questions here. Do individuals have their own intrinsic K and Linfinity? Are those individual parameters constant or do they vary with environmental factors? These are not in my list of research questions for this observational study, but I may be able to make some anecdotal observations in due course, and possibly some defensible statistical observations, after I have addressed my six primary research questions. Hampton (1991) developed a model which I will try to fit with my data. At the moment, from my ecology/biology perspective and my observations to date, I could comfortably support Wang and Thomas’ suggestion that each individual might have a different, constant, genetically determined Linfinity and a different, non-constant, environmentally determined K.
In addressing these issues, I will use several models, ranging from the Gulland-Holt and Ford-Walford linear models already introduced in earlier blogs, to various nonlinear models designed to deal with seasonal growth, individual variability and other “error” factors. The process will be iterative for the simpler models to mirror the more complex models which can simultaneously incorporate site, size, sex and season differences.
All the models seem to predict well enough, for fisheries management purposes, over median size and length increment ranges (my early observations and Haddon et al 2008). The differences are in how they perform for small and very large sizes. Length-at-tagging and time-at-liberty biases in tagging data can therefore be critical for addressing confounding factors such as seasonal growth, outlier data, and individual variability.
Being uncomfortable with non-linear models and matrix algebra, but knowing that’s where I must go, I prefer starting with graphical tools to get a feel for these issues in my data before applying more complex analyses. I am not alone. Both Wang and Thomas (2008) and Haddon (2021) used the basic plot of increment vs length-at-tagging, which is quite illustrative. My basic plots are increment/month (rather than increment) vs length at tagging and increment/month vs midLength – the length halfway between release and recapture (Gulland-Holt plot). Previous blog posts and this one have numerous examples. Both methods will have some bias because the rate is calculated over a range and is not instantaneous.
And I remind readers again that there are often “tagging disturbance” rings on the shells of recaptured snails, so all my growth rate estimates are potentially too small. I have thousands of length frequency measurements for L. torquata so I may eventually see some differences between size-transition modes/means and tagging increments.
An example of using a graph to see what is happening is Figure 1 below for female Lunella torquata which results from converting parameters obtained from a linear analysis to an age-length growth curve. This has been the dominant practice in fisheries science until recently: using growth models from tagging to produce an estimated age-length key for use in age-structured stock assessments. (More recently, size-structured models have been developed using size-transition matrices rather than parametric equations, because an age length key cannot be produced from tagging data alone. Age for at least one given length must be known, assumed [e.g. tzero is 0mm], or simulated.)
L
Figure 1. von Bertalanffy growth curves for female Lunella torquata from all sites combined, based on recapture data using times-at-liberty less than and greater than 60 days. Data points represent actual recapture lengths (L2’s) and assumed age-at-recapture calculated from theoretical age-at-release (using L1’s and VBG parameters) plus actual dT’s. The two K and Linfinity values are from Gulland-Holt plots for the two data sets (Figure 2 below).
Note that the data points at the zero-age line demonstrate one of the problems with von Bertalanffy curves from tagging data. Recaptures with a length-at-tagging (L1) greater than the estimated Linfinity cannot be plotted. These are exactly the ones which can best inform growth rates at large sizes (near the asymptote).
Note also how short-period recaptures (brown dots) do not match longer-term recaptures, reflecting the differences in Figures 2 and 3 below (for both males and females). I assume that these differences arise partly from the assumptions in the Gulland-Holt method of fitting the von Bertalanffy curve (brief discussion below and my next blog post).
Figure 2. FEMALES. Three Gulland-Holt plots for female Lunella torquata using (a) all recaptures, (b) those at liberty for less than 60 days, and (c) those at liberty for 60 or more days. Note the outlier of 6mm/month and several others in the middle plot disappear from the bottom plot, which has the highest R2 value.
Figure 3. MALES. Three Gulland-Holt plots for male Lunella torquata using (a) all recaptures, (b) those at liberty for less than 60 days, and (c) and those at liberty for 60 or more days (c). Note the unrealistic Linfinity for the middle plot.
Two possible factors in the different distributions (but not the only ones) are seasonal differences (a future blog post) and just simple extrapolation of growth rate from two or three weeks to one month. If the average measurement/estimation error is around 0.5 to 1mm, (see blog post 9) the error can be increased by scaling mm/week up to mm/month and can easily equal the overall average rate of around 1 to 2 mm/month. This “effect size” from statistical power theory is one reason for large sample sizes and longer times at liberty (but not too long?) The impact on parameter estimation can be high because the distortion will be relatively high at small sizes and relatively low at larger sizes and there are fewer data points at the extremities. The same pattern appears in both male and female data (Table 1 below).
So, one key question for me is “should very short-term recaptures be excluded from analyses?” We did this in a Bayesian analysis of my data (Kienzle et al, 2022), arbitrarily choosing 60 days after I tested outputs using 30, 60 and 90 days. (See also Fig 3 of that paper using 3,6- and 9-month subsets). I am not alone in this. Wang (1998), for example, excluded prawn tag data with dT <14 days and lobster data with dT<28 days. This seems to be a common approach.
Figure 4. Growth increment per month for Lunella torquata recaptures with short times at liberty. The data are clumped due to batches of recaptures around 15, 40 and 58 days.
In the figure above, the range of increments reduces to a more reasonable level as dT approaches 60 days.
Just as there may be no ideal minimum time to use, there seems no ideal time-at-liberty range. Using recaptures after exactly 1 year is essential for the Ford/Walford model and for fitting Haddon’s (2021) inverse-logistic model. It eliminates seasonal differences but therefore yields no information about seasonal growth. Using subsets with only annual data also greatly reduces the sample size at both ends of the time-at-liberty range. Choosing a time interval may bias the outputs due to seasonal growth if it is strong factor: a time-at-liberty of 2-3 months can only be one or two seasons whereas 9-12 months must be three or four seasons.
Consider the trend with short-, medium- and long-term data outlined below.
Figure 5. Gulland-Holt plots for short-, medium- and long-term Lunella torquata recaptures (all sites, both sexes). Note that as minimum time at liberty increases, the slope and extrapolated maximum monthly increment (Y-intercept) generally decline.(K-value declines and Linfinity increases.)
Note that there are fewer data points less than 40mm in the left in the bottom three plots simply because snails grow more over those longer periods.
The effect of filtering for time-at-liberty is illustrated in the following table of Linfinity values for Lunella.
| Males | Females | Number of male and female recaptures used. |
All Times at liberty | 92.2 | 85.7 | 837 |
dT’s < 60 days | 155.8 | 125.7 | 99 |
dT’s >= 60 days | 90.5 | 82.9 | 738 |
dT’s >= 60 days using Grotag package in R (2024 final data) | 93.4 | 87.1 | 462/375 |
dT’s>=60d Kienzle et al 2022 (early data) | 93.7 | 88.2 | 272/184 |
2 to 3 months | 87.5 | 83.9 | 74/57 |
3 to 6 months | 88.2 | 82.4 | 150/117 |
6 to 9 months | 91.8 | 88.0 | 107/69 |
9 to 12 months | 91.5 | 78.3 | 41/38 |
dT’s >12 months | 116 | 85.4 | 52/31 |
Table 1. Estimated Linfinity for tagged and recaptured male and female Lunella torquata from all sites for different times at liberty. Estimates are from Gulland-Holt plots like Figure 2 except for the two middle rows (in bold) which used Francis’ GROTAG model.
The trend for both K and Linfinity estimates is shown below.
Fig 6. Effect of time at liberty subsets (short, medium, long) on Gulland-Holt estimates of K and Linfinity for Lunella torquata (male, female, all). Note from the bottom regression that K and Linfinity values are highly correlated.
The pattern holds true with even finer scale partitioning. Gulland-Holt parameter estimates (for all sexes and all sites) using times-at-liberty of 60-69 days, 70-79, 80-89 days etc., up to 200-209 days are shown below. (See discussion about G-H plot assumptions below.)
Figure 7. K and Linfinity estimates from G-H plots using approximately equal times at liberty from 60-69 days up to 200-210 days.
Other aspects of my data also illustrate the correlations between K and Linfinity. The pattern from seasonal data (a future blog post) is shown below in Figs 8a and 8b for both of my two species.
Figure 8a. Correlation between K and Linfinity estimates for Lunella torquata using Gulland-Holt analysis of tag recapture data. Each point is a parameter pair estimated where times-at-liberty were wholly within one or two seasons (eg Summer, Spring, Summer-Autumn, Winter-Spring, etc.,)
Figure 8b Correlation between K and Linfinity estimates for Turbo militaris using Gulland-Holt analysis of tag recapture data where times-at-liberty were wholly within one or two seasons.
Gulland and Holt indicated a method of enhancing their short, equal time-at-liberty model by using a first approximation of K to calculate, for long periods at liberty, not dL/dT but dL/dT*b/(tanh b) then adding those data to the short-period dataset (b=Ka/2 where a is time at liberty.) However, I have not found studies where this has been done.
A more widely used solution to the correlation between K and Linfinity has been to adopt a “forced” G-H plot whereby the regression is forced through a fixed (known?) Linfinity to obtain a more consistent value for K. In using this model for my Lunella data, I chose values of Linfinity from the 2 GROTAG analyses shown in Table 1 above. I am more confident about those values.
Lunella torquata | FEMALES | MALES |
Linfinity used | 87.6 | 93.6 |
Mean mm/month (y) | 1.11 | 1.12 |
Mean midLength (U) | 60.24 | 62.27 |
K from above means K=y(Linfinity-U)-1 | 0.486 | 0.431 |
Mean of K values from each recapture | 0.471 | 0.430 |
Sample size | 319 | 427 |
Standard G-H plot K values from Figs 2 and 3 above | 0.58 | 0.48 |
K values from Kienzle et al 2022 | 0.41 | 0.43 |
Table 2. “Forced” Gulland-Holt K values and associated data for Lunella torquata males and females (all sites, all seasons) with dT>=60 days assuming the “true” values of Linfinity are as shown in the first row.
In selecting various models to use, I am cognisant of the assumptions for those models. A particular example is the Gulland-Holt model because it is quite widely used (eg seahorses, Harasti et al 2012; tilapia aquaculture, DeGraaf & Prein 2005; Adriatic clams, Ezgeta-Balic et al 2011; lobster, Ulmestrand & Eggert 2001). The slope of the regression equals -K *(tanh b)/b where b =Ka/2 and a is “a fixed duration of time”. If a is small, then (tanh b)/b is close to 1.0. The key assumption by Gulland and Holt (1959) was that b is usually “relatively small” and for b values up to 4, the error is 5% or less. All papers I have reviewed ignored the (tanh b)/b correction factor.
Gulland and Holt suggested that sets of data with different times at liberty would have different slopes but the same Linfinity. This seems to be the logic for “fixing” Linfinity in forced G-H plots. To explore that, since I have a lot of recaptures, I calculated K’s and Linfinity’s for subsets of recaptures in 10-day groups from 60 days to 210 days.
Mean days at liberty for the 10-day interval | K | Linfinity | Sample size |
67 | 0.49 | 102 | 68 |
74 | 0.87 | 79 | 84 |
86 | 0.37 | 98 | 65 |
92 | 0.55 | 87 | 58 |
103 | 0.59 | 84 | 81 |
116 | 0.62 | 81 | 52 |
126 | 0.45 | 89 | 22 |
136 | 0.51 | 86 | 28 |
144 | 0.61 | 84 | 54 |
154 | 0.49 | 84 | 18 |
165 | 0.63 | 86 | 25 |
172 | 0.12 | 203 | 11 |
184 | 0.49 | 90 | 28 |
196 | 0.47 | 86 | 26 |
202 | 0.51 | 89 | 83 |
Table 3. K and Linfinity estimates for Lunella torquata from Gulland-Holt plots using equal (+/- 5 days) but increasing times at liberty. (All sites, both sexes, all seasons combined).
Fig 9. K and Linfinity estimates for Lunella torquata from Gulland-Holt plots using equal but increasing times at liberty. (All sites, both sexes, all seasons combined). X axis is the mean number of days for the subset. Bars are Linfinity – left axis, line is K values – right axis.)
The results show that there is some correlation between K and Linfinity (ie Linfinity is NOT constant for these data – although there is indeed much less variation in Linfinity than in K). I am not sure if the variation shown here and above suggests Gulland and Holt’s logic is faulty. It is more likely due to different proportions in my data of sexes, sites and seasons in the various subsets. However, surely the same would have been true for species they examined, although benthic invertebrates are apparently more variable than finfish. I should stress that I am not out to critique the work of others, just understand it all better using data (and a species) that I am familiar with.
So, while the visualisations above are useful to me, parameter estimation based on graphical (linear) methods are clearly problematic. The variability highlighted by various filters reinforces the value of the Haddon et al (2008) and Francis (1988b) approaches of estimating all factors simultaneously and not blindly using the most common model.
While the larger number of sample parameters estimated simultaneously in the more sophisticated, multi-parameter models of Haddon and Francis might reduce or remove some biases, seems clear to me, and the work by Y G Wang seems to confirm it, that the possible biases in linear methods might also apply to some of the non-linear fitting methods.
Since I have plenty of data, I will exclude short-term recaptures (<60 days) from most future analyses as the default position but will also do some fits with all data included, to get a better personal understanding of the issues. Similarly, it will obviously be useful to use short, medium and long time-at-liberty datasets in those analyses as well, for males and females separately, and for different seasons. This will require many runs of each model using all data as well as subsets relevant to each research question – sex, size, season, species and site. The more models I use, with different subsets of data, the better my understanding will be.
While I will be running models for males and females separately to flesh out biologically interesting season and site differences, there is also value in fitting aggregated data for fisheries management purposes, because males and females cannot be identified separately in the field by fishers or researchers, and because a lightly fished resource requires only minimal and broad-brush management intervention, which can be based on averages or the most cautious parameters, depending on the management issue being addressed.
The next few blog posts will present the linear and non-linear model outputs.
References
De Graaf, G. and Prein, M. (2005), Fitting growth with the von Bertalanffy growth function: a comparison of three approaches of multivariate analysis of fish growth in aquaculture experiments. Aquaculture Research, 36: 100-109. https://doi.org/10.1111/j.1365-2109.2004.01191.x
Ezgeta-Balić, D., Peharda, M., Richardson, C. A., Kuzmanić, M., Vrgoč, N., & Isajlović, I. (2011). Age, growth, and population structure of the smooth clam Callista chione in the eastern Adriatic Sea. Helgoland Marine Research, 65, 457-465.
Francis R I C C 1988a. Are Growth Parameters Estimated from Tagging and Age-Length Data Comparable? Can.J.Fish.Aquat.Sci. 45(6):936-942
Francis R I C C 1988b. Maximum likelihood estimation of growth and growth variability from tagging data NZ J.Mar.Freshwater Res. 22(1):43-51 doi.org/10.1080/00288330.1988.9516276
Gulland J A and S J Holt. 1959. Estimation of growth parameters for data at unequal time intervals. J.Cons.Int.Explor.Mer 25(1):47-49
Haddon M 2021. Using R for Modelling and Quantitative Methods in Fisheries. CRC Press R Series. Taylor & Francis.
Haddon M, C Mundy and D Tarbath, 2008. Using an inverse logistic model to describe growth increments of blacklip abalone (Haliotis rubra) in Tasmania. Fish. Bull.106:58-71
Hampton J 1991. Estimation of southern bluefin tuna Thunnus maccoyii growth parameters from tagging data, using von Bertalanffy models incorporating individual variation. Fish Bull. U.S. 89:577-590
Harasti, D., Martin‐Smith, K., & Gladstone, W. (2012). Population dynamics and life history of a geographically restricted seahorse, Hippocampus whitei. Journal of Fish Biology, 81(4), 1297-1314.
Kienzle M, M Broadhurst , G Hamer. 2022. Bayesian estimates of turban snail (Lunella torquata) growth off south-eastern Australia. Fisheries Research 248, 106218. https://doi.org/10.1016/j.fishres.2021.106218
Knight W, 1968. Asymptotic Growth: An Example of Nonsense Disguised as Mathematics. Journal of the Fisheries Research Board of Canada 25(6):1303-1307 https://doi.org/10.1139/f68-114
Pilling G M, G P Kirkwood and S G Walker. 2002. An improved method for estimating individual growth variability in fish, and the correlation between von Bertalanffy growth parameters. Can.J.Fish.Aquat.Sci. 59(3) doi.org/10.1139/f02-022
Sainsbury K J 1980. Effect of Individual Variability on the von Bertalanffy Growth Equation. Can.J.Fish.Aquat.Sci. 37(2) doi.org/10.1139/f80-031
Ulmestrand, M., & Eggert, H. (2001). Growth of Norway lobster, Nephrops norvegicus (Linnaeus 1758), in the Skagerrak, estimated from tagging experiments and length frequency data. ICES Journal of Marine Science, 58(6), 1326-1334.
Wang YG. 1998. An improved Fabens method for estimation of growth parameters in the von Bertalanffy model with individual asymptotes. Can.J.Fish.Aquat.Sci. 55:397-400
Wang YG, M Thomas and IF Somers. 1995. A maximum likelihood approach for estimating growth from tag-recapture data. Can.J.Fish.Aquat.Sci. 52:252-259
Wang YG and MR Thomas. 1995. Accounting for individual variability in the von Bertalanffy growth model. Can.J.Fish.Aquat.Sci. 52:1368-1375