Geostatistical delta-generalized linear mixed models improve precision for estimated abundance indices for West Coast groundfishes

ICES Journal of Marine Science ICES Journal of Marine Science (2015), 72(5), 1297– 1310. doi:10.1093/icesjms/fsu243 Original Article Geostatistical delta-generalized linear mixed models improve precision for estimated abundance indices for West Coast groundfishes 1 Fisheries Resource Assessment and Monitoring Division (FRAM), Northwest Fisheries Science Center, National Marine Fisheries Service (NMFS), NOAA, 2725 Montlake Boulevard E, Seattle, WA 98112, USA 2 Conservation Biology Division, Northwest Fisheries Science Center, National Marine Fisheries Service (NMFS), NOAA, 2725 Montlake Boulevard E, Seattle, WA 98112, USA 3 Department of Mathematics, University of Bergen, PO Box 7800 5020 Bergen, Norway *Corresponding author: tel: +1 206 302 1772; fax: +1 206 860 6792; e-mail: Thorson, J. T., Shelton, A. O., Ward, E. J., and Skaug, H. J. Geostatistical delta-generalized linear mixed models improve precision for estimated abundance indices for West Coast groundfishes. – ICES Journal of Marine Science, 72: 1297 – 1310. Received 25 August 2014; revised 18 November 2014; accepted 5 December 2014; advance access publication 14 January 2015. Indices of abundance are the bedrock for stock assessments or empirical management procedures used to manage fishery catches for fish populations worldwide, and are generally obtained by processing catch-rate data. Recent research suggests that geostatistical models can explain a substantial portion of variability in catch rates via the location of samples (i.e. whether located in high- or low-density habitats), and thus use available catch-rate data more efficiently than conventional “design-based” or stratified estimators. However, the generality of this conclusion is currently unknown because geostatistical models are computationally challenging to simulation-test and have not previously been evaluated using multiple species. We develop a new maximum likelihood estimator for geostatistical index standardization, which uses recent improvements in estimation for Gaussian random fields. We apply the model to data for 28 groundfish species off the U.S. West Coast and compare results to a previous “stratified” index standardization model, which accounts for spatial variation using post-stratification of available data. This demonstrates that the stratified model generates a relative index with 60% larger estimation intervals than the geostatistical model. We also apply both models to simulated data and demonstrate (i) that the geostatistical model has well-calibrated confidence intervals (they include the true value at approximately the nominal rate), (ii) that neither model on average under- or overestimates changes in abundance, and (iii) that the geostatistical model has on average 20% lower estimation errors than a stratified model. We therefore conclude that the geostatistical model uses survey data more efficiently than the stratified model, and therefore provides a more cost-efficient treatment for historical and ongoing fish sampling data. Keywords: abundance index, delta-generalized linear mixed model, fishery-independent data, Gaussian random field, geostatistics, index standardization, management procedure, spatial statistics, stock assessment, template model builder. Introduction Fisheries management throughout the United States, Europe, and elsewhere is often informed by estimates of fish population abundance derived from sampling operations with predetermined sampling designs. Data derived from fishery-independent surveys are generally processed to generate an index that is intended to be proportional to population abundance (Maunder and Punt, 2004). Abundance indices are generally considered to be the most important source of information regarding fishing impacts on marine populations (Francis, 2011), and are used in data-rich and data-limited stock assessments to inform changes in management decisions, e.g. allowable catches (Methot et al., 2014). Fisheries management agencies worldwide use many different methods to estimate abundance indices and there can be wide variation in methodology even within a single fisheries management agency. Using NOAA Fisheries in the United States as an example, Published by Oxford University Press on behalf of International Council for the Exploration of the Sea 2015. This work is written by (a) US Government employee(s) and is in the public domain in the US. James T. Thorson 1*, Andrew O. Shelton 2, Eric J. Ward2, and Hans J. Skaug 3 1298 Methods Model development The goal of applying a geostatistical model to data from fisheries research surveys is to explain the catches of each species recorded in the survey, and hence to infer population density throughout the domain of the survey design. We use a delta-generalized linear mixed modelling framework (Lo et al., 1992; Stefansson, 1996; Martin et al., 2005), which separately models the probability of having non-zero catches (“encounters”) and catch rates for each encounter (“positive catch rates”). The probability that a sample encounters the target species (i.e. that catch C . 0) is approximated via a first model component: Pr[C . 0] = p, (1) where p is the probability of encounter (see Tables 1 and 2 for list of notation used). Positive catches are then approximated via a second model component: Pr[C = c|C . 0] = Gamma(c, s−2 , ls2 ), (2) where Gamma(c, x, y) is the value of a probability density function evaluated at value c given a gamma distribution with shape x and scale y, l is the expected catch given that the species is encountered, and s is the coefficient of variation of measurement errors for positive catch rates. We use the Gamma distribution here because of its flexibility, although we note that many other distributions can be explored for positive catch rates (Thorson and Ward, 2013). Each of these two components are estimated here using Gaussian Markov random fields. A random field defines the probability of a given function (e.g. densities as a function of latitude and longitude), and is analogous to a conventional random variable, which defines the probability of a given variable (Rasmussen and Williams, 2006). Specifically, a random field defines the expected value, variance and covariance of a multivariate realization from a stochastic process. In this case, the stochastic process represents the aggregate impact of environmental and biological factors that are not directly observed but still contribute to the distribution and density of the target species (Shelton et al., 2014; Thorson et al., 2015). Because the delta-GLMM framework involves two model components (probability of encounters and positive catch rates), we estimate unique random fields for each. For a Gaussian random field E, the value of the random field at a single fixed location s ¼ kx, yl (where x and y are the easting and northing for that location) follows a normal distribution, and the value of the random field at several (but a finite number of) location (...truncated)