Trophic ecology of meso- and bathypelagic predatory fishes in the Gulf of Mexico

Abstract

The trophic ecology of eight circumglobal meso- and bathypelagic fishes (Anoplogaster cornuta, Chauliodus sloani, Coccorella atlantica, Gigantura chuni, G. indica, Omosudis lowii, Photostomias guernei, and Stomias affinis) with contrasting vertical migration habits (vertical migrators vs. non-migrators) were examined using stable isotope analysis (SIA). Mean δ¹³C values of these predators were similar among species, ranging from –18.17 to –18.99 ‰, suggesting that all species are supported by a similar carbon source. This finding was supported by mixing-model analysis; all of these deep-living predators received the majority (>73%) of their carbon from epipelagic food resources. Mean δ¹⁵N values of the predators ranged from 9.18 to 11.13 ‰, resulting in trophic position estimates between the third and fourth trophic level, although significant shifts in δ¹⁵N with increasing body size suggest that some of these species undergo ontogenetic shifts in trophic position. Bayesian standard ellipses, used to estimate isotopic niche areas, differed in size among species, with those occupying the highest relative trophic positions possessing the largest isotopic niches. These results, which provide the first trophic descriptions using dietary tracers for several of these species, offer insight into the trophic structure of deep-sea ecosystems and will help inform the construction of ecosystem-based models.

bathypelagic, epipelagic, feeding ecology, mesopelagic, mixing models, stable isotopes

Introduction

The deep-pelagic zone (waters deeper than 200 m to just above the seabed) represents the largest cumulative habitat on earth and is home to a diverse array of specialized fauna adapted to its abiotic and biotic conditions (Angel, 1997; Robison, 2004, 2009). The deep sea and its inhabitants provide an array of ecosystem services that are important to humans, including carbon sequestration, nutrient regeneration, fisheries production, and waste absorption (Danovaro et al., 2008; Mengerink et al., 2014; Thurber et al., 2014). Despite its enormous volume and the economic and ecological importance of its fauna, deep-pelagic ecosystems remain chronically understudied (Webb et al., 2010) and face an increasing number of stressors including climate change, ocean acidification, overfishing, and natural resource extraction (Morato et al., 2006; Ramirez-Llodra et al., 2011; Mengerink et al., 2014). As threats to the diversity and stability of marine ecosystems increase and expand into deeper oceanic environments, there has been increasing concern regarding the status of deep-sea communities and a renewed interest in describing and understanding deep-sea ecosystem structure.

Central to our understanding of ecosystem and community structure is a thorough knowledge of foodwebs (Polis and Strong, 1996; McCann, 2000). In addition to providing important information regarding ecosystem functioning, the study of foodwebs provides understanding of how animal communities are structured and sheds light on the mechanisms underlying species coexistence and persistence. While our knowledge of deep-pelagic foodwebs has advanced considerably over the past few decades (Robison, 2009; Sutton, 2013), fundamental information in many regions, including species-specific feeding relationships, trophic position estimates, and delineations of energy pathways connecting disparate trophic levels and communities, is lacking (Mengerink et al., 2014; Drazen and Sutton, 2017).

Fishes are a dominant component of deep-pelagic ecosystems worldwide and are among the main taxa that undertake diel vertical migrations (DVM). While the standardized abundance (no. per unit volume) of meso- and bathypelagic fishes is relatively low (Angel and Baker, 1982), their global distributions have resulted in high cumulative biomass estimated at 7–10 billion tonnes (Gjøsaeter and Kawaguchi, 1980; Irigoien et al., 2014). Due to their sheer numbers and vertical migration behaviour, which can exceed 1000 m in vertical extent, it is increasingly being recognized that fishes play key ecological and biogeochemical roles in open-ocean ecosystems (Wilson et al., 2009; Drazen and Sutton, 2017). As highly abundant mid-level consumers, deep-pelagic fishes help regulate zooplankton populations (Hopkins and Gartner, 1992; Pakhomov et al., 1996). Deep-pelagic fishes also serve as trophic links between zooplankton and higher-order consumers such as epipelagic fishes (Moteki et al., 2001; Choy et al., 2013), marine mammals (Pauly et al., 1998), and seabirds (Raclot et al., 1998; Cherel et al., 2008). DVM of fishes have been shown to connect the epi-, meso-, and bathypelagic habitats with each other and with deep-benthic habitats (Porteiro and Sutton, 2007; Trueman et al., 2014).

Stable isotope analysis (SIA) has been widely used to delineate foodweb structure and provides an integrated view of an organism’s diet over time-scales relevant to tissue turnover rates rather than digestion rates (Peterson and Fry, 1987; Post, 2002). Carbon isotopes undergo relatively small amounts of fractionation during trophic transfers and are useful for determining the relative contributions of carbon sources to the production of consumers (Peterson and Fry, 1987). Stable isotopes of nitrogen undergo comparatively large levels of fractionation (∼3–5 ‰) during trophic transfer, resulting in predictable differences in the isotopic signatures of consumers and their prey (Post, 2002; Hussey et al., 2014). The relatively predictable level of enrichment of ¹⁵N during trophic transfer allows for the determination of trophic levels and can be used to identify trophic relationships within assemblages of organisms (Peterson and Fry, 1987; Post, 2002).

To date, much of the research describing the trophic ecology of deep-pelagic fishes has focused on zooplanktivorous groups (myctophids, sternoptychids, gonostomatids), while less attention has been paid to micronektonivores (stomiids, alepisauroids) that occupy higher trophic levels. The numerical importance of micronektonivores (Hopkins et al., 1996; Sutton and Hopkins, 1996a), their propensity to prey heavily on zooplanktivorous fishes (Clarke, 1982; Hopkins et al., 1996; Sutton and Hopkins, 1996b), and documented importance as prey for higher trophic level consumers (Moteki et al., 2001; Choy et al., 2013) provides the rationale for further describing their trophic dynamics. Here, we describe the trophic ecology of eight putative high-trophic-level fishes: Anoplogaster cornuta, Chauliodus sloani, Coccorella atlantica, Gigantra chuni, G. indica, Omosudis lowii, Photostomias guernei, and Stomias affinis. These species are meso- and bathypelagic fishes with circumglobal distributions, some of which have been documented as numerically important components of deep-pelagic assemblages (Sutton and Hopkins, 1996a; Moore et al., 2003; Sutton et al., 2008). Specific goals of this study are to provide estimates of trophic position, describe the isotopic niche areas and the extent of niche overlap among species, detail ontogenetic shifts in trophic position, and quantify the relative carbon contributions of particulate organic matter (POM) from the epi-, meso-, and bathypelagic zones to each of these species.

Material and methods

Sample collection and study site

Fishes were collected from the northern Gulf of Mexico (GOM) during three oceanographic cruises conducted during 2011 in spring (22 March–11 April), summer (23 June–13 July), and fall (8–27 September). All cruises were part of the Offshore Nekton Sampling and Analysis Program (ONSAP) that was implemented following the Deepwater Horizon oil spill in support of NOAA’s Natural Resource Damage Assessment (NRDA). ONSAP stations are the same as stations currently used by the long-term Southeast Area Monitoring and Assessment Program (SEAMAP) and are situated every half degree of longitude and latitude in the northern GOM (Figure 1). Specimens were collected using large midwater trawls fitted with large-mesh panels (∼80 cm) near the mouth that gradually tapered to smaller mesh (∼6 cm) sizes before the codend. Trawls were fished obliquely from the surface to depths of either 700 or 1400 m. Once the trawls were retrieved, animals were sorted, enumerated, and visually identified to species. Samples for SIA were selected haphazardly in an effort to maximize spatial and temporal coverage. All specimens for SIA were frozen whole at –20°C until processed at Texas A&M University at Galveston.

Figure 1.

Map of ONSAP sampling grid, locations of POM samples, and locations where fishes were collected for SIA (specimen number denoted by circle diameter) in the GOM.

Stable isotope analysis

SIA was conducted on 212 specimens, with sample sizes of each species ranging from 19 to 37 individuals (Table 1). White muscle tissue for SIA was dissected from the dorsal musculature of fishes and visually examined under a dissecting microscope for the presence of bones, which were subsequently removed. Cleaned samples were rinsed with deionized water, frozen, and lyophilized for 48 h. Freeze-dried samples were homogenized using a mortar and pestle, weighed, wrapped in tin capsules, and shipped to the Stable Isotope Facility at the University of California Davis for analysis. Analysis of muscle tissue δ¹³C and δ¹⁵N was carried out using an elemental analyser (PDZ Europa ANCA-GSL) coupled with an isotope ratio mass spectrometer (PDZ Europa 20-20). The long-term standard deviation of the facility at UC Davis is 0.2 ‰ for δ¹³C and 0.3 ‰ for δ¹⁵N. Stable isotope data are expressed relative to international standards of Vienna PeeDee belemnite and atmospheric N₂ for carbon and nitrogen, respectively. The C: N of fishes in this study were low (species mean C: N range 3.31–3.86; 92% of individuals C: N < 4.0) compared with C: N from similar species collected in the Atlantic and Southern oceans (C: N 3.3–12.5; Hoffman and Sutton, 2010), suggesting that lipids did not significantly confound the interpretation of δ¹³C data. Therefore, all statistical analyses were performed on uncorrected δ¹³C values.

Table 1.

Species-specific sample descriptions and bulk δ¹³C and δ¹⁵N isotope data (mean ± SD).

Species	n	Spring	Summer	Fall	Standard length range (mm)	Mean standard length (mm) ± SD	δ13C (‰) ± SD	δ15N (‰) ± SD	C: Nbulk ± SD
A. cornuta1	23	2	12	9	84–148	114.35 ± 19.09	–18.93 ± 0.67	11.14 ± 0.96	3.66 ± 0.44
C. sloani2	30	10	20	0	143–237	191.43 ± 23.32	–18.68 ± 0.43	9.51 ± 0.42	3.42 ± 0.23
C. atlantica2	19	0	19	0	44–125	89.53 ± 26.53	–18.50 ± 0.47	9.96 ± 0.83	3.67 ± 0.34
G. chuni1	24	6	9	9	34–186	134.57 ± 40.11	–18.25 ± 0.44	11.13 ± 1.08	3.31 ± 0.15
G. indica1	21	8	6	7	75–192	141.95 ± 28.83	–18.25 ± 0.94	10.70 ± 0.64	3.43 ± 0.29
O. lowii1	32	12	10	10	36–261	120.66 ± 61.60	–19.04 ± 0.32	9.79 ± 0.60	3.40 ± 0.08
P. guernei2	37	14	13	10	56–127	98.08 ± 14.03	–18.61 ± 0.40	9.18 ± 0.63	3.44 ± 0.13
S. affinis2	26	5	16	5	55–205	127.31 ± 38.66	–19.38 ± 0.84	9.98 ± 0.89	3.86 ± 0.63

¹Denotes no DVM, ²denotes asynchronous DVM (not all individuals of population migrate vertically each day). References for vertical migration patterns: A. cornuta (Clarke and Wagner, 1976), C. sloani (Sutton and Hopkins, 1996a), C. atlantica (McEachran and Fechhelm, 1998), G. chuni (McEachran and Fechhelm, 1998), G. indica (Sutton et al., 2010), O. lowii (McEachran and Fechhelm, 1998; Sutton et al., 2010), P. guernei (Sutton and Hopkins, 1996a), S. affinis (Sutton and Hopkins, 1996a).

The stable isotope data of POM used in this study are derived from the published dataset of Fernández-Carrera et al. (2016). For detailed descriptions of methodologies and sample locations, see Fernández-Carrera et al. (2016), but a brief description of methodologies follows. POM samples were collected during summer 2011 (2–21 July) in the northern GOM. In addition to samples collected in pelagic waters, the complete published dataset included samples taken from waters over the continental shelf and from waters in close proximity to the Mississippi River. In order to maximize the spatial overlap between POM samples and the collection locations of fishes, only POM data collected within close proximity to ONSAP sampling stations in waters ≥1000 m deep (Figure 1) were utilized. POM samples were collected throughout the water column using remotely fired 10-l Niskin bottles. Samples were then filtered across 47-mm glass fibre filters at low pressure and dried at 60°C for 24 h prior to isotope analysis (Fernández-Carrera et al., 2016). In order to determine if the isotopic signatures of POM samples changed with depth, we used collection depth to designate POM samples as epipelagic (0–200 m), mesopelagic (200–1000 m), or bathypelagic (>1000 m) so that statistical comparisons could be made.

Data analysis

Multivariate analysis of variance (MANOVA) was used to test for differences in δ¹³C and δ¹⁵N among species and POM depth zones. Species and season were included as factors in the linear model and tested for the presence of an interaction. If significant differences were found, univariate tests for both δ¹³C and δ¹⁵N were performed using analysis of variance among fish species and POM depth zones. A posteriori differences among means were detected using Tukey’s honestly significant difference (HSD) test. Using equation 4 from Post et al. (2007), trophic position was calculated for each species:

$${\text{TrL}}_{\mathrm{i}}=\left[\left({\delta }^{15}{\mathrm{N}}_{\mathrm{i}}–{\delta }^{15}{\mathrm{N}}_{\text{base}}\right)/{\Delta }^{15}\mathrm{N}\right]+\lambda$$

(1)

where δ¹⁵N_i is the mean species δ¹⁵N, δ¹⁵N_base is the mean δ¹⁵N of the primary producer or primary consumer being used to set the isotopic baseline, Δ¹⁵N is the trophic discrimination factor, and λ represents the trophic level of the organism being used to set the baseline. Because primary consumer data were not available for the time-period of this study, trophic position estimates made using mean δ¹⁵N values of POM collected from the epipelagic zone were compared with estimates calculated from published δ¹⁵N values of a group of primary consumers (euphausiids) collected in the pelagic northern GOM during 2007 (McClain-Counts et al., 2017). In order to explore the relationship between fish size and δ¹³C and δ¹⁵N, least-squares linear regression analysis was conducted for each species. Spatial variation in δ¹³C and δ¹⁵N of both fishes and POM was investigated using least-squares linear regression between isotopic values and longitude and latitude (0.5° intervals). Because every species was not collected at every sampling location, isotope data were pooled across species within each line of longitude and latitude (more than one site along each 0.5° of longitude). Additionally, because not all species were collected across a range of latitudes and longitudes within each season, the effect of season on the spatial relationships of the isotope data was not explored. All statistical analyses were performed in R (R Development Core Team, 2016) v. 3.3.2.

The trophic breadth of each species and trophic similarity among species were assessed by calculating standard ellipse areas (SEA) using the R package SIBER (Jackson et al., 2011) and following methods outlined by Jackson et al. (2011). Bayesian standard ellipses encompass ∼40% of the isotope data for each species are less affected by increases in sample size or statistical outliers than convex hull analysis, and represent the core isotopic niche area of a species (Jackson et al., 2011). Size-corrected SEAs (SEA_c) were calculated for each species, which adjusts for underestimation of ellipse area at small sample sizes and allows for comparison of ellipse sizes to other studies (Jackson et al., 2011). Overlap of size-corrected ellipses was used as a proxy for trophic similarity and was examined by calculating the extent of overlap between each pairwise combination of species. The percentage of overlap between species pairs was calculated by dividing the area of overlap (‰²) by the total combined ellipse area (‰²) of the two species being compared. Isotopic niche overlap was considered significant when overlap between two species was >50%. Differences in size-corrected ellipse area, a proxy for trophic breadth that assumes species with larger SEAc feed more broadly within the foodweb than those with smaller SEAc were compared among species and considered to be significantly different when 95% of posterior draws were smaller in one species compared with the other.

The Bayesian mixing model, MixSIAR (Stock and Semmens, 2015), was used to estimate the relative contribution of epi- (0–200 m), meso- (200–1000 m), and bathypelagic (>1000 m) POM to each species. Bayesian mixing models provide the most accurate estimations of source or prey contributions when tissue and species-specific discrimination factors are used (Caut et al., 2008), but discrimination factors for meso- and bathypelagic fishes are currently unknown. We chose to run mixing models using discrimination factors of 3.15 ‰ ± 1.28 ‰ and 0.97 ‰ ± 1.08 ‰ for δ¹⁵N and δ¹³C, respectively (Sweeting et al., 2007a,b), which have been previously used to study the trophic structure of meso- and bathypelagic fishes (Valls et al., 2014). Mixing models in MixSIAR estimate probability density functions using Markov chain Monte Carlo methods, and each model was run with identical parameters (number of chains = 3; chain length = 100 000; burn in = 50 000; thin = 50). Model convergence was determined using Gelman-Rubin and Geweke diagnostic tests (Stock and Semmens, 2015).

Results

Stable isotopes

Individual consumer δ¹³C values ranged from –21.49 to –16.63 ‰, while mean δ¹³C values were similar among species, with a difference of 1.13 ‰ separating the most depleted (S. affinis: –19.38 ‰ ± 0.83) and most enriched species (G. chuni: –18.25 ‰ ± 0.44 and G. indica: –18.25 ‰ ± 0.94) (Table 1; Figure 2). Individual δ¹⁵N values varied between 7.10 and 13.07 ‰, with 1.96 ‰ separating the mean δ¹⁵N values of the most enriched (A. cornuta: 11.14 ‰ ± 0.96) and depleted species (P. guernei: 9.18 ‰ ± 0.63) (Table 1; Figure 2). Species-specific differences in δ¹³C and δ¹⁵N were significant (F_{14, 382} = 17.24, p < 0.001); however, no significant seasonal effects were found (F_14 382 = 1.29, p = 0.27), and no significant interaction effect among species and season was detected (F_22 382 = 1.05, p = 0.40). Significant differences in δ¹³C values among species (one-way ANOVA; F₇₂₀₄ = 11.62, p < 0.001) were driven by G. chuni and G. indica, which were enriched in ¹³C compared with more ¹³C-depleted species such as O. lowii and S. affinis (Figure 2). Significant differences in δ¹⁵N among species (one-way ANOVA; F₇₂₀₄ = 25.55, p < 0.001) were primarily driven by A. cornuta, G. chuni, and G. indica, which were enriched in ¹⁵N compared with C. sloani and P. guernei (Figure 2). Results of all pairwise comparisons for δ¹³C and δ¹⁵N values among species are listed in Supplementary Table S1.

Figure 2.

Isotope bi-plot of δ¹³C and δ¹⁵N values from POM (squares) and fishes (circles). Data points represent means and error bars represent ± 1 SD.

The δ¹³C values of fishes were significantly correlated with latitude (r = 0.08, p < 0.01) and longitude (r = 0.04, p < 0.01), while δ¹⁵N values were not (latitude p = 0.46; longitude p = 0.19). Due to limited spatial coverage within each species, spatial trends were tested by pooling all fish species together. Because spatial variation could not be tested within each species and due to the low correlation coefficients observed among fish δ¹³C values and latitude and longitude, isotope data for each species were pooled across lines of latitude and longitude during subsequent analysis.

A total of 154 samples of POM collected from depths ranging from 1 to 2500 m were utilized (Fernández-Carrera et al. 2016). POM exhibited a wide range of δ¹³C (–27.17 to –16.41) and δ¹⁵N values (–3.58 to 11.69), with POM samples generally becoming more ¹⁵N enriched with increasing depth (Figure 2). Significant differences in POM δ¹³C and δ¹⁵N among vertical depth zones (MANOVA: F₄₃₀₂ = 14.54, p < 0.001) were observed. Significant differences in δ¹⁵N were found among depth zones (ANOVA: F₂₁₅₁ = 34.41, p < 0.001), with epipelagic POM more ¹⁵N depleted than POM collected from mesopelagic and bathypelagic depths (p < 0.001). The δ¹³C values of POM did not significantly differ across depth zones (F₂₁₅₁ = 0.42, p = 0.66). Latitudinal and longitudinal variation in POM δ¹³C and δ¹⁵N was minimal, with the only significant correlation occurring between epipelagic POM δ¹³C and longitude, although correlation coefficients were low (r = 0.04, p = 0.043). All other pairwise combinations between δ¹³C and δ¹⁵N and latitude or longitude within the epi-, meso-, and bathypelagic depth zones were non-significant (Supplementary Table S2).

Trophic position estimates

The use of primary producers or primary consumers to set the isotopic baseline had no effect on the relative trophic positions among consumers, but resulted in slight differences (0.32 TL) in calculated trophic levels. When primary producer (POM) data were used to set the baseline, consumer TPs ranged from 2.8 (P. guernei) to 3.4 (A. cornuta, G. chuni), while all species fell within the third and fourth trophic levels when primary consumers were used to set the baseline (P. guernei = 3.1, A. cornuta and G. chuni = 3.7) (Figure 3; Supplementary Table S3).

Figure 3.

Trophic level estimates calculated using δ¹⁵N data of each species. Letters represent significant differences in TL among species, with like letters being similar and non-like letters significantly different. Dashed lines represent the δ ¹⁵N threshold values of TL 3 and TL 4 when using primary consumers (euphausiids) to set the isotopic baseline; dotted lines represent the δ ¹⁵N threshold values of TL 3 and TL 4 when using primary producers (POM) to establish isotopic baseline. For species-specific TP estimates (± SD), see Supplementary Table S3.

Of the species examined, A. cornuta (r = 0.63, p < 0.001), C. atlantica (r = 0.74, p < 0.001), G. chuni (r = 0.41, p < 0.001), C. sloani (r = 0.22, p < 0.001), P. guernei (r = 0.25, p < 0.001), and S. affinis (r = 0.53, p < 0.001) exhibited significant positive relationships between δ¹⁵N and SL (Figure 4). Relationships between δ¹³C and SL were more variable than those observed with δ¹⁵N (Figure 4). Two species, G. chuni (r = 0.33, p < 0.01) and O. lowii (r = 0.44, p < 0.001) displayed significant positive relationships between δ¹³C and SL (Figure 4).

Figure 4.

Results of least-squares regression analysis between standard length (mm) and δ¹⁵N and δ¹³C values: (a) A. cornuta, (b) C. sloani, (c) C. atlantica, (d) G. chuni, (e) G. indica, (f) O. lowii, (g) P. guernei, and (h) S. affinis.

Isotopic niche breadth, calculated using SEA_c, was largest for the piscivorous S. affinis (SEA_c= 2.27), G. indica (SEA_c = 1.98), A. cornuta (SEA_c = 1.96), and G. chuni (SEA_c = 1.53), which collectively occupied the highest trophic positions within the guild of predators examined. C. atlantica (SEA_c = 0.1.19) occupied an intermediate relative trophic position and intermediate-sized trophic niche. P. guernei (SEA_c = 0.71) and O. lowii (SEA_c = 0.62), which feed primarily on crustaceans and cephalopods, respectively, occupied lower relative trophic positions and were characterized by relatively small isotopic niches (Figure 5, Table 2). Interestingly, the smallest isotopic niche also belonged to a known piscivore, C. sloani (SEA_c = 0.56), although the small calculated isotopic niche area could have been due to a limited sampled size range (Figure 4). In the 20 instances where overlap in SEAc occurred, the percentage of shared isotopic niche space ranged from 1% (A. cornuta and C. atlantica; G. indica and O. lowii) to 27% (between G. chuni and G. indica) (Table 3). Directional overlap, or the percentage of one species’ ellipse covering the ellipse of another species, varied widely from 2 to 100%. Differences in directional overlap was greatest between C. sloani and C. atlantica (81 vs. 38%), C. sloani and S. affinis (63 vs. 15%), and O. lowii and S. affinis (100 vs. 27%) (Table 3, Figure 5).

Table 2.

Metrics for estimating isotopic niche size in eight meso- and bathypelagic predators.

Species	TA	SEA	SEAc	CD
A. cornuta	6.81	1.87	1.96	0.97
C. sloani	2.68	0.53	0.56	0.46
C. atlantica	3.02	1.13	1.19	0.82
G. chuni	5.82	1.46	1.53	0.83
G. indica	5.17	1.88	1.98	0.98
O. lowii	2.33	0.60	0.62	0.58
P. guernei	2.58	0.69	0.71	0.66
S. affinis	7.01	2.18	2.27	1.04

TA, total area (expressed in ‰) encompassed by all data points of each species; SEA, standardized ellipse area for each species; SEAc, size-corrected standardized ellipse area; CD, centroid distance calculated by taking average distance of each data point from the centroid for each species.

Table 3.

Isotopic niche overlap measured in percentage of shared space (‰) between each pairwise combination of species.

	A. cornuta		C. sloani		C. atlantica		G. chuni		G. indica		O. lowii		P. guernei
C. sloani	0 (0, 0)	0
C. atlantica	1 (2, 3)	0.05	26 (81, 38)	0.99
G. chuni	9 (23, 30)	0.19	0 (0, 0)	0.99	14 (32, 25)	0.80
G. indica	16 (33, 32)	0.49	0 (0, 0)	0.99	20 (41, 25)	0.94	27 (61, 47)	0.80
O. lowii	3 (4, 12)	0	18 (39, 34)	0.66	11 (16, 32)	0.01	0 (0, 0)	0.01	1 (2, 5)	0
P. guernei	0 (0, 0)	0	30 (68, 54)	0.85	7.5 (28, 47)	0.04	0 (0, 0)	0.01	0 (0, 0)	0	11 (23, 20)	0.73
S. affinis	9 (19, 17)	0.70	12 (63, 15)	0.99	11 (33, 17)	0.98	0 (0, 0)	0.15	4 (7, 7)	0.70	21 (100, 27)	1	9 (37, 11)	1

Numbers in parentheses represent the percent overlap of species A (column) with species B (row) and vice versa. Numbers in bold represent shared overlap >50%. Second column of numbers represents the likelihood of differences in SEAc size. Numbers in bold represent statistically significant differences in SEAc size between the pair of species examined.

Figure 5.

Size-corrected SEAc plotted with mean (± s.e.) δ¹³C and δ¹⁵N values for each species.

Mean POM δ¹³C and δ¹⁵N values collected from the meso- and bathypelagic were not significantly different from each other, thus mixing models were run using epipelagic POM data and data combined from the meso- and bathypelagic zones. Mixing model results suggest that all consumers included in this study derive the bulk of their carbon from epipelagic POM (Figure 6). Relative contributions of epipelagic POM ranged from 97.87% ± 1.43 in P. guernei to 73.30% ± 3.19 in G. chuni, while contributions from meso- and bathypelagic POM were much lower, ranging from 26.70% ± 3.19 in G. chuni to 2.13% ± 1.43 in P. guernei. Diagnostic plots of posterior distributions revealed a high negative correlation between the two sources (epipelagic POM and meso-/bathypelagic POM). Considering that the producer data fully constrain consumer data when an appropriate trophic enrichment factor is applied and that model diagnostics (Gelman-Rubin Diagnostic: all variables <1.01; Gweke Diagnostic: <5% of variables outside ±1.96 for each chain) indicate that the model fully converged, the negative correlation is likely caused by the similar δ¹³C signatures of sources and not from a missing carbon source.

Figure 6.

Estimated relative contributions of POM collected from epipelagic and meso- and bathypelagic depths to (a) A. cornuta, (b) C. sloani, (c) C. atlantica, (d) G. chuni, (e) G. indica, (f) O. lowii, (g) P. guernei, and (h) S. affinis.

Discussion

Trophic structure

Trophic positions inferred through stable isotope data suggest that, within this group of fishes, the highest trophic positions are held by the largely piscivorous A. cornuta, G. chuni, G. indica, and S. affinis; intermediate trophic positions are occupied by species preying on mixtures of cephalopods and fishes (C. atlantica and O. lowii), and fishes and crustaceans (C. sloani); and the lowest trophic position is occupied by P. guernei, which eats primarily macrocrustaceans (Hopkins et al., 1996; Sutton and Hopkins, 1996b). Stomach content analysis (SCA) was performed on all samples in this study, and though sample sizes with identifiable food items were relatively small, results agree with findings from previous SCA studies and support the trophic relationships inferred through SIA (Supplementary Table S4).

For the species examined, δ¹⁵N values spanned 5.91 ‰ or 1.9 TLs, while species mean δ¹⁵N values spanned 1.96 ‰ and 0.62 TL (assuming TEF of 3.15). Using mean δ¹⁵N values and applying a TEF of 3.15, our observed range of estimated trophic levels (0.62) appears to be in line with other studies examining Mediterranean (1.1 TLs), Pacific (1.6 TLs), and GOM (1.1 TLs) fish assemblages that included both micronektonivores (stomiids, anoplogastrids) and lower trophic level zooplanktivores (myctophids, gonostomatids), which have been shown to be up to 0.6 TLs lower than micronektonivores (Valls et al., 2014; Choy et al., 2015; McClain-Counts et al., 2017).

Trophic level estimates determined using a primary consumer to set the isotopic baseline placed each species between the third and fourth trophic levels. Where species-specific comparisons of trophic positions could be made, and applying a δ¹⁵N TEF of 3.15 to reported mean δ¹⁵N values, our estimated TPs for A. cornuta (3.7) and C. sloani (3.2) were similar to estimates from the Pacific (TP = 3.5 for both species) and to the GOM (C. sloani TP = 2.8) (Choy et al., 2015; McClain-Counts et al., 2017). The observed difference in TP estimates for C. sloani in the GOM was likely caused by the inclusion of smaller C. sloani (<50 mm SL) by McClain-Counts et al. (2017). Estimates of TP for S. affinis (3.4) were similar to Stomias boa collected in the Mediterranean (TP = 3.5) (Valls et al., 2014), while estimates of the three stomiid species [C. sloani (3.2), S. affinis (3.4), and P. guernei (3.1)] were within the estimated worldwide TP range of stomiid fishes (TP = 3.0–3.5) (Choy et al., 2012). This study represents the first descriptions of trophic positions using SIA for C. atlantica, G. chuni, G. indica, O. lowii, and P. guernei, so comparisons to TP estimates in other studies were not possible.

Isotopic niche size, estimated using SEAc, was largest for fishes occupying the highest TPs within the guild (A. cornuta, C. chuni, G. indica, S. affinis) and smallest in fishes occupying intermediate and lower TPs (Figure 5). The larger SEAc of the highest TP fishes within this guild could suggest more generalized feeding compared with other species. Differences in SEAc can also be influenced by an organism’s size distribution, which was not equally comprehensive in all species. The small SEAc of C. sloani, for example, was likely affected by samples that only included the largest individuals (>140 mm SL). Isotopic niche overlap was common, although the extent of the overlap was typically non-significant (<50%) (Table 3). In species where isotopic niche overlap was significant (S. affinis and O. lowii), available SCA data suggest prey resource overlap is not as strong as isotopic niche overlap would make them appear (Hopkins et al., 1996; Sutton and Hopkins, 1996b).

Samples of POM were more ¹³C depleted at eastern longitudes that are closer in proximity to the Mississippi River, while fishes became more ¹³C enriched at southern latitudes and western longitudes. Shifts in the isotopic signatures of POM in the GOM have been observed between nearshore and offshore regions and between mesoscale cyclonic and anticyclonic oceanographic features (Wissel and Fry, 2005; Dorado et al., 2012, Wells et al., 2017). Baseline differences in POM isotopic signatures between nearshore and offshore environments of the GOM can be caused by phytoplankton in offshore regions relying more heavily on isotopically light nitrogen produced by diazatrophic cyanobacteria (Trichodesmium spp.) (Holl et al., 2007; Dorado et al., 2012), while baseline differences between cyclonic and anticyclonic regions are driven by upwelling within cyclonic features supplying ¹⁵N enriched N₂ to phytoplankton (Wells et al., 2017). The pattern of POM samples becoming ¹³C depleted at eastern longitudes and fishes becoming ¹³C enriched at lower latitudes and western longitudes is consistent with the idea that organisms collected closer to the continental shelf are more likely to be partially supported by terrestrially derived organic matter from the Mississippi River where the effects of Trichodesmium spp. and upwelling on baseline δ¹⁵N values are minimal (Dorado et al., 2012).

Ontogenetic shifts in δ¹³C and δ¹⁵N

Ontogenetic enrichment in ¹⁵N was documented in six of the eight species examined. While significant relationships between δ¹⁵N and body size suggest ontogenetic patterns in feeding ecology, observed trends in some cases were driven by a few points, and the nature of size-based relationships with δ¹⁵N could change with the inclusion of different size classes and more samples. Observed enrichment in ¹⁵N with body size could be caused by ontogenetic shifts in prey selection, as has been suggested for A. cornuta, C. sloani, and P. guernei or by ingestion of larger sized prey (Hopkins et al., 1996; Sutton and Hopkins, 1996b). For species such as S. affinis and C. sloani, which have been shown to feed on myctophid fishes across their ontogeny, observed positive relationships between δ¹⁵N and body size could be a result of ingestion of larger myctophid fishes, which have been shown to become enriched in ¹⁵N with increasing size (Sutton and Hopkins, 1996b; Cherel et al., 2010; McClain-Counts et al., 2017). The negative relationship between body size and δ¹⁵N of O. lowii contrasts with published diet data suggesting that O. lowii undergoes an ontogenetic diet shift from eating fishes as juveniles to feeding primarily on squid and fish as adults (Rofen, 1966; Hopkins et al., 1996). Omosudis lowii are known to have highly distensible stomachs and have been reported to feed on prey much larger than themselves (Rofen, 1966). However, the tendency to feed on large prey appears to occur primarily during juvenile stages, as adults have been found to feed on both large and small prey (Rofen, 1966). Thus, the lack of a relationship between SL and δ¹⁵N of O. lowii could be a function of adults and juveniles feeding on similarly sized prey or by switching to prey that occupy lower TPs. The observed relationships between δ¹⁵N and body size are not necessarily the result of ontogenetic shifts in diet and can instead reflect spatial and temporal changes in the isotopic signature of nitrogen sources at the base of the foodweb (Wells et al., 2017). Spatial variation in the isotopic signatures of primary producers has been documented in the GOM, but the increased movements and longer tissue turnover rates of fishes likely diminishes spatial variation by increasing the likelihood of an organism integrating the isotopic signatures of multiple isotopic baselines.

Relative contributions of epi-, meso-, and bathypelagic POM to deep-pelagic fishes

A paradigm of deep-sea ecology is that meso- and bathypelagic organisms feed within foodwebs largely supported by epipelagic POM and that POM suspended at deeper depths contributes little carbon to higher order consumers. Recently, through the use of compound-specific stable isotope analysis (CS-SIA) of amino acids (AAs), that paradigm was challenged by evidence which suggests that zooplankton and micronekton can partly rely on small particle (0.7–53 µm) suspended POM as a carbon source (Hannides et al., 2013; Choy et al., 2015; Gloeckler et al., 2018). Choy et al. (2015) estimated the relative contributions of epipelagic and deep-water POM to the production of four fishes (including A. cornuta) in the North Pacific and found that two meso- and bathypelagic zooplanktivores received contributions from small-particle, deep-pelagic suspended POM ranging between 39 and 81%, while contributions to the micronektonivore, A. cornuta, were far less (0–23%). Gloeckler et al. (2018) examined the δ¹⁵N values of source AAs from a micronekton assemblage and found that relative contributions of small, suspended particles to micronekton were greatest in non-migratory species with night-time distributions in the lower mesopelagic and upper bathypelagic. Species with night-time distributions within the epi- and mesopelagic, however, were found to be supported by either surface particles or large, fast-sinking particles (>53 µm) at depth (Hannides et al., 2013; Gloeckler et al. 2018).

The results from our mixing-model analyses suggest that the majority of carbon (≥73%) supporting the species examined in this study appears to be derived from epipelagic sources or from fast-sinking particles at depth which carry similar isotopic signatures to particles within the epipelagic (Hannides et al., 2013). These contribution estimates, combined with vertical distribution data which suggest the collective night-time distributions of these predatory fishes span the epi-, meso-, and upper bathypelagic (Sutton and Hopkins, 1996a; Sutton et al., 2010), are in alignment with estimations for micronekton with similar depth distributions made by Choy et al. (2015) and Gloeckler et al. (2018). It should be noted that the relative contribution of small suspended particles at depth to these species cannot be fully assessed without conducting CSIA-AA and that further investigation into the relative importance of small particles to higher trophic level consumers is warranted (Gloeckler et al., 2018). Additional support for the assertion that these species are largely supported by surface derived carbon is provided by diet studies, which suggest that these species consume migratory prey that feed within food chains supported by surface production (Hopkins et al., 1996; Sutton and Hopkins, 1996b), highlighting the extent to which spatially distinct consumers are connected in the northern GOM.

Acknowledgements

We thank Michael Novotny, Nina Pruzinski, and Matthew Woodstock for their help organizing and processing samples. We thank two anonymous reviewers for their comments, which greatly improved this paper. This research was made possible by a grant from the Gulf of Mexico Research Initiative. Data are publicly available through the Gulf of Mexico Research Initiative Information & Data Cooperative (GRIIDC) at https://data.gulfresearchinitiative.org (doi: 10.7266/N7ZS2V04).