Wide-Ranging Temporal Variation in Transoceanic Movement and Population Mixing of Bluefin Tuna in the North Atlantic Ocean

Rooker, Jay R.; Fraile, Igaratza; Liu, Hui; Abid, Noureddine; Dance, Michael A.; Itoh, Tomoyuki; Kimoto, Ai; Tsukahara, Yohei; Rodriguez-Marin, Enrique; Arrizabalaga, Haritz

doi:10.3389/fmars.2019.00398

ORIGINAL RESEARCH article

Front. Mar. Sci., 11 July 2019

Sec. Marine Megafauna

Volume 6 - 2019 | https://doi.org/10.3389/fmars.2019.00398

Wide-Ranging Temporal Variation in Transoceanic Movement and Population Mixing of Bluefin Tuna in the North Atlantic Ocean

$\r\nJay R. Rooker,*$ Jay R. Rooker^1,2*

Igaratza Fraile³

Hui Liu¹

Noureddine Abid⁴

Michael A. Dance⁵

Tomoyuki Itoh⁶

Ai Kimoto⁶

Yohei Tsukahara⁶

Enrique Rodriguez-Marin⁷

Haritz Arrizabalaga³

¹Department of Marine Biology, Texas A&M University, Galveston, TX, United States
²Department of Wildlife and Fisheries Sciences, Texas A&M University, College Station, TX, United States
³Marine Research Division, AZTI, Pasaia, Spain
⁴Regional Centre of Tangier, National Institute of Fisheries Research, Tangier, Morocco
⁵Department of Oceanography and Coastal Sciences, Louisiana State University, Baton Rouge, LA, United States
⁶National Research Institute of Far Seas Fisheries, Japan Fisheries Research and Education Agency, Shizuoka, Japan
⁷Instituto Español de Oceanografía, Centro Oceanográfico de Santander, Santander, Spain

Uncertainty regarding the movement and population exchange of Atlantic bluefin tuna (Thunnus thynnus) from the two primary spawning areas (Gulf of Mexico, Mediterranean Sea) is increasingly implicated as a major impediment for the conservation of this species. Here, two mixture methods were applied to natural chemical markers (δ¹⁸O and δ¹³C) in otoliths (ear stones) to comprehensively investigate the nature and degree of transoceanic movement and mixing of eastern and western populations in several areas of the North Atlantic Ocean that potentially represent mixing hotspots. Areas investigated occurred on both sides of the 45°W management boundary as well as waters off the coast of Africa (Morocco, Canary Islands) where both populations are known to occur. Projections of population composition (i.e., natal or nursery origin) from a multinomial logistic regression (MLR) classification method with different probability thresholds were generally in agreement with maximum likelihood estimates from the commonly used mixed-population program HISEA; however, predicted contributions for the less abundant population were occasionally higher for MLR estimates. Both MLR and HISEA clearly showed that mixing of Atlantic bluefin tuna in the Central North Atlantic Ocean was highly variable from year to year with expatriates of eastern or western origin commonly crossing into the other management area. Pronounced transoceanic movement and mixing of western migrants was also present off the coast of Africa, with the occurrence of western migrants in the Canary Islands and Morocco ranging from zero to the majority of the individuals assayed for the years examined. Results indicate highly variable rates of movement and population exchange for Atlantic bluefin tuna, highlighting the need for temporally resolved estimates of natal origin in mixing hotspots to improve population models used to evaluate the status of this threatened species.

Introduction

Highly migratory pelagic species such as billfishes and tunas are capable of transoceanic crossings that often traverse international borders or jurisdictions (Ortiz et al., 2003; Madigan et al., 2014; Rooker et al., 2016, 2019). Many of these apex predators require multiple ecosystems and/or ocean basins to complete their life cycle (Rooker et al., 2007), with individuals commonly venturing outside of their respective management area (Block et al., 2011; Wilson et al., 2015; Harrison et al., 2018). Although fisheries managers concede that knowledge of a species’ migratory history is fundamental to conserving and rebuilding populations (Kerr and Goethel, 2014), data on the movement and population structure of many large marine predators are insufficient to integrate into assessments. Unfortunately, these data gaps are common for several highly migratory species (e.g., Atlantic and Pacific bluefin tuna [Thunnus thynnus, T. orientalis], yellowfin tuna [T. albacares], blue marlin [Makaira nigricans]), limiting our ability to adequately parameterize population models and develop effective conservation plans (Collette et al., 2011).

Atlantic bluefin tuna is the prototypical highly migratory species with individuals commonly crossing the Atlantic Ocean but also returning to natal areas to spawn (Rooker et al., 2008b, 2014; Graves et al., 2015). It is widely recognized that an improved understanding of the movement and exchange of Atlantic bluefin tuna from the two primary production zones—Gulf of Mexico (western) and Mediterranean Sea (eastern)—is critically needed for management purposes (Porch et al., 2000; Kerr et al., 2017). Chemical markers in otoliths (δ¹⁸O and δ¹³C) have been used to investigate the natal origin of this species throughout the Atlantic Ocean and, previous research indicates that individuals from both production zones readily cross the 45°W management boundary and intermingle with the other population (Rooker et al., 2008a, 2014; Siskey et al., 2016); however, the degree and persistence of population mixing is unresolved in several parts of their range. The composition of Atlantic bluefin tuna on both sides of the management boundary in the Central North Atlantic Ocean (CNAO) is of particular interest to resource managers because of the strong asymmetrical production between the two populations (i.e., an order of magnitude higher in the east; Secor, 2015). Areas in the Eastern North Atlantic Ocean (ENAO) off the coast of Africa also have been identified as putative mixing areas of eastern and western populations using otolith δ¹⁸O and δ¹³C values as indicators of natal origin (Rooker et al., 2014; Fraile et al., 2015); nevertheless, the prevalence and regularity of western-migrants in this region is unresolved.

Even though it is widely recognized that population models for Atlantic bluefin tuna are highly sensitive to mixing of individuals from eastern and western production zones (Kerr et al., 2017), incomplete data on the origin of this species in presumed mixing areas of the CNAO and ENAO has delayed the integration of movement and connectivity data into population assessments. This is due in part to the lack of consensus on the primary locations and extent of population exchange, and also the absence of data regarding the temporal persistence of mixing within these zones. In addition, a variety of statistical approaches, namely mixture models and discriminant functions, have been applied to otolith chemistry data to estimate the contribution of different stocks or nurseries (e.g., Campana et al., 1999; Elsdon and Gillanders, 2003; White et al., 2008; Smith and Campana, 2010; Niklitschek and Darnaude, 2016), and the accuracy and robustness of commonly used approaches for otolith δ¹⁸O and δ¹³C data on Atlantic bluefin tuna have yet to be comprehensively investigated. The aim of the present study was to apply two different types of classification methods—multinomial logical regression (MLR) and conditional maximum likelihood estimation (HISEA)—to otolith chemistry data to assess the origin of Atlantic bluefin tuna in four areas of the CNAO and ENAO. MLR was introduced as an alternative classification method for mixed-stock analysis of Atlantic bluefin tuna based on otolith δ¹⁸O and δ¹³C because model assumptions are more flexible than HISEA. Mixture models were applied to collections made over 4–5 consecutive years in each area to comprehensively assess interannual variability in transatlantic movement and mixing of this species.

Materials and Methods

Sample Collections

Our baseline sample used to predict natal origin was comprised of yearling Atlantic bluefin tuna (approximately 12–18 months old) from eastern and western nurseries collected from 1998 to 2011 from the Mediterranean Sea/Bay of Biscay (eastern) and the Gulf of Mexico/U.S. Atlantic Ocean (western). The baseline sample was first reported in Rooker et al. (2014) and shown to be stable temporally for the primary marker (δ¹⁸O) used to determine natal origin. Eastern (Mediterranean Sea) and western (Gulf of Mexico) nurseries included areas proximal to each, including the Bay of Biscay in the east and the eastern seaboard (Mid Atlantic Bight) of the United States in the west. The baseline sample used here is comprehensive in geographic scope and includes yearling Atlantic bluefin tuna from multiple locations in both spawning/nursery areas. Moreover, samples were collected in both nurseries for the majority of years between 1998 and 2011, suggesting the baseline sample effectively accounts for temporal variability and indicates that the extended baseline does not degrade nor bias estimates of nursery origin because differences between eastern and western nurseries are retained even when samples were drawn from different years. Yearling Atlantic bluefin tuna were collected with a variety of fishing gears including baitboats, hook and line, lighted purse seines, and longlines.

The current assessment of population mixing is based on sub-adult and adult Atlantic bluefin tuna collected from 2010 to 2016 on both sides of the management boundary in the CNAO and two locations in the ENAO where Rooker et al. (2014) detected western migrants: Morocco and Canary Islands (Figure 1). All collections were made under the provision of the ICCAT Atlantic Wide Research Program for Bluefin Tuna, and samples of Atlantic bluefin tuna from all areas were comprised of medium and large size category fish (25–100 kg and >100 kg, respectively), with the majority of samples from each area and year comprised principally of individuals near or over 100 kg; hereafter both size categories are referred to collectively as adults. All areas included samples of Atlantic bluefin tuna from 4 to 5 consecutive years to evaluate interannual variability in population mixing. While the current study was based predominantly on new data for Atlantic bluefin tuna collected from 2012 to 2016, otolith chemistry data for individuals collected between 2010 and 2011 from a previous study (Rooker et al., 2014) were included in the analysis to increase the extent of our time series. Estimates of natal origin were determined by comparing δ¹⁸O and δ¹³C values in the otolith cores of adults (core material represents otolith material deposited during the yearling period) to δ¹⁸O and δ¹³C values from otoliths of yearlings (i.e., baseline sample).

FIGURE 1

Figure 1. Summary of the four sampling areas (blue shading) for Atlantic bluefin tuna (Thunnus thynnus). Collections included adult Atlantic bluefin tuna from both sides of the management boundary (45°W longitude) in the Central North Atlantic Ocean (CNAO) and two locations in the Eastern North Atlantic Ocean (ENAO): Morocco and Canary Islands. Sampling conducted from 2010 to 2016.

Otolith Processing and Chemical Analysis

For chemical analysis, one sagittal otolith from each Atlantic bluefin tuna was cleaned and rinsed following Rooker et al. (2016). Otoliths were then embedded in Struers epoxy resin (Struers A/S, Ballerup, Denmark) and a 1.5 mm thick section was cut along a transverse plane using a Buehler IsoMet saw (Beuhler, Lake Bluff, Illinois). Sections were then attached to a sample plate using Crystalbond thermoplastic glue (SPI Supplies/Structure Probe Inc., West Chester, Pennsylvania), and the region of the otolith corresponding to the first year of life was isolated following Rooker et al. (2008a). Approximately 1.0–1.5 mg of powder was collected during the milling of each otolith and a fraction of the powder was loaded into sample trays for stable isotope analysis. Otolith δ¹⁸O and δ¹³C was determined using an automated carbonate preparation device (KIEL-III; Thermo Fisher Scientific, Inc., Waltham, Mass.) coupled to a gas chromatograph – isotope ratio mass spectrometer (Finnigan MAT 252; Thermo Fisher Scientific, Inc.) at the University of Arizona. Powdered otolith samples (ca. 0.5 μg) were reacted with dehydrated phosphoric acid under vacuum at 70°C. The isotope ratio measurement was calibrated based on repeated measurements of National Bureau of Standards (NBS), NBS-19 and NBS-18, with 6 standards ran for every 40 samples; precision was ± 0.08‰ (SD) and ± 0.11‰ (SD) for δ¹⁸O and δ¹³C, respectively. Otolith δ¹⁸O and δ¹³C values are reported relative to the Vienna Pee Dee Belemnite (VPDB) scale after comparison to an in-house laboratory standard calibrated to VPDB.

Classification Methods

The stock composition program HISEA is a commonly used statistical method for estimating the proportion of populations or stocks in a sample (Millar, 1990¹), and this approach has been used extensively for assignments of Atlantic bluefin tuna to eastern and western populations (e.g., Schloesser et al., 2010; Rooker et al., 2014; Fraile et al., 2015; Secor et al., 2015; Siskey et al., 2016). Direct maximum likelihood estimation of population composition in HISEA has two key assumptions. First, it assumes that individuals used to establish the baseline are randomly sampled and representative of each population in the fishery being investigated (Millar, 1987, 1990). Second, it assumes the probability of randomly selecting an individual from a mixed sample with measured characteristics (i.e., isotopic signature) is known from the baseline. As a result, the baseline sample must be representative of all individuals of each population within a mixed sample. For highly migratory species such as Atlantic bluefin tuna, baseline data are normally sampled opportunistically and/or intermittently from a variety of nursery locations, and mixing zones comprised of individuals of unknown origin are often far removed from spawning and/or natal sites. Therefore, individuals of unknown origin in a mixed sample may not be accurately estimated from the baseline, possibly affecting the efficacy of this method for delineating the composition of Atlantic bluefin tuna in areas where populations intermingle.

The MLR is a type of classification method that is increasingly used in machine learning, psychology, and economics (Hosmer et al., 2013). This modeling framework also shows promise for estimating contribution rates in mixed-population fisheries and this study is the first to apply this approach in this manner. This MLR framework allows for estimating the probability of occurrence of each category (i.e., population) as a function of explanatory variables (i.e., δ¹⁸O, δ¹³C) in question. While model parameters in MLR are also estimated using a maximum likelihood estimator, MLR is fundamentally different than HISEA in that it does not directly estimate the proportion of each population and does not assume a classification matrix where the probability of individuals being randomly selected from a specific population in the mixed sample is known from the baseline. In fact, MLR predicts the probability for individuals randomly sampled from the mixed sample and assigns them to a category and then estimates the proportions of each population. For each population in MLR, the response variable follows a multinomial distribution and the probability of each category of the dependent variable is predicted based on information from explanatory variables included in the model. MLR assumes (1) independence among the dependent variable choices, that is the choice of one category is not related to the choice of another category (i.e., the dependent variable), (2) non-perfect separation, that is the dependent variable cannot be perfectly predicted from the independent variables, and (3) multicollinearity of independent variables is relatively low. One potential advantage of MLR is that it can be applied to a mixed sample that may not include all potential sources or nurseries by establishing probability thresholds that assign individuals with lower predicted probabilities to an undetermined source that may not be represented in the baseline. Classification in MLR is based on a vector of conditional probability, and here we applied a series of probability thresholds for positive predictions with all modeling performed in R (R Core Team, 2018).

For all MLR and HISEA mixture models, the degree of uncertainty around estimated proportions was estimated with bootstrapping (500 simulations; error term equivalent to ± 1 SD around estimated proportion). MLR-based predictions were derived for four probability thresholds: 0.5, 0.6, 0.7, and 0.8. A probability threshold of 0.5—which was the only probability also used for HISEA—indicates that when an individual’s predicted probability for the western baseline sample is above 0.5, it will be classified to the western population. As the probability threshold increases above 0.5, some individuals do not have large enough predicted probabilities to be classified to either the eastern or western population. For such cases, they are assigned to an undetermined source or nursery, and the fraction of individuals failing to meet the probability threshold of either population will increase as the probability threshold increases (highest for 0.8).

Results

Otolith δ¹⁸O and δ¹³C of yearling Atlantic bluefin tuna from eastern and western nurseries from Rooker et al. (2014) were distinct for year classes investigated (MANOVA, p < 0.01) and indicated that these markers can be used to determine an individuals’ place of origin. Otolith δ¹⁸O of eastern yearlings were consistently heavier (less negative) than western yearlings and a demarcation at −0.9%, although some degree of overlap occurred for this marker from approximately −1.1 to −0.8 (Figure 2). Despite the overlap in otolith δ¹⁸O and a slight separation of otolith δ¹³C for eastern and western yearlings, cross-validated classification success remained relatively high (QDFA = 87%), supporting the use of otolith δ¹⁸O and δ¹³C as a baseline (reference signatures) for determining the origin of Atlantic bluefin tuna from presumed mixing zones in the CNAO and ENAO.

FIGURE 2

Figure 2. Otolith δ¹⁸O and δ¹³C values for Atlantic bluefin tuna of unknown origin collected from four locations from 2010 to 2016: (A) Central North Atlantic Ocean- West of 45°W, (B) Central North Atlantic Ocean- East of 45°W, (C) Eastern North Atlantic Ocean off Morocco, (D) Eastern North Atlantic Ocean off Canary Islands. Confidence ellipses shown are based on otolith δ¹⁸O and δ¹³C values of yearling Atlantic bluefin tuna (baseline) and represent 1 SD of the western (blue) and eastern (red) sample. Otolith δ¹⁸O and δ¹³C values of yearlings used to create confidence ellipses and for mixed-stock analysis from Rooker et al. (2014).

Central North Atlantic Ocean

Otolith core δ¹⁸O and δ¹³C values of adult Atlantic bluefin tuna run in HISEA and MLR modeling frameworks indicated that individuals collected in the CNAO readily crossed the 45°W management boundary with expatriates being relatively common both east and west of the boundary (Figure 2). For both HISEA and MLR (all probability thresholds) based assignments, the majority of adults collected west of 45°W from 2011 to 2014 were assigned to the western population, although contribution rates varied widely among years and mixture models (Figure 3). For HISEA, Atlantic bluefin tuna of western origin accounted for 13% to 96% of the samples collected west of 45°W over the four years investigated (Table 1). Although findings from MLR models were generally in agreement with HISEA, MLR-derived estimates of western contribution were noticeably reduced in 2012 across all probability thresholds. Conversely, the proportion of western migrants west of 45°W in the 2014 sample was greater for all MLR models relative to HISEA.

FIGURE 3

Figure 3. Estimated proportions of western (blue), eastern (red), and undetermined (gray) Atlantic bluefin tuna in samples collected from 2011 to 2014 in the Central North Atlantic Ocean—West of 45°W. Results are based on mixed-population analysis using HISEA and multinomial logistic regression (MLR) frameworks. H.5 denotes predictions from HISEA while MLR estimates are shown as M.5, M.6, M.7, M.8 to represent the four different probability thresholds (0.5–0.8). Uncertainty around predictions shown with error bars which represent 1 ± SD around estimated proportion.

TABLE 1

Table 1. Estimated proportion of “western migrants” detected in samples of Atlantic bluefin tuna from HISEA and multinomial logistic regression (MLR) models.

Wide-ranging interannual variability in otolith core δ¹⁸O and δ¹³C values was observed for Atlantic bluefin tuna collected east of 45°W in the CNAO (Figure 2), with contributions of western migrants variable across the years investigated. The presence of western migrants in this mixing zone was generally higher for MLR (all probability thresholds) relative to HISEA (Figure 4). For samples collected from 2011 to 2015, HISEA indicated that Atlantic bluefin tuna of western origin were consistently in the minority east of 45°W, often accounting for less than 10% of the sample in certain years (2011, 2014) (Table 1); however, nearly half of the sample from 2013 was determined to be of western origin, suggesting that large numbers of western migrants may occasionally enter the eastern management zone. MLR models showed less year-to-year variability than HISEA in the population composition of Atlantic bluefin tuna but also highlight the potential for meaningful temporal shifts in the population structure within this mixing zone. MLR (0.5 probability) showed the lowest amount of variability among the 5 years investigated with western contribution ranging from 21 to 52% (mean: 35%). Increasing MLR probability led to a greater proportion of samples assigned to the undetermined category (Figure 4) and generally resulted in minor reductions in western contribution rates (ca. < 5%).

FIGURE 4

Figure 4. Estimated proportions of western (blue), eastern (red), and undetermined (gray) Atlantic bluefin tuna in samples collected from 2010 to 2014 in the Central North Atlantic Ocean—East of 45°W. Results are based on mixed-population analysis using HISEA and multinomial logistic regression (MLR) frameworks. H.5 denotes predictions from HISEA while MLR estimates are shown as M.5, M.6, M.7, M.8 to represent the four different probability thresholds (0.5–0.8). Uncertainty around predictions shown with error bars which represent 1 ± SD around estimated proportion.

Eastern North Atlantic Ocean

Atlantic bluefin tuna of western origin were clearly evident in the ENAO, and these transatlantic migrants were present in multiple years in the Moroccan trap fishery, possibly indicating that these waters occur within a common migratory pathway of the western population (Figure 2). Of particular note was highly variable mixing of western Atlantic bluefin tuna in Morocco, with all mixture models showing fluctuating rates and the strong contribution of western migrants was only observed in two of the 5 years examined (Figure 5). In both 2011 and 2014, approximately 60–75% of the sample from Morocco was predicted to be of western origin, and this finding was consistent between HISEA and among all probability thresholds in MLR models (Table 1). Interestingly, the presence of western migrants in the Moroccan sample for the other 3 years investigated was negligible (0 to 2%; 2012, 2013) or minor (15%; 2015) in HISEA. Population contributions from MLR models were similar to HISEA for the 2012 sample; however, the proportion of western migrants was consistently higher (ca. 15–30% higher in MLR models for 2013 and 2015 samples), signifying that MLR may elevate the contribution of the less abundant, western population in this mixing zone.

FIGURE 5

Figure 5. Estimated proportions of western (blue), eastern (red), and undetermined (gray) Atlantic bluefin tuna in samples collected from 2011 to 2015 off Morocco. Results are based on mixed-population analysis using HISEA and multinomial logistic regression (MLR) frameworks. H.5 denotes predictions from HISEA while MLR estimates are shown as M.5, M.6, M.7, M.8 to represent the four different probability thresholds (0.5–0.8). Uncertainty around predictions shown with error bars which represent 1 ± SD around estimated proportion.

Population mixing of Atlantic bluefin tuna in waters off the Canary Islands also varied temporally with consistent but low contributions of the western population to this fishery in three of the 4 years investigated (Figure 6). In 2013, 2015, and 2016 samples, the western contribution from HISEA ranged from 12 to 21%. Again, MLR estimates for the same years at most probability thresholds were consistently 10–20% higher than HISEA, often with approximately a third of the sample attributed to the western nursery (Table 1). No western migrants were detected in the 2014 sample from the Canary Islands for HISEA or for MLR at probability thresholds of 0.7 and 0.8, although western fish contributed marginally at probability thresholds of 0.5 and 0.6.

FIGURE 6

Figure 6. Estimated proportions of western (blue), eastern (red), and undetermined (gray) Atlantic bluefin tuna in samples collected from 2013 to 2016 in the Canary Islands. Results are based on mixed-population analysis using HISEA and multinomial logistic regression (MLR) frameworks. H.5 denotes predictions from HISEA while MLR estimates are shown as M.5, M.6, M.7, M.8 to represent the four different probability thresholds (0.5–0.8). Uncertainty around predictions shown with error bars which represent 1 ± SD around estimated proportion.

Discussion

The current study applied well-established chemical markers in otoliths to significantly expand our understanding of the movement and mixing of Atlantic bluefin tuna in areas identified as potentially important mixing hotspots. As documented in previous studies, otolith δ¹⁸O was primarily responsible for discriminating individuals from eastern and western nurseries (Rooker et al., 2008a, 2008b). Nursery-specific differences in chemical signatures of yearlings followed expected dissimilarities in seawater δ¹⁸O between eastern (Mediterranean Sea) and western (Gulf of Mexico/U.S. Atlantic) nurseries (LeGrande and Schmidt, 2006), with yearlings from the Mediterranean Sea displaying heavier otolith δ¹⁸O values than their western counterparts (Rooker et al., 2014).

Regardless of the modeling framework applied to δ¹⁸O and δ¹³C values derived from the milled otolith cores of adult Atlantic bluefin tuna, results from this study clearly showed that the 45°W management boundary that divides the two populations is readily crossed and high levels of mixing occur on both sides of this longitudinal demarcation. Our findings are in accord with earlier research using otolith δ¹⁸O and δ¹³C to demonstrate that both eastern and western migrants commonly ‘cross the line’ (Rooker et al., 2008b, 2014), and further corroborate an earlier finding that aggregations targeted by high seas longline fleets in the CNAO are consistently comprised of both populations on each side of the management boundary (Rooker et al., 2014). In response, there does not appear to be a biological basis for the construct of the current boundary used to delineate eastern and western populations, and persistent, high levels of mixing observed in both management zones in the CNAO further emphasize the importance of east-west exchange.

The presence of western migrants off Morocco and the Canary Islands in the ENAO indicates that movement extends beyond the CNAO, with individuals completing oceanic crossings and entering waters off the coast of Africa. Transoceanic migrations by Atlantic bluefin tuna have been previously reported using chemical markers in otoliths (Rooker et al., 2008b, 2014; Secor et al., 2015; Siskey et al., 2016), organochlorine tracers (Dickhut et al., 2009; Graves et al., 2015), and electronic tags (Block et al., 2005; Kurota et al., 2009; Galuardi et al., 2010). However, crossing of western migrants into the ENAO is generally assumed to be very low and expatriates from the west are even more rare inside the Mediterranean Sea. Rooker et al. (2014) evaluated the origin of Atlantic bluefin tuna from the ENAO and detected no western migrants from collections west of the Strait of Gibraltar in waters off Portugal; however, a small fraction of individuals from a trap fishery off Morocco was reported to be comprised of western migrants (6%). Similarly, Rodríguez-Ezpeleta et al. (in press) documented western migrants in fisheries from several areas of the ENAO including Morocco and the Canary Islands using single nucleotide polymorphism (SNP) markers. Still, western contribution estimates from otolith δ¹⁸O and δ¹³C were generally higher than genetic assignments, particularly for MLR-based estimates, possibly leading to overestimates of western contribution. In addition, Fraile et al. (2015) examined samples collected over 3 years from the bait boat fleet in the Bay of Biscay, and western migrants were only detected in a single year and at relatively low levels for juveniles and adults (2–5%). Here, we advance our understanding of transoceanic movements of Atlantic bluefin tuna by demonstrating that western migrants commonly enter waters off the coast of Africa year after year and mix with the eastern population. Both HISEA and MLR based estimates indicated that significant numbers of western migrants enter waters off Morocco (e.g., 2011 and 2014) and to a lesser degree the Canary Islands; nevertheless, it is important to note that observed contribution estimates may not represent the precise population composition in either fishery because samples sizes were limited in certain years and occasionally from abbreviated collection periods. In addition, early egress out of the Mediterranean Sea and into the North Atlantic Ocean by age-0 Atlantic bluefin tuna may produce chemical signatures that bear some resemblance to the western baseline (e.g., more depleted otolith δ¹⁸O), potentially overstating the contribution of western migrants in samples from Morocco and the Canary Islands.

It is also important to note that western migrants were detected in all 3 years with samples available from both Morocco and the Canary Islands regardless of the modeling framework used. Given the general proximity of both locations in the ENAO, it is possible that western migrants detected in these two locations were the result of shared migration pathway(s) across the Atlantic Ocean. Common transoceanic or inter-oceanic migration routes have been reported for other highly migratory species including southern bluefin (T. maccoyii), Pacific bluefin tuna (T. orientalis), and several other large marine vertebrates (Boustany et al., 2010; Block et al., 2011; Hobday et al., 2015). Shared population-specific migration routes across ocean basins are often correlated with dynamic ocean features, particularly biological and/or physicochemical fronts such as the Transition Zone Chlorophyll Front (TZCF). This feature is assumed to define the migration routes of highly migratory species such as tunas, sea turtles, and squid that are known to forage along this basin-wide frontal feature in the Pacific Ocean (Polovina et al., 2017) and similar dynamic features in the Atlantic Ocean (e.g., subtropical gyre front) may represent ecological bridges (sensu Briscoe et al., 2017), leading to common migration routes for Atlantic bluefin tuna to waters off Morocco and the Canary Islands.

Spatial shifts in the distribution of marine fishes occurs at multiple spatial scales ranging from small-scale habitat shifts (Moulton et al., 2017; Landry et al., 2018; Rooker et al., 2018) to large-scale journeys across ocean basins (Galuardi et al., 2010; Block et al., 2011; Kraus et al., 2011). Migratory patterns of pelagic fishes and other marine vertebrates at the basin scale (i.e., transoceanic) are often influenced by both intrinsic and external factors, including spatial shifts in prey quantity/quality, thermal tolerance/sea surface temperatures, and other incentives related to reproduction and/or energetics (Secor, 2015). In response, it is plausible to assume that movements across boundaries or oceans reported here for Atlantic bluefin tuna are physiologically and/or ecologically motivated behaviors to locate favorable physicochemical conditions and optimize foraging. Even though oceanic conditions or features are inherently dynamic in nature and known to alter the timing and extent of movements of pelagic fishes (Block et al., 2011; Sequeira et al., 2018), the common occurrence of western migrants in the mixing zones investigated suggests that this species may exploit certain areas of the CNAO and ENAO in fairly predictable ways and share common drivers that initiate these movements.

Significant temporal variation in the proportion of eastern and western migrants in the mixing zones investigated indicates that basin-scale movements of Atlantic bluefin tuna are highly variable from year to year, and our findings are in accord with earlier studies that revealed the potential for interannual changes in large-scale movements of this species. Sibert et al. (2006) contrasted movement parameters of Atlantic bluefin tuna with satellite tags and reported strong interannual variation in movement parameters of this species, albeit their assessment was based on a limited number of specimens and tagging years. Mather (1980) also reported shifts in transoceanic migration behaviors for different size classes of Atlantic bluefin tuna and hypothesized that changes in migration patterns were linked to ocean-scale events. Temporal shifts in the distribution and migration of pelagic fishes have been related to large-scale oceanic phenomena such as the Atlantic Multidecadal Oscillation (AMO) and El Niño–Southern Oscillation (ENSO), which alter the physical environment (e.g., sea surface temperatures, upwelling, circulation) and primary/secondary productivity, often times leading to changes in prey availability, the distribution of predators, and recruitment success (Otterson et al., 2001; Stenseth et al., 2002; Harford et al., 2017). While AMO is generally considered to fluctuate at the decadal scale, interannual variability is also associated with this phenomenon (Enfield et al., 2001). Alheit et al. (2014) observed that the migration and population dynamics of small pelagic fishes (i.e., forage fish of tunas) varied in synchrony with different AMO phases, and recently Faillettaz et al. (2019) demonstrated that AMO phases also influence the basin-scale distribution of bluefin tuna in the North Atlantic Ocean. Therefore, interannual changes in the AMO index during the years investigated (2011–2016²) or other large-scale oceanic/weather phenomena that alter environmental conditions may have contributed in part to changes in movement and population mixing of Atlantic bluefin tuna observed in the CNAO and ENAO.

It is also possible that year to year variation in the movement and exchange of Atlantic bluefin tuna was due in part to random exploratory and/or physiologically motivated movements to optimize feeding (e.g., dynamic optimal foraging theory) or energetics, which may or may not be linked to ocean-scale events described above (Giuggioli and Bartumeus, 2010; Hays et al., 2016; Thygesen et al., 2016). Optimal foraging theory predicts that movements by this species should occur in a manner to most efficiently use resources (i.e., prey) that are inherently spatially heterogeneous throughout the North Atlantic Ocean. Mariani et al. (2016) demonstrated the link to foraging for Atlantic bluefin tuna using a game theoretical approach and reported that migrants benefit from highly productive environments with individuals often visiting a series of foraging hotspots connected by a series of migratory pathways. Similarly, energy intake and the spatial distribution of the forage base was shown to influence migration patterns for the Pacific congener T. orientalis; however, movements of this species were not always consistent with maximizing energy and physiological constraints linked to water temperature were also important drivers of movement—both in its timing and degree (Whitlock et al., 2015).

The development and application of numerical approaches to improve our understanding of the movement and exchange of eastern and western populations of Atlantic bluefin tuna is critical to the management of this species. To date, a variety of mixture models have been applied to genetic and chemical markers to predict the contribution of eastern and western recruits to different fisheries in the Atlantic Ocean (e.g., Rooker et al., 2014; Fraile et al., 2015; Puncher et al., 2018). Here, we show that a commonly applied mixture model (HISEA) provides comparable estimates to a new MLR-based approach; however, MLR-derived proportions of western migrants east of the management boundary and into the ENAO were occasionally higher than estimates from HISEA, particularly at the 0.5 probability threshold. This suggests that MLR identified more of the less dominant population in these mixing zones, which may lead to under- (HISEA) or over-estimates (MLR) of western contribution depending upon the modeling framework being adopted. Hanke et al. (2016) applied linear discriminant, quadratic discriminant (equivalent to direct MLE in HISEA), and randomForest classifiers to assess mixed samples of this species in U.S. Atlantic fisheries, and the authors noted that predictions of natal origin were comparable between the latter two classifiers, with both classifying fewer fish from the less abundant population—in this case eastern Atlantic bluefin tuna —in a mixed sample. While results from Hanke et al. (2016) indicated that conditional MLE (i.e., HISEA) and randomForest classifiers produced comparable estimates, it again raises the question of whether mixed-population approaches commonly used accurately represent the proportion of less abundant population(s) in a mixed sample. Here, the contribution of western migrants or the less dominant population in a mixed sample was often higher for MLR-based classifications for samples from Morocco and the Canary Islands. Moreover, this result did not change meaningfully for MLR models with higher probability thresholds (i.e., less uncertainty in eastern or western assignment), indicating that HISEA may underestimate the contribution of the western population in mixed samples in the ENAO. Increasing probability thresholds in MLR also led to another interesting observation; namely, a significant fraction of the sample was characterized as undetermined and not assigned to either the eastern or western population. While this may be a function of natural overlap in otolith δ¹⁸O and δ¹³C values for yearlings from eastern and western nurseries, it may also indicate that individuals from outside sources or nurseries may be present in the mixed samples examined.

An improved understanding of the migration and mixing of western migrants into the eastern management zone has important implications because of the depressed abundance and recruitment of the western population (Standing Committee on Research and Statistics [SCRS], 2017). The sustained presence of Atlantic bluefin tuna from the smaller, western population in waters east of the management boundary intensifies the fishing pressure on this vulnerable population, and therefore highlights the need for spatially and temporally resolved estimates of population composition throughout their range. Natural chemical markers in otoliths were used here to elucidate the nature and extent of population exchange in mixing hotspots in the Atlantic Ocean, and extensive and persistent east-west exchange observed here emphasizes the need to further refine and incorporate mixing into future assessments. We also provide contribution estimates from two useful mixed-population frameworks (HISEA and MLR) with a series of probability thresholds for the latter to better understand the range of conceivable mixing scenarios. In turn, this can then be used in virtual population analysis and simulation modeling to explore sensitivity of yields and population dynamics against different mixing scenarios (Kerr et al., 2017; Morse et al., 2018). However, it is important to note that an improved understanding of population mixing and the environmental drivers responsible for trans-Atlantic movement will require more comprehensive time-space coverage of collections that accurately represent the population composition in a given area. Perceptions of stock mixing from otolith chemistry and other natural markers (genetics) are almost certainly biased to some degree by sampling limitations (e.g., collections restricted in time and space from a region) that only provide a “snapshot” of the stock composition, which may not accurately reflect east-west exchange.

Data Availability

Data supporting the conclusions of this manuscript will be made available by the authors, without undue reservation, to any qualified researcher.

Author Contributions

JR, IF, NA, TI, AK, YT, ER-M, and HA contributed to the design and collection of the biological samples. JR, HL, IF, and MD processed and analyzed the data, and performed the mixed-stock analyses. JR wrote the manuscript.

Funding

This study was funded in part by the European Community (Grant SI2/542789), Canada, Croatia, Japan, Norway, Turkey, United States (NMFS NA11NMF4720107), Chinese Taipei, and the ICCAT Secretariat. The study was also supported in part by the McDaniel Charitable Foundation and by grant CTM2011-27505 from the Spanish Ministry of Economy and Competitiveness (MINECO).

Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Acknowledgments

This work was carried out under the provision of the ICCAT Atlantic Wide Research Program for Bluefin Tuna (GBYP). The contents of this manuscript do not necessarily reflect the point of view of ICCAT or of the other funders.

Footnotes

References

Alheit, J., Licandro, P., Coombs, S., Garcia, A., Giráldez, A., Garcia Santamaría, M. T., et al. (2014). Atlantic Multidecadal Oscillation (AMO) modulates dynamics of small pelagic fishes and ecosystem regime shifts in the eastern North and Central Atlantic. J. Mar. Sys. 131, 21–35. doi: 10.1016/j.jmarsys.2013.11.002