Morphology of bar-built estuaries: empirical relation between planform shape and depth distribution

. Fluvial-tidal transitions in estuaries are used as major shipping fairways and are characterised by complex bar and channel patterns with a large biodiversity. Habitat suitability assessment and study of interactions between morphology and ecology therefore require bathymetric data. While imagery offers data of planform estuary dimensions, only for a few nat-ural estuaries bathymetries are available. Here we study the empirical relation between along-channel planform geometry, obtained as the outline from imagery, and hypsometry, which characterises the distribution of along-channel and cross-channel 5 bed-levels. We ﬁtted the original function of Strahler (1952) to bathymetric data along four natural estuaries. Comparison to planform estuary shape shows that hypsometry is concave at narrow sections with large channels, while complex bar morphology results in more convex hypsometry. We found an empirical relation between the hypsometric function shape and the degree to which the estuary width deviates from an ideal convergent estuary, which is calculated from river width and mouth width. This implies that the occurring bed level distributions depend on inherited Holocene topography and lithology. Our new 10 empirical function predicts hypsometry and along-channel variation in intertidal and subtidal width. A combination with the tidal amplitude allows an estimate of inundation duration. The validation of the results on available bathymetry shows that predictions of intertidal and subtidal area are accurate within a factor 2 for estuaries of different size and character. Locations with major human inﬂuence deviate from the general trends because dredging, dumping, land reclamation and other engineering measures cause local deviations from the expected bed-level distributions. The bathymetry predictor can be used to characterise 15 and predict estuarine subtidal and intertidal morphology in data-poor


Introduction
Estuaries develop as a result of dynamic interactions between hydrodynamic conditions, sediment supply, underlying geology and ecological environment (Townend, 2012;de Haas et al., 2017).One model for the resulting morphology is that of the 'ideal estuary' that is hypothesised to have along-channel uniform tidal range, constant depth and current velocity, and a channel width that exponentially converges in landward direction such that the loss of tidal energy by friction is balanced by the gain in tidal energy by convergence (Savenije, 2006;Townend, 2012;Savenije, 2015;Dronkers, 2017).One would expect that in this case the along-channel variation in hypsometry is also negligible.However, natural estuaries deviate from the ideal ones as result of varying degree of sediment supply, lack of time for adaptation and sea-level rise (Townend, 2012;de Haas et al., 2017) and locations wider than ideal are filled with tidal bars (Leuven et al., 2017) (Fig. 1).Differences in bed-level profiles between ideal and non-ideal estuaries are further enhanced by damming, dredging, dumping, land reclamation and other human interference (e.g.O'Connor, 1987;Wang and Winterwerp, 2001;Lesourd et al., 2001;Jeuken and Wang, 2010;Wang et al., 2015).All these natural deviations from the ideal estuary mean that there is no straightforward relation between the planform geometry of the estuary and the hypsometry or distribution of depths.
Hypsometry captures key elements of geomorphological features (Strahler, 1952;Boon and Byrne, 1981;Dieckmann et al., 1987;Kirby, 2000;Toffolon and Crosato, 2007;Townend, 2008Townend, , 2010;;de Vet et al., 2017) (Fig. 2).The hypsometric method was developed by Strahler (1952) and Boon and Byrne (1981) to relate planform area of a basin to elevation.Later the results were used to predict the influence of basin morphology on the asymmetry of the horizontal and vertical tides, to predict flood-or ebb-dominance and maturity of an estuary (Boon and Byrne, 1981;Wang et al., 2002;Moore et al., 2009;Friedrichs, 2010) and to characterise the trend of saltmarsh development (Gardiner et al., 2011;Hu et al., 2015).Furthermore, hypsometry was used as a data reduction method to characterise entire reaches spanning bars and channels in estuaries (Toffolon and Crosato, 2007) and shapes of individual tidal bar tops (de Vet et al., 2017).Hypsometry has also been used to describe dimensions of channels and tidal flats in an idealised model (Townend, 2010).Here, a parabolic shape was prescribed for the low water channel, a linear profile for the intertidal low zone and convex profile for the intertidal high zone.While these profiles are valid for perfectly converging channels, it is unknow to which extent they are applicable to estuaries with irregular planforms and whether the currently assumed profiles are valid to assess flood-or ebb-dominance.Hypsometric profiles and derived inundation duration are also relevant indicators for habitat composition and future transitions from mudflat to saltmarsh (Townend, 2008).In order to predict and characterise morphology and assess habitat area, we need along-channel and cross-channel bed-level predictions for systems without measured bathymetry (Wolanski and Elliott, 2015).
For only a few natural estuaries bathymetry is available, which leaves many alluvial estuaries with irregular planforms from all around the world underinvestigated.However, many are visible in detail on satellite imagery, which raises the question whether there is a relation between planform geometry and depth distribution.Such a relation is known to exist in rivers in the form of hydraulic geometry depending on bar pattern, and meander pool depths depending on planform channel curvature (e.g.Kleinhans and van den Berg, 2011;van de Lageweg et al., 2016).It therefore seems likely that such a relation between the horizontal and vertical dimensions exists for sandy estuaries, but this is not reported in literature.Morphological models can simulate 3D bed-levels with considerable accuracy (van der Wegen andRoelvink, 2008, 2012;van Maren et al., 2015;Braat et al., 2017), but these models are computationally intensive and need calibration and specification of initial and boundary conditions.To study unmapped systems with limited information, it would be useful to be able to estimate bed-level distributions from planform geometry.Here we investigate this relation.
Previously, we showed in Leuven et al. (2017) how locations and widths of tidal bars can be predicted from the excess width, which is the local width of the estuary minus the ideal estuary width.The summed width of bars in each cross-section was found to approximate the excess width.This theory describes bars as discrete recognisable elements truncated at low water  level on what is essentially a continuous field of bed elevation that changes in along-channel direction (Leuven et al., 2016).
However, predicting morphology in more detail requires predictions of along-channel and cross-channel bed elevations.While hypsometry can summarise bed elevation distribution as a cumulative profile, it is not clear whether the shape of the profile is predictable.Our hypothesis is that the along-channel variation in hypsometry depends on the degree to which an estuary deviates from its ideal shape.Therefore, we expect that locations with large excess width and thus large summed width of bars have a more convex hypsometry (Fig. 3c,d).In case of ideal estuary with (almost) no bars (Fig. 3a,b), we expect concave hypsometry.
The aim of this manuscript is to investigate the relation between estuary planform outline and along-channel variation in hypsometry.To do so, first hypsometric curves are used to summarise the occurring bed elevations in a cumulative profile.
Then, we use the original function of Strahler (1952) to fit the data obtained from bathymetry of four estuaries (Fig. 1).In the results, we develop an empirical function to predict hypsometry.The quality and applications of the predictor are assessed in the discussion.

Methods
In this section, first the definition of an ideal estuary is given with a description of how we derived a geometric property that characterises deviation from an ideal shape.Second, the general form of a hypsometric curve is described.Then, the available datasets that were used for curve fitting are given.Last, the methodology to fit a hypsometric function to bathymetry in systems is presented.

Deviation from 'ideal'
A useful model to describe the morphology of estuaries is that of the 'ideal estuary'.In an ideal estuary the energy per unit width remains constant along estuarine channels.This ideal state can be met when tidal range and tidal current are constant along-channel, such that the loss of tidal energy by friction is balanced by the gain in tidal energy per unit width by channel convergence (Pillsbury, 1956;Dronkers, 2017).In case the depth is constant along the channel, the ideal estuary conditions are approximately met when the width is exponentially decreasing in landward direction (Pillsbury, 1956;Langbein, 1963;Savenije, 2006;Toffolon and Lanzoni, 2010;Savenije, 2015), which also implies an along-channel converging cross-sectional area.However, when depth and friction are not constant along-channel, but for example linearly decreasing in landward direction, less convergence in width is required to maintain constant energy per unit of width.Many natural estuaries are neither in equilibrium nor in a condition of constant tidal energy per unit width.They deviate from the ideal ones as result of varying 10 degree of sediment supply, lack of time for adaptation to changing upstream conditions and sea-level rise (Townend, 2012;de Haas et al., 2017).Whether continued sedimentation would reform bar-built estuaries into proper ideal estuaries remains an open question.For our application, the concept of ideal estuaries is useful to assess the degree of deviation from ideal because of the width variations observed as bars, tidal flats and saltmarsh.While we expect a somewhat different degree of convergence such that the ideal state of constant energy per unit of width is approximately maintained, we do not study the deviation of this convergence length from that in ideal estuaries.
Ideally, we would want to assess the degree to which an estuary is in equilibrium from an aerial photograph, because this is often the only data available.However, the only indicator derivable from aerial photography is channel width and thus deviation from a converging width profile.Therefore, in Leuven et al. (2017), we defined the excess width, which is the local width of the estuary minus our approximation of the potential ideal estuary width.Here, the ideal estuary width is approximated as an exponential fit on the width of the mouth and the width of the landward river.While the empirical measure of 'ideal width' should not be confused with the 'ideal state' of an estuary, it is the only practical way to estimate deviation from an ideal estuary based on the estuary outline only.Moreover, it proved to be a good indicator of occurring bar patterns (Leuven et al., 2017) and will therefore be applied in this paper to study hypsometries.

General hypsometric curve
In the past, multiple authors have proposed empirical relations for the hypsometric shape of terrestrial landscapes (Strahler, 1952) and (partially) submerged bodies (Boon and Byrne, 1981;Wang et al., 2002;Toffolon and Crosato, 2007;Townend, 2008) (see Townend, 2008, for review).All equations, except for Wang et al. (2002), predict a fairly similar hypsometric curve based on the volume and height range of the landform (Townend, 2008).While it is of interest to use these empirical relations to predict the occurring altitude variation of a landform, the framework here is different, because in this case the hypsometric curve is applied across channel width: it is a cross-sectional hypsometric curve that we change along the system.We aim to use the general hypsometric curve to characterise the occurring cross-sectional hypsometry, which is similar to the approach of Toffolon and Crosato (2007) who fitted a power function to 15 zones along the Western Scheldt.However, the zoned approach smooths out all differences between bar complex and channel-dominated zones, which are of interest for this study.To that end, it is less relevant for which environment the hypsometric relation was proposed, as long as it is capable to describe the range of occurring hypsometries.For the case of the estuarine environment (Fig. 3), the hypsometric curve should be able to describe variations in concavity and variations in the slope of the curve at the inflection point.Here we use the original Strahler (1952) formulation, which is capable to do so, but in principle any equation that fits well could be used.Strahler (1952) formulated the general hypsometric curve as: in which h z is the value of the bed elevation, above which fraction y of the width profile occurs.In other words, h z is the proportion of total section height and y the proportion of section width.r sets the slope of the curve at the inflection point in a range of 0.01-0.50,with sharper curves for lower values of r (Fig. 2a,b).z determines the concavity of the function in a range of 0.03-2, with lower values giving a more convex profile and higher values giving a more concave profile (Strahler, 1952) (Fig. 2).Our approach changes the original definition of r and z to make them fitting parameters.It is expected that z-values depend on excess width, because the fraction of the width occupied by bars becomes larger with excess width, resulting in a more convex hypsometric profile (Fig. 3c,d), or the presence of bars generates excess width.
Excess width is defined as the local width minus the ideal width, which is given by: in which x is the distance from the mouth, W m the width of the mouth and L W is the width convergence length (Davies and Woodroffe, 2010), which can be obtained conservatively from a fit on the width of the mouth and the landward river width (Leuven et al., 2017): in which W m is the local width measured at the mouth of the estuary, W s is the width measured at the landward side of the estuary and s is the distance between these locations measured along the centreline.This practical method makes the convergence length somewhat sensitive to the selected position of the seaward and landward limit.
The landward limit was selected where the width ceases to converge on an image at the resolution of the full estuary scale and the landward width was measured between the vegetated banks.The seaward limit was selected as the location with the minimum width in case bedrock geology, human engineering or a higher elevated spit confined the mouth, because in these cases the minimum width limits the inflow of tidal prism.In other cases, the mouth was chosen at the point where the first tidal flats were observed in the estuary or where the sandy beach ends at the mouth of the estuary.However, when the mouth is chosen at a location where sand bars are present, the ideal width will be overestimated and the width of intertidal area underestimated.It is therefore recommended to either chose the mouth at a location where bars are absent or subtract the width of bars from the measured with at the mouth to obtain the ideal width profile.

Data availability and classification
Detailed bathymetries were available for four systems: the Western Scheldt (NL), Dovey (Wales), Eems (NL) and the Columbia River estuary (USA) (Fig. 1, Table 1).Data for the Western Scheldt and Eems were obtained from Rijkswaterstaat (NL), for the Dovey estuary from Natural Resources Wales and for the Columbia River estuary from the Lower Columbia Estuary Partnership.Bed elevations were extracted from these bathymetries as follows.First, the estuary outline was digitised, excluding fully developed saltmarshes, and subsequently a centreline was determined within this polygon (following the approach of Leuven et al., 2017).Bed elevations were collected on equally spaced transects perpendicular to the centreline of the estuary.The data extracted at these transects was subsequently sorted by bed level value and made dimensionless to obtain hypsometric profiles (see Fig. 3 for examples).
We classified the transects by morphological characteristics and potential susceptibility to errors.The following morphological classes were used: mouth, bar junction, bar complex, narrow bar, pointbar, channel, pioneer marsh.The mouth is the location where the estuary transitions into the sea.A bar junction is the most seaward or most landward tip of tidal bars.A bar complex is a location where a large bar is dissected by barb channels (Leuven et al., 2016) or multiple smaller bars are present.Narrow bar is used when the bars present were narrow along their entire length and often also relatively flat on their top.Pointbar is a bar in the inner bend of a large meander.Channel was assigned when bars were largely absent.Pioneer marsh was assigned when aerial photographs or bathymetry gave visual indications of initial marsh formation, such as the presence of small tidal creeks and pioneering vegetation.Fully developed marsh is excluded from the outline.
The following classes were used to indicate possible errors: presence of harbours, major dredging locations, presence of a sand spit, presence of drainage channels for agriculture, constraints by hard layers, human engineering works.Either a locally deep channel or scour occurred at one of the sides of these transects or they lacked a natural transition from channel to estuary bank, thus ending in their deepest part on one side of the transect.Major dredging locations have unnaturally deeper channels and shallower bars, resulting in a hypsometric shape that is relatively flat in the highest and lowest part and is steep in between (Fig. 3e,f).Furthermore, in a few cases side channels were perpendicular to the orientation of the main channel of the estuary.
This resulted in transects being along-channel of these side channel, which biases transect data towards larger depth and creates a flat hypsometric profile at the depth of the side channel.

Data processing
Least-squares fits resulted in optimal values of z and r in Eq. 1 (Fig. 2) for each transect, using three different approaches.
First, a regular least-squares curve fitting was used for r and z, which resulted in along-channel varying values for z, but an almost entirely constant along-channel value for r of 0.5 (Fig. 4b,c, solid lines).In the second approach we set r to a constant value of 0.5 and only fitted to obtain z (Fig. 4c, dashed line).We found that the quality of the fit was the same, as indicated by the root-mean-square error (RMSE) (Fig. 4d) and therefore apply this second approach in the remainder of this paper.
Locations where the RMSE was relatively large correspond to locations where major dredging occurred in the past century.
This possibly resulted in a hypsometry characterised by a larger fraction of the width occupied by high tidal flats, a larger fraction of the width occupied by deep channels and a smaller fraction of the width occupied by the zone between channels and bars (Fig. 3e,f).Because the hypsometric curve at these locations deviated from the original Strahler function (Eq.1), our third approach was to apply a modified function to find optimal values for z and r.To do so, the original formulation of Strahler (1952) was inverted to allow for hypsometries that describe steep transitions from bar top to channel bottom, because the original does not nearly fit as well: Applying this modified function resulted in better fits, but only at locations that were classified to be excluded because of possible errors.Therefore, results from this approach are not shown here and it is suggested to study the effect of dredging and dumping on hypsometry in more detail in future studies.
In principle, both the bed elevation (h z ) and the width fraction y in Eq. 1 are dimensionless.To compare the resulting predictions with measured values, the prediction needs to be dimensionalised.Values for y are scaled with the local estuary width.We test three options to scale h z .The first option is to scale h z between the highest bed elevation and lowest bed elevation in the given cross-section, which is sensitive to the precise cut-off of the bathymetry.The second option is to scale h z between the local high water level (HWL) and the maximum estuary depth in that cross-section, which is sensitive to bathymetric information that is usually not available in unmapped estuaries.
The third option requires a prediction of depth at the upstream or downstream boundary.Width-averaged depth profiles along estuaries are often (near-) linear (Savenije, 2015;Leuven et al., 2017), which includes horizontal profiles with constant depth.
Therefore, only the channel depth at the mouth of the estuary and at the upstream river have to be estimated and subsequently a linear regression can be made.Channel depth at upstream river (h r ) is estimated with hydraulic geometry relations (e.g.Leopold and Maddock Jr, 1953;Hey and Thorne, 1986): h r = 0.12W r 0.78 .The depth at the mouth is estimated from relations between tidal prism and cross-sectional area (e.g.O'Brien, 1969;Eysink, 1990;Friedrichs, 1995;Lanzoni and D'Alpaos, 2015;Gisen and Savenije, 2015;Leuven et al., 2017).Here we used: in which P is the tidal prism, which can be estimated by multiplying the estuary surface area with the tidal range (Leuven et al., 2017): this assumes a flat water surface elevation along the estuary and neglects portions of the estuary that might get dry during the tidal cycle (Boon, 1975).Locally, the maximum depth may be deeper or lower than predicted, due to the presence of resistant layers in the subsurface or where banks are fixed or protected.While this may affect the accuracy of the locally predicted maximum channel depth, it has a minor effect on the calculations of sub-and intertidal area.Moreover, the upper limit for dimensionalisation is chosen as the high water line, which implies that supratidal area is not included in the predictions.We will show results with all methods, but only the third can be applied when information about depth is entirely lacking.
Statistical analyses in the remainder of paper were approached as follows.In linear regressions, we minimised residuals in both the x-and y-directions.This results in regressions that are more robust than when residuals in only one direction are minimised.In case regressions are plotted, the legend will specify the multiplication factor that the confidence limits plot above or below the trend.R 2 values are given to indicate the variance around the regression.In cases where the quality of correlation of two along-channel profiles is assessed, we used the Pearson product-moment correlation coefficient (r).

Results
We found a strong relation between along-channel variation in hypsometry and the degree to which an estuary deviates from its ideal shape.Below, we will first show how the hypsometry of typical channel deviates from that of bar complexes.Subsequently, the data is presented per system, classified on their morphology and potential for errors due to human interference, method and other causes.Then we combine all data to derive an empirical formulation to predict hypsometry.Last, we apply this formulation to predict the along-channel variation in intertidal and subtidal width and validate the results with measurements from bathymetry.

Empirical relation between morphology and hypsometry
As hypothesised, it is indeed observed that channel-dominated morphology results in more concave hypsometry profiles (high z-value), while bar complex morphology results in more convex hypsometry (low z-value) (Fig. 5).Values for z in Eq. 1 range from 0.83 to 1.14 for channels, with an average value of 1.0.In contrast, z ranges from 0.36 to 0.41 for bar complexes, with an average of 0.39.
Clustering of morphological classes strongly suggests a relation between hypsometry and planform estuary shape (Fig. 6).
Mouth and channel-dominated morphology typically plot at the right-hand side of the plots in Fig. 6a,c,e,g, thus being locations close to ideal width.In case of the Western Scheldt this results in the highest values of z.In the case of the Dovey and Columbia River, the mouth region was influenced by respectively a spit and human engineering, which resulted in the formation of tidal flats on the side and thus lead to a lower z-value.Bar complexes occur at the other end of the spectrum; these locations are generally much wider than the ideal shape and are characterised with hypsometries with a z-value well below 1. Bar junctions as well as narrow bars are generally found at the transition from channel-dominated morphology to bar complex morphology and therefore also occur between these types in the plots.The point bar in the Western Scheldt (the Plaat van Ossenisse) shows hypsometry comparable to bar complexes (Fig. 6a), which reflects the complex history of formation by multiple bar amalgamations.Also the locations in the Columbia River where pioneer marsh is present, show the same trend as the locations where unvegetated bar complexes occur (Fig. 6g).
In a few cases, the transects used to extract bathymetry were not perpendicular to the main channel of the estuary.For example, landward and seaward of the pointbar in the Western Scheldt (Plaat van Ossenisse) transects were inclined, covering a larger part of the channel than perpendicular transects, resulting in higher z-values as a consequence of the apparent channeldominated morphology.Immediately landward of the spit in the Dovey, transects are almost parallel to the shallow sidechannel.Fitting hypsometry at these locations resulted in relatively low z-values, because it is a relatively shallow side-channel.
For the Western Scheldt and Eems it is known at which locations major dredging and dumping takes place (e.g.Swinkels et al., 2009;Jeuken and Wang, 2010;Bolle et al., 2010;Dam et al., 2015;Plancke and Vos, 2016).Even though the resulting z-values at these locations do not cause major outliers, the quality of fits are typically lower and the inverted Strahler function (Eq.4) fitted better.These points were therefore excluded from further analysis.
The filtered data shows quasi-cyclicity in along-channel hypsometry (Fig. 7).In general, the width at the mouth of the estuary and at the upstream estuary is close to ideal and the hypsometry is concave, except in systems with wide mouths and bars in the inlet.The part in between is characterised by variations in the local width and therefore gradual increases and decreases in the ratio between local width and ideal width.In some cases, quasi-cyclic loops are visible (e.g.Fig. 7c) caused by the asymmetry in bar complexes.In other cases, the points show more zigzag or clustered patterns, which indicate minor variation on the bar complexes or scatter in the fit applied to the bathymetry.

Hypsometry predictor
The relations between excess width and hypsometric function are similar for all estuaries, which suggest that a universal function is of value.Combining all the filtered data resulted in a regression between the extent to which an estuary deviates from the ideal shape and the predicted z-value in the hypsometry formulation (Fig. 8).Data from Columbia River, Eems in 1985, Western Scheldt in 2013 and Dovey were used to obtain this relation.Other data are shown in Fig. 8 but not used in the regression.
These results mean that we found a predictive function for hypsometry, where r is set to a constant value of 0.5 and the z(x) is calculated as: in which W ideal (x) is the ideal estuary width (Eq.2) and W (x) is the measured local width.The confidence limits of the regression plot a factor 1.9 higher and lower than the regression, which indicates that the z-value can be predicted within   a factor 2 (see Fig. 9 for an example of prediction with uncertainty).While not used in the regression, hypsometry from bathymetry in other years show similar trends and scatter as the data used in the regression.
The predictor (Eq.6) was applied to the Columbia River estuary, Western Scheldt, Dovey and Eems to check the quality of the resulting along-channel predictions of intertidal high, intertidal low and subtidal width (Fig. 10, Fig. 11).These zones can be derived after dimensionalising hypsometry and imposing a tidal range (Fig. 9b).For almost the entire along-channel profile, 5 the predictions are within a factor 2 of the measured value (Fig. 11) and best agreement was obtained when the hypsometry was dimensionalised between the minimum and maximum measured bed level for each transect (Fig. 10).

Discussion
Results from this study illustrate that bed level distributions of channel and bar patterns in estuaries are topographically forced.
The estuary outline that is observable from the surface translates into the three dimensional patterns below the water surface.
Bar built estuaries typically have a quasi-periodic planform (Leuven et al., 2017), in which major channel confluences occur at locations where the estuary is close to its ideal shape (Leuven et al., 2018).The parts between the confluences are typically filled with intertidal bar complexes.These findings are consistent with hypsometry zonations previously found for the Western Scheldt with more concave hypsometries in for channel-dominated morphology and more convex hypsometries for bar complex morphology (Toffolon and Crosato, 2007).Our cross-sectional approach additionally revealed quasi-periodic behaviour within these zones.
In contrast to an empirical description, ideally, a physics-based determination of the hypsometry would be favourable.
However, with the current state of the art of bar theory (Leuven et al., 2016) and relations for intertidal area, tidal prism, crosssectional area and flow velocities (O'Brien, 1969;Friedrichs and Aubrey, 1988) it is not yet possible to derive a theoretical prediction of hypsometry.For example, bar theory (Seminara and Tubino, 2001;Schramkowski et al., 2002) could predict occurring bar patterns on top of an (ideal) estuary shape, but current theories overpredict their dimensions (Leuven et al., 2016) and it is still impossible to scale these to bed level variations, because the theories are linear.In addition to that, the resulting predictions would need to meet the requirement that the predicted bed levels and the intertidal area together lead to hydrodynamic conditions that fit the estuary as well.2.
Previously, hypsometry was used to summarise the geometry of entire tidal basins or estuaries (Boon and Byrne, 1981;Dieckmann et al., 1987;Townend, 2008).The whole system descriptions are consistent with the original Strahler (1952) concept of a basin hypsometry based on plan area, which is a valid description in a landform context.However, these descriptions oversimplify the along-channel variability in estuaries that are relatively long.These estuaries typically have a linear bed profile varying from an along-channel constant depth to strongly linear sloping (e.g. the Mersey in UK).In the latter case, the elevation at which subtidal and intertidal area occur varies significantly along-channel (Blott et al., 2006).Additionally, friction and convergence may cause the tidal range to either dampen or amplify causing variation in tidal elevation, subtidal area and tidal prism (Savenije, 2006).Consequently, the along-channel cross-section hypsometry should be assumed to be relative to an along-channel varying high water level or mean sea level rather than an along-channel fixed vertical datum.Interpreting these along channel variations remains an open question because of the reasons outlined above.Nevertheless, if desired, alongchannel varying hypsometry predictions can be converted in one single summarising curve (Fig. 12), which shows that also the basin hypsometry can be predicted when limited data is available.
Our results show that hypsometry is not only a tool to predict morphology when limited data is available, but that hypsometry can also be used to reduce a large dataset of bathymetry and to study evolution of bathymetry over time.In the case of Strahler (1952), hypsometry fits result in along-channel profiles of zand r-values, but in practice any function or shape could be fitted.
For example, the locations along the Western Scheldt where major dredging and dumping took place showed weaker correlation with the original Strahler (1952) shape (Fig. 13).In these cases, fits with higher quality (lower RMSE) were obtained when we used the inverted hypsometric function (Eq.4) (Fig. 13b,e).So in practice, one could fit a range of different hypsometry shapes and subsequently find out which of these shapes fit best on the used dataset.It can indicate that certain parts along the estuary require a separate hypsometric description.The fitting parameters are a method to describe the along-channel variation (Fig. 13b,d).Hypsometry can be fitted to compare data from nature, physical experiments and numerical modelling and subsequently study for example the effect of vegetation, cohesive mud and the influence of management on these systems.

Implications for management of estuaries
In many estuaries from around the world subtidal channels are used as shipping fairways, while the intertidal bars (or shoals) form valuable ecological habitat (e.g.Bouma et al., 2005).For example in the Western Scheldt, the shipping fairway is now maintained at a depth of 14.5m below Lowest Astronomical Tide (de Vriend et al., 2011;Depreiter et al., 2012).With empirical hypsometry predictions we can estimate the width below a certain depth required for shipping, which gives estimates of what volume to dredge and at what locations along the estuary, which is relevant for construction of future shipping fairways in estuaries for which we may have limited data.
In contrast, low-dynamic intertidal areas are valuable ecological habitat, for example for the Western Scheldt there is an obligation to maintain a certain amount of intertidal area (Depreiter et al., 2012).Previously, Townend (2008) showed that basin hypsometry can be a tool to design breaches in managed re-alignment sites and can provide an indication of habitat composition.Hypsometry analysis per cross-section shows that estuary outline translates into intertidal area, which implies that locations where the estuary is relatively wide have a relatively wide intertidal area.The ecological value is determined by the area of low-dynamical undeep water and intertidal areas (for settling and feeding) (Depreiter et al., 2012).This means that the edges should neither become steeper nor higher (leading to permanent dry-fall) or deeper.Hypsometry fits (in case of The occurrence of vegetation species depend on bed elevation, salinity, maximum flow velocity and sediment type (de Jong, 1999;Gurnell et al., 2012).Even though predicted hypsometry only gives bed elevations, a comparison of the height interval in which Salicornia and Spartina can occur (Mckee and Patrick, 1988;Davy et al., 2001;van Braeckel et al., 2008), showed similar trends and the same order of magnitude as the measured vegetation from ecotope maps of the Western Scheldt in 2012 (Fig. 14).Some underpredictions arise in parts along the estuary where bed elevations above the high water level occur, such as at the Drowned land of Saeftinghe.However, in general the vegetation width is overpredicted because (1) hypsometry is stretched between the high water line and channel depth and (2) other constraining biotic and abiotic factors were excluded.

Conclusions
We studied the relation between along-channel planform geometry of sandy estuaries and their hypsometry, which characterises the distribution of along-channel and cross-channel bed-levels.The vertical dimensions were found to relate to the horizontal dimensions.In other words, the degree to which the estuary width deviates from an ideal converging estuary shape reflects in the occurring hypsometry.At locations where the width is much larger than ideal, convex hypsometric shapes are observed, contrary to the locations where the estuary width is close to ideal, where concave hypsometric shapes are observed.In between these extreme end members, a gradual transition with quasi-periodic variation was observed.This implies that it is possible to predict the along-channel varying hypsometry of estuaries, which is relevant for estuaries for which limited data is available.

Figure 3 .
Figure 3. Example cross-sections and hypsometry, suggesting that channel-dominated morphology (a,b) generally results in concave hypsometric functions (high z-value in Fig. 2) while bar complexes (c,d) generally result in convex hypsometry (low z-value in Fig. 2).In case of narrow bars with a flat top and relatively steep transition from bar top to channel, the original hypsometry function by Strahler (1952) is less appropriate and the inverted function fits better.

Figure 4 .
Figure 4. (a) Width along the Western Scheldt, with the maximum ideal converging width profile indicated.The green area is defined as the excess width (Leuven et al., 2017).For each along-channel transect of the estuary, the optimal fit of z and r in the Strahler (1952) function (Fig. 2) was determined.(b,c) Results for the Western Scheldt when both z and r are freely fitted as well as the results when r is fixed to a constant value of 0.5.(d) Quality of fits remains about the same when r is set to a fixed value of 0.5 as indicated by the root-mean-square error (RMSE).(e) Fitted z-values show similar trends as the ideal width divided by local width.

Figure 5 .
Figure 5. Hypsometric curves as extracted from bathymetry in the Western Scheldt, Dovey, Eems and Columbia River.Cross-sectional profiles were extracted along the centreline (Section 2.3) and were subsequently classified as channel when sand bars were (mostly) absent and as bar complex when one larger bar dissected by barb channels or multiple smaller bars were present.Channel-dominated morphology generally results in concave hypsometric profiles and bar complexes in convex profiles.

Figure 6 .
Figure 6.Results of hypsometry fitting, where clustering indicates a relation to planform geometry.(a,c,e,g) z-values were fitted for crosschannel transects in bathymetry of four estuaries and plotted by morphological classification.(b,d,f,h) Regressions for z-value as a function of ideal width divided by local width.Data points that were influenced by human interference, bedrock geology or errors in methodology or data indicated in red in panels a,c,e,g were excluded.Confidence limits are plotted at two standard deviations above and below the regression and their multiplication factor compared to the trend is given in the legend.

Figure 7 .
Figure 7. Fitted z-values as a function of deviation from the ideal width.Colours indicate the location along the estuary, with dark blue colours at the mouth transitioning into dark red colours at the landward end.Some zones show scatter in fitted z-values, while some other zones (e.g.green to orange to red in c) show quasi-periodic behaviour.

Figure 8 .
Figure 8. Fitted z-values of filtered data increase with the fraction of ideal width and local width, which indicate that hypsometric shapes become progressively more concave when the local width approaches the ideal width and become more convex when the local width becomes larger than the ideal width (i.e. the excess width increases).The data shown as asterisks is used for the regression.Confidence limits are plotted at two standard deviations above and below the regression and their multiplication factor compared to the trend is given in the legend.

Figure 9 .
Figure 9. Illustration of uncertainty of predicted hypsometry from Eq. 6 with uncertainty margins (Fig. 8).Resulting prediction for hypothethical location where (a) local width is equal to the ideal width and (b) local width is five times larger than the ideal width.Results are compared against a typical tidal range in order to show uncertainty of predicted intertidal high width, intertidal low width and subtidal width as a fraction of the total estuary width.

Figure 10 .
Figure 10.Comparison of measured and predicted values of intertidal high, intertidal low and subtidal width for the Columbia river estuary (a,c,e) and the Western Scheldt (b,d,f).Predicted nondimensional hypsometry was dimensionalised for each cross-section using three methods (explained in methods) and uncertainty margins are given for one of the predictions (solid black line).In the legend, r indicates the Pearson product-moment correlation coefficient and dev the average factor of deviation between the predicted (TP-HWL) and measured lines.

Figure 11 .
Figure 11.Comparison of measured and predicted width of intertidal and subtidal width.The (solid) line of equality indicates a perfect fit and dashed lines indicate a deviation of a factor 2. Percentage of measurements within these margins are indicated in Table2.

Figure 12 .
Figure 12.Hypsometry as summarised in a single curve for the entire estuaries.Solid lines are measured from bathymetry, dashed lines based on the predictions.

Figure 13 .Figure 14 .
Figure 13.(a) Evolution of cross-section in the Western Scheldt where significantly has been dredged and dumped (The Drempel van Hansweert), (c) Measured hypsometric profiles of the same time steps.(e) Best fit hypsometries using original Strahler equation [solid] and inverted Strahler function (Eq.4) [dashed], (b) Quality of the two types of fits shows that the shape of the best fitting hypsometric curve changes from convex to the inverted equation in the 1970s-1980s.(d) Fitting coefficients for z increase over time for both hypsometry types and both transects.(f) Intertidal high area increased over this period while intertidal low area remained constant [Transect 1] or decreased [Transect 2].
and is partly funded by the Ministry of Economic Affairs.This work is part of the PhD research of JRFWL.We acknowledge a data processing contribution by Andy Bruijns as part of his MSc guided research.Reviews by Ian Townend and one anonymous reviewers helped to improve the manuscript.

Table 1 .
Characteristics of estuaries used in this study.

Table 2 .
Percentage of points predicted within a factor 2 from the measured value.