Defining Seasonal Functional Traits of a Freshwater Zooplankton Community Using δ¹³C and δ¹⁵N Stable Isotope Analysis

Abstract

Functional-based approaches are increasingly being used to define the functional diversity of aquatic ecosystems. In this study, we proposed the use of δ¹³C and δ¹⁵N stable isotopes as a proxy of zooplankton functional traits in Lake Maggiore, a large, deep subalpine Italian lake. We analyzed the seasonal pattern of δ¹³C and δ¹⁵N signatures of different crustacean zooplankton taxa to determine food sources, preferred habitats, and trophic positions of species throughout one year. The cladocerans Daphnia longispina galeata gr., Diaphanosoma brachyurum, and Eubosmina longispina were grouped into a primary consumer functional group from their δ¹³C and δ¹⁵N isotopic signatures, but while the former two species shared the same food sources, the latter exhibited a more selective feeding strategy. Cyclopoid copepods occupied a distinct functional group from the other secondary consumers, being the most ¹⁵N enriched group in the lake. The δ¹⁵N signature of calanoid copepods showed trophic enrichment in comparison to Daphnia and Eubosmina and linear mixing model results confirmed a predator-prey relationship. In our study, we have demonstrated that the use of δ¹³C and δ¹⁵N stable isotopes represented an effective tool to define ecological roles of freshwater zooplankton species and to determine functional diversity in a lake.

1. Introduction

Functional-based approaches are increasingly being used to study aquatic ecosystems as an alternative to traditional taxonomy-based approaches. Functional diversity is a biodiversity measure based on the ecological role of the species present in a community. Species-specific functional traits, or “what they do” [1], allows species to be defined by their interactions within an ecosystem [2] in terms of their ecological roles and how they interact with the environment and with other species [3].

Many recent ecological studies [2,4,5] suggest the importance of species ecological roles, and not just the number of taxonomic species, in the relationship between biodiversity and ecosystem functioning. This is a central concept if we are to understand and predict the resilience of a community to perturbations. If two species are deemed to be functionally alike and to occupy a similar trophic niche [6,7], the loss of one of those species is not likely to have an impact on the resource pool, as the other will increase its activity accordingly. The loss of functional diversity in the ecosystem is mitigated, as the species lost does not possess unique functional traits. Thus, the sum of organism functional traits within an ecosystem can be said to represent an indirect measure of its functional diversity [1,8]. Although the importance of functional diversity is widely recognized, there is no consensus on how to quantify functional diversity within a community, as relationships between the various indices have not yet been established [9].

Freshwater zooplankton play a key role in aquatic ecosystems in the transferal of biomass and energy from phytoplankton to top predators, e.g., [10,11]. While studies on phytoplankton have increasingly highlighted the importance of functional traits and functional classification in ecological studies, only a few have attempted this approach with zooplankton. Barnett et al. [2] applied the measure of functional diversity [12] to crustacean zooplankton communities. Quantitative functional traits considered were the C:N ratio, mean body size, and preferred food size range, while qualitative traits described the preferred habitat, food selectivity, and trophic position of each species.

Analysis of δ¹³C and δ¹⁵N stable isotopes is widely used to quantify food sources, trophic positions, and the interactions of organisms, as the δ¹³C of a consumer can infer the assimilated source of dietary carbon [13,14] and δ¹⁵N the trophic role [14,15]. Considering a functional-based perspective, δ¹³C and δ¹⁵N could identify ecological relations among taxa, such as competition and predation, ecological niche, habitat preference, and taxa redundancy, and help to define functional groups in an ecosystem. In this study, we propose the use of δ¹³C and δ¹⁵N stable isotope analysis to quantify some of the “qualitative functional traits” [2] for pelagic crustacean zooplankton taxa in Lake Maggiore, a large, deep subalpine lake in Italy. As species composition, diversity, and biomass of zooplankton can change significantly seasonally, especially in temperate lakes [16], the seasonal variation of δ¹³C and δ¹⁵N was used to define seasonal patterns in habitat preference, food sources, and taxa trophic position. This allowed an interpretation of taxa redundancy and a hypothesis of bottom-up and top-down mechanisms potentially driving the observed changes.

One use of stable isotopes is to determine the proportional contributions of several sources in a mixture. An example of a source proportion calculation includes the determination of various food sources in an animal’s diet [17]. Linear mixing models are used to estimate proportions for two sources using isotopic signatures for a single element (e.g., δ¹³C), or for three sources using isotopic signatures for two elements (e.g., δ¹³C and δ¹⁵N; [18]). In this study, we have used a linear mixing model [18,19] in order to discriminate the relative contribution of different preys in the diet of zooplankton consumers such as calanoid and cyclopoid copepods and the predatory cladoceran Leptodora kindtii. This further contributed to the understanding of food preference and taxa trophic position in the zooplankton community of Lake Maggiore.

2. Materials and Methods

Lake Maggiore (45°57′ N 8°32′ E 3°47′ W) is the second deepest (d_max 370 m) and largest (area 212.5 km², volume 37.5 km³) subalpine lake in Italy. Being phosphorus-limited (TP_max ca. 10 μg L⁻¹), the lake is oligotrophic and has recovered from eutrophication of the late 1970s [20,21].

Except for September, vertical zooplankton hauls from the surface to a depth of 50 m were collected. This follows the standard routine sampling for deep subalpine lakes, in which samples are collected within the upper 50 m depth, as previous research on the vertical distribution of zooplankton showed that this is the water layer in which zooplankton live [22]. Monthly samples were collected from April to November 2009, when total zooplankton biomass was ≥3 mg·m⁻³, using a wide-mouth 450 μm mesh zooplankton net of diameter 0.58 m, filtering 13 m³ lake water from three pelagic stations (G: 45°58′30″ N 8°39′09″ E, B: 45°54′28″ N 8°31′44″ E, L: 45°49′70″ N 8°34′70″ E) [16].

Zooplankton samples for the quantification of biomass were collected with a Clarke-Bumpus plankton sampler of a 126 μm mesh size and fixed in ethanol 96% to estimate the taxa-specific population density (ind·m⁻³) and standing stock biomass (dry weight, mg·m⁻³) [22]. Organisms for isotopic analyses were kept overnight in filtered (1.2 μm GF/C filters) lake water for gut clearance, before sorting into taxa and quantities suitable for isotopic analyses. The taxa analyzed were Daphnia longispina galeata gr., Eubosmina longispina, Diaphanosoma brachyurum, Bythotrephes longimanus, Leptodora kindtii, adults of the calanoid copepods Eudiaptomus padanus and Eudiaptomus gracilis, and of the cyclopoid copepods Mesocyclops leuckarti and Cyclops abyssorum.

Samples were oven-dried for 24 h at 60 °C, before homogenizing and transferal into tin capsules of 5 × 9 mm in size. Depending on body mass, 50 to 700 individuals of each taxa were pooled to reach a minimum dry weight (DW) of 1 mg per sample. Three replicates of each taxa were run from each of the three sampling stations, as among-station differences were statistically non-significant (p > 0.05, Friedman Analysis of Variance, ANOVA test; [16]). The isotopic composition of organic carbon and nitrogen was determined from the analyses of CO₂ and N₂ by the G. G. Hatch Stable Isotope Laboratory at the University of Ottawa, Ontario, Canada, using a CE 1110 Elemental Analyser (Vario EL III manufactured by Elementar, Germany) and a DeltaPlus Advantage isotope ratio mass spectrometer (Delta XP Plus Advantage manufactured by Thermo, Bremen, Germany) coupled to a ConFlo III interface (Conflo II manufactured by Thermo, Bremen, Germany). The standard deviation of the analyses (SD) based on laboratory internal standards (C-55) was < 0.2‰ for both δ¹³C and δ¹⁵N. Isotope ratios were expressed as the parts per thousand (‰) difference from a standard reference of PeeDee Belemnite for carbon and atmospheric N₂ for nitrogen:

$$\mathsf{\delta }{}_{}{}^{13}\mathrm{C},\mathsf{\delta }{}_{}{}^{15}\mathrm{N}=\left[\left(\frac{{\mathrm{R}}_{\mathrm{sample}}}{{\mathrm{R}}_{\mathrm{standard}}}\right)-1\right]\times 1000$$

(1)

where R is the isotopic ratio: ¹³C/¹²C and ¹⁵N/¹⁴N.

Lipids can be δ¹³C -depleted as a consequence of fractionation during lipid synthesis [23], which can lead to a misrepresentation of results as differences in predator-prey δ¹³C could be greater than the expected 0.8‰ [24]. The C:N ratio was used as an indicator of lipid content. Invertebrates, including crustacean zooplankton, tend to have a C:N ratio of 4 [25], but C:N varies seasonally [26] and was as high as 7, so we used a revised version of the lipid normalizing procedure based on the C:N ratio [27], substituting the corrected parameters into Equations (2) and (3):

$$\mathrm{L}=\frac{93}{1+{\left(0.246\times \left(\mathrm{C}÷\mathrm{N}\right)-0.75\right)}^{-1}}$$

(2)

$$\mathsf{\delta }{}_{}{}^{13}\mathrm{C}\prime =\mathsf{\delta }{}_{}{}^{13}\mathrm{C}+\mathrm{D}\times \left(\mathrm{I}+\frac{390}{1+\frac{287}{\mathrm{L}}}\right)$$

(3)

where L is the proportion of lipid in the sample; C and N are the proportions of carbon and nitrogen in the sample, respectively; δ¹³C′ is the lipid normalized sample signature; δ¹³C is the measured sample signature; D is the isotopic difference between the protein and lipid (7.018 ± 0.263); and I is a constant of 0.048 ± 0.013 [28].

Zooplankton δ¹³C and δ¹⁵N isotopic signatures were referred to that of the pelagic baseline, which was expressed by the primary consumer, Daphnia longispina galeata gr. This choice came from previous stable isotope studies in Lake Maggiore [16], showing that Daphnia’s δ¹³C signature in the different seasons was closely correlated with the signature of seston (r = 0.86; p < 0.01; N = 13), confirming that Daphnia was an appropriate proxy for the pelagic baseline against which the carbon isotopic signals of other zooplankton can be compared. Δ¹³C was used to detect seasonal changes in taxa specific feeding behavior and assess the origin of carbon sources fueling the pelagic food web.

The carbon fractionation between consumer and resource (F = δ¹³C_cons − δ¹³C_diet) is ≤ 0.8‰ (±1.1‰ S.D.) [24]. The δ¹⁵N of consumers has been shown to be enriched 2.55‰ [29] for zooplankton, and was used to assess seasonal change in taxa-specific trophic position (T), as a consumer’s carbon signature is related to the baseline (F ≤ 0.8 ± 1.1) by:

T = (E/λ) + 2

(4)

where λ is the stepwise enrichment, E = 2.55‰ [29], and 2 is the value commonly assigned to the deviation of primary consumers from the pelagic isotopic baseline. A trophic level of T = 3 indicates that a consumer is feeding on a primary consumer, whereas T = 4 suggests that there is an intermediate prey.

When δ¹³C of the predator lies between that of two different prey taxa, suggesting a simultaneous use of both sources, the percent carbon contribution (p; q) of each prey to the predator’s diet was calculated by the 2-endmember linear mixing model (2-em LMM), [18,19] as:

p = (δ¹³C_predator − δ¹³C_prey2)/(δ¹³C_prey1 − δ¹³C_prey2); q = 1 − p

(5)

where p and q are the relative contributions (%) of prey1 and prey2 carbon signatures to the predator δ¹³C carbon signature (δ¹³C_predator).

When three potential prey sources were assessed, their isotopic signatures were partitioned by applying a 3-end member mixing model [18] to calculate the fractional contribution (p; q; z) of each of the three food sources to the predator’s diet as:

p = ((δ¹⁵N_prey3 − δ¹⁵N_prey2)(δ¹³C_predator − δ¹³C_prey2) − (δ¹³C_prey3-δ¹³C_prey2)(δ¹⁵N_predator − δ¹⁵N_prey2))/

((δ¹⁵N_prey3 − δ¹⁵N_prey2)(δ¹³C_prey1 − δ¹³C_prey2) − (δ¹³C_prey3 − δ¹³C_prey2)(δ¹⁵N_prey1 − δ¹⁵N_prey2));

q = ((δ¹³C_predator − δ¹³C_prey3) − (δ¹³C_prey1 − δ¹³C_prey3)p/(δ¹³C_prey2 − δ¹³C_prey3);

z = 1 − p − q

(6)

where p, q, and z are the relative carbon contribution (%) of prey (prey 1, 2, and 3). As required by the mixing model, δ¹³Cpredator and δ¹⁵Npredator were corrected for trophic fractionation, by weighting the isotopic signature of the prey against their percentage contribution to total biomass on each sampling date.

The software IsoError 04 (https://www.epa.gov/eco-research/stable-isotope-mixing-models-estimating-source-proportion) [18] was used to perform all 3-end member LMM calculations. Statistical analyses (Shapiro-Wilkinson W-test, Spearman-Rank correlation, Hierarchical Cluster Analysis) were performed using the software Statistica 12 (version 12, TIBCO Statistica Company, Palo Alto, CA, USA) and Sigmaplot 11.0 (version 11, Systat Software Inc., San Jose, CA, USA).

Cluster analysis of the seasonal variation in δ¹³C and δ¹⁵N for each taxa was performed with the software Sigmaplot 11.0. An Euclidean distance measure was used as the data were continuous [30].

3. Results

3.1. Seasonal Variation in δ¹³C and the Determination of Consumer Resources

The seasonal variation of δ¹³C for the different zooplankton species is shown in Figure 1. The variation of δ¹³C in Daphnia was most depleted in spring, with a value of −36.3‰ ± 0.6 (SD), and became most enriched during the summer in August, with a value of −26.0‰ ± 0.1 (SD). Diaphanosoma was present in the lake in October and November and its δ¹³C signature overlapped with a strong correlation (p = 0.06; R = 0.66) with the δ¹³C signature of Daphnia. Also, the δ¹³C signature of Bosmina overlapped with that of Daphnia in June, July, and October. This similarity in δ¹³C indicates that herbivorous cladocerans are sharing the same type of resources. When more δ¹³C-depleted values were recorded for Bosmina (−32.0‰ ± 0.1, SD) than for Daphnia (−29.7‰ ± 0.1, SD), this might suggest a deviation in dietary sources over the winter months (November).

The δ¹³C signature of copepods was more δ¹³C -depleted than that of the pelagic baseline (annual mean of −26.4‰ ± 4.4, SD) all year round, with a δ¹³C annual mean of −32.4‰ ± 4.6 (SD) for calanoid copepods and −34.5‰ ± 2.4 (SD) for cyclopoid copepods.

The procedure of lipid normalization for cyclopoid copepods decreased their average value of the δ¹³C signature, but the relative values between species were unaffected. The seasonal trend of δ¹³C signature for calanoid copepods generally traced the δ¹³C pelagic baseline signature (F < 0.8). The δ¹³C signature of Bythotrephes overlapped with the pelagic baseline (F < 0.8) from June to October, suggesting a tight dependence of this predatory cladoceran on pelagic food resources. However, in November, the δ¹³C signature of Bythotrephes was less δ¹³C -depleted (−26.8‰) and the least ¹³C-depleted of all zooplankton taxa. The δ¹³C signature of the other predatory cladoceran Leptodora was also related to the pelagic baseline (F < 0.8).

3.2. Seasonal Variation in δ¹⁵N and Determination of Trophic Levels

The δ¹⁵N seasonal pattern was the same for all zooplankton taxa, with more N-enriched values in early spring and late autumn, and less δ¹⁵N-enriched values during the warm months (Figure 2). The seasonal variation in δ¹⁵N was more pronounced in primary than in secondary consumers, as the predatory taxa Bythotrephes, Leptodora, and cyclopoid copepods had a small δ¹⁵N range (NR) of 1.27‰, 1.32‰, and 1.74‰, respectively, whilst Daphnia, Bosmina, Diaphanosoma, and Eudiaptomus had a wider NR range of 3.75‰, 3.32‰, 3.57‰, and 2.50‰, respectively.

Cyclopoid copepods were the most δ¹⁵N-enriched group, with δ¹⁵N ranging between 7.92‰ in October and 9.67‰ in May, and with a F_max of 6.6‰.

On average, predatory cladocerans were more δ¹⁵N-enriched than the pelagic isotopic baseline. The δ¹⁵N signature of Daphnia, Bosmina, and Diaphanosoma overlapped (≤5‰), with δ¹⁵N of the latter two taxa significantly correlated with Daphnia (p < 0.001; Spearman Rank Correlation coefficient R > 0.9; N = 8, 12, respectively).

Calanoid copepods had enriched δ¹⁵N signatures, ranging from 5.81‰ in October to 9.38‰ in April. Δ¹⁵N enrichment with respect to the pelagic baseline varied between a maximum of 4.45‰ in October and a minimum of 2.65‰ in April, suggesting a change in trophic feeding level and differential exploitation of food resources.

Using mixing models to quantify the contribution of different prey to predators’ diet, the contribution of the zooplankton preys assimilated by the consumer did not always match the zooplankton prey biomass present in the lake (Figure 3a,b). For example, in July, the estimated proportion of Daphnia and Bosmina in the diet of calanoid copepods was 53.2% and 46.8%, respectively, when Daphnia was present in the lake with a biomass of 80%. In October and November, the diet of calanoids copepods maintained a similar estimated contribution of Daphnia 57.7% and Bosmina 42.3%, when Daphnia was present in the lake with a biomass of 16.9% of total zooplankton biomass.

The diet of Leptodora in July was partitioned between Daphnia and Bosmina with 87.2% and 12.8%, respectively, while in autumn, the 3-linear mixing model indicated a higher consumption of Diaphanosoma (76.6%) than Daphnia (5.9%) and Bosmina (17.5%).

The diet of cyclopoid copepods in October and November was estimated by the mixing model to be 83.1% Bosmina and 16.9% Daphnia, when Bosmina was present in the lake, representing 75% of the total zooplankton biomass.

3.3. Determination of Functional Roles from δ¹³C and δ¹⁵N

Cluster analysis of the seasonal variation in δ¹³C for each taxon (Figure 4a) identified three major groups. The first split of the ordination separates the cyclopoid copepods from the other taxa, most likely as a result of the group utilizing deeper carbon sources in the pelagic zone. Lipid correction applied to δ¹³C did not affect the seasonal pattern observed in the taxa. The second and further split in the cluster analysis grouped together the primary herbivorous cladocera Daphnia longispina galeata gr., Diaphanosoma brachyurum, and Eubosmina longispina, and the secondary consumers Bythotrephes longimanus, Leptodora kindtii, and calanoid copepods.

Cluster analysis of the seasonal variation in δ¹⁵N for each taxon (Figure 4b) clearly grouped the taxa into two functional groups, the primary consumers, D. longispina galeata gr., D. brachyurum, and E. longispina, and the secondary consumers B. longimanus, L. kindtii, and copepods.

4. Discussion

Functional biodiversity measures are based on the functional traits of the species in a community, rather than species richness, as it is interactions that determine the response of the ecosystem to a perturbation [31,32,33]. A methodology to determine the functional diversity within an ecosystem is quantifying the functional relationships between species as the distance measure of the branch length of the connecting functional dendrogram [12]. But it is not just the distance measure or clustering algorithm used that is important, it is the choice of which functional traits to include that is crucial [34].

In this study, we propose to use δ¹³C signature as a proxy for determining habitat preference and foraging zone by inferring pelagic vs. littoral feeding preferences. We further analyzed the relationship between δ¹³C taxa signature and phytoplankton succession during the studied year of 2009.

Confounding factors in a dynamic ecosystem like Lake Maggiore include seasonality and predation pressure variation. Seasonal changes in the littoral and pelagic δ¹³C isotopic baseline were identified by a temporal shift towards less ¹³C-depleted values in the summer, which is a trend commonly observed in thermally-stratified lakes [35,36]. In lacustrine systems, δ¹³C and δ¹⁵N of suspended particulate matter (seston) varies seasonally [37] because of differences in the allochthonous input, phytoplankton species composition, and primary productivity [38,39]. In Lake Maggiore, it has been demonstrated that δ¹³C of the primary consumer Daphnia longispina galeata gr. tracks the isotopic composition of seston (50 ≤ size ≤ 126 μm) and of the pelagic baseline [16,35]. The choice of Daphnia as a useful isotopic baseline is supported by other studies [27], considering this taxa as a short-lived organism suited for fine scale temporal integration of pelagic δ¹³C or δ¹⁵N signatures.

If we consider the phytoplankton component of seston, there is variation in the fractionation of δ¹³C between phytoplankton groups, with Chrysophycee and Bacillariophycee being more depleted in δ¹³C than Cyanobacteria [37]. In Lake Maggiore, Bacillariophycee are the most dominant phytoplankton group with a high biomass throughout the year [37,40]. In the year of our study (2009), a peak in the Bacillariophycee biomass of Aulacoseira sp., Asterionella sp., and Fragilaria sp. was recorded in spring [37], when the cladoceran Daphnia longispina galeata gr. had the most δ¹³C-depleted values, while when Cyanobacteria biomass increased in the lake during the summer and autumn, Daphnia longispina galeata gr. had more δ¹³C-enriched values. This is suggestive of the opportunist nature of Daphnia longispina galeata gr., feeding on the most abundant component of the phytoplankton.

Rather than simply attributing organisms to different trophic groups per se from isolated values of δ¹⁵N, we analyzed the seasonal variation in δ¹³C and δ¹⁵N fractionation to identify changes in trophic interactions and feeding niche change. Daphnia longispina galeata gr., Diaphanosoma brachyurum, and Eubosmina longispina can be grouped into a “primary consumer” functional group from their δ¹³C and δ¹⁵N isotopic signatures and known trophic interactions. Determination of functional groups is important and consensus is being placed on the importance of how species are likely to react to a perturbation, which is particularly crucial in considering the consequences of species loss. In the cluster analysis we performed, the grouping of the δ¹³C and δ¹⁵N isotopic signature of the herbivorous cladocerans Daphnia longispina galeata gr. and Diaphanosoma brachyurum, indicates a use of the same seasonal carbon pool. This is also confirmed by the similar values of the two cladoceran taxa in their δ¹³C signature. The position of Eubosmina in the cluster analysis might show a slightly different use of food resources, a hypothesis confirmed by its δ¹³C depletion values observed in November. These results indicate a deviation in dietary sources and a possible shift towards a more selective feeding strategy for phytoplankton [29,41] over the winter, as phytoplankton tends to be more δ¹³C-depleted than detritus [28,41,42].

A redundant species can be generically defined as one that co-exists with an overlapping functional role and trophic niche to other species [5] in a community.

Bythotrephes is a cladoceran known to actively predate on Daphnia [43], which was indicated by its trophic δ¹⁵N enrichment in comparison to Daphnia. However, in our study, when Daphnia was scarce in the lake in October and November, Bythotrephes seemed to switch its feeding preference to Diaphanosoma as its prey. This diet switching from Daphnia to Diaphanosoma during periods of Daphnia scarcity was also observed for Leptodora, as the δ¹⁵N of Leptodora did not indicate that they fed on calanoid or cyclopoid copepods, rather that they preferentially exploit cladocerans.

Thus, the results in our study might indicate that Daphnia longispina galeata gr. and Diaphanosoma brachyurum are redundant species in the zooplankton community of Lake Maggiore, as their trophic roles are overlapping and interchangeable. Our hypothesis is further supported by the well-known interspecific competition of these two taxa for the same food sources in lakes [44].

The δ¹³C signature of calanoid copepods in our study was more depleted than the Daphnia δ¹³C signature. This has been observed in previous studies on freshwater zooplankton [28,38,41], and has been attributed to the calanoid copepods omnivorous diet of seston and small zooplankton like Eubosmina and rotifers [45,46,47,48].

The δ¹⁵N isotopic signature of calanoid copepods showed a trophic stepwise enrichment in comparison to Daphnia and Eubosmina, indicative of a predator-prey relationship. This was confirmed by the results of the linear mixing model [18], showing that calanoid copepods preyed on Daphnia and Eubosmina in similar percentages in spring and in autumn. Although herbivorous during some life stages, copepods are one of the main groups of invertebrate predators in limnetic and littoral inland waters [45].

In our study, cyclopoid copepods occupied a distinct functional group from the other secondary consumers, being the most ¹⁵N -enriched zooplankton group in the lake during the year. We recorded highly δ¹³C depleted values for cyclopoid copepods, and a high fractionation with respect to the pelagic baseline (F ~ 7.6), suggesting that cyclopoids may be relying upon deeper carbon sources than those exploited by Daphnia. In deep lakes, the carbon isotopic signature is strongly influenced by depth, with organisms living in deeper layers characterized by more negative values than those living closer to the surface and in the littoral zone [42,49,50]. Moreover, in oligotrophic clear lakes, especially in periods of water column stratification, it is common to have a gradient of δ¹³C POM (Particulate Organic Matter) with depth [51,52]. Matthew and Mazumder (2006) conducting a study in the oligotrophic Council Lake, recorded a δ¹³C of POM decrease with depth in the water column, and zooplankton taxa at deeper depths also had a lower δ¹³C. Our hypothesis is reinforced by the cluster analysis with δ¹³C values, which clearly show that cyclopoid copepods represented a distinct group, indicating a different utilization of carbon sources than other zooplankters, possibly found in deeper parts of the pelagic zone [49,53,54,55].

Because seasonal changes in the δ¹³C isotopic baseline were identified by a temporal shift towards less ¹³C depleted values in the littoral in the summer, this enabled the determination of the preferred foraging habitat of zooplankton taxa. The seasonal trend of cyclopoid copepods δ¹³C values could also be related to the observed presence of the plankivorous fish roach (Rutilus rutilus) in the lake during the summer, as a predation-avoidance strategy by cyclopoid copepods from the increase in predator pressure in the pelagic. The δ¹³C -enrichment observed in cyclopoid copepods in autumn may be explained by their return to the upper part of the pelagic zone when roach have migrated to the littoral zone to spawn [35].

The variation in δ¹⁵N of Bythotrephes could also be related to the presence or absence of plankivorous fish in the pelagic zone, as they exhibit δ¹⁵N depletion for the entirety of the period when plankivorous fish are feeding in the pelagic [35]. δ¹⁵N enrichment increases when the predation-pressure decreases. As the variation of δ¹⁵N between prey-predator decreases with niche occupation of increasing trophic levels, it is likely to be led by a top-down mechanism.

Cluster analysis with δ¹⁵N seasonal variation clearly split the zooplankton taxa into two functional groups, the primary consumers, Daphnia longispina galeata gr., Diaphanosoma brachyurum, and Eubosmina longispina, and the secondary consumers, Bythotrephes longimanus, Leptodora kindtii, and the copepods. The secondary consumers were more ¹⁵N -enriched than primary consumers. Cyclopoid copepods were feeding on the highest trophic level, which can be estimated to be 2–3 trophic levels higher than the primary producers, depending on whether a fractionation factor of 2.55 [29] or 3.4‰ [14] is used.

In this study, three or four trophic levels were identified in crustacean zooplankton of Lake Maggiore. In bio-magnification and energy and matter transfer in a pelagic food web, zooplankton is considered a crucial link between the primary producers and fish. In bio-energetic models, zooplankton are categorized as a “source” [56], or pooled into a singular grouping [57] of the trophic level of 2 [58,59], without taking into account species-specific traits or differences in life stage. Because the relative biomass of zooplankton taxa can significantly change during a calendar year [60], and because trophic positions of taxa and their relationships are dynamic during the year, we can conclude that freshwater zooplankton cannot be clustered together in the same ecological compartment. We propose the inclusion of seasonality and the dynamic of species trophic roles of zooplankton in the construction of models predicting bio-magnification capability.

In our study, we have demonstrated that the use of δ¹³C and δ¹⁵N stable isotope analysis represents an effective way to investigate the relationships present in the zooplankton community of a lake. In fact, the preferred habitat, food selectivity, and trophic position of each species could be defined. Through the use of δ¹³C and δ¹⁵N stable isotopes as a proxy of zooplankton functional traits, we can gain a better understanding of the ecological roles of the zooplankton species in the lake and thus define the functional diversity of the ecosystem. Moreover, combining the use of δ¹³C and δ¹⁵N stable isotopes with dietary analysis, we provided evidence that this could be an effective approach to infer functional groups, helping us understand the impact of functional differences in resource use. In particular, we have demonstrated that seasonal variations of δ¹³C and δ¹⁵N stable isotopes indicated a dynamic process of change in the relationships among zooplankton taxa, according to the different availability of food sources and of potential bottom-up and top-down (in particular fish predation) mechanisms present in the lake.