Resource partitioning may limit interspecific competition among Arctic fish species during early life

Arctic ichthyoplankton assemblages are generally composed of a high proportion of Arctic cod (Boreogadus saida) and low numbers of other common Arctic taxa, typically of the families Cottidae, Stichaeidae and Liparidae (e.g., Norcross et al., 2010; Parker-Stetter et al., 2011; Suzuki et al., 2015; Bouchard et al., 2016). Larvae of boreal species, such as Ammodytes spp., also co-occur with Arctic species in the major gateways to the Arctic Ocean (Chukchi Sea, Baffin Bay, Barents Sea). Due to climate change, the occurrence of boreal species increases in High Arctic locations, where they could enter in competition with Arctic cod (Norcross et al., 2010; Falardeau et al., 2014, 2017; Munk et al., 2015). Besides the introgression of boreal species, climate change could modify environmental conditions and favor the growth, feeding success, and survival of some Arctic ichthyoplankton species at the expense of others (e.g., Laurel et al., 2018). Changes in interspecific competition among early life stages of Arctic fishes can impact their recruitment and have cascading impacts on the entire ecosystem.

With the acceleration of climate change (IPCC, 2021) at high latitudes and the recognition of Arctic cod as a linchpin species in Arctic marine ecosystems, understanding how abiotic and biotic factors influence its recruitment has become a major goal in recent years (e.g., Bouchard et al., 2017; Gjøsæter et al., 2020; Mueter et al., 2020). Transport of the early life stages has been identified as an important factor for the feeding success, growth, and survival of young Arctic cod. Modelling studies of the Barents and the northern Chukchi seas showed high interannual variability in larval Arctic cod growth and dispersal trajectories in relation with local environment and climate indices (Huserbråten et al., 2019; Vestfals et al., 2021). East of Svalbard, a northward trend in projected spawning location was associated with anomalous drift patterns where most larvae were advected northwards along the margins of the continental shelf and into the Arctic Ocean interior instead of into their usual nursery habitats in the Barents Sea, with unknown consequences on survival (Huserbråten et al., 2019).

Another important determinant of Arctic cod recruitment is prey availability, as suggested by a positive relationship between hydroacoustic estimates of mesozooplankton density and the biomass of age-0 Arctic cod in August (LeBlanc et al., 2020). The timing of sea-ice breakup, through its effect on the phytoplankton bloom timing and subsequent secondary production, is another major determinant in the abundance of Arctic cod at the end of their first summer (LeBlanc et al., 2020). A breakup occurring too early or too late negatively affects the quantity of zooplankton prey available for age-0 Arctic cod, resulting in poor recruitment (Dezutter et al., 2019; LeBlanc et al., 2020). High quality of the available prey, in terms of size and lipid content notably, also appears crucial for Arctic cod recruitment. In the Canadian Arctic, Arctic cod larvae and juveniles actively select the lipid-rich eggs, nauplii, and copepodites of Calanus glacialis (Bouchard and Fortier, 2020). In the Chukchi Sea, areas with high potential for the accumulation of Calanus glacialis lipids were associated with increased condition of age-0 Arctic cod (Copeman et al., 2020).

In this study, we documented age-0 fish assemblages from the Greenland Sea in August–September 2017, co-dominated by gadids (mostly Arctic cod, with a low number of ice cod Arctogadus glacialis), bigeye sculpin Triglops nybelini, and gelatinous snailfish Liparis fabricii. First, we used a particle tracking method to test the hypothesis that age-0 Arctic cod have hatched in the Northeast Water Polynya. Second, we used gut content analyses to determine the main prey (in terms of both numerical and carbon contributions) and the preferred prey items (zooplankton taxa positively selected from the environment, as determined by a selectivity index) of each ichthyoplankton species to estimate their diet overlap. We then used multinet and acoustic data to examine the spatial overlap, both horizontally and vertically, of age-0 ichthyoplankton with their zooplankton prey.

The northern part of the Greenland Sea is characterized by a wide continental shelf and the presence of sea ice and icebergs throughout the year. Among the multiple recurrent polynyas occurring along the northeastern coast of Greenland, the Northeast Water Polynya (79–8l°N, 5–15°W) is the largest and northernmost one. It generally opens in April or May, closes in September and benefits from enhanced primary and secondary production compared with surrounding waters where production is generally low (Smith et al., 1990; Hirche et al., 1994; Lara et al., 1994; Boertmann et al., 2020). The Northeast Water Polynya has one of the lowest ichthyoplankton diversity recorded among Arctic seas. The latest ichthyoplankton survey of the region in the summer of 1993 reported an assemblage composed of only four taxa: Arctic cod (83%), Liparis fabricii (13%), Triglops sp. (3%), and one unidentified species (Michaud et al., 1996). The most prominent circulation feature in the Greenland Sea is the East Greenland Current (EGC), flowing southward from Fram Strait (79°N) to Cape Farewell (60°N) with a maximum strength above the continental slope, another area of increased production (e.g., Boertmann et al., 2020).

From August 22 to September 12, 2017, RV Dana surveyed the Greenland Sea in a study area covering bottom depths ranging from approximately 60 to 3000 m (Figure 1). The cruise track and sampling grid departed from the original plan as the ship needed to avoid sea ice present in the area. A CTD-rosette system (Seabird SBE 911plus®) was deployed at 82 stations. Temperature and salinity data from the CTD were post-processed using the Seabird Seasave® data acquisition and processing software. For each station, the presence or absence of different water masses was determined from CTD profiles.

Depth-stratified zooplankton samples were collected using a 0.25-m² aperture Hydrobios MultiNet®, a sampler equipped with five 50-µm mesh nets (hereafter referred to as multinet). The multinet was hauled vertically from a depth of 1000 m or 10 m above the seafloor to the surface at a speed of 0.5 m s^–1 and generally sampled in the following depth intervals: bottom–200 m, 200–100 m, 100–50 m, 50 m–surface. Samples were immediately preserved in a solution of 4% formaldehyde borax-buffered seawater. Known aliquots were taken using a Folsom splitter, and approximately 400–600 animals per sample were enumerated, measured, and identified to species or to the lowest taxonomic level attainable. For copepods, stage and sex were also determined. Calanus hyperboreus were identified to species for copepodite stages C3–C6, while C. finmarchicus and C. glacialis were identified to species for C5 and C6 based on prosome length (Swalethorp et al., 2011; Nielsen et al., 2014). Younger stages of Calanus were identified to genus only. Prosome lengths were measured for 10 individuals of each taxon and stage in each sample. Abundance of selected taxa (ind. m^–2) were integrated over the first 100 m of the water column. Station 13 was excluded from analyses as the 0–50 m sample was lost at sea.

To determine the best apparatus to sample ichthyoplankton during the study, the first two stations were sampled with bongo nets (two 0.6 m-diameter frames equipped with one 335-µm and one 500-µm mesh nets) and a Methot Isaac Kidd net (MIK net), a 2-m diameter frame carrying a 1500-µm mesh black net equipped with a General Oceanics® flowmeter and a Scanmar® real-time depth sensor. The size of the individual fish collected (all Arctic cod) ranged from 12.0 to 22.7 mm in standard length and from 1.4 to 3.1 mm in body depth at the anus. The size distribution of Arctic cod did not differ significantly between samplers (t-test; t₍₁₂₎ = 1.27; p = 0.23), suggesting no difference in size selectivity between the two nets. The large aperture of the MIK net ensures catches when ichthyoplankton abundances are low, which is often the case in High Arctic seas in August (e.g., Suzuki et al., 2015). We thus used the MIK net to collect ichthyoplankton for the remainder of the cruise and used only MIK net samples in further data analyses. At 25 stations, the MIK net was deployed from the stern at a ship speed of 3.01 ± 0.25 knots to sampling depths of 20–144 m (Figure 1). Onboard, age-0 fish (larvae and juveniles) were sorted from the macrozooplankton sample, identified to the lowest taxonomic level possible, measured fresh (standard length, SL, for 25 individuals per species per station), and preserved in 95% ethanol.

Boreogadus saida and Arctogadus glacialis are sympatric in the Greenland Sea (Mecklenburg et al., 2018; Bouchard, 2020). The early life stages of both species are impossible to distinguish using external morphology, and accurate species identification relies on either otolith microstructure or genetics (Madsen et al., 2009; Bouchard and Fortier, 2011; Bouchard et al., 2013, 2016; Nelson et al., 2020). To determine the proportion of each species encountered during sampling, we collected tissue from a subset of 119 Gadidae from 21 stations, selected to represent the spatial coverage and size range of the total sample. Tissue samples were stored in 97% ethanol and were sent to the Laboratory of Biodiversity and Evolutionary Genomics of the University of Leuven, Belgium for genetic analysis. There, DNA was extracted and amplified using primers designed to target a fragment of the cytochrome c oxidase subunit I (COI) gene in the mitochondrial genome (DNA barcoding). PCR amplicons were sequenced and species were determined following the methods described in Bouchard et al. (2020).

The hatch dates of age-0 Arctic cod collected during the study were estimated from individual standard length and capture date, assuming a larval growth rate of 0.22 mm d^–1 and a size-at-hatch of 5.0 mm typical for the circumpolar Arctic (Bouchard and Fortier, 2011; Bouchard et al., 2016). The estimated hatch dates varied between May 22 and July 28, and the median hatch date (June 23) was used as starting date in the particle tracks simulation. Daily output of surface currents (horizontal velocity components, June–September 2017) was extracted from the CMEMS-TOPAZ4 Arctic Ocean Model data set (Xie et al., 2017; https://marine.copernicus.eu). To test the hypothesis that Arctic cod larvae originated from the Northeast Water Polynya, we estimated particle tracks, termed progressive vector diagrams (PVDs; Emery and Thomson, 2001), from integrated time series of the modelled daily surface currents. PVDs were used to analyze the potential downstream displacements (tracks) of passive particles originating at different locations in the Northeast Water Polynya on June 23, 2017. PVDs were calculated at 24 initial locations close to the coast, on the inner shelf, mid-shelf and the outer shelf, respectively, until September 9, 2017, the last sampling day of the study.

The MIK net was deployed 34 times, mostly to about 100 m (mean ± standard deviation: 101.5 ± 3.2 m, n = 22), but also 10 times at shallower depths (20–51 m) and twice at deeper depths (125 and 144 m). Excluding the two deepest deployments, the filtered volume (V, in m^–3, measured with the flow meter) had a strong relationship with deployment duration (T, in min):

To account for difference in sampling depths, the filtered volume was standardized to 100-m depth with the following procedure. For each deployment, the time it would have taken to reach 100 m (T100 in minutes) was estimated using the ratio of deployment duration to sampled depth (e.g., T100 = 30 min for a 15-min deployment to 50 m). The filtered volume was then adjusted to 100 m (V100 in m^–3) using Equation 1. The two deepest deployments were discarded from the analyses. Larval fish abundances were estimated for the 32 remaining deployments using the number of fish larvae collected divided by V100, assuming that no larvae would have been collected between the sampling depth and 100 m. This assumption is justified because Arctic fish larvae are mainly found in the upper 30 m of the water column (Munk et al., 2015; Bouchard et al., 2016). Supporting this assumption, larval fish abundances (standardized to 100 m) at the seven stations where the MIK net was deployed twice, but always > 20 m, were not significantly different between deployments (paired t-test, t₍₆₎ = –0.74, p = 0.49). Larval fish abundances were averaged over the two deployments for these seven stations.

The diets of age-0 Arctogadus glacialis, Boreogadus saida, Liparis fabricii and Triglops nybelini were estimated from the gut contents of 8–25 individuals per species, selected to represent the size range and spatial distribution of the collection. The mouth gape of the individuals analyzed for gut contents were measured after preservation in ethanol. Maximum gape height was measured from the anterior-most tip of the premaxilla to the anterior-most tip of the dentary with the mouth opened to form a 90° angle between the upper and the lower jaw (Shirota, 1977). The digestive tracts were dissected in glycerol and individual prey were identified to the lowest taxonomic level possible. The most abundant prey categories were used in diet description, selectivity index calculations, and diet overlap calculations: Pseudocalanus spp. nauplii, Calanus spp. nauplii, copepodite stages of Triconia borealis, Oithona similis, Microcalanus spp., Pseudocalanus spp., and Calanus spp. The category “other/unidentified copepod nauplii” composed of Oithona spp./Triconia spp. nauplii, Microcalanus spp. nauplii, Metridia longa nauplii and unidentified copepod nauplii, was used in the diet and diet overlap calculations. The remaining taxa, representing less than 3% of the total number of prey items for all species, were grouped into the category “others” composed of Amphipoda, Isopoda, Polychaeta, Limacina helicina, Metridia longa copepodites, Mormonilloida, Harpacticoida, unidentified copepods and unidentified copepod pieces, and used in the diet analysis only. Oithona similis egg cases were identified and counted but not considered as prey items in the analyses as they were most likely ingested while attached to a female. In some cases, unidentifiable digested material was detected but was not considered in any analysis. The carbon content of each prey was estimated using published carbon-length, length-weight and carbon-weight relationships following Bouchard et al. (2016).

The preferred prey of each species were determined using the Chesson’s α selectivity index (Chesson, 1978):

where N is the number of prey taxa considered; d_j/p_j is the ratio of the relative frequencies of prey j in the diet (d_j) and in the in situ plankton assemblage (p_j); and Σ (d_i/p_i) the sum of this ratio over all prey taxa. For each station, the relative frequencies of prey were calculated from the abundance of the six considered taxa between 0 and 100 m estimated from the multinet. Prey selectivity was computed for each fish and then averaged for each species.

The diet overlaps between fish species were quantified using the Schoener’s index D (Schoener, 1970):

where p_ij and p_ik are the relative proportions of prey item i in the diets of species j versus species k. To calculate the Schoener’s index, the relative proportion of a given prey item (i) within the gut contents of each individual was determined, then these proportions were averaged across all individuals of a given species (j). The Schoener’s index varies between zero (no overlap) and one (complete overlap) and values greater than 0.6 are typically interpreted as an ecologically significant degree of diet overlap between taxa (Wallace and Ramsey, 1983).

To provide additional insight into the vertical distribution in age-0 fish and their zooplankton prey, hydroacoustic data were recorded with a hull-mounted echosounder and an acoustic probe (Simrad® Wideband Autonomous Transceiver, WBAT) using settings described in Møller et al. (2019). The EK60 and the WBAT echosounders were calibrated following the standard sphere methods using a single tungsten carbide 38.1-mm sphere (Demer et al., 2015).

A hull-mounted echosounder was used for measuring the vertical distribution of age-0 fish. The Simrad® EK60 split-beam echosounder with hull-mounted transducers at 18 kHz, 38 kHz and 120 kHz recorded data continuously during the cruise, except during the deployment of the WBAT. Acoustic data from the EK60 were analyzed following Bouchard et al. (2017). CTD profiles were used to determine sound speed in water (Mackenzie, 1981) and the coefficient of sound absorption (Francois and Garrison, 1982) for the acoustic analysis. Echograms were scrutinized to discard noise. The difference in mean volume backscattering strength (ΔMVBS) at 38 and 120 kHz (ΔMVBS_120–38) was then used to isolate the signal of fish in echo-integration cells of 5 min long by 5 m deep between 14.5 m and 100 m. Only cells with ΔMVBS_120–38 ranging from –10 to 5 dB re 1 m^–1 were kept for further analyses (Benoit et al., 2014; Geoffroy et al., 2016). Backscatter profiles of fish were then converted into density profiles (ind. m^–3) (Text S1, Table S1). We assumed that most of the fish signal originated from ichthyoplankton because, as adults, Arctic cod remain close to the bottom at this time of the year (e.g., Geoffroy et al., 2016) whereas the other fish species present are demersal.

A WBAT deployed as an acoustic probe was used to document the vertical distribution of zooplankton. The WBAT was deployed at selected stations and was mounted on the CTD-rosette system frame, which descended and ascended at a rate of 0.25 m s^–1. WBAT profiles were used to estimate vertical abundance profiles of mesozooplankton in the upper 100 m of the water column. The WBAT was equipped with a horizontally looking depth-rated 200 kHz split-beam transducer to detect mesozooplankton > 3 mm, such as the abundant Calanus spp. The S_v values at each depth were then averaged between 5 and 90 m from the transducer, and the average S_v values were converted in density profiles (ind. m^–3; Text S2).

Near-surface waters at depths < 100 m were dominated by cold (–1.73 < Θ < 2ºC) and fresh (S < 33) Polar Surface Water (PSW, Figure 2a). PSW originates in the Arctic Ocean and propagates southward on the East Greenland shelf with the EGC. The upper part of PSW (< 50 m depth) in the study area was less saline and slightly warmer because of seasonal heating and ice melt. Along the East Greenland continental slope in the depth range of 100 to 200 m, Arctic Atlantic Water (AAW) was observed at salinities 34 to 34.8 and temperatures –1.73°C to 2°C, respectively (Figure 2b). The higher salinities of AAW result from injections of more saline Atlantic waters to the southward flowing Arctic waters. The warmest and most saline water was Return Atlantic Water (RAW) in the depth range of 100 to 300 m and was mainly found along the East Greenland shelf break with salinities > 35 and temperatures of up to 4°C (Figure 2b). RAW propagates northward with the West Spitsbergen Current and becomes entrained in the southward flowing East Greenland Current system at 78°N (Håvik et al., 2017). Deep waters at depths > 400 m were mainly occupied by Upper Polar Deep Water (UPDW) at salinities < 35 and temperatures < 2°C (Figure 2b).

During our survey, average surface current speeds associated with the EGC were up to 30 cm s^–1 between 79°N and 76°N (Figure 2c). Warmer and more saline Atlantic waters entered the Greenland Sea as a distinct current with current speeds of up to 25 cm s^–1. Part of this current recirculated west of Spitsbergen and was incorporated in the ECG, thus providing the source properties of the warmer and more saline AAW and RAW on the East Greenland shelf and slope.

Gadidae dominated the larval fish assemblages, representing 53% of the total number of larvae collected (Table 1). Successful DNA extraction, COI gene amplification and sequencing allowed the identification of 112 of the 119 Gadidae selected for genetic analyses. Based on published GenBank sequences, the COI gene of eight fish corresponded to the reference sequence for Arctogadus glacialis (ice cod) whereas sequences of the remaining individuals corresponded to Boreogadus saida. Results suggest that most A. glacialis were selected for genetic identification, even if only 26% of the Gadidae collected were analyzed genetically (Table S2). A. glacialis were found at six of the 21 selected stations (Table S2). At five of these stations located on the shelf (stations 38, 41, 45, 54) and in Dove Bay (station 57), 100% of the Gadidae were genetically identified, of which 43–100% belonged to A. glacialis. At the sixth station (station 15) located on the slope, 56% of the Gadidae were genetically identified, of which 4% (1/23) belonged to A. glacialis, and 96% (22/23) to B. saida (Table S2). Based on these results, we treated all Gadidae not identified using genetics as B. saida. Potential bias resulting from this procedure would be minimal, because Arctic cod occurred in highest abundances over the continental slope, in low abundances over the continental shelf, and were absent from the two stations sampled in the warm waters of Atlantic origin (Figure 3b). Ice cod occurred mostly in the northern part of the shelf and at station 57 in Dove Bay (Figure 3a). The remaining ichthyoplankton consisted of the Cottidae Triglops nybelini (29%) and Gymnocanthus tricuspis (0.2%) and the Liparidae Liparis fabricii (17%) of similar size ranges (Table 1). These three species were almost exclusively (91% of total number) distributed over the continental shelf in waters depths < 500 m (Figure 3c and d). At stations where they occurred, the densities of age-0 A. glacialis, B. saida, L. fabricii and T. nybelini ranged from 0.02 to 0.3, from 0.09 to 12.6, from 0.04 to 3.0, and from 0.12 to 8.6 ind. 1000 m^–3, respectively (Figure 3). The smallest Arctic cod larvae occurred along the slope in the northern part of the study area, while the largest individuals were found in the southern part of the slope and the shelf (Figure S1). The abundance of L. fabricii was negatively related with the presence of upper Polar Surface Water (Mann–Whitney sum rank test, p = 0.034). No other significant relationships were detected between the abundance of any of the age-0 fish species and any of the water masses (Mann–Whitney sum rank tests, p > 0.05).

Virtual Arctic cod larvae originating on June 23, 2017, in the Northeast Water Polynya near the coast were retained close to source region or drifted toward the north (Figure 4a). The virtual larvae from an inner shelf source region only drifted a relatively short distance toward the south between the release date and September 9 (Figure 4b). The virtual larvae originating at mid-shelf around the Northeast Water Polynya drifted relatively long distances towards the south (Figure 4c). The virtual larvae from the easternmost source region on the outer shelf were displaced over the longest distances towards the southeast (Figure 4d). Hence, age-0 Arctic cod collected in August–September 2017 in our study area most likely originated from locations on the mid- and outer shelf of the Northeast Water Polynya.

The diet of age-0 gadids differed from the diet of non-gadids (Figure 5a). Arctogadus glacialis and Boreogadus saida consumed primarily small prey (< 530 µm) such as copepod nauplii and Oithona similis copepodites whereas Liparis fabricii and Triglops nybelini preyed on larger prey (> 875 µm), primarily Calanus spp. copepodites (Figure 5a). The most frequent prey of A. glacialis were, by far, O. similis copepodites (63%), with Pseudocalanus spp. copepodites as the second most frequent prey (12%). The diet of B. saida was composed mostly of copepod nauplii (44%) including Calanus spp. nauplii (25%), Pseudocalanus spp. nauplii (9%) and other or unidentified nauplii (10%), as well as O. similis copepodites (33%). The diet of L. fabricii was composed of Calanus spp. copepodites (53%), Calanus spp. nauplii (11%) and Pseudocalanus spp. copepodites (9%). The most frequent prey of T. nybelini were Calanus spp. copepodites (31%), O. similis copepodites (23%) and Calanus spp. nauplii (16%), with other taxa contributing between 1 and 10% to the diet (Figure 5a). Calanus spp. copepodites were the most important source of carbon for all species, contributing between 45% (for B. saida) and 98% (for L. fabricii) of the energy to larval diets (Figure 5b). Other significant sources of carbon for A. glacialis included O. similis copepodites (19%) and Pseudocalanus spp. copepodites (18%). The second most important source of carbon for B. saida were Calanus spp. nauplii (21%), whereas O. similis copepodites and Pseudocalanus spp. copepodites accounted for 14% and 5% of the species’ carbon intake, respectively (Figure 5b). The average mouth gape heights of the fish selected for stomach content analysis differed significantly among species (ANOVA, p < 0.001; Holm-Sidak multiple comparisons): L. fabricii and B. saida had the smallest mouth gape, A. glacialis had a gape larger than L. fabricii but not significantly different from B. saida, and T. nybelini had by far the largest mouth gape (Figure 5a).

As shown by the distinct frequencies in the fish guts and in the plankton assemblages at the stations where the fish were collected (Chesson’s index; Figure 6), Pseudocalanus spp. nauplii were neutrally selected by A. glacialis and B. saida and negatively selected by L. fabricii and T. nybelini. Calanus spp. nauplii were positively selected by A. glacialis and B. saida and negatively selected by L. fabricii and T. nybelini. Pseudocalanus spp. copepodites were positively selected by A. glacialis, B. saida and T. nybelini and neutrally selected by L. fabricii. Calanus spp. copepodites were positively selected by L. fabricii and T. nybelini, neutrally selected by B. saida and negatively selected by A. glacialis. All other taxa (Triconia borealis copepodites, Oithona similis copepodites, Microcalanus spp. copepodites and Appendicularia) were negatively selected by all species (Figure 6).

As indicated by a Schoener’s index of 0.63, the diets of both gadids overlapped significantly. The diets of both non-gadids showed an even higher overlap of 0.91 (Table 2). However, Boreogadus saida and Arctogadus glacialis diets did not significantly overlap with the diets of Liparis fabricii and Triglops nybelini (Table 2).

The prosome length of ingested Pseudocalanus spp. copepodites differed among fish species (Kruskal-Wallis analysis of variance on ranks, p < 0.001; Dunn’s multiple comparison test). Smaller individuals ingested by Boreogadus saida (median: 627 µm prosome length; Figure 7) corresponded to the length of Pseudocalanus spp. copepodite stage C2 found in the study area (Table S3), whereas the median length of the individuals ingested by Arctogadus glacialis, Liparis fabricii and Triglops nybelini (Figure 7) corresponded to Pseudocalanus spp. C3 (Table S3). The length of Calanus spp. copepodites ingested also differed among fish species (Kruskal-Wallis analysis of variance on ranks, p < 0.001; Dunn’s multiple comparison test). The median length of the Calanus spp. copepodites ingested by B. saida corresponded to the length of Calanus spp. C1 found in the study area; the median length of the individuals consumed by A. glacialis corresponded to Calanus spp. C2; the median length of the individuals consumed by L. fabricii corresponded to Calanus spp. C4; and the median length of the individuals consumed by T. nybelini corresponded to Calanus finmarchicus C5 (Figure 7, Table S3).

There was correspondence between the horizontal distribution of Boreogadus saida (Figure 3b) and that of Calanus spp. nauplii (Figure 8a), with both taxa occurring mostly along the continental shelf break. A linear relationship also demonstrated the link between the abundances of B. saida and its main prey Calanus spp. nauplii (r² = 0.45, p < 0.001, n = 22, station 9 excluded as an outlier, both variables square-root transformed). Calanus spp. C1, a zooplankton taxon and developmental stage identified as an important source of carbon for B. saida (Figures 5b and 7), also distributed mostly along the continental shelf break (Figure 8b). Oithona similis copepodites was a main prey of Arctogadus glacialis and B. saida (Figure 5a) and were omnipresent in high densities over the study area (Figure 8c). Pseudocalanus spp. copepodites, a prey selected by A. glacialis, B. saida and Triglops nybelini, were most abundant around Dove Bay and the southernmost stations (Figure 8d) which were highly influenced by water of Atlantic origin indicated by warmer surface temperatures (Figure 1). Calanus spp. copepodites, the most important source of carbon for all fish species, and likely composed of C. finmarchicus and C. glacialis, were distributed relatively uniformly in the study area, with the exception of a very high abundance of C. finmarchicus and the absence of C. glacialis copepodites at one station with water of Atlantic origin (Figure 8e and f).

For 14 of the 15 stations with both EK60 and WBAT profiles, the densities of mesozooplankton and fish larvae overlapped and peaked in the top 30 m of the water column (Figure 9). At some stations (e.g., stations 36, 38, 41 and 61), zooplankton density peaked at 10–15 m deeper than fish larvae, whereas the opposite was observed at station 74 where zooplankton density peaked at 20 m and fish larvae at 60 m. At station 66, the peak in zooplankton densities was not in the top 30 m, but at 55 m (though data are lacking in the upper 35 m). At some stations, an additional, smaller peak was observed for the zooplankton (e.g., stations 6 and 38) or the fish larvae (e.g., stations 18 and 27), below 55 m. The vertical profiles of fish larvae did not differ among the assemblages, and stations dominated by gadids had similar patterns as stations dominated by non-gadids (Figure 9).

Maximum peak density of ichthyoplankton measured from the EK60 was relatively low, varying from nearly absent at station 32 to 0.74 ind. m^–3 at station 49. Maximum peak density of mesozooplankton measured from the WBAT was higher and varied from 6.2 ind. m^–3 at station 57 to 48 ind. m^–3 at station 74 (Figure 9). Surprisingly, although there was a clear overlap between the depth of the peak density of ichthyoplankton and that of mesozooplankton, there was no relationship between the peak density values (Spearman rank correlation; p = 0.91).

Arctic cod larvae gradually increased in size and decreased in numbers from north to south along the continental slope. Based on local current patterns, we conclude that the larvae originated from the north and drifted with the EGC into the study area. As our calculated particle tracks suggest, the larvae could have hatched in the open waters of the Northeast Water Polynya (Michaud et al., 1996). There, high production of Calanus glacialis nauplii (Ringuette et al., 2002) would have provided favorable feeding conditions during early larval development, a beneficial association maintained by the subsequent co-transport of the larvae and their prey in the East Greenland Current. More precisely, the calculated particle tracks suggest an origin of the larvae on the mid- and outer shelf of the polynya. Because progressive vector diagrams are based on spatial projections of currents obtained from a single location, they can only provide a first order estimate of the true particle trajectories of Arctic cod larvae. PVDs assume spatially uniform currents and changes along the particle trajectories are not considered. Moreover, particles are assumed to drift passively throughout the entire simulation period. Despite these limitations, the EGC constitutes a continuous slope current throughout almost the entire year (Woodgate et al., 1999), and we thus estimate that PVD tracks provide a quasi-realistic approximation of the particle propagation during the survey period. A portion of the larvae probably continued drifting southward and out of the study area, while some potentially drifted inshore over the shelf due to entrainment by local eddies. In the Canadian Beaufort Sea, the Laptev Sea, and the East Siberian Sea, age-0 Arctic cod also occur in large proportion over the continental slope (Suzuki et al., 2015; Mishin et al., 2018).

In High Arctic seas, Calanus glacialis nauplii generally dominate the diet of age-0 Arctic cod < 25 mm in standard length. Calanus spp. copepodites are ingested in low number, but also account for a large proportion of the carbon intake (Bouchard and Fortier, 2020). In the Greenland Sea in August–September 2017, age-0 Arctic cod positively selected Calanus spp. nauplii, one of their main prey items, which represented 25% of abundance in gut contents. The positive relationships between the horizontal distribution of age-0 Arctic cod (Figures 3c) and Calanus spp. nauplii (Figure 8a) and between the vertical distribution of ichthyoplankton and mesozooplankton (Figure 9) support the idea of a beneficial co-transport of Arctic cod larvae and their zooplankton prey in the EGC. At 9.6-mm SL, Arctic cod have developed the mouth gape and swimming capacity to capture C. glacialis spp. copepodite stage C1 and start to positively select them from the plankton assemblages (Bouchard and Fortier, 2020). In the present study, Calanus spp. copepodites also accounted for a large fraction of the carbon intake of Arctic cod. The prosome length of the Calanus copepodites ingested by Arctic cod corresponded with the length of Calanus spp. C1, confirming that Arctic cod larvae can ingest Calanus copepodites and supporting the idea that they are important for their feeding success (Bouchard and Fortier, 2020). The copepodite stages of the cyclopoid Oithona similis constituted 33% of Arctic cod prey in the current study, but were, as in other studies, negatively selected from the mesozooplankton (Michaud et al., 1996; Falardeau et al., 2014). Different factors could explain the avoidance of age-0 Arctic cod and other fish species for O. similis copepodites despite their high densities and widespread occurrence in the circumpolar Arctic. The swimming mode of cyclopoids is characterized by pause-search movements, compared with the longer and more consistent swimming mode of calanoids. The resulting differences in prey perception have been suggested as a possible explanation for the preference of larval Atlantic cod Gadus morhua for calanoids over cyclopoids: less active prey are less likely to be seen, thus reducing their vulnerability to predators (e.g., Pepin and Penney, 1997). A second and non-mutually exclusive explanation for the avoidance of age-0 fish of O. similis copepodite stages (and of Triconia borealis, another cyclopoid) is linked to their low energetic value. The ingestion of cyclopoids has been associated with a decrease in feeding success of age-0 Arctic cod in the High Arctic (Caissy et al., 2020), concordant with a central principle of the optimal foraging theory that prey selectivity is a mechanism allowing predators to optimize energy gained over energy spent (Pyke, 1984).

Arctic cod dominated the ichthyoplankton assemblage over the continental slope whereas Arctogadus glacialis, Liparis fabricii and Triglops nybelini occurred mostly on the continental shelf. Differences in spawning areas, reproductive strategies and early life histories among the species could be the main determinants of larval drift and subsequent distribution patterns. For example, liparids, fish of the genus Triglops, and A. glacialis have demersal eggs whereas B. saida have pelagic eggs (e.g., Pietsch, 1993; Graham and Hop, 1995; Fahay, 2007; Overdick et al., 2014), which likely contributed to different drift durations and patterns among the species. However, all species co-occurred on the shelf, and Arctic cod numbers may be limited there because of the strong competitive advantages (larger size and larger mouth gape) of the other species. All species had Calanus spp. copepodites as a primary carbon source and three species (A. glacialis, B. saida and T. nybelini) positively selected Pseudocalanus spp. copepodites. This finding suggests that competitive interactions among age-0 fish likely occurred in the study area, unless food resources were either: 1) not limiting, 2) partitioned among species through vertical segregation, or 3) partitioned among species by size.

Low secondary production may have limited the feeding success and recruitment of age-0 fish in the Greenland Sea in 2017. Using a threshold of 40,000 nauplii m^–2 under which the feeding success would be limited in Hudson Bay, Michaud et al. (1996) concluded that the densities of copepod nauplii seldom limited the feeding success of early-stage Arctic cod larvae in the Northeast Water Polynya in 1993. In 2017, the densities of copepod nauplii in the top 100 m ranged from 21,480 to 525,124 nauplii m^–2 (average 273,302 nauplii m^–2) and were lower than 40,000 nauplii m^–2 at only 3 stations. However, a large part of these nauplii were of Oithona similis and Triconia borealis which, due to their small size, represent suboptimal prey for Arctic cod of about 20-mm SL. The combined densities of nauplii of Pseudocalanus and Calanus spp., which are of suitable size for later stage Arctic cod larvae, ranged from 280 to 50,160 nauplii m^–2 (average 25,220 nauplii m^–2) and were lower than 40,000 nauplii m^–2 at all but one station. These values indicate potential food limitation for Arctic cod. The densities of Calanus spp. and Pseudocalanus spp. copepodites required for age-0 B. saida, A. glacialis, L. fabricii or T. nybelini collected in this study to sustain their respective energetic requirement is unknown. The Greenland Sea is a low productivity region, with low mesozooplankton biomasses and low abundances of Calanus copepods compared to waters of similar latitude off West Greenland (Møller et al., 2018). With the exception of a few stations in Dove Bay and offshore, the abundances of Calanus finmarchicus, C. glacialis and Pseudocalanus spp. copepodites in the study area remained low. As such, it can be expected that the feeding success of age-0 fish relying on large calanoid copepods as their most energetic prey might be limited by this low level of secondary production.

Interspecific competition for zooplankton resources would be limited if the different ichthyoplankton species were vertically segregated. However, the vertical profiles of larval fish abundances obtained from acoustic data suggest that all species concentrated in the first 30 m of the water column, in correspondence with the depths of maximum zooplankton densities. Indeed, the profiles at stations dominated by Boreogadus saida (e.g., station 6) did not differ from the profiles at stations dominated by Liparis fabricii (e.g., station 45), Triglops nybelini (e.g., station 61), or co-dominated by Arctogadus glacialis and T. nybelini (e.g., station 54; Figure 9). The slight offset observed between the vertical profiles of larval fish and zooplankton densities could result from an avoidance of the zooplankton from their predators, or from a zooplankton depletion at depths where age-0 fish are highly abundant.

Interspecific competition for food likely occurred among age-0 fish in the study area, because 1) the diet of both gadids overlapped, 2) the diet of both non-gadids overlapped, 3) two species positively selected Calanus spp. copepodites and three species positively selected Pseudocalanus spp. copepodites, and 4) all species depended strongly on Calanus spp. copepodites as a main source of carbon. However, differences in the size and developmental stages of the copepodites ingested among the species suggest that resource partitioning may limit interspecific competition to some extent. The differences in prey size can be explained, at least partly, by differences in length and mouth gape. In age-0 marine fish, the size of prey consumed is generally correlated with fish length and mouth gape, but to various degrees depending on species (Pepin and Penney, 1997). Of the four co-occurring species, Arctic cod had the smallest standard length, a small mouth gape and ingested the smallest Pseudocalanus spp. and Calanus spp. copepodites (average size corresponding to developmental stages C2 and C1, respectively). On the other hand, Triglops nybelini had the largest mouth gape and preyed on the largest Calanus spp. copepodites. Liparis fabricii ingested larger Calanus spp. copepodites than Arctogadus glacialis despite a smaller mouth gape, consistent with studies showing that several factors can contribute to prey-capture performance and size-dependent prey ingestion patterns in age-0 fish besides gape size (e.g., Pepin and Penney, 1997). The size range of each species in our age-0 fish collection was relatively large. Therefore, ontogenetic changes in the diets may have contributed to some of our observations. These potential changes could not be resolved adequately due to the limited number of gut contents analyzed, but future studies documenting the diet of key Arctic fish species from the early larval to the juvenile stage would improve our understanding of resource partitioning and interspecific competition among age-0 fish in the Arctic.

In High Arctic seas, the feeding success of Arctic cod ≥ 25 mm depends strongly on the adults and final copepodite stages (C4–C6) of Calanus glacialis and C. hyperboreus (Bouchard and Fortier, 2020). These taxa make up the bulk of the summer zooplankton biomass in the surface layer of the Greenland and Canadian Arctic seas (Darnis et al., 2008, 2012; Møller and Nielsen, 2020), and likely contributed to most of the acoustic backscatter attributed to zooplankton in both LeBlanc et al. (2020) and the current study. The relationship found in the Canadian Arctic between zooplankton backscatter and the biomass of age-0 Arctic cod hold for surveys performed in August (LeBlanc et al., 2020), but not in September, as a large part of the Calanus populations are deeper than 100 m by then, having initiated or completed their seasonal vertical migration to their overwintering depths (e.g., Hirche, 1997; Darnis and Fortier, 2014). Hull-mounted echosounders thus appear as appropriate instruments to study the dependence of Arctic cod recruitment on their prey at large spatial and temporal scales, e.g., when comparing one region to the other (LeBlanc et al., 2020), but net validation is required to resolve their predator-prey interactions at smaller scales.

In the Greenland Sea from August 22 to September 12, 2017, age-0 Arctic cod fed mostly on Calanus nauplii and Oithona similis copepodites, not acoustically detectable at 120 kHz (EK60) and only detectable at 200 kHz (WBAT) when occurring in high densities. To detect single individuals, the acoustic wavelength should be of similar size or smaller than the size of the individual (Holliday et al., 1989). With a wavelength of 12.5 mm (120 kHz) or approximately 7.5 mm (200 kHz) at a sound speed of 1500 m s^–1, the frequencies used in the present study cannot detect individual nauplii and small zooplankton species (e.g., Oithona similis). However, dense aggregations of small zooplankton can be detected at 200 kHz and the acoustic data from the WBAT can be used to document the vertical distribution of these aggregations with a better vertical resolution than what can be achieved with nets (centimeters versus meters). Yet, nets are still needed for ground-truthing the species and stages of zooplankton detected by the acoustics because they provide a better taxonomic resolution of the mesozooplankton community (e.g., Wiebe et al., 1996).

Zooplankton nets likely provide more accurate density estimates of zooplankton than acoustics. In this study, the blind zone of the acoustics was 14.5 m for the EK60 and approximately 12 m for the WBAT (the latter due to noise/reflection from the surface, as the transducer was pointing horizontally). However, copepods form dense swarms near the surface at high latitudes (Basedow et al., 2019), which are sampled by the nets but not the acoustics. This difference could explain why peak densities of zooplankton and fish larvae were not statistically correlated, although both were generally distributed in the upper 30 m (Figure 9). Diel vertical migrations by the zooplankton in and out of the blind zone likely exacerbated the variation in density estimates by acoustics. In the High Arctic during late summer/autumn, diel vertical migration behavior has been observed in zooplankton species including Calanus spp. (e.g., Fortier et al., 2001; Falk-Petersen et al., 2008). Future studies should ensure that the blind zone of the surface acoustics is reduced, for example, by using autonomous surface vehicles. Another option is to have the transducers on the WBAT pointing downward measuring vertically and not sideward/horizontally as done here, as reflection from the surface resulted in lack of data in the upper approximate 12 m. This approach could suffer, however, from other issues, such as avoidance behavior of the zooplankton when the WBAT is approaching. The approach applied in the present study, where backscatter is measured horizontally up to 100 m away from the WBAT likely causes less avoidance. Nonetheless, by providing high vertical resolution information on the distribution of both zooplankton and ichthyoplankton, acoustics confirmed the vertical overlap between ichthyoplankton and their zooplankton prey, at least below 14.5 m. We thus conclude that unveiling the fine-scale interactions between larval Arctic cod, other ichthyoplankton species, and their zooplankton prey necessitates both the vertical resolution of acoustic instruments and the higher taxonomic resolution of zooplankton nets.

In Arctic and Subarctic seas, changes in sea surface temperature, sea-ice cover and salinity can have strong effects on the growth, feeding success and survival of age-0 fish, and benefit some species at the expense of others (e.g., Kono et al., 2016; Laurel et al., 2018; Spencer et al., 2020). The oceanographic conditions of the Greenland Sea are highly dynamic, spatially heterogeneous, and characterized by important interannual and interdecadal variabilities (Boertmann et al., 2020). The recent trend has been a general increase in the temperature and salinity (averaged over August–September and 0–200 m depth) and mixed layer depths compared to the 1991–2017 average, but with important variability among the East Greenland Current, the Greenland Sea and the East Greenland shelf (Boertmann et al., 2020). In 2017, the East Greenland shelf and the Greenland Sea were characterized by relatively strong positive temperature anomalies, whereas the temperature in the EGC was close to the long-term average (Boertmann et al., 2020). If the temperature in the top 30 m has the same trend as in the top 200 m, age-0 fish in the EGC, mostly Arctic cod, would have experienced relatively cold (normal) temperatures while fish on the shelf, mostly Arctogadus glacialis, Liparis fabricii and Triglops nybelini, would have experienced warmer conditions than usual. Optimal temperature ranges of these three species during early ontogeny are unknown but increasing sea surface temperature in the study area could promote their growth and survival rates, therefore increasing their competitive advantage over Arctic cod. Despite different sampling areas and months, it is tempting to speculate that such changes in water temperatures could explain why the proportion of Arctic cod in the ichthyoplankton assemblage of the Greenland Sea was 83% in 1993 (Michaud et al., 1996) and only 52% in 2017. Our results suggest that direct competition for prey from other ichthyoplankton species is likely not detrimental to Arctic cod for now. However, further increases in sea surface temperature could modify the phenology of copepods and zooplankton assemblages could become more boreal, as observed in West Greenland (Møller and Nielsen, 2020). This change could modify predator-prey relationships and potentially disrupt resource partitioning among fish species which, in turn, could exacerbate the reduction in relative abundance of Arctic cod to the benefit of other fishes, including boreal species.

DNA sequences: Genbank accessions MZ169678–MZ169789.

CTD data: Polar Data Catalogue CCIN 13252. DOI: https://doi.org/10.21963/13252.

Gut contents data: Polar Data Catalogue CCIN 13250. DOI: https://doi.org/10.21963/13250.

Mesozooplankton data: Polar Data Catalogue CCIN 13253. DOI: https://doi.org/10.21963/13253.

EK60 hydroacoustic data: Polar Data Catalogue CCIN 13251. DOI: https://doi.org/10.21963/13251.

WBAT hydroacoustic data: Polar Data Catalogue CCIN 13256. DOI: https://doi.org/10.21963/13256.

The supplemental files for this article can be found as follows:

Supplementary methods

Text S1. Calculation of age-0 fish density from EK60 backscatter profiles.

Table S1. Details MVBS_120–38 calculations for non-swimbladdered ichthyoplankton species.

Text S2. Calculation of mesozooplankton density from WBAT backscatter profiles.

Supplementary results

Table S2. Summary of barcoding results for the Gadidae collected in the Greenland Sea in August–September 2017.

Figure S1. Average standard length (mm) of age-0 Boreogadus saida in the Greenland Sea in August–September 2017.

Table S3. Median lengths of zooplankton taxa representing potential prey for age-0 fish in the Greenland Sea in August–September 2017.

We would like to thank the captain and crew of RV Dana for their professionalism during the Greenland Sea 2017 mission, Sarah Maes (UK Leuven) for the genetic analyses, and Frédéric Maps (Université Laval) for advice on data analyses.

This research project was funded by the Mineral Licence and Safety Authority (MLSA) and the Environmental Agency for Mineral Resource Activities (EAMRA) of the Government of Greenland as part of the Joint Northeast Greenland Strategic Environmental Study Program. JC and MG are financially supported by the Ocean Frontier Institute through the Canada First Research Excellence fund, ArcticNet a Network of Centres of Excellence of Canada, and the Natural Sciences and Engineering Research Council of Canada.

The authors have no competing interests to declare.

Contributed to conception and design: CB, JC, MG, EFM, MDA.

Contributed to acquisition of data: CB, AK, EFM, CM, MDA.

Contributed to analysis and interpretation of data: CB, JC, AK, MG, EFM, CM, MDA.

Drafted and/or revised the article: CB, JC, AK, MG, EFM, CM, MDA.

Approved the submitted version for publication: CB, JC, AK, MG, EFM, CM, MDA.

How to cite this article: Bouchard, C, Chawarski, J, Geoffroy, M, Klasmeier, A, Møller, EF, Mohn, C, Agersted, MD. 2022. Resource partitioning may limit interspecific competition among Arctic fish species during early life. Elementa: Science of the Anthropocene 10(1). DOI: https://doi.org/10.1525/elementa.2021.00038

Domain Editor-in-Chief: Jody W. Deming, University of Washington, Seattle, WA, USA

Guest Editor: Marcel Babin, Université Laval and CNRS, Quebec City, QC, Canada

Knowledge Domain: Ocean Science

Part of an Elementa Special Feature: Four Decades of Arctic Climate Change: A Tribute to Louis Fortier