Long-Term Patterns in the Population Dynamics of Daphnia longispina, Leptodora kindtii and Cyanobacteria in a Shallow Reservoir: A Self-Organising Map (SOM) Approach

The recognition of long-term patterns in the seasonal dynamics of Daphnia longispina, Leptodora kindtii and cyanobacteria is dependent upon their interactions, the water temperature and the hydrological conditions, which were all investigated between 1999 and 2008 in the lowland Sulejow Reservoir. The biomass of cyanobacteria, densities of D. longispina and L. kindtii, concentration of chlorophyll a and water temperature were assessed weekly from April to October at three sampling stations along the longitudinal reservoir axis. The retention time was calculated using data on the actual water inflow and reservoir volume. A self-organising map (SOM) was used due to high interannual variability in the studied parameters and their often non-linear relationships. Classification of the SOM output neurons into three clusters that grouped the sampling terms with similar biotic states allowed identification of the crucial abiotic factors responsible for the seasonal sequence of events: cluster CL-ExSp (extreme/spring) corresponded to hydrologically unstable cold periods (mostly spring) with extreme values and highly variable abiotic factors, which made abiotic control of the biota dominant; cluster CL-StSm (stable/summer) was associated with ordinary late spring and summer and was characterised by stable non-extreme abiotic conditions, which made biotic interactions more important; and the cluster CL-ExSm (extreme/summer), was associated with late spring/summer and characterised by thermal or hydrological extremes, which weakened the role of biotic factors. The significance of the differences between the SOM sub-clusters was verified by Kruskal-Wallis and post-hoc Dunn tests. The importance of the temperature and hydrological regimes as the key plankton-regulating factors in the dam reservoir, as shown by the SOM, was confirmed by the results of canonical correlation analyses (CCA) of each cluster. The demonstrated significance of hydrology in seasonal plankton dynamics complements the widely accepted pattern proposed by the plankton succession model for lakes, the PEG (Plankton Ecology Group), and may be useful for the formulation of management decisions in dam reservoirs.


Introduction
Dam reservoirs are formed and modified by human activity for specific purposes, and they therefore have different hydrological characteristics and thermal and chemical regimes compared to natural lakes. The specific features of reservoirs, especially unstable hydrological conditions, determine the species composition and dynamics of the plankton communities inhabiting them [1]. In this study, the long-term patterns in the seasonal dynamics of three types of planktonic organisms were investigated, dependent upon their interactions, the water temperature and the hydrological alterations. The study was conducted in the lowland dam of the Sulejow Reservoir, which is regarded as a model ecosystem for ecological research due to its location and characteristics [2]. We examined 1) the dominant species of daphniid in this ecosystem, Daphnia longispina (O.F. Müller); 2) the predatory cladoceran Leptodora kindtii (Focke); and 3) cyanobacteria, which are a serious problem in the reservoir due to their creation of toxic blooms.
The importance of daphniids in ecosystem functioning is related not only to their effects on the abundance, species composition and size structure of phytoplankton [3] but also to bacterial and protistan dynamics [4,5] and thus to the epilimnetic cycle of organic matter and nutrient balance [6]. The relatively well-known biology and ecology of daphniids (e.g., feeding behaviour, reproduction and both predator-prey and host-parasite interactions) underlie their particular suitability for studies on interactions between trophic levels and on the propagation of trophic effects in aquatic ecosystems [7,8]. In this context, Daphnia spp. are proper model organisms for ecological and evolutionary studies in aquatic ecosystems as well as in laboratory conditions [8]. The density of Daphnia populations varies strongly within the growing season and depends on both environmental factors and biotic interactions [9]. Because density is the essential parameter determining the functional importance of Daphnia in an ecosystem, we focused our study on this parameter and on selected factors that could potentially significantly affect density. Daphnia longispina attains a relatively large body size (0.8 to 2.1 mm) in the studied ecosystem, which indicates the potential strength of the top-down effect through their grazing. Unfortunately, during summer months, this effect is unsettled due to annual toxic cyanobacterial blooms in the Sulejow Reservoir. Microcystis aeruginosa (Kutzing), which produces microcystin-LR, microcystin-YR and microcystin-RR, is the dominant bloom-forming cyanobacteria in this reservoir [10]. Cyanobacterial filaments and colonies are inadequate and/or inedible food for grazers due to their large size, poor nutritional value and toxicity [11]; thus, Daphnia spp. usually exhibit decreased survival and reproduction potential in the presence of cyanobacteria [12]. However, our previous studies have demonstrated that D. longispina, which have coexisted with a high biomass of toxic cyanobacteria in the Sulejow Reservoir, have effective protection mechanisms against the toxic effects of microcystins [13,14]. These new facets may indicate that the Daphnia-cyanobacteria relationship is more complex than previously thought and that the presence of toxic cyanobacterial blooms is not always related to the collapse of a Daphnia population. To verify their impact on D. longispina population dynamics, we included cyanobacteria as a biotic parameter in this long-term study. We did not evaluate the grazing pressure of D. longispina on phytoplankton in our study; however, chlorophyll a was used as the indicator of the food base availability for daphniids because it provides a reasonable estimate of algal biomass [15,16].
As the next biotic factor that could potentially limit the density of zooplankton species, including D. longispina, the predatory Leptodora kindtii was taken into account. This decision resulted from our earlier studies indicating that the role of fish in the reduction of daphniid density in the pelagic zone of the Sulejow Reservoir was much smaller than that of L. kindtii. We found that the short but intense predation of L. kindtii during the summer significantly affects the dynamics of Daphnia populations by causing their seasonal decline [17]. Similar predatory pressure of L. kindtii on Diaphanosoma brachyurum (Lievin) was observed in Lake Neusiedl [18].
In this paper, the effects of selected abiotic factors on D. longispina, L. kindtii and cyanobacteria were studied. Freshwater ecosystems undergo a dual-regulation process in which abiotic conditions (e.g., hydrology and temperature) strongly influence biota and vice versa [19]. As far as the hierarchy of these factors is concerned, the abiotic conditions play a dominant role [20], mainly by limiting the growth season of organisms (primarily due to unfavourable temperature and light conditions). The dynamics of plankton populations are closely related to the thermal gradient [21]; thus, we considered temperature as an important environmental factor that generally influences the physiology, life history and demography of planktonic organisms. The seasonal response of daphniids to temperature consists of modifications of species-specific life-history traits (such as emergence from resting stages, spawning or generation times), depending on the timing of favourable environmental conditions [22]. In this regard, one of the most critical periods for cladocerans is spring, when the temperature determines the biomass of both Daphnia and L. kindtii, e.g., resting eggs of L. kindtii begin to hatch to the metanauplius stage at temperatures of 8-10°C [21]. Spring temperatures also determine the timing and the strength of predators' pressure. Wagner and Benndorf [21] found that the abundance of L. kindtii above which its predatory impact on daphniids became significant occurred at water temperature of 16°C and increased rapidly at water temperatures between 19 and 24°C. Increased temperatures in autumn as well as during mild winters may have important consequences for the overwintering strategies of daphniids [23]. The mean temperature in January-March was found to determine the future seasonal dynamics of the Daphnia population, which, consequently, determines the springtime occurrence of the clear-water phase and the duration of the daphniid midsummer decline [24]. Cyanobacteria have a relatively slow growth rate at moderate temperatures, which considerably increases in warmer conditions, often in excess of 25°C. Consequently, in temperate ecosystems, cyanobacterial dominance and bloom formation occur during hot summers [25]. These observations demonstrate that temperature is one of the major abiotic factors shaping the patterns of plankton population dynamics in water bodies of the temperate zone.
However, the seasonal events in cladoceran and cyanobacterial density dynamics may appear mostly in ecosystems with relatively stable hydrological conditions (e.g., in lakes with low inflow rates). In dam reservoirs, the abiotic environment is not as predictable as in lakes, and the dynamics of plankton populations appear to largely depend on the hydrological regime. Hydrological instability may decrease the importance of biotic interactions, including the predator-prey relationship between Daphnia and Leptodora, through the effects of flushing on zooplankton species during flash floods and through the reduction in zooplankton biomass peaks and/or shifts in the time of their appearance caused by a short retention time [26]. Moreover, during periods of frequent flushing, the plankton community is dominated by zooplankton organisms that reproduce fast enough to offset their removal by flushing, such as rotifers and small cladocerans, and replace Daphnia populations [27]. During stable hydrological conditions, biotic regulation may become dominant [20]. A long residence time and enhanced water stratification usually promote cyanobacterial growth and bloom formation [28]. Therefore, both inflow and retention time were included in our study as crucial factors that influence the plankton population dynamics in a dam reservoir.
The aforementioned abiotic and biotic factors interact with each other in a complex manner, underlying the reason behind the typical strongly non-linear relations between these factors and the population dynamics of plankton [29]. Thus, the seasonal dynamics of D. longispina, L. kindtii and cyanobacteria should be analysed with appropriate mathematical techniques. Hence, we decided to apply a Kohonen self-organising map (SOM), which is an unsupervised artificial neural network (ANN). SOMs effectively recognise patterns in complex data sets, even if 1) the relationships between variables are non-linear; the SOMs, in contrast to classical statistical methods, are applicable to data expressed at different measurement scales [30], and/or 2) the distribution of the data is non-normal. The latter problem is encountered particularly often in environmental studies because of the presence of many zeroes in data sets with an abundance of species (e.g., organism counts). Such variables, even when transformed, seldom satisfy the assumptions for normality of distribution [31]. In contrast to the classical statistical approach, analysis with SOMs requires no a priori specification about the underlying model. Moreover, SOMs reduce n-dimensional data to a two-dimensional map, which is interpreted quite easily, i.e., similar to the results of principal component analysis (PCA; similar objects are neighbouring, whereas different objects are distant) [32,33].
The aims of this study, accomplished through use of an SOM, were as follows: i) to recognise long-term patterns (1999)(2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008) in the seasonal dynamics of biotic elements (cyanobacteria, Daphnia longispina and Leptodora kindtii) despite the high variation of abiotic and biotic parameters between years and ii) to determine the influence of seasonal changes in abiotic (temperature and hydrological conditions) and biotic parameters on the population dynamics of D. longispina and L. kindtii and on the biomass dynamics of cyanobacteria in the Sulejow Reservoir.

Ethics statement
No specific permits were required for the field studies described herein. There was no activity involving endangered or protected species in this study.

Study area
The Sulejow Reservoir is a shallow reservoir that comprises 138.9 km of the Pilica River (the Vistula River catchment) in central Poland (Fig 1). The reservoir was built in 1973. Its maximum length is 15.5 km, and its maximum width is 2.1 km. At maximum capacity (75 x 10 6 m 3 ), the reservoir covers 1980 ha, with a mean depth of 3.3 m and a maximum depth of 11 m. The reservoir sub-catchment area is covered mainly with agricultural land (50% arable land, 13% meadows and pastures, 1% orchards) and forests (31%). The shoreline length is approximately 54 km. The Sulejow Reservoir is a non-stratified, eutrophic ecosystem [2]. In the initial stage of seasonal succession, the phytoplankton community consists mostly of diatoms, cryptophytes and green algae. The phytoplankton community in the reservoir has a typical composition for Central Europe, characterised by the dominance of the microcystin-producing genus of Microcystis [34]. In the summer months, blooms of mainly Microcystis aeruginosa Kütz. and Aphanizomenon flos-aquae (L.) Ralfs have been frequently observed [35].
Gillnet catches have been regularly conducted since 1993 and revealed a very stable composition of the fish stock during the 1990s, with dominance (in terms of biomass) of common bream, Abramis brama L.; pikeperch, Stizostedion lucioperca L.; white bream, Abramis bjöerkna L.; roach, Rutilus rutilus L. and perch, Perca fluviatilis L. [37].
Studies were conducted at each of the three sampling stations located along the longitudinal reservoir axis. The sampling stations (Fig 1) included the following: an upper (Tresta = TR), a middle (Bronisławów = BR) and a lower (Zarzęcin = ZA) station. The studies were conducted every week between April and October during the years 1999-2008 for a total of 412 sampling terms. The database is presented in S1 Table. Plankton sampling and analyses A total of 1236 samples (3 sampling stations × 412 sampling terms) of zooplankton and phytoplankton were collected. The samples were collected weekly in the pelagic zone of the reservoir at the three sampling stations. Water samples were taken using a 5-L sampler from each metre of the entire water column and integrated. Then, a 1-litre subsample of the integrated water was collected for the determination of cyanobacterial presence and biomass. These water samples were preserved in Lugol's solution and sedimented in the laboratory. Algae were analysed using a Fusch-Rosenthal counting cell (procedure according [38]).
Daphnia longispina and Leptodora kindtii were sampled using a 5-L Bernatowicz sampler at four depths: one, two, three and four metres (maximal depth at the sampling stations did not exceed four metres). Samples were taken with four or five replicates for each depth to properly evaluate the density of L. kindtii [39]. Because diurnal vertical migration of zooplankton has not been observed in this shallow reservoir and zooplankton were almost evenly distributed throughout the entire water column ( [36] and unpublished data based on the annual monitoring of the reservoir), samples from all four depths were integrated. Samples were then filtered with a 50 μm mesh size plankton net. Collected zooplankton were preserved in 4% Lugol's solution. The samples were identified and counted under a Nikon 115 microscope using a Sedgewick Rafter counting chamber. Morphological analyses of collected Daphnia species were performed according to Benzie [40].

Chlorophyll a analysis
Water samples were taken weekly using a 5-L sampler from each metre of the entire water column at the three sampling stations in the pelagic zone of the reservoir. Samples from different depths were integrated, and a 1-litre subsample was then obtained for the analysis of chlorophyll a. Chlorophyll a concentrations were analysed using a spectrophotometric method [41].

Hydrological parameters
The water retention time of the reservoir was calculated as the ratio of the reservoir volume to the volume of water inflow to the reservoir. These data were available from the Smardzewice branch of the Regional Water Management in Warsaw (http://warszawa.rzgw.gov.pl/theregional-water-management-authority-in-warsaw), which administers the Sulejow Reservoir and conducts permanent statutory monitoring of the hydrological parameters according to Polish and EU standards.
The hydrological data relate to the same dates as the collection of other parameters. The same values of retention time and water inflow to the reservoir were used for hydrological characteristics of the sampling station. We did not consider the pattern of water inflow propagation along the reservoir and assumed that, due to the specific reservoir shape, the biota were impacted by inflowing water in a similar way at all three sampling stations, especially during flood events. The coefficient of water inflow variation (CV) was used to classify years as hydrologically stable or instable and was calculated as the ratio of the standard deviation and the mean of the water inflows to the reservoir in June-September, expressed as a percentage.

Temperature
Water temperature (°C) was measured between 10:00 and 11:00 a.m. at each sampling term in the offshore area of the Sulejow Reservoir at 1.5 m depth using WTW 340i/SET multisensors.

SOM analysis
The patterns in the biotic components (nine variables: the biomass of cyanobacteria and the densities of D. longispina and L. kindtii at the three stations) were recognised using Kohonen's artificial neural network, i.e., a self-organizing map (SOM) [32]. An SOM is constructed from two (input and output) layers of data processing units (neurons). The SOM training was performed using the SOM Toolbox available at http://www.cis.hut.fi/projects/somtoolbox/ [42], developed in the Laboratory of Information and Computer Science at the Helsinki University of Technology. During the network training (using the sequential algorithm), the logtransformed data (nine variables × 412 sampling terms; ST) were repeatedly displayed on nine input neurons (their number = the number of variables). The input neurons, each connected to all of the output neurons, transmitted information to the output layer of the neurons. On this basis, a virtual ST (understood as a set of values of the nine analysed variables) was created in each output neuron. The virtual STs (and thus the respective output neurons) were clustered with hierarchical cluster analysis (Ward linkage method, Euclidean distance), which is a standard procedure in the SOM Toolbox [43]. The relative position of any two virtual STs on the SOM reflects their similarity; virtual STs in neighbouring neurons are similar, whereas those in distant neurons are different. Finally, each real ST was assigned to the best matching virtual ST (and the respective output neuron and cluster of output neurons); similar real STs were located in the same neuron or in adjoining neurons and were considerably different from the real STs located in distant neurons.
We decided to input STs (with values for the three stations treated as separate variables) to the SOM instead of individual samples. The latter solution was more intuitive but would cause problems in further analyses mostly because of the unbalanced statistical influence of a given ST on the values calculated for output neurons containing different numbers of dependent samples from that ST. The adopted solution allowed us to assign STs (all of equal merit) univocally to output neurons (and their (sub-)cluster) and to analyse different effects, including the differences between sampling stations, with more conventional analyses.
We finally adopted the classification of the SOM output neurons into three clusters and six sub-clusters, resulting in a fairly symmetrical division of the SOM, i.e., the particular clusters or sub-clusters of virtual units corresponded to similar parts of the total mapped variability.
The SOM Toolbox allowed visualisation of the intensity of particular variables in output neurons in the form of a greyness gradient over the SOM [33]. It was possible for both density and biomass (biomass of cyanobacteria, densities of D. longispina and L. kindtii, all input to the SOM) to explain the remaining variables (chlorophyll a, water temperature, water inflow to the Sulejow Reservoir and the retention time of water in the reservoir, not presented directly to the SOM). The greyness intensity reflected 1) the values of density/biomass variables in virtual STs (i.e., in respective output neurons) and 2) the mean values of the remaining variables calculated for real STs assigned to a given output neuron.
The number of output neurons (24 on a two-dimensional, 6 × 4 grid) was arbitrarily selected from the other tested options (from 4 × 4 to 12 × 9). Although there was no universal principle to determine the optimum map size, we considered quality measures for SOMs such as the quantisation (QE) and topographic (TE) errors [43]. The former is the average distance between each real ST and its best matching virtual ST, whereas the latter is the proportion of all real STs for which the first and second best-matching virtual STs were not in adjacent neurons [44]. In our study, similar QE and TE were observed for larger maps compared to the 6 × 4 map that was ultimately adopted (S2 Table). Some of the larger maps corresponded to the often-applied rule stating that the number of output neurons should be close to 5 p n, where n is the number of real STs [45]. Our decision on the actual number of output neurons in the map was influenced by 1) obtaining clear patterns and 2) avoiding output neurons that were either empty (i.e., without any real STs) [43] or contained 1-2 real STs assigned, such that the aforementioned greyness used to visualize the values of variables in such neurons did not assume a random intensity, which would disturb the general greyness gradients.

Canonical correlation analyses
We used canonical correlation analyses (CCA) to explore the relationship between abiotic and biotic sets of variables that were measured in the reservoir during the entire investigation period. Inflow, retention time, temperature and chlorophyll a were chosen as predictor variables, and cyanobacterial biomass, D. longispina density, and L. kindtii density were used as response variables. Data from each sampling term represented the average values of the three sampling stations. Prior to the analyses, all data were log (x+1) transformed for variance stabilisation. CCA were performed separately for the three main clusters of SOM output neurons. By analysing data subsets characterised by relatively reduced variability, we expected to reveal more subtle effects that were presumably different for particular SOM clusters.

Statistical tests
We decided to explore the data at a lower typology level by comparing the biotic and abiotic variables between 1) the six sub-clusters of SOM output neurons and 2) the three sampling stations within each SOM sub-cluster. Because the SOM does not provide a statistical verification of the observed differences, the significance of differences between the SOM sub-clusters was assessed by the Kruskal-Wallis and post-hoc Dunn tests [46] in the Statistica package (Statistica 12, StatSoft, Inc.). The significance of differences between the three sampling stations was assessed with the Friedman test for multiple groups of dependent samples [46] in the XLSTAT package. The Kruskal-Wallis and Friedman tests are non-parametric (based on ranks) equivalents of the analysis of variance that allow for testing of whether more than two groups of samples differ. The groups compared with the Kruskal-Wallis test may contain different numbers of samples, whereas groups compared with the Friedman test must be of the same size. Similar to the SOM, both tests can be applied to data of a normal or skewed distribution. Because these tests indicate the existence of a significant difference but do not indicate how many groups of samples and which of them differ, a post-hoc test should be subsequently applied for multiple comparisons between the groups [46].

of study
In all the studied years, a high water stage in the spring months (April-May) was observed. However, later on, the hydrological regime varied within the ten-year period (Fig 2A and 2B). kindtii among the three stations in all studied years (Fig 3B). The values of chlorophyll a were diverse over the 10-year period, but the dynamics of its concentration followed a similar pattern in all years; high values usually persisted from May to September (Fig 4A). Cyanobacterial blooms developed mainly from July to September and were the most intense at the TR and BR stations. In 2000, blooms were not observed (Fig 4B).

Ordination of sampling terms with the self-organising map
In the output layer of the SOM, the following three main clusters of neurons were distinguished ( Fig 5): 1. CL-ExSp (abbr. extreme/spring; including neurons A1-B3), referring to an "unstable cold season" characterized by extreme values and a high variability of abiotic factors, resulting in a dominance of the abiotic control over the biota.
2. CL-StSm (abbr. stable/summer; including neurons B4-D4), characterised by stable nonextreme abiotic conditions, making biotic interactions more important. 3. CL-ExSm (abbr. extreme/summer; including neurons E1-F4), corresponding to a "warm season with thermal or hydrological extremes" when the role of biotic control was decreased due to extreme values and/or a high variability of abiotic factors.
In the above, we intentionally used alternative terms for spring and summer because of the differences between years, which were large enough to make usage of the words "spring" and "summer" and analyses based on such divisions misleading. This is best shown by the uneven distribution of sampling terms from some years over the SOM. For example, the coldest and Long-Term Patterns in Plankton Dynamics-SOM Approach most hydrologically unstable year, 2001, was almost entirely a spring in the ecological sense (most of the sampling terms from 2001 are assigned to ESp ("early spring"), Sp ("spring") and LSp/Sm ("late spring/summer)). The year 2000 was very homogenous because the majority of its sampling terms were grouped in LSp/Sm ("late spring/summer"). In contrast, the year 2008 was almost completely devoid of the ecological spring phase because few sampling terms were assigned to ESp ("early spring"), Sp ("spring") and LSp/Sm ("late spring/summer"). In the summer of that year, similar to 2004, the most extreme conditions were recorded (sampling terms from 2004 and 2008 are the most frequent in cluster CL-ExSm "Extreme/Sumer") (S1 Fig). In this paper, we distinguished ecological states or ecological seasons (because the SOM classification was based only on density/biomass data for plankton components) and then looked for explanations of the observed patterns. The above differences between years show the advantage of our approach over others based on the calendar and/or classical seasons.
The order of sub-clusters LSp/Sm and Sm was determined on the basis of the gradients observed in abiotic factors (Figs 6 and 7). The water temperature at the three sampling stations increased for the successive five sub-clusters, from ESp to HotSm (ESp, Sp, LSp/Sm, Sm and HotSm), and the temperature of HunLSp/Sm was lower than that of HotSm. The temperature of HunLSp/Sm was similar to that of LSp/Sm or Sm. A similar upward trend was recorded for five sub-clusters, from ESp to HotSm, in the retention time of water in the reservoir; the median for HunLSp/Sm was similar to that for Sp. The inflow of water to the reservoir decreased in five successive sub-clusters ESp-HotSm and increased in HunLSp/Sm (Fig 7).  ). An asterisk indicates SOM sub-clusters for which a significant difference between sampling stations was recorded (see Table 1 for details). Abbreviations of the sub-clusters names: 1 ESp-early spring; 2 Sp-spring; 3 Sp/Sm-spring/summer; 4 Sm-summer; 5 HotSm-hot summer; 6 HunLSp/Sm -hydrologically unstable late spring/summer.
Regarding the biotic variables, the patterns that best reflected those mentioned above (for the abiotic factors) were recorded for cyanobacterial biomass, which increased in five successive sub-clusters ESp-HotSm and decreased in HunLSp/Sm (Figs 6 and 8). Increases in cyanobacterial biomass were also visible for four sub-clusters when moving downstream (from ZA to TR); the most distinct of these increases were recorded for Sm, LSp/Sm and HunLSp/Sm. Chlorophyll a increased in three successive sub-clusters, LSp/Sm, Sm and HotSm, and the differences were more pronounced downstream (Fig 8, Table 1).
Another pattern was observed in the densities of the two Cladocera species. The densities of Daphnia longispina and Leptodora kindtii generally increased in sub-clusters ESp, Sp and LSp/ Sm, followed by (after they dropped in Sm compared to LSp/Sm) HotSm and HunLSp/Sm (Figs 6 and 8). Additionally, L. kindtii decreased downstream in HotSm and HunLSp/Sm (Fig  8, Table 1).
Detailed results of the Kruskal-Wallis and post-hoc Dunn tests used for statistical verification of the significance of observed differences between the SOM sub-clusters are presented in    Table 1.

Characteristics of abiotic and biotic factors in sampling terms assigned to the SOM sub-clusters
All sampling terms were assigned to the neurons in the CL-ExSp, CL-StSm and CL-ExSm clusters according to their biotic and abiotic characteristics (more details in Table 2): 1. Sub-cluster ESp ("early spring") grouped sampling terms with very low densities of D. longispina (median of 1.19 ind dm -3 ) and L. kindtii (median close to zero). The median concentration of chlorophyll a amounted to 11.84 μg dm -3 , but the median cyanobacterial biomass was 0.09 mg dm -3 . Sampling terms grouped in ESp were characterised by high values of inflow (median > 40 m 3 s -1 ), very short retention times (median < 21 days) and low temperatures (13.86°C at three stations). Such unstable hydrology, low temperatures and low values of biotic parameters were specific to the period of early spring. Sampling terms grouped in ESp did come predominantly from May and April.
2. Sub-cluster Sp ("spring") grouped terms in which the median density of D. longispina increased to 5.79 ind dm -3 but the density of L. kindtii remained close to zero. Values for the chlorophyll a concentration and cyanobacterial biomass were similar to those in sub-cluster ESp and reached medians of 12.38 μg dm -3 and 0.14 mg dm -3 , respectively. Sp included the sampling terms of May and June with persistent unstable hydrology (median inflow value > 30 m 3 s -1 and retention time < 29 days). The median temperature was 17.98°C.

Results of canonical correlation analyses
For cluster CL-ExSp (extreme/spring), only the first canonical correlation (0.678; p<0.01), which was calculated between the first and second correlation pairs, was significant ( Table 3). The highest contribution to the formation of the U1 canonical root for predictor variables belonged to retention time, with a canonical weight of 1.252. The highest contribution to the formation of the V1 canonical root for the response variables was the D. longispina density (0.845). Considering the values of the factor loadings, which showed the overall correlation of the respective variables with the canonical roots, the highest loading for the U1 and V1 canonical roots was temperature (0.890) and D. longispina density (0.897), respectively. The redundancy of response variables calculated for the first canonical correlation indicated that 17.5% of the biotic factor variance was explained by predictor variables. When the correlation coefficients between abiotic and biotic factors were examined (Table 3), most of them were statistically significant; the correlation coefficients between biotic factors were very low, and all but one was statistically insignificant.
For cluster CL-StSm (stable/summer), both the first and second canonical correlations were statistically significant (0.666; p<0.01 and 0.391; p<0.01, respectively). The highest contributions to the U1 and V1 canonical roots had a retention time and cyanobacterial biomass with canonical weights of 1.911 and 0.656, respectively. Water inflow (1.730) and D. longispina density (-0.972) were the main contributors to the U2 and V2 canonical roots. In turn, the highest loading for the U1 and V1 canonical roots was the chlorophyll a concentration (0.840) and cyanobacterial biomass (0.919), respectively, and for the U2 and V2 canonical roots was water inflow (0.742) and L. kindtii density (0.633), respectively (Table 4). For the two significant canonical correlations, 27.4% of the biotic factor variance was explained by predictor variables. Contrary to cluster CL-ExSp, in cluster CL-StSm, most of the correlation coefficients between abiotic and biotic factors were weak, and the correlation coefficients between biotic variables were much higher and reached the highest (negative) value (-0.640) for cyanobacterial biomass and D. longispina density (Table 4).
For cluster CL-ExSm (extreme/summer), all three canonical correlations were significant (0.678; p<0.01, 0.560; p<0.01, and 0.236; p<0.05, respectively). The highest contributions to the U1 and V1 canonical roots had a water inflow and cyanobacterial biomass with canonical weights of -0.667 and 0.776, respectively. Water inflow (-0.957) and L. kindtii density (-0.901) were the main contributors to the U2 and V2 canonical roots, while temperature (-0.862) and D. longispina density (-1.060) were the main contributors to the U3 and V3 canonical roots, respectively. Similar to the canonical weights, the highest canonical factor loadings for U1 and V1 belonged to water inflow (0.715) and cyanobacterial biomass (0.925). In turn, the highest Long-Term Patterns in Plankton Dynamics-SOM Approach canonical factor loading was inflow (-0.689) for U2 and D. longispina density (-0.639) for V2. For U3 and V3, the highest factor loadings were found for temperature (-0.681) and D. longispina density (-0.752), respectively (Table 5). For the all three canonical correlations, the predictor variables explained 33% of the variance of biotic factors. All correlations between abiotic and biotic factors were significant. Among the biotic factors, the highest (negative) correlation was found for L. kindtii density and cyanobacterial biomass (-0.505) ( Table 5).

Discussion
In this study, we used an SOM to distinguish fairly homogeneous classes (clusters and sub-clusters) of sampling terms with similar biotic states, which then enabled the identification of the crucial abiotic factors responsible for the seasonal sequence of events. This allows us to now discuss the main patterns of seasonal variations in D. longispina, L. kindtii and cyanobacteria observed in our study in relation to the steps proposed by the classical model of plankton succession, the PGE model. The conceptual framework for this model concerning the seasonal dynamics of plankton communities in an idealized temperate lake was developed by an international group of limnologists, the Plankton Ecology Group (PEG), in the 1980s. The original PEG model described the roles of different abiotic and biotic factors driving the succession of zoo-and phytoplankton in 24 sequential steps [47].
The PEG model is mainly based on the temporal sequences and assumes the progressive growth of zooplankton during spring, when small, fast-growing algae are available [47].  (Table 3). In reservoirs, increases in water levels are associated with varying dynamics in the zooplankton population [48] due to both hydrological advection and an increased concentration of suspended sediments [49]. Instability in water depth may change the underwater light climate as well as the nutrient dynamics, which both affect phytoplankton biomass and species composition [50] and can indirectly modify the succession of daphniids [51]. High water inflow conditions directly remove the large-sized species of cladocerans and copepods and favour the development of rotifers [52]. Such flushing effects on macrozooplankton may be periodically stronger than the top-down control exerted by planktivores [27]. In the Sulejow Reservoir, the negative effect of the hydrological conditions (high inflow and short retention time) most likely determined the very low density of D. longispina in April and May (Table 3). Unfortunately, the hydrological regime, which is an important driver of plankton seasonality in dam reservoirs, was not included in the PEG model focused on lakes. Moreover, the model (even in the updated version [53]) does not consider the temperature as the crucial factor for zooplankton spring phenology [53,54], even though the close dependence of daphniid population dynamics on epilimnetic water temperatures during spring has been documented by many authors [21,55]. In addition, some recent analyses of longterm data have indicated that spring development of Daphnia has mainly been regulated by temperature and not by food availability [56]. The temperature directly (via changes in life processes) and indirectly (via changes in habitat properties) influences Daphnia spp. populations [56]. The response of Daphnia to temperature variability depends on species-specific thermal tolerance ranges. Previous observations of D. longispina under natural conditions in a temperate zone indicated that this species became highly abundant at temperatures above 14°C, and their favourable temperature range was 15-21°C [57]. In terms grouped in sub-cluster ESp ("early spring"), the temperature reached a median of 13.86°C. Thus, our results confirmed that the very low water temperature (< 14°C) was another critical factor in the development of the Daphnia population in early spring, in addition to hydrology (Fig 7). The terms of sub-cluster Sp ("spring") were characterised not only by a raised temperature (median range for the three stations of 17.80-18.19°C) but also by a lower inflow (median of 30.57 m 3 s -1 ) and longer retention time (median of 28 days) compared to ESp. Because the chlorophyll a concentrations in Sp were comparable to those in sub-cluster ESp, our results indicate that the increase in Daphnia density in spring (May-June) was still determined by abiotic conditions: an increase in water temperature and more stable hydrological conditions compared with early spring (Figs 7 and 8; Table 3). One of the most important sequential statements of the seasonal succession of plankton is the "clear water phase" (CWP) that usually occurs in late spring. The PEG model describes the CWP as a decline in phytoplankton biomass towards a mid-season biomass minimum due to the intensive grazing pressure of zooplankton [47]. This event is connected with a high peak of planktonic herbivores, especially daphniids and copepods [58]. A similar pattern was observed in our study: the highest density of Daphnia was observed in June [sampling terms grouped in LSp/Sm ("late spring/summer")] (Fig 8), i.e., at the same time as the appearance of the lowest concentration of chlorophyll a in stations BR and TR (median of 3.40 μg dm -3 for BR and 5.68 μg dm -3 for TR). The low concentration of chlorophyll a in late spring or early summer, i.e., when Daphnia populations were observed to be at their highest densities, could be attributed to the grazing of Daphnia on algal biomass, resulting in the spring clear-water phase (CWP). This phenomenon, characteristic of the plankton phenology of many temperate lakes [47], has been regularly observed in the Sulejow Reservoir, which is a eutrophic (with high levels of nutrients during the spring), shallow and polymictic ecosystem [2]. Thus, CWP development appears to be regulated in the reservoir mainly by selective grazing of Daphnia on smallsized algal populations, rather than by thermal stratification or seasonal changes in nutrient concentrations [59]. The duration of the CWP in the Sulejow Reservoir differed between years and was usually terminated by a collapse in the Daphnia abundance and the intensive expansion of cyanobacteria [35], as was also recorded for sampling terms in sub-cluster Sm ("summer") (Figs 6 and 8).
The terms of sub-cluster Sm ("summer") (September and July) were characterised by the lowest inflow, the longest retention time and a high median water temperature (> 20°C; Figs 6 and 7). Consequently, the analysed abiotic factors were completely different than in the sampling terms grouped in sub-clusters ESp, Sp and LSp/Sm corresponding to spring conditions. A stable hydrological regime and high temperature may increase the importance of biotic interactions [19], which is consistent with the assumptions of the PEG model [47]. In eutrophic lakes and reservoirs, such abiotic conditions are favourable for cyanobacterial growth [28]. We observed a clear relationship between cyanobacterial biomass and hydrological conditions in a spatial gradient of terms grouped in sub-cluster Sm ("summer"); blooms were more frequent at the TR and BR stations than at the ZA station (Figs 2, 4 and 8). However, in a temporal gradient, the highest cyanobacterial biomass appeared at the highest temperatures and in the most stable hydrological conditions represented by sub-cluster HotSm ("hot summer") (June, August) ( Table 2). Cyanobacteria generally grow better at higher temperatures (above 20°C) than other phytoplankton organisms, which gives them a competitive advantage in warm periods [25]. Apart from rising temperatures and nutrient concentrations, stable hydrological conditions are recognised as one of the most important factors in determining the development of bloom-forming cyanobacteria [60,61]. An increased retention time favours cyanobacterial dominance [62], which was also confirmed by our analyses (CCA, Table 4). Analyses of the correlation matrices revealed that the abundance of cyanobacteria was closely dependent on abiotic factors (temperature and hydrology) within all seasons in the studied years (Tables 3-5).
In the summer months (sub-clusters Sm and HotSm), the dynamics of the D. longispina population were dependent on biotic factors, mostly cyanobacteria. Cyanobacteria are well-known poor food sources for daphniids due to their low nutrient content, toxicity and colony and filament formation [63,64]. A negative correlation between the density of daphniids and the biomass of toxin-producing cyanobacteria has been previously reported [65,66]. The PEG model also assumed growth limitation of herbivore crustaceans resulting from the dominance of bloomforming, less edible algae or cyanobacteria during summer [47]. In our study, the sampling terms associated with sub-cluster Sm ("summer") were characterised by a high cyanobacterial biomass and a low density of D. longispina (Fig 8). Indeed, a strong negative correlation between D. longispina and cyanobacteria in the CL-StSm cluster was found (Table 4). We suppose that the negative impact of cyanobacteria on D. longispina resulted from their low nutritional value or to mechanical damage to the filtering apparatus from their colonies and filaments [67], rather than to cyanobacterial toxins. These conclusions are based on our previous studies indicating that generations of daphniids from the Sulejow Reservoir contained effective antioxidant systems to protect themselves against the accumulation of cyanobacterial toxins and their harmful effects [13,14]. In the Sulejow Reservoir, the concentrations of microcystins reached an average of 5.83 μg dm -3 , and, in the areas of the bloom, the concentrations were recorded at levels up to 22-30 μg dm -3 [2]. Despite this observation, the population of D. longispina was eventually reduced in abundance, but it did not completely decline during the blooms (as presented in this study and [14]).
In the PEG model, the effects of higher trophic levels on herbivorous zooplankton are only restricted to fish predation [47]. In addition, the updated version of the model [53] does not include the impact of invertebrate predators on filter feeders. However, our previous study on the Sulejow Reservoir revealed that the abundance of daphniids in the summer months was more dependent on L. kindtii predation than on fish pressure [17,39]. In June-July, Leptodora, due to their high densities in the pelagic zone, were able to eliminate 10-51% of Daphnia biomass, whereas fish eliminated only 0.1-5.4% [17]. The predatory rate of L. kindtii depends on temperature [21] and the density of the prey population [18]. Therefore, the density and predatory pressure of Leptodora increased during a period with a relatively high abundance of potential prey (Daphnia). This effect was particularly evident during summer months at the ZA and BR stations (Fig 3). The SOM showed that both the spatial distribution (in particular stations) and the temporal dynamics of daphniids and L. kindtii were similar in sampling terms grouped in sub-cluster HunLSp/Sm ("hydrologically unstable late spring/summer") (Fig 8). A comparable sequence of growth and decline in the densities of these animals corresponded to the classical coupled system of predator and prey oscillations observed in other lakes [68]. Thus, the SOM results confirmed that the abundance of Daphnia was regulated by the predation of L. kindtii during summer. However, because the peaks in Leptodora density were short in duration and their frequency was different within the studied years (Fig 3), the precise estimation of L. kindtii pressure on daphniid population dynamics should be measured individually for each season. Analysis of the correlation matrices revealed that Leptodora, similar to Daphnia, avoided summer blooms. Interestingly, there was also a positive influence of increased inflow on L. kindtii density ( Table 5). The positive correlation between inflow and L. kindtii density might be related to the increased turbidity of the water during periods of intensive flow. Zettler and Carter [69] observed that a reduction in transparency protected this large cladoceran from visually seeking fish predators and contributed to its increased abundance in Lake Temiskaming. Moreover, in reservoirs, the raised concentration of suspended organic matter caused by inflow is often associated with peaks in bacteria, protists and their consumers, such as Bosmina spp. [70]. These small cladocerans can be good food resources for L. kindtii [71]. Our assumptions, however, require confirmation by further detailed research.
In this paper, we analysed and discussed the role of selected mechanisms that drive temporal changes in the structure of plankton communities. However, a more comprehensive study of the Sulejow Reservoir may be necessary for expanding the protection and restoration strategies for this eutrophic ecosystem. Thus, further research projects focused on light availability and nutrient uptake as determinants of phytoplankton growth, the effect of juvenile fish on zooplankton during summer and the importance of the microbial food web for herbivorous zooplankton community dynamics have been developed.

Conclusions
In conclusion, we recommend SOMs as an appropriate method for the analysis of complex patterns in reservoir habitats. We believe that our results are important in providing a better understanding of L. kindtii-Daphnia-cyanobacteria interactions in an ecosystem that is highly dependent on the hydrological regime. Long-term studies generate comprehensive and detailed databases, which, when analysed by the appropriate mathematical approaches, can reveal ecosystem patterns, dynamics and control processes. The results of long-term research provide knowledge to the broader scientific community concerning the function of aquatic ecosystems. They are also essential in the development of management and policy decisions regarding the protection of ecosystems [72]. The application of an SOM in our research allowed the determination of a cause-and-effect understanding of the relationship between selected biotic and abiotic parameters. The importance of hydrology (flow rate, retention time, water-level fluctuations, etc.) as a plankton-regulating factor may be essential for the formulation of management decisions in dam reservoirs. This mainly refers to the regulation of retention time with respect to the control of the seasonal dynamics of large grazers, such as Daphnia spp., and/or as a tool for decreasing the occurrence of algal blooms. This method is especially recommended for anthropogenically modified ecosystems with toxic cyanobacterial blooms [28].
Supporting Information S1 Fig. The 412 sampling terms assigned to the 24 SOM output neurons. The code for each term consists of two digits for the day, two digits for the month and two digits for the year of sampling, e.g., 190799 = 19 th of July 1999. (TIF) S1