Variation in the population dynamics of small pelagic fishes has been widely documented, and progress has been made over the last few decades in understanding the factors driving long-term fluctuations (e.g., Lluch-belda et al., 1992; Chavez et al., 2003; Alheit et al., 2012). However, such changes occur across a variety of scales (e.g., seasonal, inter-annual, decadal), and encompasses fluctuations in both abiotic (temperature, precipitation, light, nutrients) and biotic processes (e.g., migration, growth, reproduction, trophic interactions).

Pelagic food webs can be highly dynamic across space and time, especially at the regional scale and in coastal areas (Hernández-Carrasco et al., 2018), where riverine inputs, upwelling processes, and anthropogenic stressors can deeply modify environmental conditions and the ecological state of the water column (Stenseth et al., 2006; DiBattista et al., 2020). These alterations can have cascading effects on fish populations and pelagic food webs (Cushing, 1995; Russo et al., 2022; Vassallo et al., 2022).

In this context, an important question is whether variability in trophic relationships at the local scale is more strongly linked to spatial or temporal factors in coastal systems. This question has fundamental implications for the ecological study and management of marine ecosystems, including the provision of ecosystem goods and services and the sustainability of economic activities such as fisheries (Lomartire et al., 2021; Nissar et al., 2023).

Small pelagic fish, due to their rapid growth, short life span, and plankton-based diet, are often considered excellent indicators for assessing changes in resource availability associated with rapid, climate-driven environmental shifts in marine systems (e.g., Drinkwater et al., 2010; Peck et al., 2013). In particular, small pelagic fishes such as anchovy (Engraulis encrasicoulus) and sardine (Sardina pilchardus) are key elements both from an economic point of view, representing important species for commercial fishing (FAO, 2022), and for marine ecosystems, linking lower and upper trophic level species (Cury et al., 2000).

Studies on the trophic ecology of these small pelagic fish have largely been conducted based on stomach content analysis (SCA), which provides insight into the most frequently consumed and abundant prey items (e.g., Borme et al., 2009; Costalago and Palomera, 2014; Costalago et al., 2014). In recent decades, scientists have increasingly used alternative techniques, including stable isotope analysis (SIA), to determine the relative contributions of prey types or trophic positions assimilated over longer periods than those captured by SCA (e.g., Pethybridge et al., 2018). This method provides information about the assimilated prey over a specific time window—ranging from a few days to several months—depending on the tissue type analysed (e.g., liver, white muscle, blood) and its isotopic turnover rate (Hesslein et al., 1993; Bearhop et al., 2004). The ability to isotopically analyse muscle tissue from anchovy and sardine at different life stages offers a time-integrated view of dietary shifts, as this tissue has an isotopic turnover rate ranging from several days to one month (Suzuki et al., 2005; Guelinckx et al., 2007; Buchheister and Latour, 2010). Because stable isotope enrichment occurs along the trophic web, nitrogen stable isotope (δ¹⁵N) have been utilised as an indicator of the trophic position of organisms along the food web, as well as an indicator of the seasonal and spatial variability in trophic habits (Tanaka et al., 2008; Malzahn and Boersma, 2009). Similarly, as different energy sources often have distinct δ¹³C signatures, carbon stable isotopes can be used to trace the original organic carbon sources at the base of the food web, assess seasonal variation, or detect geographic differences (Vander Zanden and Rasmussen, 1999; 2001; Overman and Parrish, 2001; Castro et al., 2020, and references therein).

In the present study, we compared the isotopic carbon and nitrogen (δ¹³C and δ¹⁵N) signals in white muscles of juvenile and adult individuals of both sardine and anchovy to investigate the spatiotemporal variability in their feeding behaviour across two Mediterranean areas—the Tyrrhenian Sea and the Strait of Sicily —in the years 2011, 2014, 2015, and 2017. These two areas differ markedly in terms of hydrographic conditions, morpho-bathymetry, and productivity (Bonanno et al., 2016; 2018). The SS area is characterised by highly dynamic circulation, a coastal upwelling system (Bonanno et al., 2014), low coastline complexity, minimal human impact, and limited terrigenous input from river runoff. In contrast, the TR area features a large continental shelf, strong riverine influence from several rivers (e.g., Tiber, Garigliano, Volturno, Sarno, and Sele), and high levels of urban, agricultural, and industrial activities as well as a high degree of coastline complexity with the presence of Gulfs.

To our knowledge, only a few studies have evaluated the spatiotemporal variability in the feeding habits of these two small pelagic species in the Mediterranean Sea (e.g., Bachiller et al., 2020; Navarro et al., 2020; Zorica et al., 2021; Fanelli et al., 2023), and these have generally analysed only 1 year of data and/or focused on areas located in close geographic proximity.

Taking all these factors into account, we collected fish samples from both consecutive and non-consecutive years in order to assess the spatiotemporal variability in the feeding habits of anchovy and sardine. We first analysed the inter- and intraspecific differences in the isotopic composition of both species across life stages, areas, and years, focusing on the degree of overlap and the width of the isotopic niches. In a second step, we assessed which factor had the greatest influence on their feeding behaviour (i.e., place or time).

Summary statistics for δ¹⁵N and δ¹³C values (mean, standard deviation, median, interquartile range, minimum/maximum values, and number of observations) by area, year, species, and life stage are reported in Table 1 . Overall, the differences between the mean and median values of δ¹⁵N and δ¹³C were small (maximum absolute differences were 0.6‰ and 0.2‰, respectively), highlighting the rather symmetric distributions for both δ¹⁵N and δ¹³C values in each group (i.e., area, year, species, and life stage).

Pearson correlation analyses between δ¹⁵N (and δ¹³C) fish total length (TL) showed both positive and negative correlations for anchovy and sardine across years and areas, when significant (p < 0.05) ( Supplementary Figure S2 ). For anchovy in the SS area, δ¹⁵N showed a positive correlation in 2014 and 2017, while δ¹³C showed a negative correlation in 2017 for δ¹³C; in TR it was positive in 2017 but negative in 2011 and 2014 in terms of δ¹⁵N and always negative in 2014, 2015 and 2017 for δ¹³C. For sardines, when significant, the correlations were always negative in terms of δ¹⁵N in both areas (2015 in SS; 2011 and 2017 in TR). For δ¹³C, the correlation was negative in 2014 and 2015 and positive in 2017 in SS, and positive only in 2015 in TR.

Preliminary fits of GLS models highlighted the effect of unequal variance and the lack of independence of residuals for both δ¹⁵N and δ¹³C. Consequently, the within-group heteroscedasticity structure was modelled according to area and year factor, while the within-group correlation structure accounted for the autocorrelation at catch operation level. The predictors selected for δ¹⁵N and δ¹³C best models were area and year factors as well as different interactions involving length class and other factors ( Table 2 ).

The obtained results showed that differences in isotopic signals between the two species were greater between areas than within areas across all sampled years—except for δ¹³C in adult anchovy in 2015 ( Supplementary Table S2 ). Specifically, both nitrogen and carbon stable isotopes showed higher values in the Tyrrhenian area (TR) compared with the Strait of Sicily (SS) for both species ( Figures 3 , 4 ). It is known that inter-specific variations in the feeding habits of marine fishes are consistent with inter-specific differences associated to local factors that provide variations in the baseline isotopic values in food webs and affect prey presence and abundance (e.g., Costalago and Palomera, 2014; Costalago et al., 2014; Rumolo et al., 2016).

The SS area is less impacted by continental waters, with fewer major cities, effluents, and coastal industries compared to the TR area. Coastal upwelling along the southern coast of Sicily and the presence of offshore banks (Piccioni et al., 1988) generate organic matter during upwelling events, which—after remineralisation—supports primary production (Alvarez-Salgado et al., 2002; Bode et al., 2004). In sectors as the SS area, phytoplankton generally serves as the main foundation of the food web. Productivity is primarily regulated by oceanographic processes such as upwelling and seasonal changes in stratification and mixing, which fertilise the euphotic zone (Bode et al., 2006).

Although phytoplankton in upwelling regions often becomes enriched in ¹⁵N (Montoya et al., 1990; Goering et al., 1990; Rolff, 2000), resulting in enriched zooplankton, Rumolo et al. (2016) reported low δ¹⁵N values in particulate organic matter collected in the SS area. This may be explained by the discontinuous nature of coastal upwelling along the southern coast of Sicily (Piccioni et al., 1988; Bonanno et al., 2014), where summer phytoplankton growth likely relies on regenerated nitrogen—typically depleted in ¹⁵N (Montoya et al., 1990). This mechanism may account for the lower δ¹⁵N values found in both anchovy and sardine in the SS area.

In addition, the mean isotopic signals of anchovy and sardine in the SS area were lower than those reported in the literature ( Table 3 ) for other parts of the Mediterranean Sea or similar upwelling-influenced oceanic systems (Bode et al., 2003; Sholto-Douglas et al., 1991; Costalago et al., 2012; Le Bourg et al., 2015; Bachiller et al., 2020, 2021; Lloret-Lloret et al., 2022; Fanelli et al., 2023). In this context, the absence of a consistent pattern in the correlation ( Supplementary Figures S1 , S2 ) between the δ¹⁵N and fish total length (positive, null, or negative across years for both species; Supplementary Figures S1 , S2 ), along with the narrow δ¹⁵N range, reinforces the hypothesis that anchovy and sardine in the SS area may exploit similar resources—possibly feeding directly on particulate organic matter with low isotopic signatures and/or on prey occupying lower trophic positions depending on resource availability.

However, while juvenile and adult sardines in the SS area exhibited similar feeding behaviour, as indicated by their isotopic niche overlap, anchovies showed greater feeding variability between life stages, with relatively weak isotopic overlap (never exceeding 30%; Supplementary Table S4 ). These results align with previous studies (e.g., Tudela and Palomera, 1997; Costalago et al., 2012), which describe juvenile and adult sardines as opportunistic filter-feeders with a more heterogeneous diet than anchovy. For anchovy, several authors have reported high variability in feeding between life stages, sometimes involving ontogenetic shifts (Bacha and Amara, 2009; Pepin and Penney, 1997; Scharf et al., 2000; Bachiller and Irigoien, 2013; Palomera et al., 2007; Plounevez and Champalbert, 1999, 2000), while others reported no ontogenetic shift (Bulgakova, 1993; Tudela and Palomera, 1995, 1997; Plounevez and Champalbert, 1999, 2000; Borme et al., 2009; Van der Lingen et al., 2009; Bacha and Amara, 2009; Catalan et al., 2010).

Regarding δ¹³C, which is generally considered a conservative tracer of the primary producers at the base of the food web and thus a useful indicator of foraging habitat (France, 1995; Hobson, 1999), the values presented in the SS area for both species were comparable to those reported from other Mediterranean regions (see Table 3 ). This suggests that anchovy and sardine in SS likely explored similar habitats and shared food sources.

In contrast to the SS area, the higher δ¹⁵N and δ¹³C values observed in anchovy and sardine in the TR area may be attributed either to a strong feeding preference for larger prey items or to feeding within an environment characterised by a more isotopically enriched baseline. Although particulate organic matter (POM) and zooplankton samples were not collected during the 2014, 2015, and 2017 surveys, this study focused on isotopic niche width and overlap in consumers (anchovy and sardine), which revealed greater differences between areas than among years. These area-based differences were also confirmed by Rumolo et al. (2016), who, analysing isotopic data from a single year in the Strait of Sicily and the Tyrrhenian Sea, found higher δ¹⁵N values of POM in the Tyrrhenian Sea. This was attributed to riverine input supplying high levels of inorganic nutrients, which can lead to ¹⁵N enrichment in POM (e.g., McClelland et al., 1997), and in turn to enriched δ¹⁵N in zooplankton and fish.

Moreover, Rumolo et al. (2016) reported higher isotopic values in adult anchovies and sardines from the Tyrrhenian Sea compared to other Mediterranean regions, linking these values to a higher dietary contribution of mesozooplankton species, particularly chaetognaths and euphausiids. Although the presence of chaetognaths and euphausiaceans such as Meganyctiphanes norvegica in the TR area has not been extensively studied, Wiebe and D’Abramo (1972) found Euphausia kcrohnii, Nematoscelis megalops, Meganyctiphanes norvegica, and Stylocheiron abbreviatum to be abundant in the western Tyrrhenian basin, especially in May and June. Moreover, Mussi et al. (2004) reported the presence of crustacean exoskeletons—belonging to Meganyctiphanes norvegica—in fin whale faecal material collected near Ischia Island ( Figure 1 ), close to the Cuma submarine canyon, where marine mammals have been repeatedly observed (Mussi et al., 1999; 2014).

The co-occurrence of large predators and small pelagic fish in the TR area (Leonori et al., 2021) is likely linked to the abundance of large zooplankton species, potentially triggered by the high load of nutrients provided by different mechanisms such as submarine canyon advection and river discharge. Isotopic values of anchovy and sardine in the TR area—particularly δ15N—were consistent with those reported in the northern Adriatic Sea (Table 3), a region with comparable environmental characteristics, including the presence of large rivers and densely populated coastlines (Fanelli et al., 2023; Micheli et al., 2013). In the TR area, anchovy showed no consistent correlation between δ¹⁵N and total length (TL); in some cases, higher δ¹⁵N values were found in juveniles than in adults ( Supplementary Table S4 ). For sardine, results indicated a negative correlation between δ¹⁵N and TL, though not always statistically significant. This trend is difficult to interpret due to the limited availability of isotopic and stomach content data for juvenile small pelagic fish in the Mediterranean (Bachiller et al., 2020). However, some authors (Fanelli et al., 2023; Bode et al., 2003 and references therein) found this negative trend in small pelagic fish interpreting such relationship as the consequence of the increasing consumption of phytoplankton by adults, as most clupeids are known to employ both particle-feeding and filter-feeding (James, 1988).

In the TR area, δ¹³C values were, on average, more enriched than in other areas of the Mediterranean Sea, suggesting that both species may explore more coastal areas and carbon sources ( Table 3 ). However, our results, which examined dietary variations based on assimilated food (in about one month), did not allow us to estimate individual specialisation linked to habitat use. Such an assessment would require considering the temporal and spatial variability in the δ¹³C of both prey and consumers.

Taking into account the investigated inter- and intraspecific differences in isotopic composition—by comparing area, life stage, and year—the patterns of isotopic niche overlap and width for the two species in each area and between areas ( Supplementary Tables S2–S6 ) showed the presence of different feeding systems explored by anchovy and sardine. Although the greater niche overlap observed in TR compared to SS may indicate high interspecific competition, the broader isotopic niche widths indicate dietary plasticity. This plasticity may enable the development of alternative feeding strategies (e.g. division of the niche between juveniles and adults) in response to competition (Andrades et al., 2019; Cathcart et al., 2019).

Moreover, the higher isotopic values and trophic overlap observed in sardine and anchovy in the TR area could also be related to increased consumption of larger prey (e.g., decapods and euphausiids), as previously reported by Rumolo et al. (2016), and/or driven by the presence of other high-trophic-level prey such as siphonophores and jellyfish, as noted by Bachiller et al. (2020) in the western Mediterranean Sea. In contrast, in the SS area, the low isotopic values together with smaller isotopic niche width, suggested the presence of lower trophic levels organisms and less food resources explored by the fish species. Furthermore, the variability observed in isotopic niche overlap may reflect differences in the spatial distribution of the two fish species (Barra et al., 2015).

In addition to the clear difference evidenced between the SS and TR areas, it is important to highlight that observed differences are temporally stable, indicating a certain degree of stability in feeding behaviour. This stability likely reflects underlying ecosystem characteristics rather than interannual variability in environmental conditions. Nonetheless, the high variability observed in the SS area in terms of niche overlap, along with the lower variability in isotopic values compared to TR—suggests that the SS area may be a more sensitive ecosystem, more vulnerable to environmental disturbances, including those induced by climate change.

Given the key role of small pelagic fish in energy transfer along the marine food web, careful management of these resources is should be carefully considered. On the contrary, the TR area could appear to be a more robust ecosystem, capable of handling natural or anthropogenic disturbances. To support this hypothesis, future studies should include isotopic analyses of POM and zooplankton assemblages, using a more targeted sampling design spanning multiple years.