Otolith shape variation provides a marker of stock origin for north Atlantic bluefin tuna (Thunnus thynnus)

Additional keywords: elliptical Fourier analysis, population structure, stock mixture analysis.

The stock assessment models on which the management of most fish species is based rely on the critical assumption that the unit of management is a discrete self-sustaining component with a single rate of productivity (Kerr et al. 2010). However, in reality, most exploited fish stocks exhibit some degree of spatial and biological complexity (Booth 2000; Guan et al. 2013). For highly migratory species like the Atlantic bluefin tuna (Thunnus thynnus), sustainable management requires knowledge of how movements and mixing of the targeted stocks correlate with seasonal and geographic patterns of fishing (Secor 2014). When biological processes such as migration underlie fluctuations in catches, changes in fishing mortality may be decoupled from variations in stock size, leading to erroneous estimates of stock status (Fromentin and Kell 2007; Kell and Fromentin 2007). Fisheries that coincide with overlap in the distribution of multiple stocks (mixed stock fisheries) are particularly vulnerable to mismanagement because a failure to properly account for the occurrence of multiple components within the fishery can lead to overexploitation of less-productive stocks (Fu and Fanning 2004; Ying et al. 2011). Population diversity can dampen fluctuations in abundance and confer stability to the stock complex (Schindler et al. 2010). Therefore, erosion of stock complexity may threaten the resilience of a species or metapopulation and its ability to recover from overfishing and stock collapse (Petitgas et al. 2010; Fromentin et al. 2014b).

Atlantic bluefin tuna are managed as two separate stocks by the International Commission for the Conservation of Atlantic Tunas (ICCAT): the western stock spawns in the Gulf of Mexico, whereas the eastern stock spawns in the Mediterranean Sea (Fromentin and Powers 2005). Severe overfishing in recent decades led to marked declines in abundance (Safina and Klinger 2008; MacKenzie et al. 2009), prompting fishing restrictions and stock-rebuilding plans (Fromentin et al. 2014a; ICCAT 2014). According to the latest assessment, both stocks are showing increases in spawning stock biomass in recent years that are indicative of a steady recovery. However, the assessment is subject to uncertainty, particularly with regard to stock-mixing rates (Fromentin et al. 2014a; ICCAT 2014)). Improving the accuracy of mixing rate estimates is acknowledged by ICCAT as a research priority (ICCAT 2008, 2013) The 45° meridian is the assumed boundary between the eastern and western stocks, although trans-Atlantic movements are known to occur (Block et al. 2005; Rooker et al. 2014). Although stock discreetness is maintained by high rates of natal homing, there is considerable overlap in the distribution of the two stocks outside of spawning time (Carlsson et al. 2006; Rooker et al. 2008). Consequently, commercial fisheries are likely to exploit a mixed assemblage of eastern and western origin fish, particularly in the central north Atlantic and on the east coast of the US (Rooker et al. 2007). The relative contributions of the two stocks to mixed assemblages vary spatially, temporally and ontogenetically (Taylor et al. 2011). Exploitation in mixed fisheries is likely to have a disproportionate effect on the western stock because its productivity is estimated to be one-tenth that of the eastern stock (Fromentin and Powers 2005). Stock simulation results indicate that the incorporation of mixing rates into stock assessment within a mixed-stock modelling framework is essential for sound evaluation of quotas and effective stock-rebuilding plans, particularly for west Atlantic bluefin tuna (Taylor et al. 2011).

Support for the existence of two discrete stocks of Atlantic bluefin tuna is provided by several genotypic and phenotypic markers. Analyses of DNA microsatellites (Carlsson et al. 2007), mitochondrial DNA (Boustany et al. 2008) and single nucleotide polymorphisms (Albaina et al. 2013) all show significant genetic divergence of the two populations. In addition, there is evidence of genetic heterogeneity within the Mediterranean, and the existence of separate spawning populations in the western and eastern Mediterranean has been proposed (Carlsson et al. 2004; Riccioni et al. 2010). The occurrence of mature adult Atlantic bluefin tuna in areas other than the main spawning grounds during the spawning season may also be indicative of a more complex structure than the currently accepted two-stock model (Galuardi et al. 2010). Young-of-the-year from the main spawning areas in the Gulf of Mexico and Mediterranean show differences in otolith trace elements (Rooker et al. 2003) and stable isotopes (Rooker et al. 2008), as well as concentrations of organochlorines and polychlorinated biphenyls (PCBs) in tissues (Dickhut et al. 2009). Stock-specific differences in otolith stable isotopes have facilitated the estimation of stock composition in potentially mixed assemblages (Schloesser et al. 2010; Fraile et al. 2014; Rooker et al. 2014). However, the degree of uncertainty associated with this method of assignment (Rooker et al. 2014) warrants the investigation of other stock-identification approaches to corroborate estimated mixing rates and ultimately to develop a multimarker approach to stock assignment for greater precision in mixed-stock projections. As more population markers are identified, it may also be possible to resolve finer-scale structure and complex migration pathways in Atlantic bluefin tuna.

Otolith shape is known to vary both between and within species (Lombarte and Castellon 1991) because of the combined effects of genetic and environmental factors (Vignon and Morat 2010). Fish with different life histories show variation in otolith shape (Vignon and Morat 2010), thus otolith shape measurements can be used to discriminate between stocks. Thanks to advances in image analysis, variation in otolith shape is now readily captured using geometric measurement of digitised otolith outlines (Stransky 2014). Multivariate analysis of otolith shape data can be used to characterise fish from different stocks (Begg et al. 2001; Paul et al. 2013) or to detect underlying structure in a mixed assemblage of unknown stock composition (Keating et al. 2014). Geographic variation in the shape of Atlantic bluefin tuna otoliths tuna has not been investigated previously. Given the range of environmental conditions that occurs across the wide distribution of the species, stock-specific differences in otolith shape may exist and could prove useful as markers of stock origin.

The present study examined the feasibility of using otolith shape variation as a marker of stock origin for Atlantic bluefin tuna. Spatial, temporal and length-based variability in otolith shape descriptors were investigated. After accounting for the effects of size and interannual variability, otolith shape descriptors were used to characterise the western and eastern stocks and to estimate the relative contributions of the western and eastern stocks to mixed aggregations. The accuracy of the method for discriminating between stocks and estimating mixing rates was evaluated in relation to other approaches and some recommendations for future applications of the method are presented.

Atlantic bluefin tuna were collected from 11 locations in the western Atlantic, central north Atlantic, eastern Atlantic and Mediterranean Sea in 2011, 2012 and 2013 (Fig. 1; Table 1). Samples were obtained under the provision of the ICCAT Atlantic Wide Research Program for Bluefin Tuna (GBYP) or directly from national fisheries sampling programs. In the central north Atlantic, east Atlantic and Mediterranean Sea, the fish were collected during commercial fishing operations using a combination of capture methods (traps, long-lining, bait boats). Samples from the west Atlantic were collected from commercially caught bluefin tuna (rod and reel) in the three principal fishing regions (Gulf of St Lawrence, Newfoundland and the Scotian Shelf). Sagittal otoliths were removed and fork lengths (FL) were recorded. For some fish, FLs were not available and were estimated using a monthly length–weight conversion (ICCAT 2006) and snout length conversion (Secor et al. 2014).

Fig. 1.  Map showing the sampling locations for Atlantic bluefin tuna (Thunnus thynnus) used in the otolith shape analysis: (a) Canada (CD), central north Atlantic (CA); (b) Bay of Biscay (BB), Portugal (PO), Straits of Gibraltar (GI), Morocco (MO); (c) Gulf of Lyon (GL), Ligurian Sea (LI), Adriatic Sea (AS), Malta (MA) and Levantine Sea (LS). The positions of the four regions referred to in the text (western Atlantic, central north Atlantic, eastern Atlantic and Mediterranean Sea) are shown in (a). The authors acknowledge the use of the Maptool program for the production of the sampling distribution map. Maptool is a product of SEATURTLE.ORG (information available at http://www.seaturtle.org, accessed 12 August 2015).

Table 1.  Summary details of Atlantic bluefin tuna (Thunnus thynnus)
AS, Adriatic Sea; MA, Malta; GL, Gulf of Lyon; LS, Levantine Sea; LI, Ligurian Sea; BB, Bay of Biscay; MO, Morocco; PO, Portugal; GI, Gibraltar; CA, central north Atlantic; CD, Canada; FL, fork length

Otolith images were captured using a stereomicroscope connected to a digital camera with a PC interface. Sagittal otoliths were photographed as a white object on a black background in a standard orientation, with the sulcus side uppermost. Otoliths were excluded from shape analysis when their outline was obscured by breakage or adhering dirt or tissue. To increase the availability of unbroken otoliths, both left and right otoliths were used and images were digitally rotated to produce a standard orientation. From the available material, 718 otolith images were selected in a stratified random manner that maximised the spatial, temporal coverage and the length range of the samples.

Otolith images were edited to standardise their orientation and to remove visual artefacts using Paint.NET v3.5.10 (http://www.getpaint.net, accessed 12 August 2015). Using the ImageJ software package (National Institutes of Health, Bethesda, MD, USA; http://imagej.nih.gov/ij/, accessed 12 August 2015), a set of morphological shape indices was obtained from physical measurements of each otolith image, as follows: (1) circularity, calculated as (4π × (area/perimeter²)); (2) aspect ratio, the ratio of the major and minor axes of the ellipse that binds the outline; and (3) roundness, calculated as (4 × area/(π × major axis²)).

Using the TPSdig utility (http://life.bio.sunysb.edu/morph/software.html, accessed 12 August 2015 ), images were converted to binary and otolith outlines were traced using edge detection and saved as a series of x,y coordinates. Elliptical Fourier harmonics were extracted from smoothed otolith outlines (using 200 smoothing iterations) using the momocs package in R (http://cran.r-project.org/web/packages/Momocs/Momocs.pdf, accessed 12 August 2015. The Fourier power equation (Crampton 1995) showed that over 99.9% of the variability in shape was captured by the first 12 harmonics. Each harmonic is composed of four coefficients (a_n, b_n, c_n and d_n). The first three coefficients of Harmonic 1 (a₁, b₁ and c₁) were used to standardise each outline for size, orientation and starting point. Thus, a total of 45 coefficients was included in the subsequent analysis.

The shape of the otolith is under ontogenetic control and is known to change as a fish grows (Hüssy 2008). Otolith shape could also vary from year to year within a region because of variations in the environment or the age structure of the population (Vignon 2012). These potential sources of variation could confound the interpretation of regional differences in otolith shape. Therefore, before attempting to characterise fish from the east and west Atlantic using otolith shape, the relative contributions of sampling location, collection year and fish length on otolith shape variability were investigated.

The 45 elliptical Fourier coefficients and three shape indices (henceforth collectively referred to as the shape descriptors) were tested for normality and transformed when necessary. Eleven coefficients showed significant deviation from normality (based on visual inspection of the probability distribution) that could not be corrected by transformation, leaving a total of 37 shape descriptors in the subsequent analysis. Principal component analysis (PCA) was conducted to reduce the dataset to a manageable number of orthogonal descriptor variables that summarised the variability in otolith shape. The scree plot (eigenvalues against principal component number) and cumulative variance values were used to identify the principal components that explained most of the variability in otolith shape. Correlations between these principal components and FL were tested statistically using Pearson’s correlation coefficient. A series of general linear models (GLMs) compared principal component scores between locations and years (fixed factors), with FL included as a covariate. Principal components that were significantly correlated with length and showed no spatial heterogeneity in the size–shape relationship (i.e. length × location interaction P > 0.05) were standardised using the common within-group slope, according to the following equation:

where Y_c is the corrected variable, Y is the original variable, b is the common within-group slope of the shape–size relationship and L is the measurement of fish size (FL, cm).

In order to visually represent how otolith shape changes with length, the otolith outlines were grouped into 60-cm length classes (<60, 61–120, 121–180, 181–240 and 241–300 cm) and average outlines for each length class were reconstructed from the first 12 elliptical Fourier harmonics using the mean Shapes function in the momocs package in R.

To investigate possible structuring within the east Atlantic stock, length-corrected principal components were compared between years and sampling locations in the Mediterranean Sea, Bay of Biscay and Straits of Gibraltar in a GLM. To further minimise the potentially confounding effects of differences in age classes between samples, only fish <160 cm were included in the analysis. In all, 393 fish from six sampling sites (Bay of Biscay, Gibraltar, Adriatic Sea, Gulf of Lyon, Ligurian Sea, Levantine Sea) were included in the comparison.

From the full dataset, a subset was selected for which stock origin could be assumed with a reasonable degree of certainty. Samples of adult bluefin from the western Atlantic spawning grounds in the Gulf of Mexico, which would constitute the ideal baseline for this stock, were not available. A previous study that measured stable isotopes in the otoliths of Atlantic bluefin tuna collected from the Gulf of St Lawrence and adjacent areas over three decades showed that over 99% the fish in this area originate from nurseries in the western Atlantic (Schloesser et al. 2010). Therefore, the samples collected from Canadian fishery in September 2013 were assumed to be representative of the west Atlantic stock. Adult Atlantic bluefin collected from a known spawning area in the Mediterranean Sea (off the coast of Malta) during the spawning season (May–July; Corriero et al. 2003; Rooker et al. 2007) were assumed to belong to the east Atlantic stock. Samples of Atlantic bluefin collected from Malta in 2011 were shown to comprise 100% eastern origin fish based on stable oxygen isotopes (Rooker et al. 2014). In order to minimise confounding effects of length on otolith shape, only fish >200 cm FL were used. In addition, only images of right otoliths were used. Thus, the baseline samples included 50 Atlantic bluefin from the west Atlantic (Canada) and 60 from the east Atlantic (Malta; Table 1).

Shape descriptors were tested for normality and transformed when necessary. Shape descriptors that were significantly correlated with length (with a Pearson correlation coefficient ≥0.3) and showed no spatial heterogeneity in the size–shape relationship (i.e. length × location interaction P > 0.05) were standardised using the common within-group slope as described above. A series of GLMs was used to identify which shape descriptors showed significant variation between the east Atlantic and west Atlantic baseline samples. These shape descriptors were included in a second PCA to reduce the dimensionality of the dataset. The principal components that captured the majority of the variation in otolith shape were used in a stepwise discriminant function analysis (DFA) to distinguish between fish from the western and eastern baseline samples.

A mixed sample was selected from the dataset that included Atlantic bluefin >200 cm FL from the east Atlantic (Portugal, Morocco and the Straits of Gibraltar) and the central Atlantic (Table 1). All images used were of right otoliths. The principal components that were selected in the DFA to distinguish between the two baseline groups were used in a Bayesian stock mixture analysis to estimate the proportion of fish in each sample that originated from the east and west Atlantic stocks. The analysis was conducted using the package mixFish in R as described by Smith and Campana (2010). When using the Bayesian approach, the observations from the mixed samples classified to the base populations (the western and eastern stocks in this case) are used to update the parameter estimation of the base population (unconditional estimation). Bayesian credible intervals (95% CI) were calculated as a measure of the uncertainty associated with the estimated proportions. The fit of the classification model was evaluated by comparing the means from the posterior predictive distribution with the original means for each principal component from the base populations (post-predictive check).

The PCA of the full dataset was based on 34 harmonics and three shape indices. The first 11 principal components explained 70% of the total variance in the otolith shape descriptors. Principal component (PC) 1, PC2, PC3 and PC4 were significantly correlated with length and were therefore standardised for the length effect using the common within-group slope.

For the fish <160 cm FL, there were small (R² values ranged from 1.7 to 6.5%) but significant (P < 0.05) differences between sampling locations in the east Atlantic and Mediterranean for PC2, PC4, PC6, PC8, PC10 and PC11, whereas PC3 and PC5 showed significant differences between locations and sampling years. For the principal components that showed no interannual variation, further examination of the mean and 95% confidence limits (Fig. 2) showed that the Adriatic Sea samples diverged from the Bay of Biscay samples (PC4, PC6 and PC8) and from the Ligurian Sea and Levantine Sea samples (PC2). There was also a significant difference in PC10 between fish from the Bay of Biscay and the Straits of Gibraltar.

Fig. 2.  Mean ± 95% confidence intervals for principal component scores for fish <160-cm fork length from the Mediterranean (Adriatic Sea (AS), Gulf of Lyon (GL), Ligurian Sea (LI), Levantine Sea (LS)) and East Atlantic (Bay of Biscay (BB), Straits of Gibraltar (GI)). Corpc, length-corrected principal component; PC, principal component.

The nature of the relationship between otolith shape and fish size is shown in Fig. 3. Between 60 and 120 cm FL, the shape of the otolith changes markedly, narrowing at the rostrum and post-rostrum. The antirostrum is substantially less pronounced in the larger length classes and the post-rostrum is more angular. From 120–300 cm FL, the notch between the antirostrum and the post-rostrum (excisura major) gradually flattens and the otolith becomes more slender and elongated at the rostrum and post-rostrum. In the largest size class, the post-rostrum is particularly angular and is characterised by a protrusion on the posterior end of the ventral margin.

Fig. 3.  Atlantic bluefin tuna (Thunnus thynnus) otolith outlines reconstructed using the first 12 elliptical Fourier harmonics and averaged across five length classes (60, 120, 180, 240 and 300 cm). The position of major morphological landmarks are indicated, namely the rostrum (R), post-rostrum (PR), antirostrum (AR) and excisura major (EM).

Stock-specific differences in otolith shape are apparent in the mean outlines of otoliths from Atlantic bluefin >200 cm FL from the east Atlantic (Malta) and west Atlantic (Canada) baseline samples (Fig. 4). This variation is less marked than the size-based changes described above. At the anterior end of the otolith, the width of the rostrum is similar in both types but the east Atlantic otoliths are narrower in the mid-portion of the otolith and at the dorsal margin of the posterior end of the otolith. The angle of the posterior end of the otolith is steeper in the fish from the west Atlantic and the post-rostrum is narrower on the ventral margin.

Fig. 4.  Mean otolith outlines reconstructed using the first 12 elliptical Fourier harmonics for Atlantic bluefin tuna (Thunnus thynnus) from the east Atlantic (black) and west Atlantic (grey). The mean shapes describe the baseline samples used in the stock composition analysis. AR, antirostrum; PR, post-rostrum; R, rostrum.

In all, 25 elliptical Fourier coefficients and three shape indices showed significant variation between the east and west Atlantic (GLM P < 0.05) and were not significantly correlated with length (in some cases after standardisation). The first 12 principal components from a PCA of these shape descriptors explained 84% of the variability in the dataset. The loading plot for the first two principal components (Fig. 5) highlights the shape descriptors that make major contributions to the variability in the dataset (e.g. roundness, aspect ratio, Elliptical Fourier coefficients D1, A11, C11, D9, A12 and B10).

Fig. 5.  Loading plot for the first two principal components, derived using 25 elliptical Fourier coefficients (D1, B2, A3 etc.) and three shape indices (circularity, roundness and aspect ratio) from otolith outlines of Atlantic bluefin tuna (Thunnus thynnus) >200-cm fork length from the east, west and central Atlantic and the Mediterranean Sea. Variables that were transformed for normality before the analysis are indicated in italic, whereas those that were corrected for correlations with length are indicated in bold.

Six principal components were retained in the DFA by stepwise selection (PC1, PC2, PC5, PC6, PC8, PC11) producing one canonical function that distinguished between otoliths from east Atlantic and west Atlantic fish (Wilk’s λ = 0.49, approximate F = 15.9, P < 0.0001). The canonical coefficients showed that PC1 and PC2 made the largest contribution to the separation of the groups (Table 2). The canonical scores plot showed some separation of the two baseline samples but with overlap (Fig. 6a). The distribution of the canonical scores for the mixed sample overlapped with the distributions for both the west and east Atlantic baseline samples but was most similar to the east Atlantic (Fig. 6b). The canonical function distinguished between fish of eastern and western origin with a mean jackknife classification success rate of 83% (Table 3).

Table 2.  Standardised coefficients of the canonical discriminant function, showing the relative contribution of each principal component to the classification of the east Atlantic and west Atlantic baseline samples of Atlantic bluefin tuna (Thunnus thynnus)
The principal components (PC) are ranked according to the order in which they were included in the model by stepwise selection

Fig. 6.  Density distributions of the canonical scores from the discriminant function analysis (DFA) for (a) baseline samples of east Atlantic and west Atlantic bluefin tuna (Thunnus thynnus) and (b) with the scores for the mixed sample overlaid.

Table 3.  Jackknife classification matrix from the discriminant function analysis using principal component (PC)1, PC2, PC5, PC6, PC8 and PC11 to discriminate between Atlantic bluefin tuna (Thunnus thynnus) from the east and west Atlantic baseline samples

The results of the Bayesian stock mixture analysis indicated that all the mixed samples were predominantly of eastern origin (>90% on average; Table 4). The mean estimated proportions suggested that the sample from Morocco contained the greatest proportion of western origin fish; however, the 95% CI were overlapping, indicating no significant difference in the proportions across the four mixed samples (Fig. 7). The posterior predictive check confirmed that the classification model fit the data reasonably well because the means from the posterior predictive distribution corresponded with the distribution of the original means and the probability of observing a more extreme value for the mean did not fall outside the 0.90 bounds for any of the principal components. However, the probabilities of observing an extreme value (Fig. 8) were higher for the base samples from the east Atlantic than from the west Atlantic, suggesting that the individuals in the mixed sample that were classified as being of eastern origin were less similar to the baseline samples than the individuals that were classified as being of western origin.

Table 4.  Mean predicted percentages (±1 s.d.) and 95% Bayesian credible intervals (CI) for eastern and western origin fish in samples of Atlantic bluefin tuna (Thunnus thynnus) collected from different locations in the central and west Atlantic based on conditional Bayesian estimation (mixFish program)

Fig. 7.  Credible intervals (95%) for the posterior estimates of the proportions of bluefin tuna (Thunnus thynnus) in the mixed samples assigned to the east Atlantic (E) and west Atlantic (W) base groups. Triangles (E) and circles (W) represent the position of the posterior mean and the upper and lower limits of the estimates, dotted lines (E) and solid lines (W) represent the intervals. CA, Central Atlantic; GI, Straits of Gibraltar; MO, Morocco; PO, Portugal.

Fig. 8.  Probability of observing a more extreme value for the mean of each of the principal components (PCs) of the otolith shape descriptors from the marginal posterior distribution relative to the distribution for the base samples from the east Atlantic (E; circles and solid line) and west Atlantic (W; triangle and dashed line). Bounds for the 0.90 probability limits are indicated with a dashed line at the top and bottom of the graph. Probabilities falling between the 0.05 and 0.95 limits indicate a reasonable model fit.

The results of the present study confirm that the shapes of otoliths from west Atlantic and east Atlantic bluefin tuna are sufficiently distinct to allow the discrimination of the two stocks based on elliptical Fourier descriptors and morphometric indices. Although size-based variation in shape was more pronounced than stock-specific differences, within a restricted size range fish could be classified to their parent stock with a reasonably high level of accuracy (83%). We have demonstrated the feasibility of using otolith shape descriptors to estimate the relative contributions of the western and eastern stocks to mixed aggregations. The method has the potential to complement other population markers, such as genetics (Carlsson et al. 2006), otolith stable isotopes (Rooker et al. 2014) and organochlorine markers (Dickhut et al. 2009; Graves et al. 2015), in order to improve the precision of mixed-stock estimates.

The strong effect of fish length on otolith shape variation in Atlantic bluefin tuna is not surprising and has been reported in many other species (Smith 1992; Mérigot et al. 2007; Capoccioni et al. 2011). The relationship between fish size and otolith shape reflects the combined effects of ontogeny and the environment on otolith shape (Vignon 2012). The overall species-specific shape of the otolith and the nature of its development during ontogeny from the circular larval otolith to the more complex morphology of the adult otolith is genetically determined (Hüssy 2008; Reichenbacher et al. 2009; Vignon and Morat 2010), with exogenous factors having a modulating effect through, for example, the effects of feeding and growth on the nature of crystal formation (Gauldie and Nelson 1990) and the rate of protein accretion in the otolith (Hüssy 2008). The shape of the otolith may also reflect its physiological function in hearing and balance, and such variation can have an adaptive significance. Ecomorphological studies have related intra- and inter-specific variation in otolith shape to swimming performance and feeding behaviour (Kishida et al. 2011), trophic niche (Lombarte et al. 2010) and habitat preferences (Volpedo et al. 2008; Volpedo and Fuchs 2010). For example, a long elongated otolith with a well-developed rostrum is characteristic of pelagic species, whereas benthic ecotypes typically have rounder wider otoliths (Volpedo and Echeverria 2003). Therefore, the marked change in the shape of the Atlantic bluefin otolith between the 60- and 120-cm length classes could be related to ontogenetic development and associated changes in diet and migratory behaviours (Fromentin and Powers 2005; Rooker et al. 2007). Regardless of the underlying mechanisms, the occurrence of these size-related changes in shape necessitate the careful consideration of size distribution when using otolith shape to resolve stock structure and estimate mixing rates in Atlantic bluefin tuna. In the present study it was necessary to restrict the analysis of stock composition to the size range for which baseline samples were available for both stocks. A priority for future applications of the method should be to obtain samples of known spawning origin from a wider size range that is fully representative of the size distribution in the mixing areas.

The small but significant differences in otolith shape between sampling locations within the eastern Atlantic and Mediterranean Sea may reflect the existence of groups of fish with different life histories in the samples. It is important to bear in mind that these differences do not necessarily reflect diversity in spawning origin because the otolith outline reflects the conditions experienced by the fish throughout its life. Differences between locations could be explained by variation in the environments experienced by fish from different year-classes, as has been observed in other species (Campana and Casselman 1993; Bolles and Begg 2000) or could occur if contingents within a single spawning stock pursued divergent migration pathways. It has been suggested, based on evidence from electronic tagging, that the eastern stock of bluefin tuna comprises a resident component that spends prolonged periods foraging in the Mediterranean and a more migratory component that crosses the Straits of Gibraltar to feed in the Atlantic Ocean (Fromentin and Lopuszanski 2013; Cermeño et al. 2015). Such heterogeneity in life histories could contribute to the differences in otolith shape observed between fish from the Bay of Biscay and the Adriatic Sea, for example. Further speculation regarding the origin of the subtle differences in otolith shape that were observed within the Mediterranean and East Atlantic would be beyond the limitations of currently available data. However, the results indicate that variation in otolith shape could potentially be used to distinguish between Atlantic bluefin spawning in different areas of the Mediterranean and to resolve the structuring that has been observed here (Carlsson et al. 2004; Riccioni et al. 2010). This would require the collection of baseline samples from putative spawning areas during the spawning season, which were not available for the present study. An analysis that combined information from tagging, genetics otolith chemistry and otolith shape could allow for the characterisation of contingents within the stock and provide powerful insights into variations in life histories both within and between spawning populations of Atlantic bluefin tuna.

East Atlantic and west Atlantic bluefin were discriminated with a level of accuracy (jackknife classification success of DFA = 83%) that is comparable to that achieved with other methods, namely 87% for otolith stable isotopes (Rooker et al. 2014) and 85% for otolith trace elements (Rooker et al. 2003). With all these methods there is overlap between the two stocks in the traits that are used to characterise them. This introduces uncertainty when assigning samples of unknown origin to their parent stock. Combining information from different sources in a holistic approach to stock identification can improve the accuracy of stock discrimination and provide greater insight into stock structure and the mechanisms that underlie it (Begg and Waldman 1999; Cadrin et al. 2010; Baldwin et al. 2012). Proper integration of the information provided by multiple stock markers requires careful consideration and coordination of sampling strategies (Abaunza et al. 2008). A limitation of the present study lies in the sampling design; although there is some overlap in the years and locations used in the otolith shape analysis and in previous analyses of otolith stable isotopes (Fraile et al. 2014; Rooker et al. 2014), because the otolith shape study was initiated after most of the sampling had been conducted, the majority of otoliths were not analysed using both methods. The simultaneous collection of otolith shape and chemistry data is logistically very feasible because otoliths can be photographed without interfering with subsequent chemical analysis. It is recommended that this approach be taken in future sampling programs for Atlantic bluefin tuna.

The estimates of stock composition in the mixed samples that were obtained using otolith shape descriptors are within the ranges reported by Rooker et al. (2014) for Atlantic bluefin collected from the same locations (albeit in different years and from a wider size range of fish). According to our estimates, the sample from the central north Atlantic collected in 2012 comprised 94 ± 7% (mean ± s.d.) eastern origin fish. Although this is higher than the 63.9 ± 9.6% reported by Rooker et al. (2014) for the same location in 2010, it is consistent with their estimates of 90.7 ± 5.3% for 2011. Our estimate is based on just 10 fish within a restricted size range and should therefore not be taken as representative of the entire fishery for that year. Samples of bluefin (200–297 cm FL) from the Straits of Gibraltar and Portugal were predominantly of eastern origin according to the estimates presented here (91 ± 10% and 94 ± 7% respectively) and in Rooker et al. (2014). Although our mean estimate of the proportion of eastern origin fish in catches from Morocco is lower than that reported by Rooker et al. (2014; 79 ± 13% v. 93.9 ± 4.7% respectively), the confidence limits overlap, indicating no significant difference.

When integrating information from multiple stock markers, it is important to consider the specific perspective provided by each marker (Cadrin et al. 2014). In this case, the analysis of otolith shape integrates information over the entire life of the fish, whereas the chemical composition of the otolith core reflects only the larval environment. Therefore, spatial and temporal variation shape may not always coincide with variation in otolith composition and vice versa. Our Bayesian CIs were wider than the bootstrap limits presented in Rooker et al. (2014), particularly for Portugal and the Straits of Gibraltar, which were unequivocally (100 ± 0%) assigned to the eastern stock by Rooker et al. (2014). Similar observations have been made by others (Bolker et al. 2003; Smith and Campana 2010), but it has been suggested that Bayesian CIs are more robust than bootstrap confidence intervals, particularly when the contribution of many of the base populations to some of the mixed samples is low (Bolker et al. 2003). However, the greater precision obtained by Rooker et al. (2014) may reflect greater overlap in the otolith core signatures between the baseline and mixed samples than we obtained with the shape descriptors. In the present study, the posterior predictive check revealed that agreement between the baseline and mixed sample was better for the fish that were predicted to be of western origin than those that were predicted to be of eastern origin. This suggests that there is greater variability in the eastern population than is captured by the eastern baseline samples. This could arise if fish of eastern origin from the mixed sample represent multiple spawning groups or if they contain contingents with diverse environmental histories, as discussed above. The latter would not affect assignments based on otolith core composition and may explain the lower precision of the estimates obtained using otolith shape descriptors. It is recommended that future sampling efforts focus on obtaining samples from other putative spawning grounds in the Mediterranean to enable the eastern baseline to be more accurately characterised.

It is important to consider the possibility that the baseline samples do not all belong to the spawning populations from which they are assumed to originate. Ideally, samples of spawning adults collected from known spawning grounds should be used to characterise the stocks. For the western stocks, it was not possible to obtain samples from the spawning area in the Gulf of Mexico. Previously published analyses of otolith stable isotope signatures provided indirect verification of the likely spawning origin of these fish (Schloesser et al. 2010; Rooker et al. 2014). If some of the assumed western origin fish are actually of eastern origin (and vice versa), this would most likely increase the similarity of the two baselines in terms of otolith shape and reduce the discriminatory power of the otolith shape descriptors. The fact that the two baseline samples can be distinguished with a reasonably high rate of accuracy provides some reassurance that most of the individuals do originate from distinct stocks. Although possible errors in the composition of the baseline do not invalidate the main conclusions of the study, a priority for future applications of this approach should be to obtain baseline samples of indisputable origin. This could be achieved by collecting spawning adults from all known spawning areas in multiple years and by obtaining information from several stock markers.

To conclude, the results of the present study confirm that otolith shape analysis can be used to distinguish between east Atlantic and west Atlantic bluefin tuna with a high rate of success that is comparable to that achieved with other methods. The potential exists to use otolith shape descriptors in combination with other population markers to provide greater insight into complexity and to improve the accuracy of mixing-rate projections. This would greatly assist the incorporation of mixing rates into stock assessment and the development of spatially explicit stock assessment models to improve the reliability of mortality and stock size estimates and allow management scenarios to be properly evaluated.