Skip Navigation


ICES Journal of Marine Science: Journal du Conseil Advance Access originally published online on January 5, 2007
ICES Journal of Marine Science: Journal du Conseil 2007 64(2):234-247; doi:10.1093/icesjms/fsl025
This Article
Right arrow Abstract Freely available
Right arrow FREE Full Text (PDF) Freely available
Right arrow All Versions of this Article:
64/2/234    most recent
fsl025v1
Right arrow Alert me when this article is cited
Right arrow Alert me if a correction is posted
Services
Right arrow Email this article to a friend
Right arrow Similar articles in this journal
Right arrow Alert me to new issues of the journal
Right arrow Add to My Personal Archive
Right arrow Download to citation manager
Right arrowRequest Permissions
Google Scholar
Right arrow Articles by Cotter, A. J. R.
Right arrow Articles by Piet, G. J.
Right arrow Search for Related Content
PubMed
Right arrow Articles by Cotter, A. J. R.
Right arrow Articles by Piet, G. J.
Social Bookmarking
 Add to CiteULike   Add to Connotea   Add to Del.icio.us  
What's this?

Crown Copyright © 2007. Published by Oxford Journals on behalf of the International Council for the Exploration of the Sea. All rights reserved

Estimating stock parameters from trawl cpue-at-age series using year-class curves

A. J. R. Cotter1,, B. Mesnil2 and G. J. Piet3

1 Cefas Laboratory, Pakefield Road, Lowestoft, Suffolk NR33 0HT, UK
2 Département EMH, IFREMER, BP 21105, F44311 Nantes Cedex 03, France
3 Institute for Marine Resources and Ecosystem Studies (Wageningen IMARES), Haringkade 1, IJmuiden, Netherlands 1970 AB

Correspondence to A. J. R. Cotter: tel: +44 1502 524323; fax: +44 1502 513865; e-mail: john.cotter{at}cefas.co.uk

Cotter, A. J. R., Mesnil, B., and Piet, G. J. 2007. Estimating stock parameters from trawl cpue-at-age series using year-class curves. – ICES Journal of Marine Science, 64: 234–247.

A year-class curve is a plot of log cpue (catch per unit effort) over age for a single year class of a species (in contrast to the better known catch curve, fitted to multiple year classes at one time). When linear, the intercept and slope estimate the log cpue at age 0 and the average rate of total mortality, Z, respectively. Here, we suggest methodological refinements within a linear least squares framework. Candidate models may include a selectivity term, fleet-specific parameters, and polynomials in year to allow for gradual variations of Z. An iterative weighting method allows for differing precisions among the different fleets, and a forward (one-step ahead) validation procedure tests predicted cpue against observed values. Choice of the best approximating model(s) is made by ranking the biological credibility of each candidate model, then by comparing graphic plots, precision of prediction, and the Akaike Information Criterion. Two example analyses are (i) a comparison of estimated and true results for five stock simulations carried out by the US National Research Council, and (ii) modelling three beam trawl surveys for plaice (Pleuronectes platessa) in the North Sea. Results were consistent with known, age-related, offshore migrations by plaice. Year-class curves are commended as a widely applicable, statistically based, visual, and robust method.

Keywords: analysis of relative residual variance, cpue, fish stock assessment, forward validation, iteratively weighted least squares, North Sea, plaice, year-class curve

Received 5 May 2006; accepted 26 October 2006; advance access publication 5 January 2007.


    Introduction
 Top
 Introduction
 Theory
 Examples
 Discussion
 Appendix
 References
 
A more or less linear decline in the logarithm of catch per unit effort (cpue) indices over age is typically observed for many commercially important marine fish caught in trawls (Jensen, 1939; Silliman, 1943; Beverton and Holt, 1957, section 13), and has even been reported for the larval stages of Atlantic mackerel caught in plankton nets (Sette, 1943). It has also been seen for several species of European trawl-caught fish, at least for the older, fully selected ages. A noteworthy feature, also reported by Beverton and Holt (1957, p. 182), was that the slopes, which estimate a time-averaged value of Z, the coefficient of total mortality, showed little variation either across different year classes or over a period of more than a decade (Cotter, 2001). It should be pointed out, because there is often confusion, that these year-class curves, as they are named here, are not the same as catch curves (Ricker, 1975; Jensen, 1985). A year-class curve is fitted to the successive ages of one year class, but a catch curve is fitted to successive year classes observed at their respective ages in one catch, i.e. when no time-series is available. Year-class curves allow estimation of relative annual recruitment, whereas catch curves allow estimation of average recruitment only.

As Z is the sum of fishing and natural mortality (Z = F + M), year-class curves, with Z described by a constant average value, stand in contrast to methods of fish stock assessment that estimate variations in F, and perhaps M as well, over A age classes and Y years. The models used for these methods replace the single parameter Z with up to (A + Y) or even (A * Y) parameters for F, and perhaps more again for M (Deriso et al., 1985; Gavaris, 1988; Megrey, 1989; Shepherd and Nicholson, 1991; Shepherd, 1999; Grønnevik and Evensen, 2001). Cotter et al. (2004) argue that assessment models requiring estimation of large numbers of parameters can signal spurious changes in F or Z because of weaknesses in the data, the assumptions, or the model itself. It is concluded that, if there are major doubts about any of these aspects, a simpler modelling method is advisable and year-class curves have much to offer. Other uses for year-class curves could be as visual screening for cpue data preparatory to some more elaborate modelling method, or when limited computing resources are available.

Year-class curves are easily estimated using ordinary least squares linear regression. Here, we propose the following developments with the aim of improving estimation and prediction without loss of the computational and statistical benefits of the linear least squares method:

  • use of an iterative weighting method for data sets from different sources (Cotter and Buckland, 2004);
  • a set of nested, candidate models for the curves, and a general, although somewhat subjective, protocol for selecting the best;
  • use of polynomials in year to permit Z to vary gradually over time with the minimum number of additional parameters;
  • use of an analysis of relative residual variance (Cotter, 2001) to estimate the precision of different fleets, given the chosen model; and
  • a forward (one-step ahead) validation procedure for testing the predictive abilities of different models and for estimating prediction errors.
Two examples of using the year-class curve method are presented briefly. The first tests precision of estimation, using simulated stocks (National Research Council, 1998) for which true results were available. The second uses cpue data for plaice (Pleuronectes platessa) from Dutch beam trawl surveys carried out by the Institute for Marine Resources and Ecosystem Studies (IMARES) in the Netherlands. We suggest how an analysis of year-class curves might be carried out, and demonstrate one of the strengths of the year-class curve method, namely that alternative models can easily be tested and compared. Selecting the "best approximating model" is an important component of statistical modelling (Burnham and Anderson, 2002), but few existing stock assessment methods allow much flexibility for that purpose. The analyses were carried out with a package written in the R programming language and referred to as YCC (Year-Class Curve). It is freely available from the first author, and is being built into the FLR software suite for fish stock assessments (http://flr-project.org).


    Theory
 Top
 Introduction
 Theory
 Examples
 Discussion
 Appendix
 References
 
Derivation of the year-class curve model
The following derivation aims to be pragmatic and simple, rather than conformable with the advancing mathematics of catch-at-age equations (Xiao, 2006). The usual model of mortality over time t, assuming no net migration to or from the stock, is


Formula 025UM1

where Z, the instantaneous rate of total mortality, is expected to have a negative value. [The absence of a minus sign before Z is unconventional in fisheries work, but it leads to Equation (2) having all terms positive, as is conventional for regression models.] Solving this gives


Formula 025M1

(1)
We now assume that cpue (denoted U) is a constant proportion of N, i.e. U = qN for all ages included in the analysis, and that Z represents a constant, average value over time. Then, taking natural logarithms of Equation (1), restricting attention to one year class, c, substituting age for t, and adding a random error term, e, gives the basic model for a year-class curve:


Formula 025M2

(2)
where U0,c is the cpue index for age zero, a is the age class, i.e. the age in years as an integer index, and age is age in years as a real number. The term e is assumed to be normally distributed around zero with residual variance, {sigma}e2.

Allowing for non-linearity with age
Frequently, the observations appear to be better fitted by a curve over age than by the straight line represented by Equation (2). This could often, but not exclusively, be explained by the inclusion of young age groups that are not fully selected by the trawl. They can either be removed from the analysis, or a term can be added to Equation (2) to allow for the curvature, a convenient possibility being the term log (age + 1). It increases gradually from zero with age in a curve that seems to simulate trawl cpue data well when young fish are less well caught than older fish, whereas a negative coefficient simulates the opposite effect. For lack of a better name, it is here referred to as the selectivity term, but there is no clear relationship with the well-known logistic selectivity function commonly used for trawls. That function is avoided here because its parameters cannot be estimated within a linear least squares framework. The regression coefficient for selectivity is denoted S. Equation (2) thus becomes


Formula 025M3

(3)
Note that estimates of Z' and S will tend to be highly correlated because both relate to age and, for this reason, Z' estimated from fitting Equation (3) is likely to be a biased estimate of total fishing mortality, Z. Z can be estimated numerically as the slope of the fitted log Ua,c for the age at which full selectivity is thought to occur, often the oldest age for trawl fisheries. The estimation is over a small interval of age, i.e. as {Delta} ln U/{Delta}age. The validity of the estimate depends on having an adequate range of ages in the data to find the age of full selectivity, and on the number of fish observed around that age.

Alternative models for different fleets
The term fleet is used here to refer either to a single vessel (such as a research vessel carrying out a survey) or to a fleet of many vessels (such as a commercial fishing fleet). Different fleets are likely to exhibit different cpues for given stock densities because of different fishing powers. This can be allowed for by adding a fleet factor, Vf, to Equation (3):


Formula 025M4

(4)
where f = 1, ... , NF-1 and NF is the number of fleets. Vf brings all fitted log Ua,c,f to the level of one fleet treated as the standard (Cotter, 2001).

Equation (4) implies that Z' and S are common to all fleets. This may be appropriate if all fleets fish the same subset of the stock with gears having similar selectivity functions at age. On the other hand, nesting of one or more of these parameters within fleet may be more appropriate. This is denoted by subscripting with f. Therefore, if each fleet fishes the same subset of the stock, but with a different selectivity function, we would estimate NF parameters Sf, or if, say, catches of some fleets are affected by age-dependent migrations into or out of the survey area, we would estimate NF parameters Zf. The model with all terms nested within fleets, including the log cpue index at age 0, is referred to here as the global model:


Formula 025M5

(5)
It is equivalent to fitting Equation (3) separately to each fleet, except that only one common residual variance is estimated.

Allowing for changing slopes over time
Z may vary over time because of changing fishing practices, rates of natural mortality, or migrations. The Z or Z' slope of a set of year-class curves can be made linearly variable over time by adding an age*year variable to one of the models aforementioned. Curvature can be added with the polynomial terms age*year2 and age*year3. This is more parsimonious with parameters than using Za,y for each year class and age class, and leaves more degrees of freedom for estimating precision. Polynomials also have the advantage that they are easily estimated within a linear least squares framework. Random walks or autocorrelated processes might be used instead, but they would require a more elaborate estimation procedure.

Polynomials can result in huge numbers (e.g. year3) in the predictor matrix that can cause serious rounding errors during least squares matrix inversion. To reduce the problem, the year variables should be transformed to yi = yeari–mean(year), where i indexes rows in the predictor matrix. Next, the polynomials should be orthogonalized to reduce changes in Z' from values previously estimated without polynomial terms. This is achieved with yi' = (yi–mean(yi)), y'i2 = (yi2–mean(yi2)), etc. (This transform may not be important for the odd powers because the mean is then expected to be zero for balanced data.) The transformed variables are orthogonal because, if the same ages are present in every year, {Sigma} yi' (0 + 1 + ··· + A) = 0, as required, and the same for the other powers of yi'. Therefore, instead of Z age in models (2)–(5), we can use Z0 age + Z1 age y' + Z2 age y'2 + Z3 age y'3, or a subset of these terms.

Weighting different fleets
Different series of cpue data are likely to estimate stock parameters and year-class curves with different precision depending on the season and area covered by the fleet, on the precision of age determination and other practical aspects, and on how well the chosen model fits the data. Weighting of different data sets to reflect their precision with respect to the chosen model is therefore desirable. Cotter and Buckland (2004) suggest that the weighting estimated for each fleet's data set should be balanced with the reciprocal of the estimated residual variance specific to that fleet, computed after the model is fitted, i.e. wf {propto} Formulaf–2. They describe how the method can be implemented using iteratively weighted least squares (IWLS), taking into account the degrees of freedom contributed by each fleet to the estimates of each parameter. Usually, two or three iterations produce stable values. Additionally, using the fleet-specific residual variances, the relative precision of the different fleets can be compared with F-tests (Cotter, 2001). Note that biased cpue series will produce biased weights (Quinn and Deriso, 1999, p. 353). Fleets that appear exceptionally precise should be scrutinized to see whether biased sampling may be the cause, e.g. attributable to clustering of observations in restricted times or places (Cotter and Buckland, 2004).

Finding the best model: AIC and forward validation
The possible models for year-class curves have widely different numbers of parameters, and too few or too many in a model can both lead to biased estimation. Burnham and Anderson (2002) make a good case for using the Akaike Information Criterion (AIC) to select the best approximating model or models, rather than a sequence of F- and t-tests for the statistical significance of parameters. A summary of their approach is given below.

  1. Create a list of candidate models based on "thoughtful, science-based, a priori" reasoning. One of these models should be the global model that includes all feasible and important parameters so far as is consistent with the principle of parsimony. Biologically infeasible models should be excluded from the list.
  2. Fit all the candidate models and estimate the AIC using exactly the same set of data for all cases. The small sample AIC, denoted AICc, should be used when the number, n, of independent observations is less than approximately 40 times the number, K, of estimated parameters, Formula:


    Formula 025UM2


  3. Compute the AIC differences for each model relative to the best (minimum AIC), i.e. {Delta}i = AICi–AICmin, and use these to compute the odds against each of the models relative to the best, i.e. 1/exp (–0.5{Delta}i).
  4. Select the best model if it is clearly more likely than any other, or alternatively, select the set of R best supported models and use the methods of multi-model inference (MMI) based on Akaike weights.

When modelling fish stocks, the a priori reasoning referred to concerns the biology and fishery for the species. For example: are there good reasons to expect that selectivity, apparent mortality over time taking into account possible migrations, and apparent recruitments by year class will be similar or different among different fleets? Equation (4), with or without nested parameters and polynomials in year, is intended to provide candidate models for cpue data, with Equation (5) being the global model. However, it may be possible to eliminate some of these models on prior scientific grounds and, if so, this should be done. A clearly best model may not appear at step (4) which, as Burnham and Anderson (2002) point out, is perfectly reasonable because of the complexity of biological systems.

Here, we partially adopt Burnham and Anderson's (2002) approach to model selection and estimation, and we have not attempted MMI. One reservation is that a model that fits all available data well is not necessarily good for predicting next year's stock, the primary task of a stock assessment. A second is that the AIC statistic for year-class curves may be invalidated because of dependence among observations of cpue across different ages within a year and within a fleet. Estimation of random year effects would remove some within-year dependence, e.g. due to weather or cruise-leader effects, but it appears to be a complicated task for the range of nested models suggested here, and may not improve the predictive ability of fitted models (as a random effect is not predictable). For these reasons, we did not attempt it.

The method of forward validation was developed to supplement the AIC method because of these reservations. Starting from an early year and proceeding forward in the time-series, it finds the differences between the predicted log cpue and the observed log cpue for one year after the time domain of the data used to fit the model. The preferred model is that whose mean difference is closest to zero, and for which the mean square of the differences is lowest. This is merely a simulation of a fish stock assessment working group making predictions each year for the coming year, then checking them when the outcome is known, and, for this reason, we are proposing forward-validation in preference to the more widely used cross-validation.

Forward validation starts with selection of a starting year, v, near the beginning of the data series and for which there are sufficient observations to estimate all parameters of the selected model. All available data up to and including year v are fitted, and predictions formed for year v + 1. For linear models, the fitted model may be used directly for prediction, but for polynomial models, extrapolation beyond v may produce erratic results. A Taylor series prediction method based on differential coefficients estimated close to, but behind the v th observation is then preferable. This is described in the Appendix. The predicted log cpues are compared with observed log cpues for each fleet and each age in year v + 1 to form prediction errors, {delta}v + 1 log U. Next, the same model is fitted in the same way to v + 1 observations, and predictions and errors are prepared for year v + 2, and so on, until the penultimate data year is reached and prediction errors are formed with the final data year. Note that successive prediction errors from forward validation are not independent. An outlier at the beginning of a series could cause serial correlation of prediction errors for that fleet for several years afterwards.

A weighted estimator of the mean square prediction error for next year's log cpue for age-class a and fleet f is


Formula 025M6

(6)
where ni is the number of years of observations fitted for the ith forward validation step. This estimator assumes a linear relationship between the number of years of observations and the reliability of the model fitted to them, while also giving more weight to recent observations. The MSPE would often be greater than the residual variance, {sigma}e2. Similarly, a mean prediction bias factor can be estimated from


Formula 025M7

(7)
This estimator is anti-logged, so has value 1 when prediction is perfect.

A complication with forward validation arises when Z is allowed to vary over time. Here, we use low degree polynomials in year for this purpose. It can be expected that, when the index i is small, the polynomial terms will fit short period fluctuations of Z over time, but when i approaches the penultimate data year, the same polynomial terms will fit longer period trends in Z leaving the short period fluctuations to appear as noise. For this reason it is suggested that MSPE for polynomial functions should be estimated using only the last half of the data series, depending on the length of the series available and the amount of short period noise that is thought to exist in the data.

A further choice to be made is whether to estimate prediction errors separately or jointly for each age and/or fleet. Initial screening of candidate models is much simpler if interest is restricted to one measure for each fleet. For this purpose, the average MSPE per age class can be calculated:


Formula 025M8

(8)
where Af is the number of age classes in fleet f for which prediction errors are available. Another approach that might be more suitable for skewed prediction errors would be to examine quantiles of the prediction errors for each fleet pooled across all ages. Later evaluation of the shortlist of preferred models could involve examination of MSPE and MPB for age classes individually, particularly those of most importance for predicting the future size of the stock.


    Examples
 Top
 Introduction
 Theory
 Examples
 Discussion
 Appendix
 References
 
NRC simulated stocks
The Committee on Fish Stock Assessment Methods of the US National Research Council (NRC) published details of catch-at-age simulations of five stocks that were used as part of a review of various assessment methods (National Research Council, 1998). Summaries of the principal varying features of the NRC data sets are given in Table 1. Here, we use the 30-year series of simulated survey abundance indices of numbers-at-ages 2–14, and weight-at-age matrices. Candidate models were fitted with the YCC package to each set of data and compared, then the preferred model was used to estimate numbers-at-age and recruitment annually without knowledge (at the time) of the true values underlying the simulations. As both biomass and recruitment were estimated on a per unit effort basis by YCC while the NRC values were absolute, it was necessary to standardize both the NRC true values and the YCC estimates to units of standard deviations from their respective mean values in order to compare the two time-series on the same scales. (Note that this also standardized the variances of the series.)


View this table:
[in this window]
[in a new window]

 
Table 1. Summary of the principal varying features of five stock simulations used to test the YCC method.

 
Details of the preferred models, their estimated MPB factor, and, very briefly, the reasons for preferring them over other candidate models plus any reservations about them are given in Table 2. Some of the comments relate to voluminous, unreported diagnostic results. True and estimated trajectories of biomass and recruitment obtained with the preferred model for each data set are shown in Figures 1a and 1b, respectively. For biomass, the trajectories of the simulated survey observations are also shown in Figure 1a, because they formed the input to YCC and their variance around the truth could have caused some of the imprecision of the YCC estimates. Inspection of Figure 1a and b indicates that biomass and recruitment were estimated with reasonable precision for most of the five simulations, especially considering that only the survey data were used. Set 3 gave most difficulties, as expected because of the step change in survey catchability after year 15 (Table 1). That was handled by splitting the survey into two fleets. Forward validation was extended back to 25 years (instead of to 10 years for the other sets), so that the MSPE/Age class included results for the first survey. Catchability of the survey over years 16–30 was estimated to be 2.4x that during years 1–15. The true value was 2.0. Figure 1a shows a large overshoot of estimated biomass in year 15 apparently through overestimation of 6–9 year-olds in that year, which in turn was linked to overestimation of recruits in years 6–9 (Figure 1b). Bearing in mind the irregularities in set 3, the analysis with YCC was considered acceptable.


Figure 1
View larger version (33K):
[in this window]
[in a new window]
[Download PowerPoint slide]
 
Figure 1. Comparison between NRC simulations for stocks 1–5 (National Research Council, 1998) and blind estimates made with preferred models of year-class curves. Ordinates are normalized to standard deviations from the mean for each series. For set 3, data were analysed as two surveys changing at year 16. (a) Total biomass estimated for ages 2–14; (b) recruits aged 1.

 


View this table:
[in this window]
[in a new window]

 
Table 2. Brief details of YCC models preferred for fitting the five stock simulations prepared by NRC (1998).

 
Plaice in the North Sea
This example uses survey cpue (tuning) data for North Sea plaice taken from an ICES stock assessment report (ICES, 2005). It illustrates how an analysis of year-class curves might be applied in a stock assessment, suggests uses for some of the outputs offered by YCC, and gives results consistent with known migrations of plaice. Data were for three beam trawl surveys carried out by IMARES for the stock assessments of plaice and sole. The short names for the surveys and details of the data they provided are:
  • BTS Isis, for fish aged 1–9 years, for 1985–2003;
  • BTS Tridens for fish aged 2–9 years, for 1996–2003; and
  • SNS for fish aged 1–3 years, for 1982–2002.

BTS Isis covered the southeastern North Sea with twin 8-m beam trawls with eight tickler chains. BTS Tridens was an expansion into the central and northern part of the North Sea. The same gear was used, but with the addition of a flip-up rope to cope with rough ground. SNS covered transects within coastal waters of the southeastern North Sea using twin 6-m beam trawls and aiming at younger age groups than the two BTS surveys. Figure 2 shows the coverage of the three surveys. More information about them is available in two project reports (Beare et al., 2002; Piet, 2004). In all, 291 observations of cpue at all ages were analysed. Forward validation was started 10 years behind the final data year using the Taylor series method. Note that only six or seven forward predictions were available for BTS Tridens, and only two age classes for SNS, suggesting that MSPE for those surveys was estimated poorly relative to those for BTS Isis.


Figure 2
View larger version (51K):
[in this window]
[in a new window]
[Download PowerPoint slide]
 
Figure 2. Coverage of three beam trawl surveys carried out by IMARES (the Netherlands) as used for the YCC analysis of cpue abundance indices-at-age for plaice.

 
Table 3 contains the shortlist of six candidate models considered, notes on the assumed biological meaning of each model, and in each case, comments on its scientific implications, plus a three-level subjective ranking of its overall credibility prior to any statistical analysis. Models including polynomials in year were omitted at this stage because the biological credibility of the various terms is little affected by whether or not Z is allowed to vary gradually over time.


View this table:
[in this window]
[in a new window]

 
Table 3. Plaice in the North Sea. Preliminary, subjective assessment of candidate models for fitting year-class curves to cpue data from three beam trawl surveys carried out by IMARES.

 
After fitting the different models in Table 3, models (4) and (10) were preferred. Table 4 shows selected diagnostic results, plus reasons for preferring them over other candidate models, together with any reservations. Models (4) and (10) showed the lowest AICc, good MSPE/age class, and similar levels of MPB factor/age class. These models were therefore re-fitted with the addition of polynomial terms in year to see whether better fits or predictive abilities could be found. For model (4), AICc was improved by polynomial terms, but MSPE/age class was either unchanged or became worse, and the unsteady behaviour of residuals on age observed with the linear models was not improved. Giving priority to MSPE/age class and the principle of parsimony, the preferred model was the linear model (4). For model (10), first and second degree polynomials improved neither AICc nor predictive abilities. The third degree polynomial did improve AICc but not MSPE/age class, nor the unsteady behaviour of residuals. For the same reasons as for model (4), the preferred model was the linear model, (10).


View this table:
[in this window]
[in a new window]

 
Table 4. Plaice in the North Sea. Selected results from fitting the preferred pair of models, (4) and (10), see Table 3, to cpue data from three beam trawl surveys (BTS Isis, BTS Tridens, SNS) carried out by IMARES (The Netherlands). Some comments refer to results not presented.

 
Model (10) might be preferred to model (4) at this stage because of its greater prior credibility (Table 3), its lower AICc, and its comparable MSPE/age class for each of the three fleets (Table 4). Taking model (10) as AICmin, {Delta}4 = 555.699–554.848, and the odds against model (4) are 1/exp (–0.426) = 1.53 to 1, making model (10) only slightly preferable to model (4).

The values of Z estimated numerically from models (4) and (10) at the maximum observed ages for each survey are compared in Table 5. BTS Tridens showed the shallowest slopes (least negative values). This is consistent with BTS Tridens covering the largest area of sea consisting of mainly offshore stations, because estimated Z is, in that case, likely to be less affected by age-dependent migrations of plaice out of the survey area. SNS showed the steepest Z (most negative values), even though the maximum observed age was only 3. Model (4) estimated Z = –1.84 year–1, substantially different from the estimates for the other fleets. This would imply very high losses of older fish from the SNS transects to offshore waters. Model (10) estimated Z for SNS as –1.34 year–1, also quite distinct from the values for the other two fleets. As SNS was designed as a survey of young fish, detection of offshore migrations is not surprising and lends support for use of the YCC method. Beverton and Holt (1957, pp. 182–183) found average Z for North Sea plaice aged 5–10 between 1929 and 1938 to be –0.83 year–1, a low value that they attributed to immigration to the fishing area off Lowestoft, UK, and to discarding. Concerning estimates of log cpue indices for plaice aged 0, models (4) and (10) gave series with different elevations but identical patterns (not illustrated), implying that the choice of model was unimportant for estimation of relative recruitment to year classes in this case.


View this table:
[in this window]
[in a new window]

 
Table 5. Plaice in the North Sea. Coefficients of total mortality, Z, as estimated numerically at the ages shown for models (4) and (10) fitted to cpue data from three beam trawl surveys (BTS Isis, BTS Tridens, SNS) carried out by IMARES (The Netherlands).

 
Figure 3 shows grey-scale representations of the correlations of log residuals across ages obtained with model (10) for each survey. The patterns for models (4) and (1), not illustrated, were nearly identical, implying that the observed correlation structure was a property of the data rather than a consequence of the model. Patterns were different for each survey. Correlations were positive and high for SNS, whereas BTS Isis showed more independence, i.e. more medium grey in Figure 3, and BTS Tridens was intermediate. High correlations imply that parameters estimated by a survey should have larger standard errors than those obtained by least squares regression, based on the assumption of independent residuals. Table 6, here called the D-table, shows the relative quantities of information contributed by each fleet to the estimate of each parameter in model (10), based on the number of observations and the fleet weightings, i.e. D in Equation (4) of Cotter and Buckland (2004). Also shown in Table 6 are the parameters and standard errors on the log scale estimated after three iterative re-weightings, the absolute fleet weights (i.e. not constrained to add to 1, as they are in Table 4), and the numbers of observations. Table 6 may be used to indicate how much confidence to place on the standard errors of estimates. As examples, the 1976–1978 and 2002 year classes and a selectivity factor were all estimated solely by BTS Isis and, because residuals appeared to be reasonably independent for this survey (Figure 3), the standard errors could be considered reasonable. The coefficient of age (i.e. log Z' in Table 6) was estimated with substantial contributions of information from BTS Tridens and SNS, both of which showed correlated residuals (Figure 3); the standard error of ±0.090 is therefore likely to be underestimated. The standard errors for the fleet and selectivity factors estimated solely by and for SNS are likely to be substantially underestimated because of the correlated residuals for that fleet (Figure 3). Unfortunately, a link between the degree of underestimation of the standard error and the degree of correlation among the residuals-at-age is not available owing to the lack of statistical theory. The D values for model (4), not shown, did not differ substantially from those for model (10) in Table 6, except for log Z', for which each fleet contributed a separate estimate.


Figure 3
View larger version (17K):
[in this window]
[in a new window]
[Download PowerPoint slide]
 
Figure 3. Plaice in the North Sea. Grey scale representations of correlations of log residual errors across ages with model (10) for three beam trawl surveys carried out by IMARES (the Netherlands). (a) BTS Isis; (b) BTS Tridens; (c) SNS.

 


View this table:
[in this window]
[in a new window]

 
Table 6. Plaice in the North Sea. D-table for model (10), i.e. relative contributions of information on a scale of 0 to 1 by each of three IMARES survey fleets to estimates of each parameter based on the numbers of observations (n) and the estimated fleet weightings (w). Three iterative reweightings. First column: ‘1976 Yclass’ = log (U0,1976) etc. Right 2 columns, log-transformed parameter estimates and corresponding standard errors (s.e.).

 
Figure 4 is a graphic presentation of the correlations between estimated parameters for model (4). That for model (10) was similar. In both cases, estimated year-class parameters were highly correlated positively. This is consistent with the two models estimating similar relative recruitments, but with different elevations, because if one estimate is high, all are high when positively correlated. With model (4) (Figure 4), all year-class estimates were positively correlated with the coefficient of age estimated for BTS Isis. This is consistent with BTS Isis contributing the longest data series and the most age groups to the analysis. With model (10) (not illustrated), year-class estimates were most positively correlated with the coefficient of age that was common to all fleets. The year-class estimates in Figure 4 were negatively correlated with the coefficient of selectivity estimated for BTS Isis because a large value depresses the estimated numbers of young fish in each year class. It may be helpful to note that the correlations among parameters, {theta}, are computed from cov ({theta}) = Formulae2 (X'X)–1, where X is the predictor matrix of the model. The patterns in Figure 4 therefore depend on the choice of parameters in the model, and of years, ages, and fleets over which cpue were observed, and not on the observed values themselves. Other strong relationships, both positive and negative, among estimated parameters are also evident in Figure 4, and although they may be worthy of further investigations, will not be commented upon here.


Figure 4
View larger version (41K):
[in this window]
[in a new window]
[Download PowerPoint slide]
 
Figure 4. Plaice in the North Sea. Grey scale representations of correlations of parameters estimated with model (4) from data for three beam trawl surveys carried out by IMARES (the Netherlands). Refer to Table 2 for parameters corresponding to model terms. Year classes are labelled by birth year.

 
Year-class curves obtained by fitting model (10) are shown in Figure 5a–c. Those for BTS Isis are concave upwards, whereas those for BTS Tridens show the opposite effect, reflecting the different signs of S (Table 6). The negative value for BTS Isis implies that old plaice were less well caught than younger fish, consistent with migrations by older fish out of the survey area. The SNS curves showed a tendency to overestimate older fish in the more recent years but, considering the shortness of the age series, this is not surprising. Plots of residuals over year, year-class, and age were also examined for all surveys (not shown). Patterns over year were evident, but were not cured by adding polynomial terms. BTS Isis and, to a lesser extent, BTS Tridens showed waves in residuals over age. The presence of patterns was not ideal, but the situation was not generally improved by other candidate models. The patterns were accepted because of other favourable properties of model (10), and a reluctance to further complicate the set of candidate models on an ad hoc basis.


Figure 5
Figure 5
View larger version (40K):
[in this window]
[in a new window]
[Download PowerPoint slide]
 
Figure 5. Plaice in the North Sea. Year-class curves fitted with the YCC package to log cpue obtained from three beam trawl surveys using model (10); see Tables 3 and 4. (a) BTS Isis, 1985–2003; (b) BTS Tridens, 1996–2003; (c) SNS, 1982–2002. Panel labels show year classes. Surveys carried out by IMARES (the Netherlands).

 
An analysis of relative residual variance (Cotter, 2001) was also carried out to see whether the residual variance exhibited by any one fleet was significantly higher than the residual variance of all the others, as judged by an F-test. BTS Tridens showed the highest fleet weightings with either model (4) or (10) (Table 4), and correspondingly the lowest fleet-specific residual variances (Table 7). The comparison of residual variances is somewhat different from that foreseen when this method was presented (Cotter, 2001), because the interest now is in whether one survey is significantly better, not worse, than two others. Here, we tested ratios of fleet-specific residual variances as F for Isis/Tridens and SNS/Tridens. The data and calculations for model (10) are shown in Table 7. The residual variances of BTS Isis and SNS were significantly higher than for BTS Tridens (p < 0.05), suggesting that BTS Tridens was achieving better precision. This result is expected given that BTS Tridens covers a much larger survey area (Figure 2) but, bearing in mind dependence among some of the residuals and the short time-series for BTS Tridens, the result should be treated with caution. Results for model (4) (data not shown) were similar.


View this table:
[in this window]
[in a new window]

 
Table 7. Plaice in the North Sea. F-tests of which of three IMARES beam trawl survey was most precise with one model.

 
To summarize, models (4) and (10) were selected on prior biological grounds, competitive predictive abilities, and low AICc. Exclusion of polynomials in year was preferred for both models because the extra terms did not notably improve prediction variance, mean bias, or homogeneity of residuals. This exclusion implies that Z did not change detectably over the survey period, from 1982 to 2003. A final choice between the two models was not possible, but it may not be particularly important for an assessment of the stock because both gave comparable estimates of relative recruitments, and the trend in –Z was clearly SNS < BTS Isis < BTS Tridens with both models (Table 5). This is consistent with the increasing geographic spreads of the three surveys and the loss of plaice offshore with age. Heincke (1905, p. 17) notes that for the North Sea: "the occurrence of the young plaice, from the coast to the open sea, is arranged like the steps of a ladder ... the smallest and youngest quite close to land, the largest and oldest the furthest out". The reader is also referred to Metcalfe et al. (2006) and Rijnsdorp et al. (2006).


    Discussion
 Top
 Introduction
 Theory
 Examples
 Discussion
 Appendix
 References
 
Year-class curves have several desirable features not all of which are available with other, more elaborate methods of stock assessment.

  • They are based on statistical theory for least squares and AIC.
  • They are equally applicable to cpue results from surveys or commercial fishing.
  • They use continuous rather than categorical variables for all variables except year-class and fleet, so reducing the numbers of parameters to be estimated and allowing best precision with the available observations.
  • They generate illustrated output that can be discussed with fishers and other, non-mathematical stakeholders, e.g. when selecting the best approximating model.
  • They allow fleets to be ranked for precision using their fitted weights and relative residual variances. Results are model-dependent but may, nonetheless, be useful for assigning priorities to different sources of data.
  • They are fitted quickly.
  • They are compatible with a range of candidate models from which one may be selected so as best to capture the biological and fishery-related factors affecting cpue.

In addition, the inclusion of polynomials in year allows year-class curves to vary gradually over time, a feature that was not suggested when they were previously described as a method for intercalibrating groundfish surveys (Cotter, 2001). Allowing only a conservative estimate of variability in Z gives fewer opportunities for erroneously interpreting sampling variance or random year effects as trends in time. The possibility of model-related bias arising from including or excluding polynomials can be checked by seeing whether predictive abilities are improved as seen by forward validation. A weakness of year-class curves is that they are totally dependent on catchabilities remaining constant over time. However, few, if any, cpue-based stock assessment methods are immune to this problem.

The large range of candidate models available for year-class curves opens up the problem of how to select the "best approximating" model. Some other methods for estimating stock parameters offer a more limited range of candidate models. For example, the XSA method offers only two, concerning whether or not abundances of the youngest age classes are proportional to eventual year-class strength (Darby and Flatman, 1994; Shepherd, 1999). The small number of options simplifies the analyst's task, but may mean that model-related bias is overlooked. The several candidate models together with the visual means of examining them offered by year-class curves encourage thorough consideration of the biological relevance of each model to the species and locations of interest, as in Table 3. This would surely represent valuable use of a stock working group's time, particularly if commercial fishers were invited to take part in the discussions and perhaps help to select the best model. Smith and Punt (1998) report involving fishers in a (Bayesian) assessment of gemfish (Rexea solandri) in eastern Australia, leading to better acceptance by the fishing industry of the science and management process.

Selecting the best model is, however, not easy. The results of the analysis of plaice in the North Sea showed that the various indicators of fit may favour different models. For example, the AICc implied that polynomials in year would be beneficial in some of the models even though the MSPE was not improved, or only slightly, and none of the models listed in Table 3 provided neat horizontal bands of residuals over year, year-class, and age, as is usually considered important for a good fit. Priority was given to the MSPE criterion here because of the importance of prediction for stock assessments, and because forward validation is based on repeated use of the model with successively lengthening data series, as it would be applied in reality from year to year. The AICc, on the other hand, relates only to the model fitted once to the complete set of data. It is also compromised by dependence of data across ages within years and fleets, as well as by possible inadequacy of the constant variance, normal distribution model of residual errors. The conflict among indicators of fit underlines the importance of biological considerations and the principle of parsimony of parameters when choosing the model with which to go forward. Elaboration of the candidate models with extra terms or different functions designed to deal with patterned residuals is always an option, but it could amount to little more than an ad hoc exercise with no lasting explanatory power unless there are clear prior justifications based on "thoughtful" science.

Year-class curves do not allow absolute stock numbers or fishing mortality, F, to be estimated but, as has been pointed out (Rivard, 1989; Cotter et al., 2004), other methods can only do this if the coefficient of natural mortality, M, is accurately known, which is not usually the case (Vetter, 1988; Hewitt and Hoenig, 2005). Statistical inferences from year-class curves are compromised somewhat by dependence among observed data, but so too are other assessment methods unless a covariance matrix of errors is an estimated component of the method. The practical relevance of dependence is that precision is generally overestimated. Analysts are consequently encouraged to add too many terms to models, and to have false confidence in the estimated parameters. This lends further support to the use of parsimonious models.

Year-class curve models could be augmented by inclusion of commercial fishing effort, E, for a stock (Beverton and Holt, 1957, section 14.3). This could be achieved by disaggregating the Z age term in the model equations:


Formula 025UM3

where F and M are the coefficients of fishing and natural mortality, respectively, and q' is a catchability coefficient. E is total commercial fishing effort, or an index of it. Adding an annual estimate (multiplied by age) to the predictor matrix for each cpue-at-age could remove the need to use polynomials in year to allow Z to vary over time. Alternatively, if polynomials were still found necessary in the model, they would indicate changes in M over time. However, if effort does not vary much from year to year (as in the North Sea, for example; Jennings et al., 1999), the two factors E age and age would be nearly collinear, meaning that estimates of q' and M would be highly correlated. Another problem is that a reliable index of total fishing effort is often not readily obtainable, particularly when different fishing gears are in use commercially and their effort measures cannot be readily added. For these reasons, we did not attempt to include effort in our analyses.


    Appendix
 Top
 Introduction
 Theory
 Examples
 Discussion
 Appendix
 References
 
Prediction using a Taylor series
Taylor's theorem states that a function of x may be predicted at a point (x + h) with


Formula 025UM4

The primes indicate successive differential coefficients. Denote the date of the final observation of any forward validation step as Y0, and let a small fraction of a year be {delta}. Let Y1 = Y0{delta}, Y2 = Y0–2{delta}, etc. Backward differences are used to keep estimation within the fitted domain, so that the differential coefficients are estimated as close as possible to the value to be predicted for forward validation. Then the coefficients are approximated by


Formula 025UM5

{delta} is set to a value, say 0.01, which is not too large for estimating the differentials accurately, and not so small as to cause significant rounding errors in the calculations. Finally, the Taylor series is evaluated with x = Y0 and h = 1 year.


    Acknowledgements
 
This paper and the YCC software were developed as a contribution to the FISBOAT project, number 502572, http://www.ifremer.fr/drvecohal/fisboat/, funded by the European Commission. We are grateful to all referees for beneficial comments. The NRC simulation data were kindly provided by Terry Quinn. AJRC is grateful to Laurie Kell for encouraging development of YCC using R. No statement in this paper should be interpreted as official policy of the EC or of the authors' employers.


    References
 Top
 Introduction
 Theory
 Examples
 Discussion
 Appendix
 References
 

    Beare D., Castro J., Cotter J., van Keeken O., Kell L., Laurec A., Mahé J-C., et al. (2002) Evaluation of research surveys in relation to management advice (EVARES). Final report. DGXIV Fisheries, European Commission, Brussels. FISH/2001/02 – Lot 1. Available from john.cotter{at}cefas.co.uk.

    Beverton R. J. H. and Holt S. J. (1957) On the dynamics of exploited fish populations. Fishery Investigations 19:533 Series II.

    Burnham K. P. and Anderson D. R. (2002) Model Selection and Multimodel Inference. (Springer, New York)488.

    Cotter A. J. R. (2001) Intercalibration of North Sea International Bottom Trawl Surveys by fitting year-class curves. ICES Journal of Marine Science 58:622–632 [Erratum, Ibid, 58: 1340].[Abstract/Free Full Text]

    Cotter A. J. R. and Buckland S. T. (2004) Using the EM algorithm to weight data sets of unknown precision when modeling fish stocks. Mathematical Biosciences 190:1–7.[CrossRef][Web of Science][Medline]

    Cotter A. J. R., Burt L., Paxton C. G. M., Fernandez C., Buckland S. T., Pan J-X. (2004) Are stock assessment methods too complicated? Fish and Fisheries. 5:235–254.

    Darby C. D. and Flatman S. (1994) Virtual Population Analysis: version 3.1 (windows/DOS) user guide. (CEFAS, Lowestoft, UK).

    Deriso R. B., Quinn T. J., Neal P. R. (1985) Catch-age analysis with auxiliary information. Canadian Journal of Fisheries and Aquatic Sciences 42:815–824.

    Gavaris S. (1988) An adaptive framework for the estimation of population size. Canadian Atlantic Fisheries Scientific Advisory Committee 88/29:12.

    Grønnevik R. and Evensen G. (2001) Application of ensemble-based techniques in fish stock assessment. Sarsia 86:517–526.

    Heincke F. R. (1905) The occurrence and distribution of the eggs, larvae and various age-groups of the food-fishes in the North Sea. 39 Rapports et Procès-Verbaux des Réunions du Conseil Permanent International pour l'Exploration de la Mer, 3, Appendix E.

    Hewitt D. A. and Hoenig J. M. (2005) Comparison of two approaches for estimating natural mortality based on longevity. Fishery Bulletin US 103:433–437.

    ICES. (2005) Report on the assessment of demersal stocks in the North Sea and Skagerrak. ICES Document CM 2005/ACFM: 07 772 http://www.ices.dk/iceswork/wgdetailacfm.asp?wg=WGNSSK.

    Jennings S., Alsvaag J., Cotter A. J. R., Ehrich S., Greenstreet S. P. R., Jarre-Teichmann A., Mergardt N., Rijnsdorp A. D., Smedstad O. (1999) Fishing effects in northeast Atlantic shelf seas: patterns in fishing effort, diversity and community structure. 3. International trawling effort in the North Sea: an analysis of spatial and temporal trends. Fisheries Research 40:125–134.[CrossRef][Web of Science]

    Jensen A. J. C. (1939) On the laws of decrease in fish stocks. Rapports et Procès-Verbaux des Réunions du Conseil Permanent International pour l'Exploration de la Mer 110:85–96.

    Jensen A. L. (1985) Comparison of catch-curve methods for estimation of mortality. Transactions of the American Fisheries Society 114:743–747.[CrossRef]

    Megrey B. A. (1989) Review and comparison of age-structured stock assessment models from theoretical and applied points of view. American Fisheries Society Symposium 6:8–48.

    Metcalfe J. D., Hunter E., Buckley A. A. (2006) The migratory behaviour of North Sea plaice: currents, clocks and clues. Marine and Freshwater Behaviour and Physiology 39:25–36.[CrossRef]

    National Research Council. (1998) Improving fish stock assessments. National Academy of Sciences, Washington, DC pp. 178 http://books.nap.edu/catalog/5951.html.

    Piet G. (2004) Development of a central database for European trawl survey data (DATRAS). Final report. DGXIV(European Commission, Brussels) EC project QLRT-2001-00025. http://www.ices.dk/datacentre/datras/Final%20report%20to%20EU.pdf.

    Quinn T. J. and Deriso R. B. (1999) Quantitative Fish Dynamics. (Oxford University Press, Oxford, UK) pp. 542.

    Ricker W. E. (1975) Computation and interpretation of biological statistics of fish populations. Bulletin of the Fisheries Research Board of Canada 191:382.

    Rijnsdorp A. D., Daan N., Dekker W. (2006) Partial fishing mortality per fishing trip: a useful indicator of effective fishing effort in mixed demersal fisheries. ICES Journal of Marine Science 63:556–566.

    Rivard D. (1989) Overview of the systematic, structural, and sampling errors in cohort analysis. American Fisheries Society Symposium 6:49–65.

    Sette O. E. (1943) Biology of the Atlantic mackerel (Scomber scombrus) of North America. 1. Early life history, including the growth, drift, and mortality of the egg and larval populations. Fishery Bulletin of the US Fish and Wildlife Service 50:149–237.

    Shepherd J. G. (1999) Extended survivors analysis: an improved method for the analysis of catch-at-age data and abundance indices. ICES Journal of Marine Science 56:584–591.[Abstract/Free Full Text]

    Shepherd J. G. and Nicholson M. D. (1991) Multiplicative modelling of catch-at-age data, and its application to catch forecasts. Journal du Conseil International pour l'Exploration de la Mer 47:284–294.

    Silliman R. P. (1943) Studies on the Pacific pilchard or sardine (Sardinops caerulea). 5. A method of computing mortalities and replacements. United States Department of the Interior, Fish and Wildlife Service 24:10.

    Smith A. D. M. and Punt A. E. (1998) Stock assessment of gemfish (Rexea solandri) in eastern Australia using maximum likelihood and Bayesian methods. In. In Quinn T. J., Funk F., Heifetz J., Ianelli J. N., Powers J. E., Schweigert J. F., Sullivan P. J., Zhang C-I. (Eds.). Fisheries Stock Assessment Models(Alaska Sea Grant College Program, AK-SG-98-01. University of Alaska, Fairbanks) pp. 245–286.

    Vetter E. F. (1988) Estimation of natural mortality in fish stocks: a review. Fishery Bulletin US 86:25–43.

    Xiao Y. (2006) Catch equations: calculating the instantaneous rate of fishing mortality from catch and back. Ecological Modelling 193:225–252.[CrossRef]


Add to CiteULike CiteULike   Add to Connotea Connotea   Add to Del.icio.us Del.icio.us    What's this?


This article has been cited by other articles:


Home page
ICES J. Mar. Sci.Home page
A. Silva, D. W. Skagen, A. Uriarte, J. Masse, M. B. Santos, V. Marques, P. Carrera, P. Beillois, G. Pestana, C. Porteiro, et al.
Geographic variability of sardine dynamics in the Iberian Biscay region
ICES J. Mar. Sci., April 1, 2009; 66(3): 495 - 508.
[Abstract] [Full Text] [PDF]


This Article
Right arrow Abstract Freely available
Right arrow FREE Full Text (PDF) Freely available
Right arrow All Versions of this Article:
64/2/234    most recent
fsl025v1
Right arrow Alert me when this article is cited
Right arrow Alert me if a correction is posted
Services
Right arrow Email this article to a friend
Right arrow Similar articles in this journal
Right arrow Alert me to new issues of the journal
Right arrow Add to My Personal Archive
Right arrow Download to citation manager
Right arrowRequest Permissions
Google Scholar
Right arrow Articles by Cotter, A. J. R.
Right arrow Articles by Piet, G. J.
Right arrow Search for Related Content
PubMed
Right arrow Articles by Cotter, A. J. R.
Right arrow Articles by Piet, G. J.
Social Bookmarking
 Add to CiteULike   Add to Connotea   Add to Del.icio.us  
What's this?