Introduction

Plankton communities are the fundamental basis of marine food webs and drive the dynamics of entire ecosystems1,2. They are a complex group of organisms represented by different taxonomic categories, from Bacteria to fish larvae, that respond rapidly to both external influences and internal dynamics3,4. Environmental factors such as temperature, salinity, and pH affect taxa composition and productivity, and thus have strong impacts on plankton biodiversity5. In general, changes in physicochemical factors can strongly influence planktonic organisms and their turnover rates6, making their communities useful for monitoring ecosystem health7,8,9. For example, resulting changes in biomass can lead to changes in the trophic structure of the plankton community10. Overall, the plankton community is able to rapidly cope with new conditions11,12, in part due to a variety of processes and functions that can be performed by the community at the right time13,14. Some basal processes, such as mixotrophy, heterotrophy and detritivory, are more or less pronounced and may be expressed in response to changing environmental conditions, for example, to maintain system resilience14 and ecologically meaningful processes13.

Especially in transitional ecosystems such as a coastal lagoon, plankton populations are subject to fairly frequent and significant habitat disturbances, such as those caused by freshwater or seawater inputs15,16. In addition, the interplay of various forces, such as wave energy, fishing activities, atmospheric disturbances, and climate change, are among the most important factors influencing and determining the exchange of matter and energy between the system components, consequently affecting the presence of resident species or the exclusion and arrival of other species17,18,19.

Observing plankton communities can therefore inform on how and to what extent the aquatic ecosystem is able to cope with sources of variability, and it is crucial to understand what the key processes involved are1,2,13. In this context, ecological network models are useful tools since they provide estimates of flows that are difficult to disentangle and measure14. The analysis of plankton networks allows for a holistic understanding of changes in the functioning of the marine system, as they represent a wide range of different taxa involved in the basic processes of microbial loop and therefore link to fundamental trophic and ecological processes13. Applying such ecological network studies to coastal areas has a number of advantages, including the increased availability of knowledge and data2 and known high variability in responses of organisms in these systems due to multiple external pressures20.

The present study compares the trophic status of the Venice Lagoon in July 2005 and 2007, by developing plankton trophic networks. Using a novel approach based on iterative random samplings of parameters within specific realistic ranges, we reconstructed the planktonic food webs in July for three main reasons: (I) the unimodal annual peak of phytoplankton biomass occurs in this month21; (II) it is the period with the best biological data coverage from microzooplankton (size between 20 and 200 µm22) to mesozooplankton23 (size between 200 and 20,000 µm22); and (III) July has been historically characterised by some economically important ecological processes, such as recruitment of small pelagic fishes24. Therefore, our study aims to understand the planktonic community structure of the Venice Lagoon, a well-studied coastal system25,26,27, and to determine the role of specific functional groups within these communities, with the resulting potential implications at the food web level as an emerging effect of factors affecting the system at the two different times.

Materials and methods

Study area

The Venice Lagoon is one of the largest lagoons in the Mediterranean and is part of the Italian Long-Term Ecological Research Network (http://www.lteritalia.it/)28. In its history, this important transitional ecosystem has undergone relevant changes over time and has been extensively monitored over decades from different aspects due to its socio-ecological importance29,30. Plankton communities, especially phytoplankton and zooplankton, have also been monitored for many years, using different methods and approaches depending on different objectives or needs related to impact assessments31,32,33.

The study focused on a limited shallow water area of the northern Venice Lagoon, i.e., the Palude della Rosa (Fig. 1), a typical lagoon area influenced by both freshwater inflow, from Canale Silone34, and saltwater input, due to incoming tides via Canale di Torcello15. Palude della Rosa covers an area of about 3.5 km2 with an average depth of about 0.5 m34 and is located in an intermediate position between the mainland coast and the sand barriers separating the Venice Lagoon from the Adriatic Sea. The plankton community at Palude della Rosa is therefore alternately influenced by river discharges and seawater intrusions34.

Figure 1
figure 1

Location of Palude della Rosa (black circle) in the Lagoon of Venice (black square), Italy. This image was created using the sf35, ggplot236, and ggsn37 packages (versions 1.0.9, 3.4.0, and 0.5.0 respectively) for R38.

Data

The available dataset covers July 2005 and 2007 and includes: taxonomic composition (where possible at species level or higher taxonomic levels), biomass (in mg C m−3) for organisms ranging in size from 1.150 to 28,000.000 µm. Plankton sampling was conducted during the neap tide to minimize variability associated with direct marine influence. Sampling and laboratory methods for biomasses of Bacteria, phytoplankton, mixoplankton, nanozooplankton, and microzooplankton are described in Pugnetti et al.39, while those of mesozooplankton, macrozooplankton, and non-living nodes in D’Alelio et al.10. Additional detailed description of methods used to collect plankton data, as well as ranges used for parameters and biomasses is reported in Supplementary Material.

Structure of the plankton networks

The plankton was classified at the lowest possible taxonomic level. For a few organisms, species-level classification was possible, while for others only larger taxonomic groups were available. Data at the lowest taxonomic level were grouped into species with similar ecological functions, known interactions, and similar biological rates to simplify the model40, resulting in a set of ecologically meaningful functional nodes (FNs) (Table 1). Each FN is characterised by its size and trophic role. Although ecological preferences of plankton in the Venice Lagoon vary widely41, two macro-preferences (pelagic and benthic) were considered for all FNs. These categories were then used for the analyses.

Table 1 List of parameter values for simulations.

Four metabolic parameters were assigned to each FN: production rate per biomass unit (μ, as d−1), consumption rate per biomass unit (α, as d−1), unassimilated fraction of biomass consumed (ε, dimensionless), and the phototrophy proportion in individual metabolism (ph, dimensionless). The latter has a value of 0 for heterotrophs, 1 for autotrophs, and a value between 0 and 1 for mixotrophs. The metabolic parameters μ, α and ε have a range with a maximum and a minimum value as extreme values for each FN, depending on the specific metabolism of each FN, which in turn is influenced by water temperature40. The proportions of flows to non-living nodes (γ) describing the fate of faeces, mortalities, and excreta also have a range.

The ordinal qualitative trophic links between FNs are ranked with four different values: 0, 1, 2, 3, representing no interaction, weak interaction, moderate interaction, and strong interaction, respectively (Table 2). These values were determined based on expert knowledge of plankton trophic ecology.

Table 2 Ordinal qualitative interactions between functional nodes (FNs).

Modelling approach

Plankton food webs were based on biomasses (B, as mg C m−3) of plankton FNs and flows between them as daily flows (mg C m−3 d−1).

Weighted plankton food webs have been developed that assume a balance between production, natural mortality, and consumption by predators for each living node k:

$${\upmu }_{\mathrm{k}}\cdot {\mathrm{B}}_{\mathrm{k}}-{\sum }_{\mathrm{j}=1}^{\mathrm{n}}\left({\mathrm{\alpha }}_{\mathrm{j}}\cdot {\mathrm{B}}_{\mathrm{j}}\cdot {\mathrm{DC}}_{\mathrm{k},\mathrm{j}}\right)-{\mathrm{m}}_{\mathrm{k}}=0$$
(1)

where μk is the production rate per biomass unit of FN k and Bk is its biomass. The first negative term is the sum of the consumptions of predator j as the product of the predator’s biomass Bj, its consumption rate per unit of biomass αj, and the proportion of living prey in the predator’s diet (DCk,j). The total number of FNs in the network is n and mk is the natural mortality of node k.

And for each non-living node d, a balance is established between flows to non-living nodes, consumption by detritivores, exports and imports:

$${\sum }_{\mathrm{i}=1}^{\mathrm{n}}\left[{\upgamma }_{\mathrm{i},\mathrm{d}}\cdot \left({\upvarepsilon }_{\mathrm{i}}\cdot {\mathrm{\alpha }}_{\mathrm{i}}\cdot {\mathrm{B}}_{\mathrm{i}}+{\mathrm{m}}_{\mathrm{i}}\right)\right]-{\sum }_{\mathrm{j}=1}^{\mathrm{n}}\left({\mathrm{\alpha }}_{\mathrm{j}}\cdot {\mathrm{B}}_{\mathrm{j}}\cdot {\mathrm{DC}}_{\mathrm{d},\mathrm{j}}\right)-{\mathrm{ex}}_{\mathrm{d}}+{\mathrm{im}}_{\mathrm{d}}=0$$
(2)

where γi,d is the proportion of flows from any node i to the non-living node d and εi is its unassimilated fraction of biomass consumed. The first negative term is the sum of the consumptions of detritivores j as the product of the predator’s biomass Bj, their consumption rate per biomass unit αj, and the proportion of non-living nodes in the diet of the predator (DCd,j). The amount of export and import of node d are exd and imd, respectively.

For each FN (i), production (P = μ ∙ B), consumption (Q = α ∙ B), and unassimilated (UN = ε ∙ α ∙ B) were related to estimate respiration (R):

$${\mathrm{R}}_{\mathrm{i}}={\mathrm{Q}}_{\mathrm{i}}-{\mathrm{P}}_{\mathrm{i}}\cdot (1-{\mathrm{ph}}_{\mathrm{i}})-{\mathrm{UN}}_{\mathrm{i}}$$
(3)

where phi is the phototrophy proportion of node i. Each time a trophic network met all conditions and constraints, it was accepted and the process began again until the ensemble of 1000 networks was reached. This procedure was applied to both 2005 and 2007.

Randomly generated networks with a-posteriori validity check

The system of equations was applied using an iterative approach in which μ, α, ε, and γ of each FN were randomly sampled from their range. The values obtained for the proportions of flows to non-living nodes were transformed so that the sum of proportions for each FN was equal to 1. Trophic links were transformed from ordinal qualitative values to quantitative values by randomly drawing two boundaries between 0 and 1 for each consumer. We constructed the matrix of the proportion of the diet (DCij, with flows from prey i to predator j) as follows: ordinal qualitative values equal to 1 were replaced by random values sampled from 0 and the first boundary; ordinal qualitative values equal to 2 were replaced by random values sampled from the first and second boundaries; and ordinal qualitative values equal to 3 were replaced by random values sampled from the second boundary and 1. Finally, the values obtained were transformed so that the sum of proportions of trophic links for each consumer was equal to 1. All samplings to determine the proportions of flows to non-living nodes, metabolic parameters, and links were performed using a uniform distribution.

Equations (1) and (2) were used for all FNs of the food web to establish a system of algebraic linear equations in which several parameters had a range (μ, α, ε, and γ) or were defined only in ordinal qualitative terms (DCij). Other parameters such as natural mortality, exports, imports, respirations, and gross food conversion efficiency were estimated using the previous parameters. Considering the range of parameters, the system of equations does not have a unique solution. To examine all potential parameter combinations when no relevant information about parameter distributions is available, we randomly sampled them from a uniform distribution over the specified ranges and analysed a posteriori the distribution of parameters for valid networks. The resulting networks were tested for ecological and thermodynamic realism based on a set of simple constraints to limit respiration and gross food conversion efficiency.

Thus for each FN (i):

$${\mathrm{R}}_{\mathrm{i}}\ge 0$$
(4)

and for each consumer (j):

$$0.15\le {\mathrm{GE}}_{\mathrm{j}}\le 0.5$$
(5)

where Ri and GEj are respectively the respiration flow and gross food conversion efficiency42,43,44, the ratio of heterotrophic production to consumption \({\mathrm{P}}_{\mathrm{j}}\cdot (1-{\mathrm{ph}}_{\mathrm{j}})/{\mathrm{Q}}_{\mathrm{j}}\), of nodes i and j, respectively.

Indicators

A set of whole system indicators and the omnivory index, which provide information on the ecological characteristics of food webs, were calculated (Table 3). Each indicator was calculated for the two studied plankton networks to obtain a distribution of values with which statistical tests were performed for comparison.

Table 3 Indicators for the study of the networks of the Venice Lagoon.

For each of these indicators, the abbreviation, the formula and a brief description of its meaning are given.

Statistical analysis

Principal component analysis (PCA) was carried out on the consumption flows. Only the non-constant flows shared by the food webs of the July of two years were selected, since they had a different number of nodes. Because the value ranges of the variables were very different, these were shifted to be zero-centered and scaled to have unit variance before analysis to avoid larger value variables dominating the PCA results. In analysing the PCA results, the first 15 loadings, ranked by their relative importance to the first two principal components, were considered to determine the FNs with the greatest contribution to the split of the simulations of the July of the two years tested in terms of consumption flows.

Comparisons between 2005 and 2007 food webs and between pelagic and benthic primary production within July of each year were made using Mood's median test. The one-sample sign test was used to verify if the values of relative ascendency were statistically different from a reference value of 0.459651. Non-parametric tests were used because the assumptions for parametric ones were not met.

The entire modelling approach, including the calculation of indicators and statistical analyses, was developed in R version 4.2.238 with RStudio version 2022.07.2 + 57652, using EnvStats53, DescTools54, and NetIndices55 packages versions 2.7.0, 0.99.47, and 1.4.4.1, respectively.

Results

Principal component analysis

PCA on the consumption flows shows that the first two principal components together account for 34.251% of the total variance (Fig. 2a). The first 15 loadings, ranked by their relative importance for the first two principal components, are shown in Fig. 2b. The consumption flows involving FN 28 (Strombididae) contribute most to the first principal component, while the consumption flows involving FN 39 (Evadne spp. and Podon spp.) contribute most to the second principal component.

Figure 2
figure 2

Principal component analysis. The plot was created using the ggplot236 and ggpubr56 packages, versions 3.4.0 and 0.6.0, respectively, for R38. (a) The first two principal components are given and together account for 34.251% of the total variance. (b) The first 15 loadings ranked by their relative importance to the first two principal components, are given. The consumption flows involving FN 28 (Strombididae) contribute most to the first principal component, while those involving FN 39 (Evadne spp. and Podon spp.) contribute most to the second principal component. "Q" indicates the consumption flow moving from the prey (the first number in brackets) to the predator (the second number).

Shifting topological roles in key functional nodes

The results show that consumers with the phototrophy proportion greater than or equal to 0.5 have higher OI values, namely mixo-dinoflagellates and Mesodinium cf. rubrum (FNs 18, 22, and 23). On the other hand, exclusive detritivores, namely Harpacticoida and Bacteria (FNs 48 and 49, respectively), have the lowest value, i.e., zero (Fig. 3).

Figure 3
figure 3

Violin plot of the omnivory index (OI) for each consumer in July in each year (2005 and 2007). Each violin was plotted to have the same maximum width, but if there is only one violin in a year, it has twice the maximum width. The plot was created using the ggplot2 package36 version 3.4.0 for R38.

Output structure of the models

The iterative process to develop ecological networks using the input data sets (and ranges) on plankton communities in the Venice Lagoon for the July of 2005 and 2007 resulted in 1000 meaningful networks for each year. These valid networks result from a random selection of parameters within the ranges and do not lead to unrealistic ecological processes (such as unrealistic respiration, mortality, etc.). The first, second, and third quartiles, listed in Table 1, are used to describe the results of the parameters since they do not have a normal distribution.

The median primary production of pelagic and benthic nodes was respectively 20.890 and 27.144 mg C m−3 d−1 for 2005, and 6.373 and 28.347 mg C m−3 d−1 for 2007. Mood's median test between pelagic and benthic primary production was highly significant in July of each year (the p-value was 0 in July of both years and the value of the z-statistic was − 31.029 for 2005 and − 44.710 for 2007). The results of Mood's median test for total imports to undissolved detritus and for parameters related to FNs 28 (Strombididae) and 39 (Evadne spp. and Podon spp.) between 2005 and 2007 are shown in Table 4. Compared to 2005, total imports to undissolved detritus (FNs 50, 51, 52) increased in 2007, OI and percent consumption of FN 28 prey that can provide "domesticable" plastids, i.e., phyto-nanoflagellates, mixo-dinoflagellates, and Mesodinium cf. rubrum smaller than 20 µm57, which are FN 16, 17, 18, and 22, decreased in the diet of FN 28, while the production rate per biomass unit and OI of FN 39 remained the same.

Table 4 The results of all comparisons between 2005 and 2007 are given.

Whole system indicators

The results of the comparison of each whole system indicator between the July of 2005 and 2007, calculated using Mood's median test, are shown in Table 4. For each whole system indicator, the test revealed significant statistical differences between July of each year. Median values of Finn's cycling index and relative ascendency increased over time, while median values of relative internal ascendency, detritivory to herbivory ratio, and primary production to community respiration ratio decreased. The distributions of some whole system indicators are shown in Fig. 4.

Figure 4
figure 4

Histograms of some whole system indicators for July of the two years (2005 and 2007): relative internal ascendency (Ai/Ci), ratio of detritivory to herbivory (D/H), Finn's cycling index (FCI), and ratio of primary production to community respiration (PP/R). Graphs were generated using the ggplot2 package36 version 3.4.0 for R38.

Discussion

This study focused on modelling summer plankton food webs in one of the most important transitional ecosystems of the Mediterranean, the Venice Lagoon. The structure of the plankton food web was numerically derived in July 2005 and 2007 based on experimental data. In this ecosystem, the plankton community is highly influenced by a mixed control mechanism that depends mainly on tidal conditions affecting salinity and nutrient gradients, but also on anthropogenic influences29,30, that may limit the approach used. Nonetheless, the choice of the same month and the sampling carried out during the neap tide seem to have minimized the environmental variability, which is supported by the fact that the whole system indicators associated with the networks of July of the two years largely converge despite the differences in the composition and abundance of the FNs.

The results show that the consumption flows involving FNs 28 and 39 were the most important in differentiating the food webs of the two years. These FNs are mixotrophic ciliates (Strombididae, unicellular facultative mixotrophs) and Cladocera (i.e., Evadne spp. and Podon spp., fast-growing metazoans), respectively. Strombididae feed on small unicellular organisms belonging to the pico- and nanoplankton, of which some photosynthetic ones are retained as plastids for photosynthesis, while the heterotrophic ones are digested for energy production57. Hence, Strombididae are called generalist non-constitutive mixotrophs58. In 2007, there was a greater diversity of autotrophic prey of FN 28 (Strombididae), resulting in a decrease in OI compared to 2005. However, the percent consumption of its prey that provide "domesticable" plastids57, i.e., phyto-nanoflagellates (FNs 16 and 17), mixo-dinoflagellates (FN 18), and Mesodinium cf. rubrum smaller than 20 µm (FN 22), decreased and its production rate per biomass unit did not change in July of the two years (Table 4). The combination of these factors may indicate a more heterotrophic behaviour of Strombididae in 2007 by increasing predation on other, more abundant, unicellular prey such as Bacteria (FN 49) or hetero-nanoflagellates (FN 34) (Table 1).

In 2007, even hetero-nanoflagellates were an important food for some metazoans (Table 2), such as Evadne spp. and Podon spp. (FN 39), whose biomass increased by an order of magnitude in 2007 (Table 1). Cladocerans are organisms with an affinity for marine and coastal waters and hardly reside in the interior of the Venice Lagoon23,59, but under certain tidal conditions, incoming marine waters cause them to extend into more interior areas of the lagoon59, such as Palude della Rosa. Cladocerans are strongly influenced by the seasonality and spatial variability of environmental conditions60, thus parthenogenesis allows them to respond very quickly to environmental changes61. Their higher growth rates, compared to other planktonic crustaceans, must be sustained by higher consumption rates of prey whose abundance is more stable over time61, such as picoplankton61,62 (size between 0.2 and 2 µm22). This ability to respond rapidly to environmental variability was also evident in our study. In fact, a greater number of FNs, with a lower abundance (Table 1), fulfilled the function of primary producers in 2007. Thus Evadne spp. and Podon spp. (FN 39) doubled the number of prey, from 16 in 2005 to 32 in 2007, to maintain their high growth rate, while their OI remained constant (Table 4). Although the abundance of these crustaceans was very low (Table 1), their presence under certain conditions may reveal structural and functional changes already observed in other coastal planktonic food webs10.

Despite differences in community composition, due to interannual population fluctuations typical of transitional environments such as coastal lagoons63,64, whole system indicators show very similar results for the planktonic food web for July of both years. In our models, an inverse trend is observed between FCI and Ai/Ci, with median FCI higher and median Ai/Ci lower in 2007 than in 2005 (Table 4). Although FCI is considered an indicator of system maturity49, it can be used in conjunction with other whole system indicators, such as Ai/Ci, to assess whether or not the system is stressed46. Although FCI is a partial measure of recycling in networks65, it is used here to compare results with other studies. However, the use of a whole system indicator that provides more accurate information about recycling within the system, such as the comprehensive cycling index65, is strongly recommended for future studies.

Our results indicate a significant increase in total import to the undissolved detritus (FNs 50, 51, and 52) from 2005 to 2007 (Table 4), which is confirmed by sampling of their total biomass, which increased from 318.362 in 2005 to 334.327 mg C m−3 in 2007 (Table 1). This led to a progressive increase of stress in the system, which increased matter/energy recycling and, consequently, FCI66. At the same time, Ai/Ci decreased so that resilience increased and the system could cope with the perturbation66. The efficiency associated with internal flows decreased so that the system became more dependent on external flows46. Systems with lower Ai/Ci ratios are more resilient because they have higher redundancy in trophic pathways, which allows them to recover disrupted ones46. However, systems with lower Ai/Ci ratios are not resistant because they have low internal stability, which makes them more susceptible to external influences that can alter their configuration46.

In our models, A/C increased even if Ai/Ci decreased over time, so overall efficiency increased even if internal efficiency decreased. The A/C median values below 0.5 in July of both years indicate that the overall resilience of the system under study was greater than its overall efficiency45. In particular, for ecological networks, an A/C value of 0.4596 has been suggested as optimal51 to represent two opposing trends in a system development, efficiency and resilience45. The A/C medians are statistically lower than the optimum, suggesting that the system likely increased its resilience at the expense of its efficiency to cope with the source of stress. This tendency of A/C is similar to eutrophication67 but the enrichment comes from organic matter rather than nutrients. The system should have increased its efficiency to improve its sustainability (robustness) in the long term51.

The Venice Lagoon is an ecosystem exposed to various natural and anthropogenic influences68. In particular, during the first decade of the 2000s, several factors led to the resuspension of sediments: dredging of new large channels, increasing number and speed of boats69, lower stabilisation of sediments due to the decline of seagrass70,71, and widespread use of mechanical fishing gears for harvesting Manila clams (Ruditapes philippinarum)71,72.

Sediment resuspension led to an increase in water turbidity and consequently to a decrease in phytoplankton31. In this regard, chlorophyll a, which is a proxy for phytoplankton biomass, showed a decreasing trend from 2001 to 200721. Our work also showed that the PP/R values were well below 1 in July of both years, suggesting that sediment resuspension simultaneously inhibited primary producers due to turbidity and favoured detritivores due to the increase in organic carbon content in the water column72. As a result, community respiration increased (respiration of Bacteria, FN 49, has medians of 84.187 and 90.905% of community respiration in 2005 and 2007, respectively). Similarly, hetero-nanoflagellates (FN 34), which fed primarily on Bacteria (Table 2), seemed to benefit from this condition, contributing to a significant increase in community respiration (median values of 13.367 and 4.537% of community respiration in 2005 and 2007, respectively). Among the mesozooplanktonic FNs, which are between 200 and 20,000 µm in size22, a similar reasoning could be applied to Acartia spp. (FN 41), as in both years their biomass alone was greater than the median of the sum of the other mesozooplanktonic FNs (Table 1).

In the Venice Lagoon, intense sediment resuspension from Manila clam harvesting was a source of stress73 that acted like an external press perturbation74 on the system. In particular, our results confirm that detrital resuspension was a necessary component for the Venice Lagoon ecosystem75,76. In fact, our models estimate a mandatory import to non-living nodes to maintain detritivore consumption. The very high values of D/H make it clear that non-living nodes were the most utilised resource at the lowest trophic level. Nevertheless, only 4 FNs performed the detritivore function (Table 2): Acartia spp. (FN 41), Oithona spp. (FN 43), Harpacticoida (FN 48), and Bacteria (FN 49). Of these, Bacteria benefited the most from this condition17, because in both years their biomass alone was higher than the sum of the other FNs, making the decomposition processes significant29. Thanks to sediment resuspension, heterotrophic Bacteria were able to maintain high densities even when the carbon source of dead phytoplankton was insufficient to sustain them17. In addition, Bacteria that thrived on detritus could strongly influence the food web, as Bacteria are intensively consumed by Protozoa, which in turn were eaten by higher trophic levels13. Since the average depth of Palude della Rosa is about 0.5 m34, there is a close coupling between the benthic and pelagic environments29, so sediment resuspension has profoundly altered the ecosystem, not only because of greater resources for detritivores75, but also because of the resuspension of benthic FNs11,31. These accounted for most of the total primary production, especially in 2007 when the primary production of benthic FNs was more than three times that of pelagic FNs (Table 1).

A hypothetical persistence of conditions, such as those highlighted in our work, would inevitably have implications for the population dynamics of consumers at higher trophic levels, and thus for the structure and functioning of entire food webs77,78,79. In fact, a sharp decline in the abundance of fish feeding on plankton was observed in the Venice Lagoon landings from 1995 to 2001, with a decrease in the ratio between pelagic and demersal fishes that showed a value of less than 1, as well as an effect of the rapid decline of higher level consumers80. This phenomenon might appear to be overfishing, but instead it was due to the direct and indirect impact of Manila clam harvesting on the entire ecosystem80.

The approach used in this study proved to be a valid tool for capturing and interpreting the major forces affecting aquatic food web dynamics at two different time points even in highly dynamic environments. However, the application of this network approach to consecutive plankton samples in time and space is necessary in the future to link them to possible interpretations or predictions of future scenarios. With the latter in mind, and measures to identify tools to mitigate and possibly prevent the various sources of impact, the idea was to capture the driving forces of change in lagoon planktonic communities by placing them in the context of the various pressures to which they may be exposed. The results confirm that the plankton community can serve to assess the health of the whole ecosystem7,8,9, as it provided results comparable to those of other studies on high trophic level networks. Indeed, we demonstrated that sediment resuspension was a source of stress73 on which the system was highly dependent75,76, and that the system increased its resilience at the expense of its efficiency in coping with the perturbation66 moving away from its optimum of robustness51.