State-space models have been promoted as the next-generation of fisheries stock assessment and evaluation of their reliability is needed. We simulated operating models that varied fishing pressure, magnitude of observation error, and sources of process error. For each operating model, we fit a range of estimating models with correct and incorrect configurations. We measured reliability of estimating models by convergence rate, accuracy of Akaike information criterion (AIC)-based model selection, estimation bias, and magnitude of retrospective patterns. All reliability measures were generally better with lower observation error, contrast in fishing pressure over time, and when median natural mortality rate is known. The magnitude of the log-likelihood gradients was not a reliable indicator of convergence. AIC can generally distinguish process error source with lower observation error and higher true process error variability. Distinguishing the stock–recruit relationship with AIC required large contrast in spawning biomass and low recruitment variation, but bias in stock–recruit parameter estimation was prevalent. Retrospective patterns were not large for mis-specified models. These findings improve our understanding of when results from state-space models will be reliable.
The study of animal diets and the proportional contribution that different foods make to their diets is an important task in ecology. Stable Isotope Mixing Models (SIMMs) are an important tool for studying an animal's diet and understanding how the animal interacts with its environment. We present cosimmr, a new R package designed to include covariates when estimating diet proportions in SIMMs, with simple functions to produce plots and summary statistics. The inclusion of covariates allows for users to perform a more in-depth analysis of their system and to gain new insights into the diets of the organisms being studied. A common problem with the previous generation of SIMMs is that they are very slow to produce a posterior distribution of dietary estimates, especially for more complex model structures, such as when covariates are included. The widely-used Markov chain Monte Carlo (MCMC) algorithm used by many traditional SIMMs often requires a very large number of iterations to reach convergence. In contrast, cosimmr uses Fixed Form Variational Bayes (FFVB), which we demonstrate gives up to an order of magnitude speed improvement with no discernible loss of accuracy. We provide a full mathematical description of the model, which includes corrections for trophic discrimination and concentration dependence, and evaluate its performance against the state of the art MixSIAR model. Whilst MCMC is guaranteed to converge to the posterior distribution in the long term, FFVB converges to an approximation of the posterior distribution, which may lead to sub-optimal performance. However we show that the package produces equivalent results in a fraction of the time for all the examples on which we test. The package is designed to be user-friendly and is based on the existing simmr framework.
Tiger Grouper ( Mycteroperca tigris ) form fish spawning aggregations (FSAs) around the winter full moons (typically January through April) in the Caribbean. Males defend territories to attract mates in a lek-like reproductive strategy. Prior studies have documented rapid declines in populations with FSA-associated fisheries. This study examines the migratory behavior of adult male Tiger Grouper in Little Cayman, Cayman Islands, to better understand the impacts of aggregation fishing. As part of the Grouper Moon Project, we acoustically tagged ten spawning male Tiger Grouper at the western end of Little Cayman in February 2015. Using a hydrophone array surrounding the island, we tracked the movements of the tagged fish for 13 months. We observed 3 migratory strategies: resident fish ( n = 2) that live at the FSA site, neighboring fish ( n = 5) that live within 4 km of the site, and commuter fish ( n = 3) that travel over 4 km for spawning. Fish began aggregating 2 days before the full moon and left 10–12 days after the full moon, from January to May. Regardless of migratory strategy, all tagged fish that aggregated after February 2015 returned to the west end FSA. However, in January 2016, one fish appeared to attend a different FSA closer to its presumed home territory. Tiger Grouper may establish multiple FSAs around Little Cayman, and males appear to attend FSAs near their home territories. Protracted spawning seasons, FSA site infidelity, and putative FSA catchments should all be considered to ensure sustainable fisheries management for this important species.
Dispersal of eggs and larvae from spawning sites is critical to the population dynamics and conservation of marine fishes. For overfished species like critically endangered Nassau grouper ( Epinephelus striatus ), recovery depends on the fate of eggs spawned at the few remaining aggregation sites. Biophysical models can predict larval dispersal, yet these rely on assumed values of key parameters, such as diffusion and mortality rates, which have historically been difficult or impossible to estimate. We used in situ imaging to record three-dimensional positions of individual eggs and larvae in proximity to oceanographic drifters released into egg plumes from the largest known Nassau grouper spawning aggregation. We then estimated a diffusion–mortality model and applied it to previous years' drifter tracks to evaluate the possibility of retention versus export to nearby sites within 5 days of spawning. Results indicate that larvae were retained locally in 2011 and 2017, with 2011 recruitment being a substantial driver of population recovery on Little Cayman. Export to a nearby island with a depleted population occurred in 2016. After two decades of protection, the population appears to be self-replenishing but also capable of seeding recruitment in the region, supporting calls to incorporate spawning aggregation protections into fisheries management.
Weighting data appropriately in stock assessment models is necessary to diagnose model mis-specification, estimate uncertainty, and when combining data sets. Age- and length-composition data are often fitted using a multinomial distribution and then reweighted iteratively, and the Dirichlet-multinomial ("DM") likelihood provides a model-based alternative that estimates an additional parameter and thereby "self-weights" data. However, the DM likelihood requires specifying an input sample size (n(input)), which is often unavailable and results are sensitive to n(input). We therefore introduce the multivariate-Tweedie (MVTW) as alternative with three benefits: (1) it can identify both overdispersion (downweighting) or underdispersion (upweighting) relative to the n(input); (2) proportional changes in n(input) are exactly offset by parameters; and (3) it arises naturally when expanding data arising from a hierarchical sampling design. We use an age-structured simulation to show that the MVTW (1) can be more precise than the DM in estimating data weights, and (2) can appropriately upweight data when needed. We then use a real-world state-space assessment to show that the MVTW can easily be adapted to other software. We recommend that stock assessments explore the sensitivity to specifying DM, MVTW, and logistic-normal likelihoods, particularly when the DM estimates an effective sample size approaching n(input).
The productivity of many fish populations is influenced by the environment, but developing environment-linked stock assessments remain challenging and current management of most commercial species assumes that stock productivity is time-invariant. In the Northeast United States, previous studies suggest that the recruitment of Southern New England-Mid Atlantic yellowtail flounder is closely related to the strength of the Cold Pool, a seasonally formed cold water mass on the continental shelf. Here, we developed three new indices that enhance the characterization of Cold Pool interannual variations using bottom temperature from a regional hindcast ocean model and a global ocean data assimilated hindcast. We associated these new indices to yellowtail flounder recruitment in a state-space, age-structured stock assessment framework using the Woods Hole Assessment Model. We demonstrate that incorporating Cold Pool effects on yellowtail flounder recruitment reduces the retrospective patterns and may improve the predictive skill of recruitment and, to a lesser extent, spawning stock biomass. We also show that the performance of the assessment models that incorporated ocean model-based indices is improved compared to the model using only the observation-based index. Instead of relying on limited subsurface observations, using validated ocean model products as environmental covariates in stock assessments may both improve predictions and facilitate operationalization.
Survival is an important population process in fisheries stock assessment models and is typically treated as deterministic. Recently developed state-space assessment models can estimate stochastic deviations in survival, which represent variability in some ambiguous combination of natural mortality (M), fishing mortality (F), and migration. These survival deviations are generally treated as independent by age and year, despite our understanding that many population processes can be autocorrelated and that not accounting for autocorrelation can result in notable bias. We address these concerns, as well as the strong retrospective pattern found in the last assessment of Southern New England yellowtail flounder (Limanda ferruginea), by incorporating two-dimensional (2D, age and year) first-order autocorrelation in survival and M. We found that deviations were autocorrelated among both years (0.53 +/- 0.09, 0.63 +/- 0.16) and ages (0.33 +/- 0.12, 0.40 +/- 0.16) when estimated for survival or M, respectively. Models with 2D autocorrelation on survival or M fit the data better and had reduced retrospective pattern than models without autocorrelation. The best fit model included 2D autocorrelated deviations in survival as well as independent deviations in M and altered estimates of spawning stock biomass by 18 % and F by 21 % in model years. In short-term projections with F = 0, including 2D autocorrelation in survival or M reduced spawning stock biomass by 48 %. We conclude that incorporating 2D autocorrelated variation in survival or M could improve the assessment of Southern New England yellowtail flounder in terms of model fit and consistency of biomass projections.
The Gulf Stream Index is calculated as the standardized first principal component time series from the empirical orthogonal function (EOF) analysis of the 200 m temperature time series from the EN4.2.1 dataset at the 20 base points, as detailed in Chen et al. [2021].
Fish spawning aggregations (FSAs) are vulnerable to overexploitation, yet quantitative assessments of FSA populations are rare. We document an approach for how to conduct such an assessment, evaluating the response of Critically Endangered Nassau Grouper (Epinephelus striatus) to protections in the Cayman Islands. We assessed pre-protection status on all islands using length data from fishery catch. We then used 17 years of noninvasive length-frequency data, collected via diver-operated laser calipers, to estimate recruitment and spawning biomass of Nassau Grouper on Little Cayman following protection. Bimodal length distributions in 2017-2019 indicated a large recruitment pulse (4-8x average) derived from spawning in 2011. Biomass recovered to 90-106% of the pre-exploitation level after 16 years, largely driven by the strong 2011 year class. Length distributions were also bimodal in 2017-2019 on nearby Cayman Brac, implying a synchronous recruitment pulse occurred on both islands. Our results demonstrate that: (i) in situ length data can be used to monitor protected FSAs; (ii) spatiotemporal FSA closures can be effective, but success takes time if population recovery depends upon sporadic recruitment; and (iii) FSA fishery management targets may need to be higher than commonly recommended (i.e. spawning potential ratio >0.6 instead of 0.4).
The rapid changes observed in many marine ecosystems that support fisheries pose a challenge to stock assessment and management predicated on time-invariant productivity and considering species in isolation. In single-species assessments, two main approaches have been used to account for productivity changes: allowing biological parameters to vary stochastically over time (empirical), or explicitly linking population processes such as recruitment (R) or natural mortality (M) to environmental covariates (mechanistic). Here, we describe the Woods Hole Assessment Model (WHAM) framework and software package, which combines these two approaches. WHAM can estimate time-and age-varying random effects on annual transitions in numbers at age (NAA), M, and selectivity, as well as fit environmental time-series with process and observation errors, missing data, and nonlinear links to R and M. WHAM can also be configured as a traditional statistical catch-at-age (SCAA) model in order to easily bridge from status quo models and test them against models with state-space and environmental effects, all within a single framework. We fit models with and without (independent or autocorrelated) random effects on NAA, M, and selectivity to data from five stocks with a broad range of life history, fishing pressure, number of ages, and time-series length. Models that included random effects performed well across stocks and processes, especially random effects models with a two dimensional (2D) first-order autoregressive, AR(1), covariance structure over age and year. We conducted simulation tests and found negligible or no bias in estimation of important assessment outputs (SSB, F, stock status, and catch) when the operating and estimation models matched. However, bias in SSB and F was often non-trivial when the estimation model was less complex than the operating model, especially when models without random effects were fit to data simulated from models with random effects. Bias of the variance and correlation parameters controlling random effects was also negligible or slightly negative as expected. Our results suggest that WHAM can be a useful tool for stock assessment when environmental effects on R or M, or stochastic variation in NAA transitions, M, or selectivity are of interest. In the U.S. Northeast, where the productivity of several groundfish stocks has declined, conducting assessments in WHAM with time-varying processes via random effects or environment-productivity links may account for these trends and potentially reduce retrospective bias.
The Northeast U.S. shelf (NES) is an oceanographically dynamic marine ecosystem and supports some of the most valuable demersal fisheries in the world. A reliable prediction of NES environmental variables, particularly ocean bottom temperature, could lead to a significant improvement in demersal fisheries management. However, the current generation of climate model-based seasonal-to-interannual predictions exhibits limited prediction skill in this continental shelf environment. Here, we have developed a hierarchy of statistical seasonal predictions for NES bottom temperatures using an eddy-resolving ocean reanalysis data set. A simple, damped local persistence prediction model produces significant skill for lead times up to months in the Mid-Atlantic Bight and up to in the Gulf of Maine, although the prediction skill varies notably by season. Considering temperature from a nearby or upstream (i.e., more poleward) region as an additional predictor generally improves prediction skill, presumably as a result of advective processes. Large-scale atmospheric and oceanic indices, such as Gulf Stream path indices (GSIs) and the North Atlantic Oscillation Index, are also tested as predictors for NES bottom temperatures. Only the GSI constructed from temperature observed at 200 m depth significantly improves the prediction skill relative to local persistence. However, the prediction skill from this GSI is not larger than that gained using models incorporating nearby or upstream shelf/slope temperatures. Based on these results, a simplified statistical model has been developed, which can be tailored to fisheries management for the NES. this study, we have developed a collection of statistical models that produce seasonal predictions of NES bottom temperature with 1–12 months lead time. Variables considered in these prediction models include local persistence of bottom temperature from prior months, bottom temperature from an upstream or nearby region, and large-scale atmospheric and oceanic indices representing the North Atlantic Oscillation or position of the Gulf Stream (GS). Only considering local persistence provides significant skill for lead times up to ∼ 5 months in the Mid-Atlantic Bight and up to ∼ 10 months in the Gulf of Maine, although the skill varies by season. Using upstream or nearby bottom temperature and the GS index both generally improve the prediction skill. However, the GS index does not provide higher prediction skill than those upstream or nearby bottom temperatures. A simplified statistical model has been developed, which can be tailored to fisheries management on the NES.
Spatiotemporal predictions of bycatch (i.e., catch of nontargeted species) have shown promise as dynamic ocean management tools for reducing bycatch. However, which spatiotemporal model framework to use for generating these predictions is unclear. We evaluated a relatively new method, Gaussian Markov random fields (GMRFs), with two other frameworks, generalized additive models (GAMs) and random forests. We fit geostatistical delta-models to fisheries observer bycatch data for six species with a broad range of movement patterns (e.g., highly migratory sea turtles versus sedentary rockfish) and bycatch rates (percentage of observations with nonzero catch, 0.3%–96.2%). Random forests had better interpolation performance than the GMRF and GAM models for all six species, but random forests performance was more sensitive when predicting data at the edge of the fishery (i.e., spatial extrapolation). Using random forests to identify and remove the 5% highest bycatch risk fishing events reduced the bycatch-to-target species catch ratio by 34% on average. All models considerably reduced the bycatch-to-target ratio, demonstrating the clear potential of species distribution models to support spatial fishery management.
The ongoing evolution of tracer mixing models has resulted in a confusing array of software tools that differ in terms of data inputs, model assumptions, and associated analytic products. Here we introduce MixSIAR, an inclusive, rich, and flexible Bayesian tracer (e.g., stable isotope) mixing model framework implemented as an open-source R package. Using MixSIAR as a foundation, we provide guidance for the implementation of mixing model analyses. We begin by outlining the practical differences between mixture data error structure formulations and relate these error structures to common mixing model study designs in ecology. Because Bayesian mixing models afford the option to specify informative priors on source proportion contributions, we outline methods for establishing prior distributions and discuss the influence of prior specification on model outputs. We also discuss the options available for source data inputs (raw data versus summary statistics) and provide guidance for combining sources. We then describe a key advantage of MixSIAR over previous mixing model software-the ability to include fixed and random effects as covariates explaining variability in mixture proportions and calculate relative support for multiple models via information criteria. We present a case study of Alligator mississippiensis diet partitioning to demonstrate the power of this approach. Finally, we conclude with a discussion of limitations to mixing model applications. Through MixSIAR, we have consolidated the disparate array of mixing model tools into a single platform, diversified the set of available parameterizations, and provided developers a platform upon which to continue improving mixing model analyses in the future.
Quantifying effects of fishing on non-targeted (bycatch) species is an important management and conservation issue. Bycatch estimates are typically calculated using data collected by on-board observers, but observer programmes are costly and therefore often only cover a small percentage of the fishery. The challenge is then to estimate bycatch for the unobserved fishing activity. The status quo for most fisheries is to assume the ratio of bycatch to effort is constant and multiply this ratio by the effort in the unobserved activity (ratio estimator). We used a dataset with 100% observer coverage, 35440 hauls from the US west coast groundfish trawl fishery, to evaluate the ratio estimator against methods that utilize fine-scale spatial information: generalized additive models (GAMs) and random forests. Applied to 15 species representing a range of bycatch rates, including spatial locations improved model predictive ability, whereas including effort-associated covariates generally did not. Random forests performed best for all species (lower root mean square error), but were slightly biased (overpredicting total bycatch). Thus, the choice of bycatch estimation method involves a tradeoff between bias and precision, and which method is optimal may depend on the species bycatch rate and how the estimates are to be used.
Increasing complexity in human-environment interactions at multiple watershed scales presents major challenges to sediment source apportionment data acquisition and analysis. Herein, we present a step-change in the application of Bayesian mixing models: Deconvolutional-MixSIAR (D-MIXSIAR) to underpin sustainable management of soil and sediment. This new mixing model approach allows users to directly account for the 'structural hierarchy' of a river basin in terms of sub-watershed distribution. It works by deconvoluting apportionment data derived for multiple nodes along the stream-river network where sources are stratified by sub-watershed. Source and mixture samples were collected from two watersheds that represented (i) a longitudinal mixed agricultural watershed in the south west of England which had a distinct upper and lower zone related to topography and (ii) a distributed mixed agricultural and forested watershed in the mid-hills of Nepal with two distinct sub-watersheds. In the former, geochemical fingerprints were based upon weathering profiles and anthropogenic soil amendments. In the latter compound-specific stable isotope markers based on soil vegetation cover were applied. Mixing model posterior distributions of proportional sediment source contributions differed when sources were pooled across the watersheds (pooled-MixSIAR) compared to those where source terms were stratified by sub-watershed and the outputs deconvoluted (D-MixSIAR). In the first example, the stratified source data and the deconvolutional approach provided greater distinction between pasture and cultivated topsoil source signatures resulting in a different posterior distribution to non-deconvolutional model (conventional approaches over-estimated the contribution of cultivated land to downstream sediment by 2 to 5 times). In the second example, the deconvolutional model elucidated a large input of sediment delivered from a small tributary resulting in differences in the reported contribution of a discrete mixed forest source. Overall D-MixSIAR model posterior distributions had lower (by ca 25-50%) uncertainty and quicker model run times. In both cases, the structured, deconvoluted output cohered more closely with field observations and local knowledge underpinning the need for closer attention to hierarchy in source and mixture terms in river basin source apportionment. Soil erosion and siltation challenge the energy-food-water-environment nexus. This new tool for source apportionment offers wider application across complex environmental systems affected by natural and human-induced change and the lessons learned are relevant to source apportionment applications in other disciplines.
Compound-specific stable isotope (CSSI) fingerprinting of sediment sources is a recently introduced tool to overcome some limitations of conventional approaches for sediment source apportionment. The technique uses the 13C CSSI signature of plant-derived fatty acids (δ13C-fatty acids) associated with soil minerals as a tracer. This paper provides methodological perspectives to advance the use of CSSI fingerprinting in combination with stable isotope mixing models (SIMMs) to apportion the relative contributions of different sediment sources (i.e. land uses) to sediments.
Many animals are considered to be specialists because they have feeding structures that are fine-tuned for consuming specific prey. For example, “smasher” mantis shrimp have highly specialized predatory appendages that generate forceful strikes to break apart hard-shelled prey. Anecdotal observations suggest, however, that the diet of smashers may include soft-bodied prey as well. Our goal was to examine the diet breadth of the smasher mantis shrimp, Neogonodactylus bredini, to determine whether it has a narrow diet of hard-shelled prey. We combined studies of prey abundance, feeding behavior, and stable isotope analyses of diet in both seagrass and coral rubble to determine if N. bredini’s diet was consistent across different habitat types. The abundances of hard-shelled and soft-bodied prey varied between habitats. In feeding experiments, N. bredini consumed both prey types. N. bredini consumed a range of different prey in the field as well and, unexpectedly, the stable isotope analysis demonstrated that soft-bodied prey comprised a large proportion (29–53 %) of the diet in both habitats. Using a Bayesian mixing model framework (MixSIAR), we found that this result held even when we used uninformative, or generalist, priors and informative priors reflecting a specialist diet on hard-shelled prey and prey abundances in the field. Thus, contrary to expectation, the specialized feeding morphology of N. bredini corresponds to a broad diet of both hard-shelled and soft-bodied prey. Using multiple lines of study to describe the natural diets of other presumed specialists may demonstrate that specialized morphology often broadens rather than narrows diet breadth.
Mixing models are statistical tools that use biotracers to probabilistically estimate the contribution of multiple sources to a mixture. These biotracers may include contaminants, fatty acids, or stable isotopes, the latter of which are widely used in trophic ecology to estimate the mixed diet of consumers. Bayesian implementations of mixing models using stable isotopes (e.g., MixSIR, SIAR) are regularly used by ecologists for this purpose, but basic questions remain about when each is most appropriate. In this study, we describe the structural differences between common mixing model error formulations in terms of their assumptions about the predation process. We then introduce a new parameterization that unifies these mixing model error structures, as well as implicitly estimates the rate at which consumers sample from source populations (i.e., consumption rate). Using simulations and previously published mixing model datasets, we demonstrate that the new error parameterization outperforms existing models and provides an estimate of consumption. Our results suggest that the error structure introduced here will improve future mixing model estimates of animal diet.