Depth diversity gradients of macrophytes: Shape, drivers, and recent shifts

Abstract Investigating diversity gradients helps to understand biodiversity drivers and threats. However, one diversity gradient is rarely assessed, namely how plant species distribute along the depth gradient of lakes. Here, we provide the first comprehensive characterization of depth diversity gradient (DDG) of alpha, beta, and gamma species richness of submerged macrophytes across multiple lakes. We characterize the DDG for additive richness components (alpha, beta, gamma), assess environmental drivers, and address temporal change over recent years. We take advantage of yet the largest dataset of macrophyte occurrence along lake depth (274 depth transects across 28 deep lakes) as well as of physiochemical measurements (12 deep lakes from 2006 to 2017 across Bavaria), provided publicly online by the Bavarian State Office for the Environment. We found a high variability in DDG shapes across the study lakes. The DDGs for alpha and gamma richness are predominantly hump‐shaped, while beta richness shows a decreasing DDG. Generalized additive mixed‐effect models indicate that the depth of the maximum richness (D max) is influenced by light quality, light quantity, and layering depth, whereas the respective maximum alpha richness within the depth gradient (R max) is significantly influenced by lake area only. Most observed DDGs seem generally stable over recent years. However, for single lakes we found significant linear trends for R max and D max going into different directions. The observed hump‐shaped DDGs agree with three competing hypotheses: the mid‐domain effect, the mean–disturbance hypothesis, and the mean–productivity hypothesis. The DDG amplitude seems driven by lake area (thus following known species–area relationships), whereas skewness depends on physiochemical factors, mainly water transparency and layering depth. Our results provide insights for conservation strategies and for mechanistic frameworks to disentangle competing explanatory hypotheses for the DDG.


| INTRODUC TI ON
Describing and explaining biodiversity gradients have been central goals of biogeography and ecology since the beginning of the respective fields (Gaston, 2000). Improving our understanding of the biodiversity gradients and their drivers is still an important requirement to deal with impending species loss. Therefore, many studies have explored environmental gradients as explanatory variables for biodiversity patterns along different geographic scales (Rahbek, 2004;Whittaker et al., 2007) such as (a) latitude (Stehli et al., 1969;Rohde, 1992;Pontarp et al., 2019;Etienne et al., 2019), (b) elevation (Colwell & Rangel, 2010;Graham et al., 2014;Hutchinson, 1953;Lomolino, 2001;Nogués-Bravo et al., 2008;Rahbek, 1995;Rahbek et al., 2019;Sanders & Rahbek, 2012), (c) tree height in forests (Petter et al., 2016), (d) depth in soils (Jakšová et al., 2019;Rendoš et al., 2016), or (e) depth in water (Gong et al., 2015;Rex & Etter, 1998;Smith & Brown, 2002). These geographic gradients share some environmental gradients, which are expected to influence spatial structuring of diversity gradients, for example, temperature, light, or seasonality. However, the shorter the spatiotemporal scales are, the less confounding biogeographical contingencies there are, such as the legacy of the glacial cycles on latitudinal gradients, and dispersal/connectivity limitations. Hence, studying gradients expressed at short spatiotemporal extents may provide valuable insights on drivers of biodiversity. Still, the short spatiotemporal gradients, like depth in freshwaters, are often overlooked.
Freshwater ecosystems have a high biodiversity with a high rate of species loss (Gatti, 2016;He et al., 2017;Strayer & Dudgeon, 2010), exceeding those of terrestrial systems (Dudgeon et al., 2006). Nonetheless, studies that focus on the short spatiotemporal diversity gradient in freshwater are surprisingly scarce, although light gradients in freshwater must represent a very strong driver. The few existing studies seem to show predominantly a general decrease of biodiversity along the depth gradient, for example, for bacteria (Cantonati et al., 2014;Zhao et al., 2019), chironomids (Zhao et al., 2019) or diatoms (Kingsbury et al., 2012;Stoof-Leichsenring et al., 2020), or hump-shaped patterns along depth, for example, for diatoms (Zhao et al., 2019) or submerged macrophytes (Ye et al., 2018). Macrophytes are, however, comparatively less studied.
Macrophytes play a pivotal role in lakes by reducing nutrient concentrations (Song et al., 2019), by providing food for a lot of other species (Bakker et al., 2016), and by giving shelter to a large number of other aquatic organisms like zooplankton, juvenile fish, and amphibians (Jeppesen et al., 1998). However, there are several knowledge gaps on macroecology of freshwater plants (Alahuhta et al., 2020). Some aspects of depth gradients in macrophytes can also be found in the literature, with previous studies mainly focus on depth limits of species (Domin et al., 2004;Middelboe & Markager, 1997;Søndergaard et al., 2013), growth of single species dependent on depth (Fu et al., 2018;Li et al., 2020;Xu et al., 2020), functional diversity along depth (Fu et al., 2017), or biomass of macrophytes along depths (Chambers & Kaiff, 1985;Collins et al., 1987). However, the depth pattern of submerged macrophytes species richness is sparsely studied (Bolpagni et al., 2016;Fu et al., 2015;Ye et al., 2018). The few studies that have looked at depth distribution of macrophytes in lakes mainly focused on Lake Erhai in Yunnan Province, China. They report a hump-shaped pattern along the water depth gradient for species richness and community biomass of submerged macrophyte species (Ye et al., 2018). Looking at all functional types including emergent species, Lake Erhai shows a significant decrease in taxonomic and functional diversity along the water depth gradient and its niche differentiation (Fu et al., 2014b(Fu et al., , 2015. Also, hump-shaped and decreasing patterns of species numbers along depth were found in four Italian lakes, changing patterns with time (Bolpagni et al., 2016). Still, it remains unclear if the described pattern of macrophytes is generalizable across multiple lakes and whether it stays robust over time.
The lack of macrophyte diversity gradient studies is intriguing because the environmental gradients along lake depth represent one of the sharpest found in nature, with strong variation over just few meters. With increasing lake depth, multiple abiotic factors that influence the growth of macrophytes (light, temperature, nutrients, water quality, disturbances/hydrologic variability) drastically change (Bornette & Puijalon, 2011). Light is gradually attenuated with increasing depth due to absorption and scattering, resulting in a specific reduction of light quality and quantity depending on depth and on the water turbidity. Water temperature in deep lakes does not decrease gradually, but rather abruptly with depth (Bornette & Puijalon, 2011). The formation of thermally stratified lakes results in an abrupt thermocline, especially during growing season. The thermocline influences the within-lake fluid dynamics in each thermal layer, further leading to stratified gradients in nutrients and water quality components during stratification (Bornette & Puijalon, 2011).
Moreover, mechanical disturbances, like wind or waves, lose their influence gradually with depth (Van Zuidam & Peeters, 2015). The probability that water-level fluctuations result in drying up the soil (Evtimova & Donohue, 2016) is also reduced, and browsing pressure by waterfowl decreases with depth (Bakker et al., 2016). How these different environmental gradients influence the species richness of macrophytes remains unclear, although knowing the processes shaping species diversity might help to predict how global change will affect biodiversity and how management strategies might mitigate potential negative diversity responses.
This study aims to describe the depth distribution of macrophyte diversity, to assess the relative importance of its drivers, and to search for recent shifts. Specifically, we ask the following questions:

1.
1.1. What is the general shape of the depth diversity gradient (DDG) of submerged macrophytes in deep lakes?
To address these questions, we use macrophytes occurrence data of 274 transects taken over 13 years across 28 natural deep lakes in Bavaria that were mapped for monitoring in relation to the European Water Framework Directive. We expect a hump-shaped DDG (question 1.1) corresponding to previous scattered empirical evidence and following the typical patterns found along elevation (Nogués-Bravo et al., 2008). We assume no strong differences between lakes and diversity components (question 1.2), as the pattern along elevation is shown to be generalizable for alpha and gamma richness (Bhatta et al., 2018) and as we consider beta richness as additive partitioning of alpha and gamma richness. To tackle question 2, we test whether the shape of the DDG can be explained by geographic and physical-chemical characteristics of the lakes. We expect water quality to have a high degree of influence, since water quality affects resource availability (light, temperature) (Bornette & Puijalon, 2011). Finally, we assess, with regard to questions 3.1-3.2, whether there have been detectable temporal changes in the DDG.
We suppose that the DDG is a quite stable pattern over time, as macrophytes react slowly to changes (Bakker et al., 2013). However, due to the overall warming in annual average water temperature during the last decades we expect that species richness increases, as invasive species are expected, and warm-adapted species might expand. Our results provide the most refined and extensive assessment of macrophyte biodiversity patterns in freshwater lakes up to date, giving insights for the development of long-term conservation strategies for freshwater systems in general.

| Data and study area
Bavaria has a wide variety of lakes, which vary in size, depth, altitude, and physiochemical parameters. Information about lake surface area, maximal lake depth, and mixing regime was provided by For example, both the concentration of total nitrogen and the concentration of nitrate or even ammonium can influence biodiversity.
The total nitrogen indicates the basic nutrient situation, while the ammonium or nitrate concentrations can affect species with certain nitrogen preferences (Nelson et al., 1981). We did not include Chl-a measurements as they were not consistently available and phosphorus concentration is supposed to be the better predictor for trophic state and lake productivity in central European lakes. Physicalchemical data include monthly measurements at the deepest point of the lake.
Macrophyte data were also extracted from the Bavarian State Office for the Environment (www.lfu.bayern.de). The macrophyte data were recorded for the EU-Water Framework Directive Monitoring following an official sampling strategy (Schaumburg et al., 2015) and include vegetation surveys for all large lakes of Bavaria (>0.5 km² surface area) for at least one and maximum five different years. At each lake, macrophyte data for several transects perpendicular to the shoreline at characteristic sections are available (see sampling strategy in Figure 1a). Each transect is considered from shoreline to the lowest macrophyte occurrence and is subdivided along the depth gradient into four depths (0 to −1 m; −1 to −2 m; −2 to −4 m; >−4 m). At each depth, the frequency of all species is sampled in five steps following the scale of Kohler (1978), an estimate of abundance.

| Data preparation
Data preparation and analyses were done in R 3.5.3 (R Core Team, 2019

| Species richness components and depth diversity gradient measures
As depth-independent component of species richness, we calculated gamma richness as the total number of species per lake and year.
As depth-dependent component of species richness, we determined for every lake and year an additive alpha richness as the number of species per depth range averaged across transects (Figure 1b).
The gamma richness for every lake and year was defined as the total number of species per depth range ( Figure 1c). We then calculated an additive beta richness as gamma richness minus alpha richness (Tuomisto, 2010) (Figure 1d).
To further characterize the diversity depth gradient, we identified the peak of the richness depth curve (Figure 1c). For each transect, we filtered the depth with the maximal species number. Thereafter, we averaged this valued across transects per lake and year, from now on termed the depth with maximal alpha richness (D (α,max) ). The corresponding maximal species number averaged across transects is termed the maximal alpha richness (R (α,max) ). Similarly, the depth with maximal gamma richness (D (γ,max) ) is the depth of the maximal gamma richness (R (γ,max) ) along depth, and maximal beta depth (D (β,max) ) describes maximal beta richness (R (β,max) ) along depth.

| Statistical analysis
We addressed the study questions with several analyses, focusing on different dataset levels dependent on data availability. The biodiversity dataset contains all macrophyte recordings (274 mapped transects in 100 field campaigns, mapping of lake in one year is called field campaign) of the selected 28 lakes. As no complete information is available for all mapped lakes and years, we compiled two subsets of the biodiversity dataset: The environmental & biodiversity dataset is a subset dataset with all macrophyte recordings for which all abiotic data (see Table 1 (GAMMs), computed with the R package "gamm4" (Wood, 2011). D (α,β,γ,max) and R (α,β,γ,max) were used as response variables, the lake as random effect.
To reduce the high correlations between abiotic factors (Pearson correlation test), we performed a principal component analysis (PCA) and named the main axis (>80% variance) after the corresponding abiotic factor, whenever an axis encompassed more than 40% of the variation of a variable. We used the loadings of the main PCA axes (>80% variance) as explanatory variables for the GAMM. We constructed a full model with all PCA axes; then, we stepwise excluded the least significant terms until obtaining a minimal model (Wood, 2008).
To answer questions 3.1 and 3.2, concerning the temporal change of the depth diversity gradient, we used the biodiversity time series dataset. First, we calculated the invariability coefficient (IC) as inverse of the coefficient of variation (CV): The IC is a statistical tool to evaluate the degree of invariability also for datasets with different means (Question 3.1). To check for temporal trends, we built simple linear regression models for depthindependent gamma richness and the DDG measures, D (α,β,γ,max) and R (α,β,γ,max) , as response variables and time as explanatory variable for (a) the complete dataset and (b) each individual lake. We identified all models that showed significant linear trends (p < .1) and characterized the direction of their slopes (Question 3.2).

| RE SULTS
For the dataset of all macrophyte recordings (biodiversity dataset), a total of 75 submerged species is documented in Bavaria. The available taxonomic groups are mainly Spermatophytes (51 species), Charophytes (20 species), Bryophytes (two species), and Pteridophytes (two species). The complete abiotic and biotic data (environmental and biodiversity dataset) cover 57 different species, whereas the biodiversity time series dataset included 66 species. The total (depth-independent) gamma richness per lake ranges from 5 to 34 species of submerged macrophytes. The mean gamma richness averaged across lakes is 15.36 species with a standard deviation of 6.27 species.

| The depth diversity gradient (DDG) patterns of macrophytes
The mean depth pattern of submerged macrophytes' alpha and gamma richness is hump-shaped, showing a peak between −1 and −2 m, respectively (Figure 2a,c). The mean alpha richness at the hump's peak is 4.5 species (SD = 2.2), whereas the mean gamma richness peak is 11.4 species (SD = 5.1). In contrast, beta richness shows a decreasing curve with its highest richness being 7.0 species

| Drivers of the depth diversity gradients
The DDG measures correlate with some of the abiotic variables (Supporting information). R (α,max) correlates highly significantly The GAMM analysis for D (β,max), R (β,max), D (γ,max), and R (γ,max) had all R 2 < 2.1% (see results in Supporting information).

| Temporal dynamics of the depth diversity gradients
The DDG measures of the different richness components show different degrees of invariability coefficient (IC) as measure for stability. The IC of Gamma richness (mean = 7.14, SD = 3.69), of D (α,max) (mean = 6.62, SD = 4.11), of R (α,max) (mean = 6.36, SD = 2.87), of R (β,max) F I G U R E 2 Depth diversity gradient (richness along depth) total mean (black line) and standard deviation (gray area) for alpha (a), beta (b), and gamma richness (c

F I G U R E 3
The PCA loadings (a-d) of the four first axes (81% variation) of a PCA for all abiotic variables were used as explanatory variables of the GAMM. The names of the axes are given by the variables with >40% of the loading (highlighted in blue). The full variable names are given in Table 1. Panels a-d are ordered corresponding to the order of a-d. PCA biplot can be found in the Supporting information. Minimal best models for generalized additive mixed models explaining D (α,max) by PC2 (e), PC4 (f), PC3 (g), and PC1 (h) and R (α,max) by PC1 (i). GAMMs without significant explanatory variables are not shown (R (β,max) , R (γ,max) , D (β,max) , D (γ,max) ). Panels e-h are ordered by decreasing drop contribution (dc). The relative contribution of each variable to the full model (% drop contribution) is measured as the drop in deviance explained by the model when the predictor is removed. Significance levels of p-values: *p < .05; **p < .01; ***p < .001. R 2 adj for D (α,max) is 73.5%, for R (α,max)   shows two significant positive trends (peak shifting toward water surface-lake Starnberg and lake Tegernsee) and two significant negative trends (peak shifting to deeper waters-lakes Großer Alpsee and Woerthsee). The D (β,max) increases significantly at Lake Riegsee, while D (γ,max) increases for lakes Riegsee, Staffelsee Nord, and Tegernsee.

| The DDG patterns of macrophytes
We showed that submerged macrophytes in deep lakes have in general a hump-shaped depth diversity gradient (DDG) for alpha richness, a prevailing decreasing pattern for beta richness and a dominantly hump-shaped pattern for gamma richness (Figure 2a-c) (question 1.1). As we had only significant differences between mid and greater depths for all richness components, an even coarser species mapping resolution might be helpful. Our results are congruent to the few existing studies, which also show a hump-shaped pattern (Ye et al., 2018) for alpha richness of submerged macrophytes.
A simple explanation for the predominantly hump-shaped pattern of alpha and gamma richness might be the mid-domain effect (Colwell et al., 2004): Niches along environmental gradients overlap and build a peak of richness following geometric constraints. Furthermore, the generally decreasing beta DDG might be explained by a change in local species between transects in shallower depths. In shallow water, disturbances resulting from the surface might be more diverse, which may increase coexistence with spatial partitioning of occurrences. As disturbances are coming from the surface, we propose the hypothesis that shallow water has a higher environmental heterogeneity which might be the reason for an increased beta diversity.
We see a high variety of DDG shapes between lakes, as all these DDGs vary in their maximum richness (R max ) and the corresponding depth (D max ), but a robust hump-shaped pattern can be seen for alpha richness (Figure 2d-f) Gr. Alpsee ; significance levels of p-values: 0 "***" 0.001 "**" 0.01 "*" 0.05 "." 0.1.

Lake Niedersonthofen
2 and 35 species per lake. This wide range of species richness and environmental conditions broadens our understanding of the DDG, previously limited to one single lake in China . Although for alpha and gamma richness, hump-shaped curves along depth are predominant ( Figure 2g-i), we also see increasing and decreasing patterns at single lakes. Increasing curves must be hump-shaped, as we can safely assume that plant richness should decrease to zero further down in deep lakes. We detected more decreasing DDGs for gamma than for alpha richness, which reflects predominantly decreasing beta richness curves. Nevertheless, besides geometry and disturbances, there must be further variables affecting the DDG, as DDG shape varies between lakes, which themselves have different properties.

| The DDG drivers
The drivers of the macrophyte DDGs strongly differed for DDG measures and richness components. Whereas pairwise correlations detected many strong relationships across richness components, multiple models revealed significant variables only for DDG measures of alpha richness.
The R max correlate of the different richness components with a very similar set of abiotic parameters. All R max correlate with area. This reflects species-area relationships (SAR) (Connor & McCoy, 1979;Lomolino, 2000;Patiño et al., 2014;Qian et al., 2007) for macrophytes, which is also shown by high correlation of area with total gamma richness. Looking at the GAMM results for nonlinear responses, R (β,max) and R (γ,max) did not show any significant results, but R (α,max) is exclusively influenced by the "SiO 2 & Conductivity axis" (PC1) (Figure 3). Area, SAC, and WLF also have a high contribution to PC1 and thus, area might be the key driving force again. The positive SAR of macrophytes was already shown in several studies (see Alahuhta et al., 2020 for a review). Moreover, a study comparing SAR of macrophytes with terrestrial plants would be very informative.
For this purpose it would be interesting to add information about lake bathymetry, as lake area is just a proxy for the colonizable area per depth. Still, very generalizing indices like volume development index are also not suitable to determine the colonizable area in this case because the lake's morphology is very diverse (from kettles with several deep funnel-shaped basins to v-shaped glacial lakes and lake basins created by glacial erosion). However, it was not shown yet that the size of lakes also influences the shape of DDG.
The D (β,max) and D (γ,max) could not be explained with abiotic variables, neither by correlations nor by a GAMM. Unlike D (α,max) , the gamma, and consequently beta, values along DDG are more variable, indicating spatial heterogeneity and possibly unsaturation (Karger et al., 2014). Still, D (α,max) correlates positively with P tot and Tempsd (higher P tot and Tempsd evoke a D (α,max) in shallower water) and negatively with O 2diss and Transp. Furthermore, looking at nonlinear influences, the D (α,max) is affected by all four PCA axes.
The PC2 (Temperature & P tot axis) shows the highest influence ( Figure 3). This means that in lakes with high P tot and/or high Temp the DDG peaks in shallower waters. Phosphorus is the limiting factor for phytoplankton growth and phytoplankton reduces the light availability for macrophytes. In contrast, it is still debatable whether the phosphorus concentration in the water is a limiting growth factor for macrophytes (Carr et al., 1997). One important point to consider in this debate is that rooted submerged macrophytes can also take up nutrients from the sediments (Lacoul & Freedman, 2006). Hence, phosphorus might affect macrophytes by promoting phytoplankton growth, which then reduces light availability and shifts DDG to shallower depths. Besides phosphorus, temperature is a major factor influencing metabolic processes as photosynthesis and respiration.
Additionally in lakes, higher and prolonged high temperature result in higher nutrient levels due to increased mineralization and internal fertilization processes (Moss, 2012). Internal fertilization processes occur when longer high water temperatures lead to increased layering stability, prolonged oxygen consumption, anoxia in deep waters, resulting in anoxic resuspension of phosphorus from the lake sediments. These resuspended nutrients promote phytoplankton growth, thus reducing light for macrophytes. Therefore, the PC2  (Huston, 2014;Rajaniemi, 2003;VanderMeulen et al., 2001).
Besides light, temperature (Temp) also seems to influence D (α,max) , via surface water temperature and its influence on light availability (as explained above) and via the lake's layering depth. This second mechanism by which temperature layering affects DDG is demonstrated by the fact that along PC3 (Tempsd-Chloride axis) D (α,max) decreases. A high Tempsd (shallow epilimnion-the upper temperature layer in a stratified lake) promotes a shallow D (α,max), while a low Tempsd (broad epilimnion) allows deeper D (α,max) . Tempsd is positively correlated to Temp demonstrating that higher temperatures can lead to a shallower upper warm layer in water bodies as the stratification is more stable (Adrian et al., 2009).
The weakest influencing effect (lowest drop contribution) is provided by PC1 (NH4 + -SiO2 & Conductivity axis). Just at very high values of PC1 D (α,max) becomes shallower. As Cond is negatively correlated with Transp (cor = −0.71, p < .001), we speculate that also here Transp is the actual mechanism that influences D (α,max) .
In summary, the main influences on D (α,max) seem to be, as expected, factors of water quality that influence light quantity (transparency, influenced by phosphorus and temperature), light quality (CDOM), and layering depth (temperature). The main influence on R (α,max) is the lake surface area.

| DDG temporal change
We showed that the stability of the pattern depends on the DDG measure (question 3.1). D (β,max) and D (γ,max) were quite variable measures over years, while D (α,max) , R (α,max) , R (β,max), and R (γ,max) are comparatively stable measures. This may be related to the fact that there is neither pairwise correlation between nor an explaining model for D (β,max) and D (γ,max) .
Contrary to our expectations, we see no general trend of increasing species richness or decreasing D max (question 3.2). Although we observe high variety in DDG temporal change between lakes, the DDG temporal change for single lakes, especially for D (α,max) , is low and develops into different directions for different lakes. Still, we see linear trends that are consistent over time within lakes. These patterns suggest that global change effects will be more complex than anticipated. In fact, climate and land use change influence all the highly connected chemical and physical gradients known to significantly affect DDG (Hossain et al., 2017). Thus, the following hypotheses can be formulated ( Figure 4): (1) As temperatures rise, so do lake surface water temperatures (O'Reilly et al., 2015;Pilla et al., 2020).
This seems to result in shallower epilimnion  and generally shallower D max and a lower R max .
(2) Furthermore, rising temperatures entail higher phosphorus content, as they promote internal fertilization. But extreme weather events combined with enriched fertilization in agriculture can also cause fertilization events (Rose et al., 2016), which might result in shallower light depth and consequently in shallower DDG.
(3) Browning, which is generally increasing due to temperature-induced decomposition rates and changes in precipitation events (Guarch-Ribot & Butturini, 2016;Sobek et al., 2007;Weyhenmeyer & Karlsson, 2009), leads to a shallower D max . (4) However, water management reduced the external nutrient loading of European lakes enormously during the last decades (Eigemann et al., 2016;Murphy et al., 2018). This trend is still ongoing and might still lead toward lower nutrient contents and thus to deeper D max . All these opposing environmental trends make it hard to draw a general trend for multiple lakes for short timespans.
However, for long timespans it seems to be a race between climate change impacts (Hypothesis 1-3 in Figure 4) that might lead to a shallower D max and thus generally less macrophytes and water management impacts that might deepen the D max via improved water quality (Hypothesis 4 in Figure 4). In summary, this study sets a good comparison for future studies once longer time series become available.

| Implications for diversity gradients and hypotheses in general
Comparing different diversity gradients might provide deeper insights into mechanisms shaping species richness. Here, the DDG of macrophytes brings advantages compared with other gradients. The DDG assembles over shorter spatial scales than the latitudinal diversity gradient (LDG) or the elevational diversity gradient (EDG), which implies a lower importance of dispersal or connectivity and an easier replicability. For LDG, the options for replicates are restricted to two hemispheres , whereas EDG comparative studies require a high logistic sampling effort (Kessler et al., 2011;Nogués-Bravo et al., 2008). Because of these advantages, other small-scale diversity gradients may also be insightful. One example is the vertical diversity gradient (VDG) from forest floor to tree crowns, which involves sharp gradients of light intensity, temperature, and humidity. Therefore, we discuss below the potential explanatory hypotheses with all mentioned diversity gradients.
One explanation of the observed hump-shaped DDG might be the intermediate productivity hypothesis (IPH). The IPH states that at low productivity level (deep waters with low light quantity and low temperature) only few specialized species survive, whereas at high productivity level (shallow waters) only few competitive species F I G U R E 4 Summary figure showing the submerged macrophytes depth diversity gradient for alpha, beta, and gamma richness as well as the main drivers (black arrows) of the alpha richness peak, characterized by R (α,max) and D (α,max) . Additionally, hypotheses for global change influences on the alpha richness peak for Bavaria are given in the gray box. The hypotheses are: (1) layering depth might become shallower due to rising water temperature.
(2) Light quantity might be lowered due to lowered transparency.
(3) Light quality is said to decrease. (4) Light quantity might be increased if water management gets adapted survive. Previous LDG study of freshwater plants revealed its peak at subtropical to low tropical latitudes (Murphy et al., 2019), thus peaking at intermediate level of solar productivity and reflecting our analysis of DDG. Intermediate light intensity and temperature would also match the mid-canopy VDG peak for vascular epiphytes (Acebey et al., 2017;Krömer et al., 2007;. Although quantification of productivity along depth should be attempted, our findings and the evidence from other diversity gradients already indicate a key role of light quantity and temperature in shaping DDG. Another hypothesis is the mid-domain effect (MDE), which proposes mid-gradient peaks due to geometric constraints (Colwell et al., 2004). If depth ranges from shallow-water species overlap with ranges of deep-water species, a species richness peak in the middle of the gradient can be expected. The MDE is used to explain hump-shaped patterns of the LDG  and the EDG (Colwell & Lees, 2000). Indeed, the overlap of light and temperature preferences may explain the subtropical peak (Murphy et al., 2019) in LDG of macrophytes and our reported mid-depth peak in DDG as well as the mid-canopy VDG peak in vascular epiphytes . Still, an adequate evaluation of the MDE requires quantification of environmental preferences for each species-an important direction for future empirical research. In this regard, the MDE evaluation may be more feasible to perform for DDG, as it considers a smaller regional species pool than the LDG and a better experimental feasibility than VDG given the faster life cycles of macrophytes compared with vascular epiphytes.
Another explanation is the intermediate disturbance hypothesis (IDH). It suggests species richness peaking at mid-levels of disturbance as species of early and late successional phases coexist (Connell, 1978). Whereas the disturbances along EDG are caused by human activities at lower elevation (Nogués-Bravo et al., 2008) and the disturbances along VDG can be associated with higher branchfall toward the outer crown of a tree (Cabral et al., 2015;, depth-dependent disturbances in water can be caused by anthropogenic use, waves, herbivory, ice cover, and water-level fluctuations (Evtimova & Donohue, 2016). Water-level fluctuation was already integrated in our study in a very simple way, but showed no strong effect on richness, thus did not explain the DDG. Nevertheless, considering that several disturbances in shallow waters should happen in lakes, future monitoring schemes should quantify more types of disturbances.

| Limitations and perspectives
The main limitation is that, in some lakes, the deeper end of the DDG up the lowest depth level to have a finer resolved depth gradient and to quantify a metric termed "the lower macrophyte limit" (Søndergaard et al., 2013). This metric is often used as indicator for water quality and might be useful to further characterize the DDG as it defines the lower limit and the occupied space.

| CON CLUS ION
Our study makes a step toward a cross-lake generalizable understanding of the depth diversity gradient (DDG) of submerged macrophytes, their regional and temporal heterogeneity as well as the drivers of the DDG shape. Submerged macrophyte richness peaks predominantly at intermediate depth forming a hump-shaped pattern for alpha and gamma richness, but a decreasing pattern for beta richness ( Figure 4). Well-known hypotheses of biogeography shape diversity gradients in general, such as mid-domain effect and meanproductivity hypothesis. The latter is already supported by our findings on the role of light and temperature as DDG drivers. The key advantage of DDG in contrasting these hypotheses is the logistic feasibility of short-distance scales and the exclusion of confounding effects associated with dispersal constraints. The key drivers of DDG we determined were area influencing the species richness peak height (R (α,max) ) and light quality, light quantity, and layering depth influencing the species richness peak depth (D (α,max) ). However, as there are many other possible factors for which we did not have data but which could play a role, further research is needed before general conclusions can be drawn from this study. Although we found that the DDG in general remained stable over the past few years for most lakes, we still found shifting trends for richness metrices for some lakes. However, these trends were shown to be diverse across lakes. Whereas climate change might be more ubiquitous, land use change may be lake-specific. This suggests that water management strategies should also consider, besides global warming, lake characteristics, and change in the surrounding land use. The interaction of these aspects also means that although higher temperatures lead to a reduction in the quantity of light available to aquatic plants in lakes, land use measures can be taken to counteract this. Nevertheless, our findings already indicate that warmer water temperatures may still lead to a shift in species along depth towards shallower waters dependent on further efforts to hold or increase water quality of lakes.

ACK N OWLED G M ENTS
We thank the LfU, in particular Christine Schranz and Anette Maetze, for data provision. A special thanks go to Ludwig Leidinger,

CO N FLI C T O F I NTE R E S T
The authors declare that they have no conflict of interest.

DATA AVA I L A B I L I T Y S TAT E M E N T
Original raw data are publicly available, provided by Bayerisches