Hypoxia in deep waters of moderately eutrophic marine lakes, Island of Mljet, eastern Adriatic Sea

Summary : In this study, we explored the impact of eutrophication and stratification on hypoxia in deep waters of moder- ately warm Croatian marine lakes. Although the Mljet Lakes (MLs) are predominantly oligotrophic, mesotrophic conditions are present at depths below 20 m in the Small Lake (SL) and below 30 m in the Big Lake (BL), along with higher apparent oxygen utilization (AOU). Hypoxia at depths ≥ 25 m in SL and and ≥ 40 m in BL was observed between October 2009 and January 2010, and in SL in summer (July and September 2010). Significant differences (p<0.05) in several physical, bio- logical and chemical parameters were detected between the lakes, while AOU, derived oxygen utilization rate (OUR) and organic carbon remineralization rate (OCRR) were not significantly different (p>0.05) between the lakes. An intense and persistent pycnocline throughout the year, comparatively high water temperature, extended water renewal time and summer phytoplankton bloom were identified as physical and biological parameters which might have significantly contributed to increased frequency of hypoxic events in a shallow SL. Significantly (p<0.05) higher ammonium concentration in SL, especially in its deep water, seems to be a long-term chemical feature related to the poor ventilation and higher sediment oxygen demand. At the current level of eutrophication and the present climate change trends, the MLs and similar systems may ex-perience more persistent and intense stratification, which could further prevent mixing between upper and deep waters, likely leading to increasing duration of hypoxia and its negative impacts on the biodiversity of benthic communities.


INTRODUCTION
Eutrophication, of either natural or anthropogenic origin, is closely connected to hypoxia in coastal seas (Nixon 1990).The cause of deficiency or absence of oxygen in aquatic systems can be of a physical origin.For example, a) sills and barriers affect circulation of water; b) excessive precipitation over evaporation intensifies the stratification that prevents mixing; and c) weak cooling of surface waters hampers water sinking and decreases water column mixing and aeration of deep waters.In addition, if a pronounced pycnocline develops in such a system, diffusion of oxygen through the density gradient can be further prevented.There are large portions of marine systems that are hypoxic, such as subsurface waters of the East China Sea (Chen et al. 2007, Wang et al. 2016) and deep waters of the Baltic Sea (Conley et al. 2009a, Rolff andElfwing 2015), the Black Sea (Konovalov et al. 1999) and the Gulf of Mexico (Justić et al. 1996(Justić et al. , 2005)).Layers of the open ocean can be hypoxic due to a significantly larger supply of organic matter than oxygen (Kamykowski and Zentara 1990), but it is hypoxia/anoxia in coastal seas (Diaz 2001) that causes greater concern regarding damage to the biodiversity of ecosystems (Ciglenečki et al. 2015), which many human activities depend on, such as fisheries and tourism.
Lakes behave like natural laboratories because of their comparatively prompt responses to environmental changes (Ciglenečki et al. 2015).Those which contain sea water originated a few thousand years ago by processes related to late postglaciation, i.e. isostasia and elevation of sea level.The marine lakes in Scandinavia were transformed into freshwater lakes due to isostasia, and hence some of them contain the trapped sea water at the bottom (Bøyum 1973).In comparison, the marine lakes of the eastern Adriatic coast are incorporated into the system of karst, containing sea water in the entire water column (Žic et al. 2012, Hrustić et al. 2013, Ciglenečki et al. 2015).The Mljet Lakes (MLs) were formed during the Mesozoic Era as karstic depressions filled with fresh water (Govorčin et al. 2001).The postglacial sea level rose around 4 to 5 thousand years ago and consequently the ingression of sea turned them into a marine environment (Bognar and Curić 1995).Unlike the well-oxygenated Jezero Mir (Lake Mir) (43.8864°N 15.1667°E), which is up to 10 m deep (Žic et al. 2012), and the meromictic lake Zmajevo oko (Dragon's Eye) (43.5309°N 15.9588°E), which is up to 15 m deep, permanently hypoxic below 13 m and occasionally completely anoxic (Ciglenečki et al. 2015), episodes of hypoxia are expected in deep waters of the MLs (Buljan and Špan 1976) within the system of coastal karst in Croatia.
The strait (depth ~0.8 m, width ~2 m) between SL (~241000 m 2 , ~3349000 m 3 ) and BL (~1450000 m 2 , ~58000000 m 3 ) and the strait (depth ~2.5 m, width ~10 m) between BL and Soline Channel (SC) (Fig. 1) hinder the exchange of deep lake waters with the adjacent coastal sea.BL is up to 49 m deep, with an average depth of 13.6 m, while SL is up to 29 m deep, with an average depth of 7.6 m (Cuculić et al. 2012).The advection from the open sea occurrs strictly in the surface layer, while circulation is dominantly driven by tides <25 cm and tidal phase shifts through the straits acting as filters at semidiurnal frequencies (Peharda and Vilibić 2008).The absence of upwelling, unlike from coastal waters (Batistić et al. 2013), constrains the mixing and supports the development of hypoxia in deep waters of the MLs.Some instabilities, i.e. internal waves in deep water of BL, were explained by severe cooling and sinking of the surface water to the bottom (Buljan and Špan 1976).Such an event has not been recorded in SL so far, but some water column instability can come from lateral and bottom subaquatic karstic springs (Wunsam et al. 1999).Cuculić et al. (2012) estimated that water renewal times (WRTs) in SL and BL during the mixing period (October-March) were 81 and 70 days and during the period of stratification (April-September) they were 44 and 36 days, respectively.WRT for the stratification period refers only to the upper waters above the depth of 10 m in SL and 17 m in BL.When SL and BL were considered as one system exchanging the water with the open sea, WRTs of the MLs during the mixing and stratifed periods were 82 and 43 days, respectively.These values, being closer to WRT in SL, indicated that the exchange between BL and SL is the limiting factor for the total exchange of water in the system (Cuculić et al. 2012).
In this study, we intended to investigate the link between natural eutrophication, water column stratification and episodes of hypoxia in deep waters of moderately warm marine lakes.The arrangement of data concerning the mixing and stratification periods and water column division was adjusted to be comparable to the recent report on WRT in the MLs (Cuculić et al. 2012), while offering an extension to parts of the research that revealed an anoxic event in deep water of BL more than a decade ago (Benović et al. 2000).Anoxia in the water of the MLs was not detected after August 1996 and vertical homogenization has been noted to occur regularly since 2003 (Vilibić et al. 2010).Having apparent oxygen utilization (AOU) as an important parameter of our study, including the report on the WRT (Cuculić et al. 2012), we were able to determine oxygen utilization rate (OUR) and organic carbon remineralization rate (OCRR) (Feely et al. 2004) in the MLs for the first time, and to provide information and conclusions that may be valuable for future studies of similar environments worldwide.

Sampling
The sampling locations and dates are presented in Table 1 and Figure 1.The stations in SL (at its deepest point), BL (Vrbovačka) and Gonoturska Bay (G) are 29, 46 and ~60 m deep, respectively.The sampling was performed between 9 am and 3 pm, repeating the sequence SL, BL and G.

Thermohaline properties, precipitation, compensation depth
Salinity (S), temperature (t) and density (σ T ) were determined by an SBE 19plus CTD probe at a depth resolution of 1 m to 27, 45 and 50 m in SL, BL and G respectively.All sensors were regularly calibrated by Sea-Bird Electronics Inc.The difference in t, S and σ T for every metre of the water column depth was calculated and named ΔS (m -1 ), Δt (°C m -1 ) and Δσ T (kg m -4 ).The depth visibility of the white Secchi disc (Sec, diameter 30 cm) was multiplied by 3 (Berman et al. 1985) to determine the depth of the photic zone (Pz).Data on precipitation at the nearest meteorological station to the MLs in Goveđari (Fig. 1) was ensured by Državni hidrometeorološki zavod (DHMZ, Croatian Meteorological and Hydrological Service, meteo.hr).

Dissolved oxygen concentration, apparent oxygen utilization (AOU), oxygen utilization rate (OUR), organic carbon remineralization rate (OCRR) and nutrients
The samples were collected by 5-L Niskin bottles at the surface (0-1), 5, 10, 15, 20, 25, 30, 35, 40 and 43 m in BL and at the same depths up to 25 m in SL.The samples in G were taken at the same depths up to 20 m plus at 30, 40 and 45 m.Dissolved O 2 (DO) was precipitated in Winkler bottles and determined by iodometric titration (Grasshoff et al. 1983).The results (ml O 2 l -1 ) were converted to mg O 2 l -1 using the conversion factor 1.42903 (Owens andMillard 1985, Garcia andGordon 1992, Sea-Bird Electronics Inc.).The threshold for hypoxia was set at <2.8 mg O 2 l -1 , as proposed by Diaz and Rosenberg (1995), which is equal to <1.96 ml O 2 l -1 and <87.5 μmol O 2 l -1 .Anoxia stands for the environment without DO.AOU was calculated as the difference between DO concentration at saturation of 100% (O 2 ') at a given S, t and depth (Weiss 1970) and the associated measured DO concentration, expressed as μmol O 2 kg -1 (Feely et al. 2004).The unit for DO i.e.O 2 (mg l -1 ) was adjusted to the criterion for hypoxia (Diaz and Rosenberg 1995).We also calculated OUR (μmol kg -1 year -1 ) as AOU/water age (Feely et al. 2004) using corresponding individual WRTs on the basis of "e-folding time" by Cuculić et al. (2012).OCRR (μmol kg -1 year -1 ) was calculated as OUR×R C:O (Feely et al. 2004) , where R C:O represents the molar "Redfield ratio" of C:O.We employed 0.69 for R C:O as it represents an R C:O of 117:170 (Anderson and Sarmiento 1994) and the average range of R C:O (0.53-0.85) reported by Li and Peng (2002).The calculation of OCRR represents only the organic carbon mineralized when oxygen is the acceptor of the electron transport system.If another acceptor such as nitrate is also used, there will be an extra OCRR not accounted for by this calculation.
Samples (50 ml) for analysis of ammonium (NH 4 ) were stabilized by addition of 2 ml of phenol (1 mol l -1 ) dissolved in 95% vol/vol ethanol (Ivančić and Degobbis 1984) and stored in the dark at 4°C.The samples (500 ml) for other nutrients were stored at -22°C.The concentrations of nitrate (NO 3 ), nitrite (NO 2 ), NH 4 , reactive silicates (SiO 4 ), orthophosphate (PO 4 ), total phosphorus (TP) and total nitrogen (TN) were determined according to Strickland and Parsons (1972) without filtration.Difference between TP and PO 4 was designated as "other P" i.e.P oth .Difference between TN and dissolved inorganic nitrogen (DIN) (NO 3 +NO 2 +NH 4 ) was designated as "other N" i.e.N oth .The persulfate oxidation of unfiltered samples in autoclave (127°C, 25 min) enabled us to measure TP (Menzel and Corwin 1965) and TN (D'Elia and Steudler 1977), which meant that P oth and N oth consisted of dissolved organic P (DOP) and N (DON), and non-dissolved (colloidal and particulate) organic and non-dissolved inorganic P and N species, respectively.

Chlorophyll a and trophic index
DIN (μg l -1 ), TP (μg l -1 ), chlorophyll a (Chl a) (μg l -1 ) detected by the CTD probe and O 2 absolute deviation (%) from saturation (O 2 ') were used to calculate TRIX: TRIX is the trophic index first developed for the Italian northwestern Adriatic coastal waters, and it has a range of values from 0 to 10, where 0-4, 4-5, 5-6 and 6-10 indicate the high, good, moderate and poor quality of the water, respectively (Vollenweider et al. 1998).The detector of fluorescence was calibrated by the Chl a measurements from the culture of a diatom Thalassiosira weissflogii at the Sea-Bird Electronics Inc. laboratories.The same concept was employed by Viličić et al. (2009).The data on Chl a in October and November 2009 were lost, so TRIX could not be calculated for October and November 2009.

Picoplankton
The samples (50 ml) were preserved in neutralized solution of formaldehyde (~4% vol/vol) on board and stored in the dark at 4°C until analysis.Abundance and sizes of identified picoplankton groups were determined by epifluorescent microscopy (Zeiss Jenalumar, magnification 1500).Heterotrophic bacteria (HBa) were detected by direct epifluorescent microscopy (Hobbie et al. 1977).Heterotrophic picoflagellates (HPF), autotrophic picoflagellates (APF) and picocyanobacteria (Cyano) were detected using proflavine (Haas 1982).Biomass of HBa was calculated by the conversion factor 20 fg C cell -1 (Lee and Fuhrman 1987).Biomass of Cyano was calculated by the conversion factor 250 fg C cell -1 (Kana and Gilbert 1987).The volume (μm 3 ) of picoplankton cell categories was calculated using the equations according to morphological classes: ellipsoid (V=πLw 2 /6) and sphere (V=4πr 3 /3), where L is the length (μm), w is the width (μm) and r is the radius (μm) of the cell.Biomass of HPF was calculated as B=nVFC, where n is abundance, V is average cell volume (μm 3 ) and FC is the conversion factor 0.22 pg C μm -3 (Borsheim and Bratbak 1987).The equation 0.433(μm 3 ) 0.863 =pg C cell -1 (Verity et al. 1992) was employed as a conversion factor for APF by taking into account average cell volumes (μm 3 ).

Statistics
Statistica for Windows (Statsoft Inc.) was used to determine the distribution of variables by Kolmogorov-Smirnov and Lilliefors criteria and to calculate Spearman Rank (non-normal distribution) and Pearson's correlation coefficients (r) at p<0.05 after log(x+1) transformation of the data.One-way analysis of variance and Student-Newman-Keuls tests were applied to determine the significant differences (p<0.05) between SL, BL and G, representing the coastal sea.T tests were employed to determine significant differences (p<0.05) between upper and deep waters at each location, using the depths 10 and 17 m (Cuculić et al. 2012) as distributors for SL and BL, respectively, while 20 m was employed as a distributor for G. Furthermore, t tests were employed to detect significant differences (p<0.05) between periods of mixing (October-March) and stratification (April-September) at each location.
The thermocline in SL appeared in March, but a deeper and more pronounced thermocline developed around the depth of 10 m between April and June (Figs 2A and 3A).Further, the thermocline advanced deeper during the summer, being the strongest (-2.716°C m -1 ) between 10 and 11 m (24.685-21.969°C,respectively) in August and the deepest in September (14-15 m, -2.22°C m -1 ).These results supported the division of the water column of SL at a depth of 10 m, as proposed by Cuculić et al. (2012).The thermocline in BL and G developed with a small time lag after that in SL (Fig. 2).The most intense Δt in BL stretched from 13-16 m in April to 22-25 m in September (Figs 2B and 3B).Δt max (-2.107°C m -1 ) in BL developed in the 23-24 m layer (18.476-16.369°C, respectively) in September (Figs 2B and 3B).Δt max did not ideally match the division of the water column in BL at a depth of 17 m, but the most intense Δσ T were detected in the 16-19 m layer, except in September (23-24 m, Δσ T 0.29 kg m -4 ) (Figs 2B and 3B), which supported the division at a depth of 17 m, as proposed by Cuculić et al. (2012).Δt max in G (-0.843°C m -1 ) developed in the 3-4 m layer (21.210-20.367°C, respectively) in August (Figs 2C and 3C).The water column of G had a notable Δσ T in May-September, being occasionally pronounced both above and below 20 m (Figs 2C and 3C).
Deep water of BL was persistently colder than coastal sea at comparable depths (Fig. 2B, C).A similar finding was reported by Benović et al. (2000).The only water that is denser than the near-bottom water of BL in October (σ T 28.81 kg m -3 , depth ≥43 m) was situated below the depth of 44 m in G (σ Tmax 28.83 kg m -3 , depth of 49 m) in March (Fig. 2C).Δσ Tmax in G (0.34 kg m -4 ) was found in May (16-17 m) (Fig. 2C).Δσ Tmax in BL (0.32 kg m -4 , 18-19 m, July) (Fig. 2B) was similar to that in G, while Δσ Tmax in SL (1.44 kg m -4 , 2-3 m, December) exceeded the values of the highest Δσ T in BL and G more than fourfold.In addition, all samplings in SL had greater Δσ Tmax (the lowest monthly Δσ Tmax 0.36 kg m -4 , 14-15 m, September) than overall Δσ Tmax in G and BL.It was noted that Δσ Tmax in SL in July was almost double (0.62 kg m -4 , 8-9 m) the overall values, 10 m deeper Δσ Tmax in BL in July.The pycnocline in the MLs plunged deeper towards the end of summer (Fig. 2A, B).
Because of the extraordinary low S in the surface water of SL in December 2009, data on precipitation before and during samplings were analysed (Fig. 4).Intense precipitation from 14 to 17 December (average, ~30 mm m -2 day -1 ) significantly impacted on the lowest S in the surface water of SL on 18 December.By applying the partial differential equation while keeping the S and t constant (http://www.csgnetwork.com/h2odenscalc.html) in two separate tests between the depths of 2 and 3 m, we concluded that ΔS (2.411 m -1 ) and Δt (2.319°C m -1 ) contributed circa 130% and -30%, respectively, to the detected Δσ Tmax in the 2-3 m layer in SL in December.

The depth of the photic zone
Mutually different Secchi disc (Sec) values decreased from the outer station towards SL (Table 2).Minimum Sec in SL (June-July) indicated a compara-

DO, AOU, OUR and OCRR
The maximum (193.90 μmol kg -1 ) and minimum (-71.89μmol kg -1 ) AOU were detected at 43 m in BL (October) and 10 m in SL (August), respectively (Fig. 5A, B).There were only two samplings in SL (April, August), when O 2 at 25 m exceeded 5 mg l -1 (5.16, 5.63 mg l -1 ) (Fig. 5A) accompanied by AOU 93.42 and 75.60 μmol kg -1 (Fig. 5A), respectively.In February, March, May and June, O 2 ranged from 3.13 to 4.16 mg l -1 (Fig. 5A), with AOU from 124.00 to 155.22 μmol kg -1 at 25 m in SL (Fig. 5A).The near-bottom waters of SL at 25 m were hypoxic from October to January and in July and September (Fig. 5A).In comparison, the near-bottom waters of BL (depths of 40 and 43 m) were hypoxic only in the period October-January, having a range of 2.06-2.76mg O 2 l -1 (Fig. 5B).In spite of having higher DO than samples from 25 m depth in SL (October-January, AOU 182.62-189.83μmol kg -1 , Fig. 5A), the deepest samples from BL (43 m) had a similar or higher AOU (184.04-193.90μmol kg -1 ) (Fig. 5B).We found that AOU in the water column of G was mostly negative, being negligibly positive only occasionally (three samples in December and two samples in February) (Fig. 5C), whereas AOU was positive in the whole water columns of SL and BL in December (Fig. 5A, B).Significantly greater consumption than production of O 2 was detected below 20 m in SL and 30 m in BL throughout the whole study period (Fig. 5A, B).The minimum of O 2 (1.70 mg l -1 ) was detected at 25 m in SL in December (Fig. 5A), when Δσ Tmax occurred in the subsurface water (2-3 m).Using individual WRTs (Cuculić et al. 2012) as the water age in SL and BL, we calculated OUR and OCRR (Table S3, S6, S9, S10).

Nutrients and estimate of phytoplankton growth
The negative correlations (p<0.05) between O 2 and inorganic nutrients in the MLs (Table 3) were most likely associated with regeneration after decomposition of the sinking organic matter excreted by Aurelia sp.(Alvarez-Colombo et al. 2009, Tinta et al. 2016) and other zooplankton and phytoplankton.The positive correlations between S and NO 3 in the MLs are in accordance with oxidation (i.e.remineralization) of organic matter during its sinking towards deep waters of greater S (Fig. 2A, B).SL had significantly (p<0.05)higher NH 4 than BL and G, especially in deep waters (Fig. 5A, B, C), while G had significantly (p<0.05)lower SiO 4 and SiO 4 /DIN than the MLs (Table S1, S2).In April, SiO 4 and NO 3 were notably depleted in the water columns and upper waters of the MLs, respectively (Fig. 6A, B).SiO 4 was significantly (p<0.05)correlated with O 2 , Chl a and t in BL.The nutricline (NO 3 and SiO 4 ) in the MLs (Fig. 6A, B), except for SiO 4 in April and NO 3 in SL in the stratification period, was accompanied by the rise of AOU under ~20 m in SL and 30 m in BL (Fig. 5A, B).Concentration of PO 4 had the maxima of 0.39, 0.50 and 0.71 μmol l -1 (not presented) in G, BL and SL at 0-1 m in February, at 0-1 m and at 5 m in January, respectively.Maximum NO 2 concentrations of 0.28, 0.24 and 3.17 μmol l -1 (not presented) in G, BL and SL were detected at 45 m in March, 30 m in March and 25 m in October, respectively.Maximum NH 4 , NO 3 , and SiO 4 concentrations (4.10, 6.55 and 52.50 μmol l -1 ) were detected at 25 m in SL in November (Fig. 5A) and 43 m in BL in October and in December (Fig. 6B), respectively.

Chl a and TRIX
Chl a had the highest (p<0.05)concentration in SL (Table S1, S2).The trend of rising Chl a concentration from G towards the more enclosed waters in SL was observed.Chl a was positively correlated with depth at deeper locations (BL, G) (Table 3).Chl a max (1.76 μg l -1 ) was detected in SL at 27 m in March (Fig. 6A).Chl a max in BL (1.09 μg l -1 ) at 40 m in July (Fig. 6B) might have been connected to the rich pool of nutrients (Fig. 6B) and photoadaptation of phytoplankton (Geider et al. 1997) in the colder (11.78°C) deep water  of BL (Fig. 2B).Seasonal development of Chl a was similar for the MLs.While SL had relatively high Chl a in the whole water column in April and July, BL had a notable concentration of Chl a (Fig. 6A, B) only below the thermocline in the same months.Chl a in BL had also positive correlations (p<0.05) with SiO 4 , P oth , NO 3 and NO 2 and a negative correlation (p<0.05) with t.A deep chlorophyll maximum (Zavatarelli et al. 1998, Viličić et al. 2008) of 0.67 μg l -1 was detected at 49 m in G in July (Fig. 6C).The spring (April) and the summer (July) blooms were evident because of the higher and more homogenous Chl a concentrations (1.04±0.05 and 0.99±0.03μg l -1 , respectively) throughout the water column of SL.The bloom in July was accompanied by weaker hypoxia than in September, although the largest difference in σ T between upper and deep waters was detected in SL in July (Fig. 2A).Hypoxia was not detected in the less eutrophicated deep water of SL in August (Figs 5A and 6A).The mesotrophic conditions (TRIX>4) dominantly developed below 20 m in SL and below 30 m in BL (Fig. 6A, B), while the investigated locations were prevailingly oligotrophic (TRIX<4).TRIX max (4.83) was found at 35 m in BL in September (Fig. 6B).

Differences between the mixing and stratification period
Although October-March was considered as the mixing period, the upper water of SL up to 7 m was highly stratified except in October-November because of the pronounced halocline (Fig. 2A).In the mixing period, SL had 5 and 10-fold higher ΔS than BL and G, respectively (Supplementary material Table S3).Unlike SL, where ΔS and Δt exchanged the leading role in defining Δσ T during the mixing and stratification period, respectively, Δt dominantly impacted on Δσ T in BL and G in April-September (Tables S3, S4).The mixing period for SL corresponds to the fact that, except in October and November, the upper waters of SL (≤7 m) were notably stratified on the sampling days in the mixing period (October-March).There was no significant difference in parameters related to O 2 between the mixing and stratification period in SL (Tables S3, S4, S5).AOU was significantly higher (p<0.05) in BL and G in the mixing period than in the stratification period (Tables S3, S4, S5).
HPF had significantly larger biomass (p<0.05) in SL than in G (Tables S1, S2).There was no significant difference in picoplankton between the specified periods in SL (Table S3, S4, S5).APF and HPF reacted positively to t increase in SL (Table 3), but other impacts (probably grazing was the most important) did not support stronger development of the picoplankton in the warmer period (Tables S3, S4, S5).An increase in t was significantly (p<0.05)correlated with HPF in BL (Table 3), resulting in larger HPF in BL in the stratification period (Tables S3, S4, S5).Also, an increase in t resulted in significantly (p<0.05)higher HBa and lower Cyano in G in the stratification period.Cyano had almost unchanged averages in the MLs between the specified periods (Tables S3, S4).
Potentially limiting nutrients, according to Redfield et al. (1963), had a uniform succession: PO 4 , DIN, SiO 4 at all locations in the contrasting periods (Table S3, S4).DIN at G in the stratification period was within the range (0.1-0.6 μmol l -1 ) of half saturation constant (K m ) for NO 3 and NH 4 uptake by the phytoplankton community in the oligotrophic environment (MacIsaac and Dugdale 1969).Taking into account the compara-tively low PO 4 (Table S4), it seems that phytoplankton growth at G during stratification was potentially N and P co-limited.PO 4 concentrations at all stations were significantly (p<0.05)higher in the mixing period than in the stratification period (Table S5).
In the stratification period, Chl a was significantly (p<0.05)higher in SL, not significantly different in BL and significantly (p<0.05)lower in G than in the mixing period.The contrast was evident between notable Chl a and O 2 concentrations in the deep water of the MLs in early spring and poor Chl a concentration in the same, but hypoxic water in late autumn (Figs 5A, B and 6A, B).TRIX was significantly (p<0.05)lower at all locations in the stratification period than in the mixing period.HPF in BL, Cyano and HBa in G had significantly different (p<0.05)values between the mixing and stratification periods (Tables S3, S4, S5).

Upper and deep waters in different periods
There was no significant difference in parameters related to S between specified layers in each lake within the stratification period (Supplementary material Tables S6, S7, S8).All investigated locations above the thermocline had similar O 2 and, according to AOU, significantly greater production than consumption of O 2 (Table S6).SL had persistently higher averages of O 2 sat.than BL, even at lower O 2 (Table S4, S6, S7), due to significantly higher t than in BL (Table S2).All locations, above the thermocline, had low PO 4 accompanied by NO 3 and NH 4 means within the range (0.1-0.6 μmol l -1 ) of K m of their uptake by the phytoplankton community in the oligotrophic environment (MacIsaac and Dugdale 1969), indicating a potential N and P co-limitation (Table S6).The most pronounced oligotrophy was detected in G above the thermocline.The significant difference in SiO 4 between specified layers during stratification in the MLs (Tables S6, S7, S8) might have been caused by the sinking detritus, which was notably composed of diatoms.PO 4 was significantly (p<0.05)higher only in deep water of BL during the stratification period (Tables S6, S7, S8).
Chl a had three and six-fold higher concentrations below the thermocline at deeper locations, BL and G, respectively, while Chl a rise in deep water of SL was not dramatic (Table S7).NO 3 was three-and nine-fold higher below the thermocline in SL and BL, respectively (Table S7).All locations had low DIN above the thermocline, while G had low DIN even below the thermocline.NO 2 , especially in SL, and TRIX were significantly (p<0.05)higher below the thermocline and accompanied by higher AOU (Tables S6, S7, S8).DIN/PO 4 below the thermocline in the MLs was notably higher than Redfield's 16.Two specified layers at G did not show significant differences in the ratios of analysed nutrients (Table S8).APF and Cyano had significantly (p<0.05)higher values above the thermocline than in deep waters in SL (Tables S6, S7, S8).Since division on mixing and stratification period did not match the reality in SL, at least for our samplings, it was preferred to determine averages and standard deviations of the study parameters in the mixing period (Tables S9, S10) associated with the same layers in the stratification period, including the test of the potentially significant (p<0.05)differences between the layers (Table S11) and within the same layers in different periods (Table S12).SiO 4 in deep water of G was an order of magnitude lower than in the MLs (Table S10), perhaps as a consequence of the decomposition of detritus in the stagnant deep waters of the more productive MLs.The upper waters at all locations had DIN/PO 4 close to Redfield's 16 (Table S9), but only deep water of G retained similar DIN/PO 4 (Table S10).Although DIN in deep waters of the MLs did not differ notably, NH 4 and NO 2 were more important constituents of DIN in deep water of SL, whereas NO 3 was more important in deep waters of BL and G (Table S10).There were no significant changes in picoplankton constituents at all stations among upper and deep waters during the mixing period (Table S11).The largest number of significant changes between the layers occurred in SL (Table S11).This is in line with stratification of SL even in the mixing period.

Same layers in different periods
SL did not have significantly different Δσ T between the mixing and stratification periods (Supplementary material Table S12), indicating a notable pycnocline throughout the year.Deep water of SL did not show significant differences in parameters related to O 2 in different periods, while in BL only O 2 was significantly (p<0.05)lower in deep water in the mixing period (Tables S7, S10, S12).The stratification period, at least in its upper waters, seemed to be more productive according to differences (p<0.05) in O 2 sat.and parameters derived from O 2 sat.(Tables S6, S9, S12).Deep waters of the MLs were moderately eutrophic in both periods (Tables S7, S10, S12).Upper water of SL did not differ significantly in PO 4 between the periods (Table S12), despite an extremely high concentration of PO 4 at 5 m in January (not presented) that lifted the average of PO 4 (Table S9) notably above the average in that layer in the stratification period (Table S6).Unlike in the MLs, the Chl a concentration in deep water of G was not significantly different between the periods (Table S12).Also, there were no significant differenc-es in parameters related to picoplankton and nutrients other than inorganic ones in water layers between the specified periods in SL (Table S12).Significant differences (p<0.05) in constituents of picoplankton (HPF in BL, Cyano and HBa in G) and N oth (BL and G) were detected in the upper waters, while concentration of analysed organic matter differed significantly (p<0.05) in deep waters of BL (N oth and P oth ) and G (N oth ) between the specified periods (Table S12).

DISCUSSION
The lack of significant differences in any of the parameters related to O 2 between the mixing and stratification periods in SL was most likely related to persistent stratification (no significant difference in Δσ T between the mixing and stratification period), which constrained the mixing between its upper and deep water.In addition, we found the lowest numbers of parameters with significant (p<0.05)changes between the periods in upper and deep water of SL (Supplementary material Table S12), which are most likely related to the same reason.In comparison, the significant (p<0.05)differences in the upper water of BL (HPF) and G (Cyano and HBa) between the periods (Table S12) are most likely a consequence of changes in the food web, which is a topic beyond the scope of this study.Of note, N oth in BL and G was significantly (p<0.05)higher in the mixing period, but we cannot offer a plausible explanation for these observations.
The thermocline was crucial for stratification at all explored stations in the stratification period, while a relatively shallow (up to 7 m) and persistent halocline was crucial for stratification of SL in the mixing period, to the point that difference in t in the 2-3 m layer in December, which impacted negatively on the pycnocline by ~30%, was compensated by an 130% positive influence of ΔS on Δσ T in that layer.High range (>10°C) between low and high t in the upper water of SL during less and more productive periods (Table S9, S6), respectively, played a role in balancing between different levels of O 2 production and solubility.These include an exponentially positive relation between the rate of metabolic processes and O 2 requirements of marine organisms at increasing t (Brown et al. 2004), along with decreasing O 2 solubility in seawater with the rise of t (Garcia and Gordon 1992).Ciglenečki et al. (2015) detected an increase in t of 1.1°C in the surface layer (0-2 m) and 0.64°C at 12 m in Dragon's Eye for 2001-2012 in comparison with 1996-2000.Therefore, we speculate that the presented trend will have greater influence on SL than on BL concerning the diminished solubility of O 2 with rising t.Decreased O 2 solubility in the upper water of SL in both periods due to higher t is mitigated by lower S in comparison with BL (Table S6, S9).An important difference between the MLs was the notably warmer (2.12 and 2°C) and somewhat saltier (0.32 and 0.38; negative effect on solubility of O 2 ) deep water of SL during periods of stratification and mixing (Tables S7, S10), respectively.Although phytoplankton blooms produced organic matter and DO throughout the water column of SL in April and July, saltier (on average 0.198 and 0.167) and warmer (on average 2.242°C and 2.933°C) deep water of SL in April and July, respectively, caused reduced solubility of O 2 in comparison with deep water of BL, while having higher O 2 requirements due to organisms at higher t at the same time (Brown et al. 2004).It was predicted by Vaquer-Sunyer and Duarte (2011) that survival times of benthic macrofauna would decrease by a mean of 35.6% under hypoxia with the average 25.5% increase of DO concentration threshold for high mortality in a 4°C warmer ocean, representing the average increase in t expected at the end of the 21 st century under moderate CO 2 emission scenarios.Along these lines, Lojen et al. (2010) have detected 4-5°C higher t (12°C in May and 13°C in October) in the pore water of SL than in the pore water of BL (8°C).Therefore, the difference of ≥2°C between deep waters of SL and BL should be monitored in the future.
Detailed discussion on the impact of p(CO 2 ) and potential synergistic effects of hypoxia and increased p(CO 2 ) on benthic invertebrates and fish can be found in Steckbauer et al. (2015), where measurement of calcification rates are stressed as an important additional tool to provide more prediction on the future outcome of marine invertebrates.Influence of hypoxia on benthic communities of the MLs (Kružić 2002, Mikac 2003) has not been investigated.BL is a home for the scleractinian coral Cladocora caespitosa, which forms a reef (650 m 2 ) between depths of 6 and 18 m near the Vrbovačka sampling station (Hrustić et al. 2013).While experimenting on another scleractinian coral reef, Reynaud et al. (2003) found out that calficification was decreased by 50% when t and p(CO 2 ) were both elevated, with no significant change in response to increased p(CO 2 ) under normal t.Their findings indicated a need to re-evaluate the projected decrease of marine calcification by the year 2100.Taking into account their findings, McNeil et al. (2004) reported that the annual average coral reef calcification rate will increase with future ocean warming and eventually exceed pre-industrial rates by about 35% by the year 2100.
Although notably stratified throughout the year, anoxia did not develop in deep waters of SL.Overturn of the water column under influence of strong cooling of the surface water has not been documented in SL so far, but such events were documented in the meromictic lake Dragon's Eye a few hundred kilometres farther north, causing the mass mortality of benthic and planktonic organisms in the early autumns of 1997 (Kršinić et al. 2000, Barić et al. 2003) and 2011 (Ciglenečki et al. 2015).These events revealed that the copepod Acartia italica has been well adapted to the extreme conditions of anoxia, while some benthic Pennatae diatoms also survived anoxia in Dragon's Eye in 1997 (Barić et al. 2003).In addition, the research in the northern Adriatic Sea at 24 m conducted by inducing anoxia revealed that meiobenthic harpacticoid copepods from the family Cletodidae were well adapted to anoxia (Grego et al. 2014).
The absence of upwelling in the MLs indicates that diffusion is the most important mechanism of transport of dissolved matter from deep to upper waters.If up-welling were a normal feature of the MLs, PO 4 released from iron-bound P in oxygen-poor sediments would reach the upper water very fast, alleviating potential P limitation according to DIN:PO 4 >16 (Table S1) and acting as a positive feedback to increase hypoxia by stimulating phytoplankton production (Conley et al. 2009b).Deep water of SL was >100% and ~150% richer in total Fe than deep water of BL and G (Buljan and Špan 1976), respectively, while iron-phosphorus still remains unexplored in sediments of this system.
The experiment of fertilization of the MLs was started via addition of 21.5 tonnes of "superphosphate" in eight instalments in BL in 1954 (Buljan and Špan 1976).In total, an average of 1.15 μmol PO 4 l -1 was introduced in SC and the MLs.The system reacted by multiplying the primary productivity more than sixfold in comparison with 1953 (Buljan and Špan 1976).Hydrogen sulfide (H 2 S) was almost persistently present at 25 m in SL between 1951 and 1956 (Buljan and Špan 1976).It was occasionally present at notable concentrations at 20 m (e.g.4.84 mg l -1 ; overall maximum 4.95 mg l -1 was detected in SL at 25 m on the same day, 11 December 1952).The fertilization experiment caused formation of H 2 S (2.93 mg l -1 ) at 13 m at Pospile Donje in BL on 20 July 1954 (Buljan and Špan 1976).These authors noted "pink and red water" at 20 m in SL in May and June 1953, which was most likely connected to the presence of anoxygenic phototrophic purple sulfur bacteria, commonly found at the chemocline and in the hypolimnion (Pjevac et al. 2015).Consequently, the MLs were recognized as a highly vulnerable ecosystem and became the protected area of the Mljet National Park (NP) in 1960.It has been shown that the widening and deepening of SC and its entrance to BL in 1960 mitigated anoxic conditions in the near-bottom water of SL, according to historical profiles of the redox-sensitive elements in sediment cores of the MLs (Sondi et al. 2017).WRT in the MLs before 1960 was most likely longer than the WRT estimated recently by Cuculić et al. (2012).Buljan and Špan (1976) did not detect H 2 S in deep water of SL in 1961, which is in line with the statement of Sondi et al. (2017), indicating poor ventilation as one of the most important parameters linked to estimated persistent anoxia in the near-bottom water of SL before 1960.
The negative correlation (p<0.05) between PO 4 and O 2 in the water columns of the MLs is common in hypoxic waters, e.g.bottom waters of the Baltic Sea (Conley et al. 2002).The periods (Table S5) rather than water layers (Tables S8, S11), except for BL in the period of stratification, differed significantly (p<0.05) in concentration of PO 4 at all stations, having higher values during the mixing period (Tables S3, S4, S5).Because of the relatively high DIN/PO 4 at high productivity in a P-enriched environment, Buljan and Špan (1976) assumed N 2 fixation in the MLs.A similar hypothesis, derived from the excess of NO 3 in SL in early spring, was proposed by Hrustić et al. (2013).N 2 fixation can increase eutrophication, especially if coupled with a flux of PO 4 from deep water or sediments (Jensen et al. 1995) to the photic zone (Conley et al. 2009a).In oligotrophic P-starved regions, N 2 fixers can increase their competitiveness by using excessive N and energy to produce exoenzymes involved in utilization of DOP (Landolfi et al. 2015).The problem of eutrophication in the Baltic Sea was a good reason for the development of a new bioassay to test the distance between N and P limitation in mesotrophic surface waters (Hrustić et al. 2017).The integrated (depth of 0-10 m) samples from the ocean acidification experiment in the Tvärminne field station (Baltic Sea) in August 2012, without a bloom of diazotrophic cyanobacteria, showed the surplus of P, P + =0.30±0.10μmol l -1 , but the effect of ocean acidification was not detected by this method (Hrustić et al. 2017).It would be interesting to use this new bioassay in the future to explore the distance between N and P limitation in the protected mesotrophic MLs, where the food web is controlled by Aurelia sp., and the particulate feeding of Cladocora caespitosa on phytoplankton from the surface currents in BL (Hrustić et al. 2013) may supply missing nutrients in the coral reef (Wyatt et al. 2010).
Anoxia below 39 m in BL in the period 26-28 August 1996 (Benović et al. 2000) did not last for long and wintertime cooling homogenized the entire water column.Quiescent or slowly moving speciemens within a swarm of Aurelia sp. were noted just above the anoxic layer.Aurelia sp. were repeatedly observed in swarms with very large numbers of individuals in BL except in September 1997.The quiescent individuals above the layer of anoxia were most likely affected by hypoxia since decreasing O 2 sat.from the maximum of 101% at 18.5 m to depleted O 2 at 39 m (Benović et al. 2000) probably caused death to many individuals of Aurelia sp., which might have been a significant source of NH 4 in BL (Tinta et al. 2010(Tinta et al. , 2016)).However, deep water of SL had significantly (p<0.01)higher NH 4 than deep water of BL, which is, together with the gradient of NH 4 concentration below the depth of 20 m in SL (Benović et al. 2000), in line with our results (Fig. 5A, B; Tables S7, S10).
These findings point out other reasons why NH 4 in deep water of SL is significantly higher than NH 4 in deep water of BL.Higher t and Chl a are positively related to phytoplankton growth (Agawin et al. 2000).Increased NO 2 in the photic zone might be notably related to excretion during phytoplankton growth (Al-Qutob et al. 2002).This fits with the detection of maximum NO 2 concentrations in the photic zone of BL and G in March, but is less likely to be relevant for the maximum value of NO 2 in deep water of SL in October.According to the highest yearly means of t and Chl a (p<0.05)(Tables S1, S2), a comparatively intense phytoplankton growth was estimated for SL.However, the difference in NH 4 between the MLs was markedly higher during lower production, i.e. in the mixing period when the difference in Chl a was lower (Tables S3, S4).Similar findings can be seen throughout the layers in different periods, epecially in upper waters during the mixing period when BL had higher Chl a but lower NH 4 than SL (Table S9).This leads us to conclude that phytoplankton is probably not responsible for significantly higher NH 4 in SL.
Because of its smaller size, SL receives more terrestrial organic and mineral detritus relative to the amount of authigenic sediment than BL (Lojen et al. 2010).Particulate organic matter (POM) undergoes remineralization through the water column.A negligible portion of POM (~1%) escapes remineralization, being passively buried in oceanic sediments (Delaney 1998), but the percentage of POM received by sediments in the shallower SL is expected to be notably higher.Sondi et al. (2017) determined oxic-anoxic boundaries in the sediment cores of SL, BL and G at depths of 2, 5 and 15 cm, respectively, followed by the profile of dissolved iron.These results led us to conclude that surface sediments of SL have significantly higher oxygen demand than sediments of BL and G, producing notably more NH 4 from the increased amount of POM.That would be an explanation for the negative correlation (p<0.05) between NH 4 and O 2 in SL (Table 3).Wang et al. (2017) reported that re-aeration and sediment oxygen demand were more important than photosynthesis and water column respiration as counterbalanced coupled processes affecting the near-bottom DO concentration in the shallow estuary.Re-aeration of deep waters in the MLs is not expected in summer (Vilibić et al. 2010).Intermittent summer (July and September) hypoxia in the warmer deep water of SL can be attributed to presumably higher organic matter supply to the surface sediments, which have an increased oxygen demand.The plunge of the intense pycnocline might have contributed to hypoxia in the deep water of SL in summer.
Positive AOU in the whole water columns of SL and BL in December (not presented), as well as higher AOU at all investigated locations in autumn-winter, were most likely related to the lower photosynthetic activity.The critical depths of dominantly higher consumption than production of O 2 , i.e.where AOU>0 μmol kg -1 , involving a notable regeneration of nutrients by organic carbon remineralization, were situated ~10 m below the water columns' division points of 10 m for SL and 17 m for BL (Fig. 5A, B).There were no significant differences (p>0.05) between MLs concerning O 2 , AOU, OUR and OCRR.The larger volume of water in BL was characterized by a net utilization of O 2 (AOU>0 μmol kg -1 ) (Fig. 5 A, B) at a higher average OUR than in SL.A larger amount of organic carbon was remineralized in deep water of BL in the mixing period, due to a higher OCRR and a larger volume of water with positive values of OCRR than in deep water of SL in the same period (results not presented, derived from AOU).This finding is in line with a presumably higher % of POM that escapes remineralization during sinking in SL.Since we did not have data on WRT for deep waters of the MLs in the stratification period or any data on WRT in G, we performed statistical analyses for OUR and OCRR to a lesser extent.
According to the vertical profiles of SiO 4 and NO 3 in the MLs in April (Fig. 6A, B), we assume that diatoms were a major player in a spring bloom.It is interesting that Chl a had high concentrations only at depths ≥30 m in BL during the whole water column blooms in SL in April and July.Since a major player in a mesozooplankton community of BL, the calanoid copepod Calanus helgolandicus, cannot reproduce in SL (Miloslavić 2012), we speculate that the higher grazing on diatoms in the upper water of BL, at least in April, may explain this result.In the warmer and more oligotrophic July, when C. helgolandicus is expected to reside below the thermocline in BL, it is possible that Aurelia sp.feeding on copepods had a positive effect on autotrophic biomass below the thermocline (Turk et al. 2008).In our study, maximum concentrations of nutrients (Figs 5A,B and 6A,B) in the MLs were higher than those in SL in 1997and BL in 1997-1998, except for PO 4 in SL and TN in BL, which were similar to the values detected in 1997-1998 (Benović et al. 2000).We could not compare TP as it was not reported by Benović et al. (2000).Average nutrient concentrations in our study seem to be similar to those reported in Benović et al. (2000) and we cannot infer more about any change in the trophic status of the study area.Under the circumstances, we could not calculate the trophic status of the MLs from the data collected two decades ago, but we present TRIX for the period October 2009-September 2010.Steckbauer et al. (2011) proposed 3.5 mg O 2 l -1 as a new operational threshold to designate hypoxic waters, representing 88% of the distribution of mean lethal DO concentrations in warm coastal waters, except for those of the 12% most sensitive species.That threshold (~52% DO saturation) is less conservative than 4.6 mg O 2 l -1 (~70% DO saturation) for warm coastal waters, representing 90% of the distribution of mean lethal DO concentrations (Vaquer-Sunyer and Duarte 2008).Therefore, it is a scientific effort of the highest priority to re-evaluate DO thresholds that define hypoxic waters, taking into account climate change trends.In this context, the MLs seem to be an excellent laboratory for studying differences in benthic communities related to a future potential t increase.
Detailed research on the thermohaline properties of SL in the mixing period should be conducted to review whether the upper water tends to be markedly stratified year after year.Persistent stratification along with the present trends in climate change, e.g. higher precipitation (Justić et al. 2005), have a potential to extend the duration and impact of hypoxia.Since small systems like marine lakes respond fast to environmental changes (Ciglenečki et al. 2015, Miloslavić andLučić 2015), research on these lakes and similar systems should be encouraged in order to gain more insights into the impacts of driving factors of hypoxia on biodiversity in moderately warm marine systems worldwide.

CONCLUSIONS
In the present study, we have demonstrated that two moderately eutrophic Mediterranean marine lakes that communicate and depend on each other, with only a few significant differences, were affected by hypoxia in different frequencies throughout the year.A long-term chemical feature of significantly higher ammonium concentration in a smaller, shallower lake more isolated from the coastal sea has been shown by comparison with data collected from the period more than a decade ago.A notable difference in temperatures of deep waters in these lakes might be one of the key parameters to focus on in future studies, for example, on the impact of driving factors of hypoxia on the benthic communities in the Mljet NP.For the first time, an apparent oxygen utilization and the OUR as well as OCRR have been presented for this study area.Although not significantly (p>0.05)different between the lakes, AOU, OUR and OCRR indicated that a larger volume of deep water in BL in the mixing period was affected by net oxygen utilization and organic carbon remineralization.These findings are in line with the expectation that more POM in the shallower water escapes remineralization, and is then subsequently oxidized in the surface sediments of SL.

Parameter
Small Lake Big Lake  S2. -Significant (p<0.05)differences among the explored locations in Mljet waters (SL, Small Lake; BL, Big Lake; G, Gonoturska Bay).See the list of abbreviations in Table S1.

Table 1 .
-Sampling dates with corresponding seasons and periods.The first three letters of each month (except for Sept.), season and period are used in the abbreviations.

Table 2 .
-The estimation of the photic zone depth in Small Lake (SL), Big Lake (BL) and Gonoturska Bay (G).Pz, the depth of the photic zone; Sec, the depth of the Secchi disc visibility.

Table S3 .
-Mean±sd of the analysed parameters in the mixing period(October 2009-March 2010); n/a, not analysed; OUR, oxygen utilization rate; OCRR, organic carbon remineralization rate.See the list of abbreviations in TableS1.

Table S4 .
-Mean±sd of the analysed parameters in the stratification period (April 2010-September 2010).See the list of abbreviations in

Table S7 .
-Mean±sd of the analysed parameters in the deep waters of the study locations within the stratification period (April 2010-September 2010).See the list of abbreviations in TableS1.