Annual, seasonal and spatial differences in the growth rate of Baltic cod larvae and early juveniles in relation to zooplankton biomass fluctuations in 2006-2014

Because of its intra-continental location and bottom topography, the Baltic Sea has very specific hydrological conditions that fluctuate under the influence of both local and global processes. Local hydrological conditions in the Baltic Sea vary depending on the season, the intensity and frequency of water inflows from the North Sea, and the intensity of precipitation and freshwater inflows from rivers. In individual sub-areas, intensive local hydrological processes lead, inter alia, to stratification that may obstruct free circulation in the water column. This may lead to increases in near-bottom anaerobic zones and decreases in the water layer that offers favourable conditions for the development of fish species that spawn pelagic eggs (e.g. cod) and also the lower survival rate of new generations (Hüssy et al. 2012). On the other hand, long-term hydrological effects may be triggered by global changes associated with climate warming (Houde 2008; Pinsky and Byler 2015). This variability of environmental conditions in the Baltic Sea has a significant impact on the inter-seasonal and spatial differentiation of the GR of fish in their early stages of development, especially eastern that of the Baltic cod. To investigate that phenomen, we distinguished three study areas located in different sub-basins of the Baltic Sea, i.e. the Bornholm Basin (BB), the Słupsk Furrow (SF) and the Gdańsk Basin (GB) (Fig. 1).

Figure 1

Larval and juvenile cod specimens (n = 835; age range 2 – 136 days; standard length (SL) 4.1 – 39.2 mm) were caught during Polish ichthyoplankton cruises of the RV Baltica conducted in spring (between April and June) and summer (between July and August) between 2006 and 2014 (Table 1). Catches at each station were made with Bongo nets (mesh sizes of 305 and 505 μm) and MIK plankton nets (Methot-Isaac-Kidd; mesh size 1 mm). At the same station, however, the larvae could be present in different nets in different years. During the entire research period 243 samples were collected from approximately 90 different stations in total, which were located in three sub-basins in the Southern Baltic (Fig. 1 and S1). Larval and juvenile specimens of cod were sorted from the ichthyoplankton samples on board the vessel immediately after the nets were hauled in and they were preserved in 96% ethyl alcohol. The remainder of the collected material was preserved in a 4% formaldehyde solution. The alcohol was changed in all samples within a few months to avoid the risk of otolith destruction, which may occur even if the samples are preserved in alcohol (Fey 2018). Since no differences in larval and juvenile growth rates were noted in the same sampling areas with BONGO 335 and 505, and MIK nets during previous studies (Fey 2015), the research material collected with these nets was treated equally in the research.

Zooplankton samples were collected during the same RV Baltica cruises during which the ichthyoplankton samples were collected in spring (between April and June) and summer (between July and August) between 2006 and 2014. 362 samples of zooplankton were collected and analysed from 19 different stations in the same areas as the ichthyoplankton stations (Fig. 1; S1). The samples were collected according to the HELCOM Guidelines for the monitoring of mesozooplankton (2017) with a WP-2 net (mesh size 100 μm) that was deployed vertically only when distinct stratification in the water column was observed. Samples were collected independently within the three layers: from 5 m above the bottom to the halocline, from the halocline to the thermocline and from the thermocline to the surface. The material collected was preserved in a 4% formaldehyde solution.

Hydrological data (i.e. water temperature) was collected at each ichthyoplankton sampling station from 2006 to 2014. Measurements were performed from the surface to the bottom with a CTD Neil-Brown sensor and a CTD Ocean Seven 316 Plus probe (Idronaut). The temperature data were analysed in two water column layers. The first layer was from the surface to the thermocline, while the second was from the beginning of the thermocline to the bottom.

After several months, the otoliths (lapillus) were extracted from each specimen and the SL (mm) of each specimen was measured to the nearest 0.1 mm with an electronic caliper. Length measurements were corrected with a correction value (+0.7 mm) for cod that permits converting length after preservation into ‘live’ length prior to preservation (Fey 2012), so to avoid the possibly negative effect of shrinkage during preservation on the results of growth rate analysis (Greszkiewicz and Fey 2018).

Once extracted, the otoliths were placed on basic microscope slides and embedded in DEPEX (Electron Microscopy Sciences Hatfield, PA). The ages of cod larvae and juveniles were determined based on the number of micro-increments between the centre of the otolith and the edge (left or right, depending on which had more legible microstructures) (Spich and Fey 2022). Assumptions regarding the daily frequency of growth increment deposition in cod had previously been verified (Geffen 1995; Steffensen 1980; Hüssy et al. 2003). Since the first growth increment is deposited in cod on the lapillus on the day of hatching (Radtke 1989), the number of micro-increments is a direct indication of age. The otoliths were read under transmitted light at magnifications of 500× or 1000× using an Eclipse 80i microscope (Nikon Corp., Tokyo, Japan). The lapillus was used in the analysis because Spich and Fey (2020) reported that this structure ensures higher precision and accuracy of age estimates than does the sagitta, and it is less prone to producing age estimate variations among readers. They also reported that the lapillus in most cases permits reading the age of larvae and early juvenile cod without the necessity of polishing, which, in the case of sagitta, is impossible. The number of micro-increments were read by the same technician twice at an interval of a few months. If the difference exceeded 15%, the sample was read once again. The values obtained were averaged. Some of the otoliths were illegible (e.g. they were damaged during extraction) or were determined to be insufficiently legible and were excluded from analysis without attempting an age reading. In the end, age was determined for 491 of the 835 specimens obtained from the ichthyoplankton samples.

The zooplankton samples were split into analysable sub-samples using the Motoda method (Motoda 1959). Taxonomic analysis was performed in a zooplankton counting chamber (Hydrobios 80 × 100 mm) with a stereo microscope under a total magnification of up to 600×. Each mesozooplankton specimen was identified to the lowest possible taxonomic level. The identified zooplankton organisms were measured (length and width) to the nearest 0.05 mm with a measuring slide.

The average GR (mm/d) was calculated using the age determined from the otoliths of individual specimens of cod larvae and juveniles with the formula: GR=SLcatch−SLhatchAge $$GR = {{\left( {S{L_{catch}} - S{L_{hatch}}} \right)} \over {Age}}$$

Based on data from the literature (Hardy 1978; Horbowa and Fey 2013), it was assumed that the value of SL_hatch was 4.9 mm. Since the obtained mean individual GR values were dependent on the size of individual specimens (S2), and were consequently incomparable among different periods, years or areas, the growth rate index (GR index) was also determined in order to account for the GR differences among larval fish sizes. For this purpose, the length-at-age data for all years were pooled and fitted with a best-fit polynomial function using a generalised regression model (GLM) (Fig. 2B): SL=5.021+0.3077×Age−0.0049×Age2+6.8519E−5×Age3−2.6686E−7×Age4n=491,r2=0.902,p<0.001 $$\matrix{ {SL = 5.021 + 0.3077 \times Age - 0.0049 \times Ag{e^2} + } \cr {6.8519E - 5 \times Ag{e^3} - 2.6686E - 7 \times Ag{e^4}} \cr {\left( {{\rm{n}} = 491,{{\rm{r}}^2} = 0.902,p < 0.001} \right)} \cr } $$

Figure 2

The polynomial function of the 4^th order was used because the 3^rd order function seemed not to provide sound data description for the smallest specimens (visual observation of the residuals). The R² was also noticeably higher for the 4^th order function (0.902) than the 3^rd order function (0.854). On the other hand, the 5^th and 6^th order functions were not considered because neither changed shape nor increased noticeably the R² compared to the 4^th function. Residuals of the age-length-relationship vs. length for the chosen 4^th order polynomial function has been presented in S3. No other method of choosing the best model, such as AIC, were used in this case because only one independent variable (age) was considered, so there was no need to reduce the number of variables. It is important to note that the intercept of the fitted function (5.0 mm) corresponds very well to the size of cod larvae at hatch. We viewed the residuals from Fig. 2 as measures of the difference of individual growth from population growth trends. The absolute values of the residuals (mm) were then converted to a percentage of the average length-at-age to allow for the pooling of residuals across different ages. Since the GR index represented relative values, which took into consideration the size of individual specimens, it was possible to make comparisons among specimens of different sizes. The GR index was therefore then analysed in terms of its variability among years (2006–2014), seasons of the year (spring and summer), and sub-areas of the Baltic Sea (BB, SF and GB).

The GR index estimated for cod larvae and early juveniles has been compared between seasons and among geographical areas with the nested-ANOVA (separate analysis for area and season effect). The seasons and geographical areas were nested within years. The differences in GR index of larval and juvenile cod among years (data pooled for seasons and geographical regions) were compared with ANOVA for all nine years together. The aggregated GR data for two groups (years 2006–2011 and 2012–2014) were compared with a t-test. The assumptions for the parametric statistics used were verified and satisfied using the Kolmogorov-Smirnov test (normal distribution, p > 0.05) and Levene’s test (homogeneity of variance, p > 0.05). The assumption of independence of cases was ensured from the design of the study. All the statistical analyses were conducted in Statistica ver. 12.0 (TIBCO Software, Inc., Palo Alto, CA). Differences were considered statistically significant at p < 0.05 (α = 0.05).

Analysis of zooplankton was performed in order to determine quantitative (biomass) and qualitative shares at different sampling stations from 2006 to 2014. In order to calculate zooplankton biomass, zooplankton abundance (n/m³) was calculated first, based on the number of specimens at a given taxonomic level at a given sampling station and the volume of water that was filtered through the WP-2 net. The analysis considered the zooplankton collected from the first two layers of the water column, i.e. from the surface to the thermocline, and from the thermocline to the halocline. This is because the majority of selected sampling stations were characterised by anoxia below the halocline. Zooplankton biomass was calculated with the Standard Size Classes method (Witek et al. 1996). Based on the measurements of the length of cylindric shape of Copepods body, each organism was classified into certain size classes defined in the specially design worksheet developed by Witek et al. (1996). The body volume of Copepods was calculated with formula presented below. The wet weight is a geometric mean of boundaries in a fixed size classes. V=π×R1×(0.5R2)^2 $${\rm{V}} = {\rm{\pi }} \times {\rm{R1}} \times {\rm{(0}}{\rm{.5}}\,{\rm{R}}2)^\cap {\rm{2}}$$

R1 – body length

It was decided to limit the analyses of zooplankton data to zooplankton of the Copepoda class, since these are considered as the main prey organisms of larval and juvenile Eastern Baltic cod (Economou 1991, Zuzarte et al. 1996, Hüssy et al. 1997, Voss et al. 2003). Six taxonomic groups were identified: Acartia spp., Centropages hamatus, Eurythemora affinis, Oithona similis, Pseudocalanus spp. and Temora longicornis. Two size fractions were also identified: the developmental stages from nauplius to copepodite CIII (NP+C3) and from copepodite CIV to adult forms (C4+AD). The biomass of zooplankton from different taxonomic and size groups were analysed with regard to season (spring– April to June, and summer – July to August) and research area (BB, SF and GB).

The biomass of the four dominant zooplankton taxa in different years (means from all stations for Acartia spp., Pseudocalanus spp., T. longicornis and C. hamatus) was compared, and the intra-annual differences among geographical regions and between seasons were analysed with nested-ANOVA, where the geographical areas and the seasons were nested within a year. The relationships between the GR index and environmental variables (zooplankton biomass and water surface temperature) and between zooplankton biomass and surface temperature, for the entire analysed period from 2006 to 2014, were estimated with linear regression analysis, where mean annual values were used for each of the relationships. The significance of the regression lines obtained was evaluated, assuming significance at the p < 0.05. Zooplankton biomass for Pseudocalanus spp., Temora longicornis, Acartia spp., Centropages hamatus, and all four zooplankton groups together. The assumptions for the parametric statistics used were verified and satisfied using the Kolmogorov-Smirnov test (normal distribution, p > 0.05) and Levene’s test (homogeneity of variance, p > 0.05). The assumption of independence of cases was ensured from the design of the study. All the statistical analyses were conducted in Statistica ver. 12.0 (TIBCO Software, Inc., Palo Alto, CA). Differences were considered statistically significant at p < 0.05 (α = 0.05).

The mean growth rates of individuals from hatching to being caught (n = 491; age 2 – 136 d; SL: 4.1 – 39.2 mm) was 0.25 mm/d ± 0.14 (SD) (S2). The growth rate depended on larval age, and the greatest variation in growth rates was noted in the age group of up to approximately 10 d.

By fitting linear regressions to the SL-at-age dependence and comparing the slopes of the regressions, we obtained the population growth rate. Two age groups were distinguished that were characterised by different growth rates (Fig. 2A). Age ≤ 60 d, GR = 0.18 mm/d (y = 0.184x + 5.4344, n = 430, R² = 0.724, p < 0.001); age > 60 d, GR = 0.26 mm/d (y = 0.2588x + 1.5274, n = 74, R² = 0.618, p < 0.001). The splitting point was estimated by description of the SL-at-age data with 2^nd order polynomial function and with exponential function first. The inflection point of those two functions, which corresponded in both cases to the age of ca. 60 days, was used as the point to split the data for the description with two separate linear functions.

The growth rates of larvae for individuals and the population were age-dependent, so the GR index (%) was calculated for each specimen from the residuals of the polynomial function fitted to the SL-at-age data (SL = 5.021 + 0.3077 × Age - 0.0049 × Age² + 6.8519E -5 × Age³ - 2.6686E-7 × Age⁴, n = 491, r² = 0.902, p < 0.001) (Fig. 2B).

The within-year comparison of the GR index between seasons (i.e. spring and summer) was possible in 2008, 2012, 2013, and 2014 (Fig.3). Although the season effect within a given year (n = 238) was not statistically significant (nested-ANOVA, season nested within years, p < 0.01 for year effect, p = 0.073 for season effect, p = 0.172 for the interaction term), faster growth in summer than in spring was observed when the data were pooled among years: GR index oscillating at the level of ca. from -30 to 0 for spring and from -20 to 20 for summer (S7).

Figure 3

The effect of geographical area (BB, SF and GB) on the GR index (Fig. 4) was compared (data availability for all the three compared areas) in 2006, 2008, 2010, 2012 and 2013. The differences for the geographical area effect (n = 360) were not statistically significant (nested-ANOVA, geographical area nested within years, p < 0.001 for year effect, p = 0.131 for geographical area effect, p < 0.05 for interaction term effect). Those results were confirmed when the GR index was plotted against geographical position (latitude), without grouping the data into three areas, i.e. after the data were pooled among years (S4). Some tendency towards faster growth was observed in such way as to present data for larvae occurring closer to the shore (lower latitude), but only in spring (S4).

Figure 4

The GR index of larval and juvenile cod (data pooled for seasons and geographical regions) (n = 491) (Fig. 5) was statistically different among years (ANOVA, p < 0.001), and the last years analysed (2012-2014) were characterised by slower growth (mean GR index = -12.5) than the previous years (2006-2011) (mean GR index = 5.2). The differences between the two time periods were statistically significant (t-test, p < 0.001).

Figure 5

Between 2006 and 2014 the dominant Calanoida representatives in the biomass were the zooplankton identified as Acartia spp., Centropages hamatus, Pseudocalanus spp. and Temora longicornis. Species such as Oithona similis and Eurythemora affinis contributed only negligibly to the zooplankton biomass (Fig. 6).

Figure 6

In spring the NP+C3 fraction of small zooplankton were dominated by Pseudocalanus spp. and T. longicornis (50 and 32% of the mean biomass for 2006-2014). Considering the general pattern, from 2007 the contribution of Pseudocalanus spp. increased while that of Acartia spp. and C. hamatus decreased. In 2014, however, the contribution of Pseudocalanus spp. decreased significantly and that of T. longicornis increased.

In spring the C4+AD fraction of large zooplankton was most abundantly represented by Pseudocalanus spp. (29%), T. longicornis (32%) and Acartia spp. (21%). From 2011, however, the contribution of Pseudocalanus spp. declined significantly. This pattern was associated with the increased contribution of T. longicornis, C. hamatus and Acartia spp. (Fig. 6).

In summers throughout the study period the NP+C3 size fraction of small zooplankton was mostly dominated by Pseudocalanus spp. (62%), followed by Acartia spp. (10%) and T. longicornis (20%). It was difficult to identify any visible trends in the changes in the contributions of zooplankton species in this group.

In summer the C4+AD size fraction of large zooplankton was dominated by T. longicornis (50%). Considering the general pattern, the contribution of T. longicornis increased from 2006.

Changes in the biomass of the four dominant zooplankton groups (Pseudocalanus spp., Acartia spp., T. longicornis and C. hamatus) among years are presented in Fig. 7, in which geographical area and zooplankton size group (NP+C3 and C4+AD) are denoted.

Figure 7

Between 2010 and 2012 there was a visible decrease in the biomass of younger development stages (NP+C3) of Calanoid in all three analysed areas, especially in the spring (Fig. 6, 7, and 8).

Figure 8

In the case of the older developmental stages of Calanoid, a decrease in biomass was recorded as early as in 2008 and continued until the end of the research period in 2014. The only area where an increase in biomass could be observed was BB (Acartia spp., T. longicornis, C. hamatus) (Fig. 6, 7, and 8).

Among all 4 groups of Calanoida, Pseudocalanus spp. NP+C3 had the highest average biomass for the entire period (2006-2014) in all 3 areas (5.90E+03 μg m^-3 in BB, 6.10E+03 μg m^-3 in GB and 5.60 E+03 μg m^-3 in SF). C.hamatus NP+C3, however, showed the lowest average biomass for the period between 2006 and 2014 from all 4 Calanoid groups (1.03E+03 μg m^-3 in BB, 6.09E+02 μg m^-3 in GB, 4.81E+02 μg m^-3 in SF) (Fig. 7).

The largest fluctuations in biomass values among individual years were recorded in T. longicornis NP+C3 (max 1.21E+04 μg m^-3 in GB 2007 and min 2.25E+02 μg m^-3 in BB 2006). Similarly, sharp biomass peaks were recorded for C. hamatus NP+C3 (max 4.74E+03 μg m^-3 in BB, min 5.84E+01 μg m-3 in SF) (Fig. 7). The biomass of the NP+C3 fraction for all four species observed throughout the study period was much lower than was observed for the C4+AD fraction. There was also a general trend of biomass reduction between 2006 and 2014 for all four species.

The biomass of the four dominant zooplankton taxa in different years (means from all stations for Acartia spp., Pseudocalanus spp., T. longicornis, and C. hamatus) was compared with the long-term mean by season (spring and summer), zooplankton size fraction (NP+C3 and C4+AD), and geographical region (Fig. 8). Intra-annual differences among geographical regions were not statistically significant for the small or large fractions (nested-ANOVA, geographical area nested within a year, p = 0.251 for NP + C1, p = 264 for C4+AD). Considering intra-annual variances among seasons, the biomass of both NP+C3 and C4+AD was higher in spring than in summer. The differences were statistically significant (nested-ANOVA, season nested within a year, p < 0.001 for NP+C3 and C4+AD).

After the 2006-2008 period a decline in zooplankton biomass was noted, especially with regard to large C4+AD zooplankton organisms. The reductions noted in the biomasses of both zooplankton size fractions after 2008 were, however, more distinct in spring than in summer.

The mean water temperature in the surface layer during the period between 2006 and 2014 was 9°C in spring and 17.9°C in summer. In the deep water layer the mean water temperature was 5.3°C in spring and 6.1°C in summer. The greatest inter-annual temperature fluctuations (Fig. 9) were noted in spring in the surface layer with the highest values in 2006, 2007 and 2014 and the lowest in 2010 and 2011. The smallest inter-annual fluctuation was noted in the surface layer in summer.

Figure 9

The statistically significant correlation in spring was found for GR index and biomass of whole fraction of NP+C3 zooplankton, specifically for T. longicornis and C. hamatus. Among four distinguished groups of C4+AD zooplankton fraction only Pseudocalanus spp. had significant (p < 0.05) correlation with GR index of cod larvae and juveniles in the area studied (Table 2; S5). No statistically significant correlation was determined for the GR index of the early developmental stages of cod and the summer zooplankton biomass.

In terms of the surface water temperature, the statistically significant correlation was determined for both size fractions of zooplankton biomasses in spring as in summer (Table 2, S6). Analysing the correlation of each of the zooplankton group separately, only Pseudocalanus spp. in NP+C3 fraction and T. longicornis, Acartia spp., and C. hamatus in C4+AD size fraction proved to be statistically significantly correlated in spring with the surface water temperature. During the summer the biomass of Pseudocalanus spp., T. longicornis and C. hamatus in NP+C3 fraction and T. longicornis and C. hamatus in C4+AD had a significant correlation (p < 0.05) with the surface water temperature. No statistically significant correlation was determined for the GR index of early developmental stages of cod and the surface water temperature during both seasons studied (Table 2; S7).

In this study the mean growth rate of larval and early juvenile cod (n = 491; age 2 – 136 d; SL 4.1 – 39.2 mm) was 0.25 mm/d ± 0.14. Very few articles address the somatic growth rate of larval and juvenile cod from the Baltic Sea, especially the growth rate represented as a SL-at-age relationship. There is least information available about the growth rate of Baltic cod larvae. Oeberst and Böttcher (1998) studied cod juveniles of body lengths between 20 mm and 180 mm TL, and average growth rates oscillated around 0.6 mm/d. A similar value of 0.78 mm/d was presented by Fey and Linkowski (2006) for cod juveniles between 40 mm -153 mm. Huwer et al. (2014) reported the most comparable growth rate data for juvenile cod with body lengths that ranged between 15 mm and 57 mm SL and an age range of between 52 and 101 days (0.86 mm/d). Although it may be concluded that the growth rates of Baltic cod larvae and early juveniles observed in our study were relatively slow, the small size (mostly between 4 mm and 30 mm) of the specimens analysed should be borne in mind.

In this study the growth rates of larval and juvenile cod were not significantly different among the designated sub-areas of the Baltic Sea (Bornholm Basin – BB; Słupsk Furrow – SF; Gdańsk Basin – GB). Differences in larval fish growth rates among these geographical areas are not surprising, but no other data for larval or early juvenile cod are available to compare with these results.

The lack of differences in growth rates among regions correspondeds to the lack of differences in zooplankton biomass. The biomass of Calanoida in the Baltic Sea is regulated by the prevailing abiotic and biotic environmental conditions including changing seasons of the year and hydrological processes in the water column and also predation by higher trophic levels. Differences in the number of Calanoida species and in the proportions among them in zooplankton assemblages occur in the different regions of the Baltic Sea (Diekmann et al. 2012; Otto et al. 2014; Daewel et al. 2014).

Although some tendency was observed towards faster growth in summer than in spring, this pattern was not confirmed when the data were analysed separately within years. The absence of a clear seasonal effect could result from too few spring-time data available for separate within-year comparisons.

The literature lacks information on the topic of differences in the larval and juvenile growth rates of cod caught in the Baltic in spring and summer seasons. Fey and Linkowski (2006) only noted differences in the residuals of size at age among sampling months in juvenile cod specimens. Samples from April were significantly lower than those collected in other months (between December and March). Steffensen (1980), however, reported higher numbers of increments on otoliths in juvenile cod of the same length class caught in spring than in specimens caught in autumn (i.e., faster growth rates in autumn than in spring).

Considering the effect of environmental conditions in this study, the temperature in spring was lower than in summer, which would explain the lower larval growth rates in that season. The beneficial effect of temperature on larval fish growth is well-known (Otterlei et al. 1999). On the other hand, zooplankton biomass was higher in spring than in summer, which should have supported faster larval growth in spring, which it failed to do. The within-year effect of temperature seems to be more important for larval growth than the abundance of zooplankton, especially if zooplankton biomass abundance is sufficient to protect larvae from starvation.

Significant differences in growth rates have been noted in this study for larval and juvenile cod among years, with lower growth rates between 2012 and 2014 than those between 2009 and 2011. Fey and Linkowski (2006) analysed inter-year differences in the size and age of young cod caught in the Gulf of Gdańsk between 1991 and 1995 and found no significant differences. The differences in growth rates observed in this study were not associated with temperature differences among years. There was, however, a statistically significant correlation with among-years differences in biomass of zooplankton in spring. Above all, this was the effect of the small zooplankton fraction of all species (NP+C3) and P. elongatus from the large zooplankton size fraction (C4+AD). The importance of the availability of zooplankton for larval fish growth is clear, and this includes the availability of appropriate zooplankton sizes since, as they grow, cod larvae feed on increasingly larger sizes of food (Jacobsen et al. 2020; Voss et al. 2003; Zuzarte et al. 1996). The condition of juvenile cod in spring is therefore regulated by the availability of the appropriate sizes of zooplankton, the bottom-up process, while in summer juvenile specimens regulate zooplankton biomass, the top-down process (Jacobsen et al. 2019).

In the 1980s climate change caused a rapid decrease in salinity and an increase in mean water temperature in the Baltic Sea environment, and this so-called regime shift triggered changes in the biomass of two main zooplankton groups. Rapid declines in the biomass of the then dominant, Pseudocalanus acuspes were accompanied by increased Acartia spp. biomass (Möllmann et al. 2005, 2008; Otto et al. 2020). This phenomenon was also observed in our data for the period 2006-2014.

Larval cod are not uniformly dispersed in the water column (Huwer et al. 2011). Cod larvae hatch at depths of > 50 m and from here move to shallower zones to begin feeding on nauplii (Huwer et al. 2011). It is important at the moment when exogenous feeding begins that the larvae have food of the appropriate size. The highest density of larval cod in the northern part of the BB in August 2007 was at 25 m to 40 m (Huwer et al. 2011), which corresponded with the depths at which they shift to exogenous feeding (Schmidt and Hinrichsen 2008). Strategically, larval cod begin foraging before the yolk-sac is fully resorbed (Biernaczyk et al. 2016; Jacobsen et al. 2020). In some parts of the Atlantic, such as the Faroe Plateau, larval cod feed on phytoplankton and Calanoida eggs until the yolk-sac is fully resorbed (Jacobsen et al. 2020), while in the BB phytoplankton is not a appropriate food for larval cod, which is likely the result of the lack of spatial and temporal co-occurrence (Zuzarte et al. 1996; Voss et al. 2003). Calanoida copepods are the preferred food of larval and juvenile cod (Economou, 1991; Zuzarte et al. 1996; Voss et al. 2003). In the early stage, larval cod prefer Calanoida nauplii, while older larval and juvenile cod forage on progressively old developmental copepod stages (Jacobsen et al. 2020). Hüssy et al. (1997) reported that copepods were an element of the diets of juvenile Baltic cod up to 70 mm in length. As cod body size increases, the share of copepods in its diet also decreases as the fish shift from pelagic food to inhabiting the bottom. Baltic cod transitioned from pelagic to benthic habitats at a body length of approximately 50 mm (Hüssy et al. 1997).

Köster et al. (2017) observed a decline in the abundancy of larvae in BB in between 2012 and 2014, especially in August. That situation could have an influence on the decline in recruitment success. In the following years, there was a drastic drop in value of recruitment to roughly 1.8 billion in 2015 from approximately 5 billion in 2011 and 2012 (ICES, 2022). The decrease in the GR index of larvae and juvenile cod observed in this study between 2012 and 2014 may not be the direct cause of the reduction in recruitment, but its influence cannot be excluded.