Spatio-Temporal Dynamics in Physico-Chemical Properties, Phytoplankton and Bacterial Diversity as an Indication of the Bovan Reservoir Water Quality

The study aimed to investigate the physico-chemical properties as well as phytoplankton and bacterial community diversity of Bovan Lake reservoir in Serbia to gain insight into the seasonal dynamic of water quality. All analyses were performed at three localities and water depths in spring, summer, autumn, and winter 2019. Seven phytoplankton phyla comprising 139 taxa were detected at all three localities (Chlorophyta 58%, Bacillariophyta 14%, Cyanobacteria 9%, Chrysophyta 5%, Dinophyta 5%, Euglenophyta 5%, and Cryptophyta 4%). Winter 2019 was characterized by the presence of 36 unique species in all phyla except Euglenophyta. Bacterial diversity analyses showed that Proteobacteria, Actinobacteria, Bacteroidetes, Cyanobacteria, and Verrucomicrobia dominated the water intake locality at all three water depths (0.5, 10, and 20 m below the surface). In general, the physico-chemical parameters, phytoplankton, and bacterial community composition depended on the season and the water depth and showed that Bovan Lake was of satisfactory ecological status and water quality at all localities, and meets the needs for all intended purposes.


Introduction
Freshwater represents an invaluable yet limited and unequally distributed resource. Surface water bodies are the most important sources of freshwater [1], among which lakes contain around 90% of all of Earth's aquatic ecosystems [2]. Lakes are classified as either natural or artificial and serve as an important water supply source for human utilization [1]. Many organisms depend on freshwater for survival, and humans frequently depend on lakes for many 'goods and services', such as drinking water, industrial activity, waste removal, fisheries, agricultural irrigation, recreation [3], commerce, transportation, tourism, and energy production [2,4].
Most freshwater aquatic ecosystems are rapidly deteriorating in quality due to changes in climatic conditions and anthropogenic activities [5,6], leading to eutrophication and habitat pollution [7,8]. Urbanization and modern, intensive agricultural practices have increased the release of untreated wastewater from households, industries, and farms into lakes [9][10][11]. This severely affects the water's physical, chemical, and biological properties and causes an imbalance in aquatic organisms' structure and functioning, and makes the water unfit for the intended purposes.
The lake ecosystems' physical, chemical, and biological components are greatly interdependent and in constant interaction [1]. It is well known that water parameters depend highly on environmental conditions. For instance, the influx of nutrients into lakes leads to the increasing productivity of aquatic ecosystems, which can cause excessive algal biomass increase. As a major producer of organic matter in aquatic ecosystems, any disruption in the community structure or changes in the number of certain phytoplankton species could affect the water quality, especially in water supply reservoirs. The most harmful are the toxin-producing cyanobacteria, which can cause acute or chronic health issues in plants, animals, and humans [12,13]. Since cyanobacteria coexist and interact with heterotrophic bacteria in freshwater ecosystems, it is not surprising that the composition of cyanobacterial and phytoplankton communities greatly impacts the structure of the associated heterotrophic bacterial communities [14,15]. Moreover, due to the diverse metabolic features and the great adaptability to environmental fluctuations, bacteria represent good and practical indicators of surface water pollution [5]. Another great concern for water quality is the presence of pathogenic microorganisms and heavy metals, which can have a devastating impact on human and animal health [6].
Furthermore, water properties are prone to changes on a day-to-day basis, with fluctuations evident on both temporal and spatial scales [1]. Regular and exhaustive monitoring of freshwater sources, which includes the analysis of physical, chemical, and microbiological parameters, is necessary to maintain the proper and desired water quality needed for all the intended uses.
Serbia has relatively extensive groundwater resources, but they are not uniformly distributed [16]. The national strategy for the water supply of the Republic of Serbia encompasses a plan for regional systems of the water supply by the construction of 29 reservoirs [17].
This study aimed to gain insight into the spatio-temporal changes in the water quality of Bovan Lake by assessing the physico-chemical parameters as well as phytoplankton and bacterial community diversity. The set goal was achieved by conducting the aforementioned assessment at different localities (at the inflow of the river Moravica, the main bathing area, and the water intake area), water depths (0.5 m, 10 m, and 20 m below the surface), and time points (spring, summer, autumn and winter 2019). In addition to conventional cultivation methods, metagenomic 16S amplicon sequencing was employed to obtain reliable information on bacterial diversity and, consequently, the water quality at the water intake locality of Bovan Lake.

Description of the Study Area
Lake Bovan, one of the 29 reservoirs mentioned above, is an artificial accumulation formed between 1978 and 1984 following the dam's construction in the middle flow of the river Sokobanjska Moravica in Serbia. The primary aim of the lake building was the regulation of the Morava River basin and protection of the Djerdap I reservoir [18]. The Bovan accumulation arrests aqueous deposit sediment and is an irrigation and electricity source. It is also pivotal for supplying drinking water to the municipality of Aleksinac with its 72 settlements [16,19]. The entire population of the Aleksinac municipality is estimated at 55,000, of which about 65% live in rural settlements. The water supply for 20 settlements and nearly 35,000 inhabitants is regulated from a unitary system, the water treatment plant Bresje, at Lake Bovan [16]. Lake Bovan is also an important tourist attraction, since it offers recreational fishing, sailing, and several bathing areas [20]. Lake Bovan, therefore, represents a multifaceted accumulation, and its recreational and water-supplying purposes make it an important object for careful and thorough water quality monitoring.

Sampling Sites
Sampling was carried out at the inflow of Moravica River (N 43 •  Seven water samples were collected in total at the three locations at Bovan Lake, at three water depths (0.5 m, 10 m, and 20 m below the surface) in spring, summer, autumn, and winter of 2019 (Table S1, Supporting Information). The sampling site 'the inflow of Moravica River' had a low water column, so the samples were taken from a depth of just 0.5 m.

Sampling Methodology
Three depths were chosen for sampling: one near the surface of the water (0.5 m), another at the middle of the water column (10 m), and the third near the bottom of the lake (20 m). The samples for quantitative phytoplankton analyses, determination of chlorophyll a concentration, and microbiological analyses, as well as those for analyses of chemical and physical parameters of the water, were collected by Van Donr's Water Sampler (volume 2 L) by lowering it to the depths mentioned above according to the Rulebook on the method of sampling and methods for laboratory analysis of drinking water [21]. The collected samples were placed in sterile plastic or glass bottles (volume 1 L) and transported to laboratories in portable refrigerators. The average time from sampling to the beginning of the analyses was eight hours.

Physico-Chemical Analyses
The physico-chemical analyses were carried out at the Institute of Chemistry, Technology and Metallurgy (University of Belgrade) and Anahem Laboratory in Belgrade. The general properties of the water samples, salinity, hardness, oxygen regime, organic load, nutrient presence, and metal concentrations were determined using the following methods: pH (

Phytoplankton Analyses
The phytoplankton net collected water samples for the qualitative phytoplankton analysis (mesh size 22-23 mm, net frame 25 cm ø), neglecting the volume of the water. The collected samples were placed in sterile plastic bottles (volume 100 mL) and preserved by adding Lugol's solution (at ratio 1:100), except those for chlorophyll a analysis, according to European standard EN 15204:2008. All samples were transported to the laboratory in a hand-held refrigerator at 4 • C in the dark. The transparency of the water was measured in situ using a Secchi disk.
The chlorophyll a concentration analysis (ISO 10260:1992(E)) was based on the extraction of filtered phytoplankton biomass by the hot ethanol and extract absorbance measurements at 665 and 750 nm before and after acidification with 3M HCl. The final concentration of chlorophyll a was calculated and expressed in µg l −1 . The trophic status of Bovan Lake was determined according to [22].
The quantitative phytoplankton analyses were done according to the Utermöhl method [33] and European standard (EN 15204:2008) using the Leica inverted microscope and Hydro-Bios sedimentation chambers. A final abundance of phytoplankton organisms was expressed as a number of individuals per liter and a number of cells per liter. After transforming the collected data into percentages, the information was used to generate a heatmap using the R package ggplot2 [34] via RStudio version 1.2.1335-1.
The Shannon-Weaver index of diversity was determined according to [35].

Cultivation Methods
To determine the number of facultative heterotrophic (aerobic organotrophic) and oligotrophic bacteria, decimal dilutions of water samples were inoculated on Tryptone Soy Agar and R2A agar plates, respectively, incubated five to seven days at 22-26 • C, and the number of cells per mL was calculated. Additionally, FO/H index and water quality class were determined by [36]. The concentration of aerobic mesophiles and psychrophiles, amylolytic, lipolytic, proteolytic and saccharolytic bacteria was determined according to [37].

Metagenomic Analysis
Metagenomic 16S amplicon sequencing was performed for the water intake locality in winter 2019 (samples S5, S6, and S7). Exactly 100 mL of each sample was filtered through 0.22-µm mixed cellulose ester membrane filters (Merck Millipore Ltd., Darmstadt, Germany). The extraction of total DNA from filters was performed by ZymoBIOMICS™ DNA Miniprep Kit (Zymo Research, Irvine, CA, USA) following the manufacturer's instructions, and the DNA concentration was determined by Qubit™ 4 fluorometer (Invitrogen, Waltham, MA, USA).
Amplification and sequencing of the 16S rRNA gene V3-V4 region was performed at Novogene Company Ltd. in Hong Kong, China. All PCR reactions were carried out with Phusion ® High-Fidelity PCR Master Mix (New England Biolabs) containing 341F (5 CCTAYGGGRBGCASCAG 3 ) and 806R (5 GGACTACNNGGGTATCTAAT 3 ) primers. The PCR products were checked for size by electrophoresis on 2% agarose gel, and purified with Qiagen Gel Extraction Kit (Qiagen, Hilden, Germany). The amplicon libraries generated with NEBNext ® Ultra™ DNA Library Prep Kit for Illumina ® were quantified by Qubit and sequenced by the IlluminaNovaSeq 600 platform using a paired-end 2 × 250 bp strategy.
The paired-end reads were trimmed to remove barcode and primer sequences and merged using FLASH (v1.2.7) [38]. The raw sequences were demultiplexed and qualityfiltered to obtain the high-quality clean tags [39], according to the Qiime (v1.7.0) [40]. The tags were compared with the reference database (Gold database, v11) using the UCHIME algorithm [41,42], which allowed the detection and removal of chimera sequences, thus obtaining the effective tags.
Sequences analysis was performed by Uparse software v7.0.1001 [43] using all of the effective tags. To the same OTUs, all sequences with ≥97% similarity were assigned. A Mothur software was performed against the SSU rRNA database of SILVA Database for each representative sequence [44] for species annotation at each taxonomic rank-kingdom, phylum, class, order, family, genus, and species (Threshold:0.8~1) [45]. To obtain the phylogenetic relationship of all OTUs representative sequences, the MUSCLE (v3.8.31) [46] was used to compare multiple sequences rapidly. OTUs abundance information was normalized using a standard sequence number corresponding to the sample with the least sequences. The data were subjected to alpha and beta diversity analyses.
Beta diversity analysis was used to evaluate bacterial diversity and species complexity differences between samples. Beta diversity on both weighted and unweighted UniFrac was calculated by QIIME software (v1.7.0). Principal component analysis (PCA) was applied to reduce the dimension of the original variables using the FactoMineR package and ggplot2 package in R software (v2.15.3). Principal Coordinate Analysis (PCoA) was based on weighted or unweighted UniFrac distance matrices to compare the differences between samples [47]. The PCoA analysis was displayed by WGCNA package, stat packages, and the ggplot2 package in R software (v2.15.3). Unweighted Pair-group Method with Arithmetic Means (UPGMA) clustering was conducted by QIIME software (v1.7.0) to interpret the distance matrix using average linkage. The data were made publicly available at the GenBank database (NCBI) under BioProject Acc. No. PRJNA719709.

Physico-Chemical Analyses
The physico-chemical parameters were measured at three different localities of Bovan Lake, three water depths, and four sampling seasons. In 2019 (Table S2a-c, Supporting Information), the water temperature increased from spring to autumn at all localities, while it decreased with depth (except for samples S5, S6, and S7 in spring) and was appropriate to the season of sampling. The lowest pH was detected in spring and the highest in summer at 0.5 m below the surface. The electrical conductivity was highest in spring (443-454 µS cm −1 ) for all localities and all water depths, showing the trend of decline in summer (371-437 µS cm −1 ) and autumn (312-445 µS cm −1 ), except for samples S4 and S7. When comparing different localities, suspended solids in spring were the highest at the inflow of Moravica River (S1) and the water intake (S6 and S7, respectively) in summer and autumn 2019. Dissolved oxygen decreased from spring to autumn (except in S1, where it was the highest in summer). The water bathing area (S2, S3, and S4) was the locality with the highest dissolved oxygen values, except in autumn when the highest dissolved oxygen was measured at 0.5 m below the surface at the water intake site.
In winter 2019 (Table S2d, Supporting Information), the electrical conductivity of water was 427 µS cm −1 in S5, 415 µS cm −1 in S6, and 430 µS cm −1 in S7, while the pH was 8.2, 8.1, and 8.3 in S5, S6, and S7, respectively. Chemical and biochemical oxygen demands were the highest at 0.5 m below the surface (47 and 15 mg L −1 , respectively), with an evident decrease with depth. Total nitrogen was 0.12 mg L −1 for S5, 0.31 mg L −1 for S6, and 0.32 mg L −1 for S7. Ammonium ion concentration ranged from 0.63 mg L −1 in S5 to 1.50 mg L −1 in S6. Nitrite concentration was less than 0.03 mg L −1 in all three samples, while the lowest concentration of nitrates was detected in S6 (0.95 mg L −1 ) and the highest in S5 (1.30 mg L −1 ). Total phosphorous concentration was 0.02 mg L −1 in S5 and 0.05 mg L −1 in S6 and S7. In addition to the common physical and chemical parameters, the concentration of 24 metals was determined in this survey for samples collected in winter 2019 (Table S3, Supporting Information).

Phytoplankton Analyses
The overall qualitative analysis of phytoplankton comprising all sampling seasons and sites showed the presence of 139 taxa belonging to seven phyla (Table S4, Supporting Information). Species from the Chlorophyta phylum constituted 58% of the total phytoplankton, followed by Bacillariophyta (14%), Cyanobacteria (9%), Chrysophyta (5%), Dinophyta (5%), Euglenophyta (5%), and Cryptophyta (4%) ( Table S4, Supporting Information). The number of phytoplankton taxa varied with the season of sampling. The highest number of different taxa was recorded in autumn at 75, followed by 68 in winter, 47 in spring, and 43 in summer 2019 (Table S4, Supporting Information). When comparing the sampling seasons, Figure 1 shows that Chlorophyta were the most represented in all seasons, ranging from 36% in spring 2019 to 71% in winter 2019. Following Chlorophyta, Bacillariophyta taxa were the highest in number in all seasons except autumn 2019, when Cyanobacteria outnumbered them. spring, and 43 in summer 2019 (Table S4, Supporting Information). When comparing the sampling seasons, Figure 1 shows that Chlorophyta were the most represented in all seasons, ranging from 36% in spring 2019 to 71% in winter 2019. Following Chlorophyta, Bacillariophyta taxa were the highest in number in all seasons except autumn 2019, when Cyanobacteria outnumbered them. The abundance of phytoplankton taxa is shown in Figure 2. The most abundant were species belonging to Chlorophyta and Bacillariophyta. Cyclotellaocellata was detected in all localities, water depths, and sampling seasons except winter 2019, and was represented mainly by the greatest number of individuals per l, apart from S1, S3, and S7 sampled in spring 2019. Except for the water intake at a depth of 20 m in winter 2019, where it was not recorded, Cyclotellameneghiniana was the second most abundant species among the Bacillariophyta. However, Navicula sp. and Nitzschiapalea stood out as highly abundant at the inflow of river Moravica in spring and autumn, respectively.
Concerning Cyanobacteria, most species were detected in autumn at the main bathing area at a depth of 0.5 m, and Aphanocapsaholsatica could only be detected in autumn at all localities and depths, and in winter at 10 and 20 m of the depth of the water intake locality.
Other phyla were characterized by an uneven distribution of species among the samples and no regularities regarding the concentration of individual species.
In general, winter 2019 could be easily distinguished among other sampling seasons by the presence of 36 unique species in all phyla but Euglenophyta (Table S4, Supporting Information; Figure 2). The most abundant phytoplankton species at the water intake locality in winter was Cryptomonaserosa at all three sampling depths ( Figure 2).
Based on the chlorophyll a concentration, the trophic levels of the water samples ranged from oligotrophic to meso-eutrophic, depending on the season of sampling, sampling site, and depth. The Shannon-Weaver index showed the diversity was the highest The abundance of phytoplankton taxa is shown in Figure 2. The most abundant were species belonging to Chlorophyta and Bacillariophyta. Cyclotellaocellata was detected in all localities, water depths, and sampling seasons except winter 2019, and was represented mainly by the greatest number of individuals per l, apart from S1, S3, and S7 sampled in spring 2019. Except for the water intake at a depth of 20 m in winter 2019, where it was not recorded, Cyclotellameneghiniana was the second most abundant species among the Bacillariophyta. However, Navicula sp. and Nitzschiapalea stood out as highly abundant at the inflow of river Moravica in spring and autumn, respectively.
Concerning Cyanobacteria, most species were detected in autumn at the main bathing area at a depth of 0.5 m, and Aphanocapsaholsatica could only be detected in autumn at all localities and depths, and in winter at 10 and 20 m of the depth of the water intake locality.
Other phyla were characterized by an uneven distribution of species among the samples and no regularities regarding the concentration of individual species.
In general, winter 2019 could be easily distinguished among other sampling seasons by the presence of 36 unique species in all phyla but Euglenophyta (Table S4, Supporting Information; Figure 2). The most abundant phytoplankton species at the water intake locality in winter was Cryptomonaserosa at all three sampling depths ( Figure 2).
Based on the chlorophyll a concentration, the trophic levels of the water samples ranged from oligotrophic to meso-eutrophic, depending on the season of sampling, sampling site, and depth. The Shannon-Weaver index showed the diversity was the highest in autumn for all samples except S4 and S6, where the highest diversity was observed in summer (Table S5, Supporting Information). in autumn for all samples except S4 and S6, where the highest diversity was observed in summer (Table S5, Supporting Information).

Figure 2.
Heatmap showing the abundance of phytoplankton species from Bovan Lake in each of the 24 individual samples. The x-axis represents samples at three localities, three water depths, and four sampling seasons. The y-axis represents recorded phytoplankton taxa. Values are presented as percentages. The color code indicating abundance ranges from brown (low abundance) to red (high abundance).

Microbiological Analyses
The number of facultative heterotrophs (aerobic organotrophs) and oligotrophic microorganisms varied among the different localities, water depth, and sampling season ( Figure 3, Table S6, Supporting Information).

Microbiological Analyses
The number of facultative heterotrophs (aerobic organotrophs) and oligotrophic microorganisms varied among the different localities, water depth, and sampling season ( Figure 3, Table S6, Supporting Information).
The concentration of organotrophs at the inflow of the River Moravica (S1) was the highest in summer (7417 CFU mL −1 ), and the lowest in autumn (1634 CFU mL −1 ), and a similar trend was observed concerning the oligotrophic microorganisms (14,955 CFU ml −1 and 2864 CFU mL −1 in summer and autumn, respectively). At the main bathing area (S2, S3, and S4), the concentration of oligotrophic microorganisms was higher than that of organotrophs, except in autumn at the 20 m depth. The same was noticed at the water intake locality (S5, S6, and S7), except for the water sample at the 10 m depth in autumn. Seasonal The concentration of organotrophs at the inflow of the River Moravica (S1) was the highest in summer (7417 CFU mL −1 ), and the lowest in autumn (1634 CFU mL −1 ), and a similar trend was observed concerning the oligotrophic microorganisms (14,955 CFU ml −1 and 2864 CFU mL −1 in summer and autumn, respectively). At the main bathing area (S2, S3, and S4), the concentration of oligotrophic microorganisms was higher than that of organotrophs, except in autumn at the 20 m depth. The same was noticed at the water intake locality (S5, S6, and S7), except for the water sample at the 10 m depth in autumn. Seasonal and water-depth variations were also recorded for aerobic mesophiles and psychrophiles, amylolytic, lipolytic, proteolytic, and saccharolytic bacteria at all localities.
According to the number of facultative heterotrophs and oligotrophs [36], Bovan Lake water belonged to the quality class II, except for samples S4 in spring, S2 and S6 in autumn (I-II class), and S2 and S5 in winter (I class). The self-purifying ability of the water, as shown by the FO/H index, was low in S6 in autumn, and satisfactory for the rest of the samples examined at different time points (Table S7, Supporting Information).
The bacterial community composition and diversity were determined based on the metagenomic amplicon sequencing of the V3-V4 region of the 16S rRNA gene for water intake samples collected at three water depths in winter 2019 only, given the importance of this locality as a source of drinking water. The average sequence length was 411 bp for S5 and S6, and 416 bp for S7, while the total number of paired-end sequence reads obtained after chimera filtering was 461,333 (Table S8, Supporting Information).
Alpha diversity analysis was used to determine the bacterial diversity within the three water samples. Species richness was estimated by Observed OTUs, Chao1, and ACE, while the bacterial diversity was determined based on Shannon and Simpson indices (Table 1). The observed OTU numbers for S5, S6, and S7 were 1480, 1533, and 1532, respectively. Accordingly, Chao1 and ACE indicated that the species richness was the highest at a depth of 20 m and the lowest at 0.5 m. Both Shannon and Simpson indices showed that the microbial diversity was the greatest at a depth of 20 m and the lowest at 10 m.
The number of shared and unique OTUs between samples is shown in the Venn diagram ( Figure S1, Supporting Information). All three depths of water intake had 1127 According to the number of facultative heterotrophs and oligotrophs [36], Bovan Lake water belonged to the quality class II, except for samples S4 in spring, S2 and S6 in autumn (I-II class), and S2 and S5 in winter (I class). The self-purifying ability of the water, as shown by the FO/H index, was low in S6 in autumn, and satisfactory for the rest of the samples examined at different time points (Table S7, Supporting Information).
The bacterial community composition and diversity were determined based on the metagenomic amplicon sequencing of the V3-V4 region of the 16S rRNA gene for water intake samples collected at three water depths in winter 2019 only, given the importance of this locality as a source of drinking water. The average sequence length was 411 bp for S5 and S6, and 416 bp for S7, while the total number of paired-end sequence reads obtained after chimera filtering was 461,333 (Table S8, Supporting Information).
Alpha diversity analysis was used to determine the bacterial diversity within the three water samples. Species richness was estimated by Observed OTUs, Chao1, and ACE, while the bacterial diversity was determined based on Shannon and Simpson indices (Table 1). The observed OTU numbers for S5, S6, and S7 were 1480, 1533, and 1532, respectively. Accordingly, Chao1 and ACE indicated that the species richness was the highest at a depth of 20 m and the lowest at 0.5 m. Both Shannon and Simpson indices showed that the microbial diversity was the greatest at a depth of 20 m and the lowest at 10 m.
The number of shared and unique OTUs between samples is shown in the Venn diagram ( Figure S1, Supporting Information). All three depths of water intake had 1127 common OTUs. The highest number of distinctive OTUs was discovered at 10 m of depth (149), followed by 130 OTUs at 20 m and the lowest at 0.5 m of depth (121). Samples S5 and S6 shared 107 OTUs, samples S6 and S7 150, and samples S5 and S7 had 125 OTUs in common.
According to the taxonomic annotation results, each sample's taxa distribution was determined at phylum, class, order, family, and genus taxonomic rank (Figure 4; Figures S2-S5, Supporting Information). 4. Actinobacteria were the most abundant phylum, with the relative abundance ranging between 32.3% in S7 and 44.7% in S5. Proteobacteria were predominant in S7 (45.8%) and almost equally represented in S5 and S6 samples (26.9% and 26.7%, respectively). Bacteroidetes, Cyanobacteria, Firmicutes, and Verrucomicrobia were slightly more represented at 0.5 and 10 m than 20 m below the surface. Seven bacterial genera that harbor species belonging to different physiological groups responsible for the decomposition of organic matter were singled out. Their relative abundances at three different depths of the water intake locality of Lake Bovan are presented in Table 2. Table 2. Relative abundance of selected bacterial genera at three water depths (S5-0.5 m, S6-10 m, S7-20 m) of the water intake of Lake Bovan.  The relative abundance of the top 10 taxa in the phylum rank is illustrated in Figure 4. Actinobacteria were the most abundant phylum, with the relative abundance ranging between 32.3% in S7 and 44.7% in S5. Proteobacteria were predominant in S7 (45.8%) and almost equally represented in S5 and S6 samples (26.9% and 26.7%, respectively). Bacteroidetes, Cyanobacteria, Firmicutes, and Verrucomicrobia were slightly more represented at 0.5 and 10 m than 20 m below the surface.

Genus
Seven bacterial genera that harbor species belonging to different physiological groups responsible for the decomposition of organic matter were singled out. Their relative abundances at three different depths of the water intake locality of Lake Bovan are presented in Table 2. Table 2. Relative abundance of selected bacterial genera at three water depths (S5-0.5 m, S6-10 m, S7-20 m) of the water intake of Lake Bovan. Beta diversity analysis was conducted to determine the differences between samples. Weighted UniFrac and unweighted UniFrac distances were used to measure the dissimilarity coefficient between pairwise samples. The heatmap of weighted UniFrac and unweighted UniFrac distances showed that the diversity of bacterial OTUs of Lake Bovan at 20 m below the surface differed from two other samples (Figure 5a). Additionally, PCoA and PCA plots based on weighted UniFrac distances confirmed the clear separation of sample S7 from S5 and S6 (Figure 5b,c).
Beta diversity analysis was conducted to determine the differences between samples. Weighted UniFrac and unweighted UniFrac distances were used to measure the dissimilarity coefficient between pairwise samples. The heatmap of weighted UniFrac and unweighted UniFrac distances showed that the diversity of bacterial OTUs of Lake Bovan at 20 m below the surface differed from two other samples (Figure 5a). Additionally, PCoA and PCA plots based on weighted UniFrac distances confirmed the clear separation of sample S7 from S5 and S6 (Figure 5b,c).

Discussion
Changes in the physico-chemical properties of the water are relatively easy to conduct, but are not enough for their quality assessment. Checking the phytoplankton and bacterial diversity as an indication of water quality is not only necessary, but very informative. While physico-chemical properties only depict the current situation, the diversity of phytoplankton and bacteria is a window into the past and a powerful tool to assess future water quality. For several decades, the Republic Hydrometeorological Service of Serbia carried out physico-chemical and biological tests of the water quality of the Bovan reservoir under the systematic surface water quality examination program, at three localities (at the dam, in the central part of the reservoir, and at the entrance of the reservoir) and three depths (0.5 m beneath the surface, in the middle part of the water column, and near the bottom of the reservoir) [48].
Several surveys that investigated the physico-chemical parameters, phytoplankton community structure, and microbiological quality of Bovan Lake water have been con-ducted so far [18,[48][49][50]. However, the obtained data are difficult to compare and interpret, since the methodology employed and the choice of sampling locations were inconsistent. Therefore, a need for a standardized and comprehensive water quality monitoring strategy is evident.

Physico-Chemical Parameters
The physical and chemical parameters investigated in the current study showed variations depending on the location, season, and depth of sampling, which was also found in many previous studies [5,48,[51][52][53][54]. The obtained results also suggested that the water of Bovan Lake is of satisfactory quality to be exploited for drinking or recreational purposes.
The average pH of Bovan Lake at all three localities and water depths was 7.0 in spring 2019, 8.1 in summer 2019, 7.8 in autumn 2019, and 8.2 in winter 2019, which suggested a neutral to a slightly alkaline environment. The pH of water depends on the water composition and is usually in the range of 6.5-8.5. However, it has no direct impact on water consumers, so no health-associated guideline values are established by the World Health Organization [55].
Variations in electric conductivity have strong relations with pollution in lakes. In Bovan Lake, the electric conductivity was lower than 1500 µS cm −1 in all four months of sampling, which is a requirement for drinking water [55]. In most cases, higher temperatures and increased evaporation during dry seasons will equate to higher electrical conductivity. On the other hand, a large amount of precipitation during wet seasons increases overall groundwater volumes, which results in decreased electrical conductivity.
Based on the concentration of dissolved oxygen, the quality of Bovan Lake water in 2019 was the best in spring, and it decreased through summer to autumn, when the lowest values were recorded. The decrease of dissolved oxygen during the summer is mainly due to increased temperature and microbial activity, as previous research shows [5,48,51]. As a chemical parameter of water quality, the extreme decrease of dissolved oxygen (below 5 mg L −1 ) can indicate increased algal activities in the aquatic ecosystem, high concentrations of phosphate and ammonia, or high amounts of untreated wastewater inputs [5].
According to the concentrations of antimony, arsenic, barium, copper, boron, zinc, iron, chrome, cadmium, cobalt, molybdenum, nickel, lead, selenium, and mercury recorded in winter 2019, the Bovan Lake water at the water intake locality belonged to the I class [56]. The low concentrations of different metals in Bovan Lake show that the water fulfills the drinking water usage requirements/could be safely used as a drinking source.
An investigation of the water quality inĆelije, Bovan, Grliška, Gruža, Bukulja, Garaši, Borkovac, and Pavlovci reservoirs, and 20 aquatic ecosystems in Vojvodina showed accelerated eutrophication [57]. A previous study found that the highest concentration of chlorophyll a was detected in October 2004, which suggested the polytrophic status of Lake Bovan [18]. However, the trophic status of the Bovan Lake showed improvement from spring 2019, when it was mainly meso-eutrophic at all localities, to oligo-mesotrophic in winter 2019 at the water intake locality. An exception was the inflow of the Moravica River, which is not surprising given that rivers bring high amounts of nutrients.

Phytoplankton Analysis
The previous research conducted in July 2003 at the same three localities of Bovan Lake showed that the taxa from Chlorophyta were dominant, followed by Bacillariophyta, while all other phyla were present with a lower number of taxa [48], which is in accordance with the results obtained in 2019. However, the phyla distribution was slightly different in winter 2019, with Cryptophyta being the second most abundant phylum. This is not surprising given that the phytoplankton community composition depends greatly on the season. Other studies investigating phytoplankton's temporal and spatial variability in temperate region lakes showed a very similar proportion of phytoplankton phyla [51,58]. However, qualitative phytoplankton analyses of Bovan Lake in November 2005 showed the presence of species from divisions Bacillariophyta (54.16%), Chlorophyta (29.1%), Cyanophyta (12.5%), and Pyrrophyta (4.16%) [49], which differed substantially from the results obtained in the current research.
In contrast to our findings, in July 2003, the composition of Cyanobacteria and the number of species differed, with Anabaena, Aphanizomenon, Lyngbya, Microcystis, and Oscillatoria species being present at each tested locality of Bovan Lake [48].
Cyanobacteria are an important phytoplankton phylum due to their ability to cause "blooming of the water" when the environmental conditions are suitable, usually in late summer and early autumn [52]. The phenomenon threatens the health of all aquatic organisms and humans, especially when it occurs in drinking water reservoirs. The production of different toxic substances by certain members of the phylum can have serious acute, but primarily chronic, consequences to human health.
According to the qualitative and quantitative phytoplankton analyses, the diversity and quantity of Cyanobacteria were the highest in autumn 2019 ( Figure 2, Table S4). Among Cyanobacteria, Aphanocapsaholsatica stood out by its abundance (464,000 ind L −1 and 51,520,000 cells L −1 20 m below the surface at the main bathing area), thereby fulfilling the criteria for induction of "blooming", which requires ≥10,000 cells mL −1 [59]. The species was not recorded in spring or summer 2019 and was detected in significantly lower numbers in winter 2019 at the water intake locality. Given that Aphanocapsa spp. are potentially harmful [60], the World Health Organization considers Aphanocapsaholsatica a candidate for surveillance [61]; this species should be closely monitored in the future.
Benthos cyanobacteria Cylindrospermumstagnale isolated from a lake in Brazil produced neosaxitoxin, decarbamoylsaxitoxin, and saxitoxin [62]. In contrast, in Bovan Lake, it was only detected at the surface of the main bathing area in spring 2019. Additionally, Pseudanabaenacatenata, known for producing microcystins [63], was only detected at the inflow of Moravica River in spring and autumn and at the water intake locality 0.5 m below the surface in summer 2019. Both species were absent from samples collected in winter 2019. None of the three potentially harmful cyanobacteria were detected in Bovan Lake before 2014, which could be attributed to severe floods in central Serbia in the summer of that same year.

Microbial Community Analysis
Since lakes represent closed aquatic bodies and isolated ecosystems, any environmental change or increase in organic and inorganic substances can strongly affect the composition of microbial communities, which implies the need for continuous monitoring of microbial communities. Even though considered a "gold standard" in water quality monitoring, commonly employed culture-dependent methods have a significant drawback, as it is estimated that less than 1% of bacteria can currently be successfully cultivated [64]. As opposed to that, metagenomic 16S amplicon sequencing can discover the most present microbial taxa, and the composition of bacterial communities can be compared [65]. However, the culture-dependent methods remain useful as a fast and straightforward way of determining the concentrations of important bacterial groups, such as fecal coliforms, as indicators of health-endangering water pollution.
Variations in the concentration of heterotrophic and oligotrophic bacteria were detected in the current study at different localities of Bovan Lake. The values obtained varied with both the season of sampling and water depth, as already shown [18,50]. Comparison of the vertical distribution of heterotrophs showed the highest concentration of bacteria at 20 m below the surface of the main bathing area and the water intake localities in autumn 2019, while another study found the greatest number of heterotrophs at the Bovan Lake dam was in August 2006 at a depth of 10 m [50]. The highest number of heterotrophs and oligotrophs was recorded in summer 2019, which is expected, since higher bacterial counts usually characterize the summer due to the increased anthropogenic influence. Another study investigated the monthly fluctuations of different bacterial groups, including heterotrophs [5]. The researchers found that minimum numbers were recorded in January and maximum in May or June, depending on the locality, but they did not proceed with the investigations into the late summer.
The enumeration of different physiological bacterial groups involved in the decomposition of organic matter in Lake Bovan was investigated in the current study, but also in February, May, October, and December in 2004 at the inflow of Moravica River [18]. Given the 15-year period between the two investigations, the obtained data are incomparable, and more recent studies are lacking. In the current study performed in 2019 (Figure 3), the number of aerobic mesophilic bacteria at all three localities was the highest in August. Amylolytic bacteria were most abundant in autumn when their numbers were recorded at all localities and all sampling depths. No amylolytic bacteria were detected at the main bathing area 10 and 20 m below the surface and at 0.5 m and 10 m of depth at the water intake in spring, as well as in sample S6 in summer 2019. Lipolytic bacteria were recorded at the highest concentration at the inflow of the Moravica River in spring, summer, and winter. In autumn, the highest abundance of these bacteria was detected 20 m below the surface at the main bathing area and the water intake locality. In spring 2019, saccharolytic bacteria were not detected at the lowest sampling depth at the main bathing area. Their presence was also not established in summer and autumn at the main bathing area and the water intake of Lake Bovan. Saccharolytic bacteria were most abundant in spring at the water intake locality, while the concentration of proteolytic bacteria was inconsistent throughout the seasons and sampling depths.
In winter 2019, it was observed that the water level of the Bovan reservoir was drastically lower than usual. Since any significant fluctuation to the water level affects its physical and biological properties [66], this occurrence suggested changes in the bacterial community composition and the distribution of various physiological bacterial groups. Given that only a small fraction of bacteria can be cultivated, and that the commonly used cultivation methods for determining the presence and concentration of specific bacterial groups vastly underestimate the diversity of bacteria, metagenomic analysis of the water intake locality was carried out in the winter season. The obtained results provided a deeper understanding of the bacterial community composition at three different depths, which is especially significant given the importance of the locality as a source of drinking water for many human settlements in the vicinity of Bovan Lake. The metagenomic 16S amplicon sequencing showed the ten most abundant phyla detected in all three samples of the water intake locality in winter 2019 ( Figure 4) are typical bacterial phyla in lakes and are often recorded in various studies, but at different relative abundances [6,[52][53][54]. Proteobacteria, Actinobacteria, Bacteroidetes, Cyanobacteria, and Verrucomicrobia dominated the water intake locality at all three depths, which is in accordance with a similar study that showed that these five phyla accounted for 91.7% of all assigned sequences of the 16S rRNA gene [53]. In the same research, Chloroflexi, Acidobacteria, and Gemmatimonadetes were also among the top 10 most abundant phyla, which were also found in the current study. Besides the eight mentioned phyla, Firmicutes were detected in all three water samples with notable abundance, whereas in another study, it represented only 1.1% of all detected phyla [52]. Among Proteobacteria, Gammaproteobacteria were dominant at 20 m below the surface of the water intake locality, while Alphaproteobacteria were represented more at 0.5 and 10 m below the surface in this study. In contrast to these findings, Alpha-and Gammaproteobacteria rarely dominate freshwater bacterial communities, except under raised organic and inorganic loads during phytoplankton blooms [15], which, however, occur mainly during spring.
When analyzing the phyla distribution at different depths of the water column, no significant differences in the diversity of bacterial communities were obvious, only slight differences in the relative abundance of detected phyla, where sample S7 stood out as shown by the beta diversity analysis ( Figure 5). These differences could be attributed to the physical and chemical properties of the water, since it is well known that factors such as season, pH, conductivity, and dissolved oxygen can influence the composition of bacterial communities [5,6,53,67].
The results obtained from the physico-chemical, phytoplankton, and microbiological analyses indicated that the water of Bovan Lake is of satisfactory ecological status and water quality, and meets the needs for all intended exploitations, especially for drinking water purposes. However, the water properties are highly dependent on the site, season, and sampling depth, which should be considered when conducting future investigations. To our knowledge, this study is the first to use metagenomic 16S rRNA amplicon sequencing to characterize the bacterial communities in the Bovan reservoir, providing deeper insight into the composition and diversity of such communities at the water intake locality. Due to the importance of water bodies such as lakes and reservoirs to the surrounding human population, the continuous and regular monitoring of microbial ecology and water quality throughout the year should be considered along with standardization of the employed methodology.
Supplementary Materials: The following are available online at https://www.mdpi.com/article/ 10.3390/w14030391/s1, Figure S1: Venn diagram showing the common and unique OTUs at the water intake locality of Bovan Lake at three different water depths (samples S5, S6, and S7); Figure S2: Relative abundance of the top ten bacterial classes at three water depths of the water intake locality of Bovan Lake in winter 2019; Figure S3: Relative abundance of the top ten bacterial orders at three water depths of the water intake locality of Bovan Lake in winter 2019; Figure S4: Relative abundance of the top ten bacterial families at three water depths of the water intake locality of Bovan Lake in winter 2019, Figure S5: Relative abundance of the top ten bacterial genera at three water depths of the water intake locality of Bovan Lake in winter 2019; Table S1: Sample information-sampling site, depth, sample names, and season of sampling; Table S2a: Physical and chemical parameters of Bovan Lake water in spring 2019 at all localities; Table S2b: Physical and chemical parameters of Bovan Lake water in summer 2019 at all localities; Table S2c: Physical and chemical parameters of Bovan Lake in autumn 2019 at all localities; Table S2d: Physical and chemical parameters of the water intake locality of Bovan Lake in winter 2019; Table S3: Heavy metal concentrations at the water intake locality of Bovan Lake in winter 2019; Table S4: Phytoplankton taxa detected in Bovan Lake in 2019 at different sampling seasons; Table S5: Chlorophyll a concentration, trophic level, and Shannon-Weaver diversity index of Bovan Lake; Table S6: The concentration of bacteria (CFU mL −1 ) at three localities of Bovan Lake, three water depths, and four sampling seasons in 2019; Table S7: Self-purifying ability of Bovan Lake water at three different localities and four sampling seasons.

Data Availability Statement:
The data supporting the findings of this study are available on request from the corresponding author.

Conflicts of Interest:
The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, or in the decision to publish the results.