How human activities affect the fine sediment distribution in the Dutch Coastal Zone seabed

Abstract The fine sediment distribution in the seabed is an important indicator for the ecological functioning of shallow coastal seas. In this paper, we investigate the processes and conditions that determine the fine sediment distribution in the Dutch coastal zone surficial seabed, while also assessing the response of the system to human interventions. An extensive sediment dataset, collected in the Dutch coastal zone from 2006 to 2014, is presented. These data are used to map the distribution of fines in the seabed of the DCZ at unique spatiotemporal scales. For the entire Dutch coastal zone, the distribution of fines generally agrees well with previous studies. The recent extension of the Port of Rotterdam, the Maasvlakte 2 reclamation, was found to locally change the distribution of fines. In the sand mining pit and directly south of the reclamation, fines percentages in the seabed increased by more than 10%. We developed a conceptual framework to analyse the distribution of fines and how it is affected by human interventions. Three components are distinguished within this framework: (1) sources of fines; (2) transport pathways; and (3) accumulation potential. These components are determined both qualitatively and quantitatively, based on high-resolution bathymetric and hydrodynamic model data. The distinction between the three components makes it possible to unravel the contributions of different human interventions to the changes in the fines distribution. In the case of Maasvlakte 2, the local increase of fines percentage in the seabed could thus be attributed to a temporary additional source of fines and enhanced accumulation potential. The high spatiotemporal resolution of the new sediment dataset proved crucial to enable development and testing of the framework to evaluate the impact of (large) engineering works on the spatial distribution of fines.


Introduction
The seabed of the southern North Sea has a significant ecological and economic value. It accommodates a substantial amount of living and non-living resources and fulfils vital ecosystem services, by providing habitat for a wide range of benthic organisms (Rees et al., 2007;Stephens and Diesing, 2015). Over the past years, human activities affecting the seabed have increased (Degraer et al., 2019;Emeis et al., 2015). These activities include sand mining, operation of offshore platforms (Stolk and Dijkshoorn, 2009), beam trawling (Rijnsdorp et al., 2008), accommodation for pipelines and cables buried in the sea bed (de Groot, 1982;Rouse et al., 2017), and the construction and operation of offshore wind farms (Breton and Moe, 2009). These human activities lead to increasing pressure on the southern North Sea ecosystem (Burdon et al., 2018;Piet et al., 2019). To make sure these human activities are carried out in a sustainable way, a balance between economic and ecological interests must be maintained. An important parameter affecting the local ecosystem is the amount of fines in the water column and seabed, which are related. This requires knowledge about the dynamics and composition of the seabed sediment (Degraer et al., 2019;Reed et al., 2012).
In a large part of the southern North Sea, the seabed mainly consists of sand (Eisma et al., 1987;Huthnance, 1991;Irion and Zollmer, 1999), containing a relatively small fraction of fines, i.e., sediment with a grain size smaller than 63 μm. However, many studies have shown that small fractions of fines can readily exert a profound influence on the behaviour of the seabed and the benthic ecosystem (e.g., Degraer et al., 2008;Heip et al., 1992;van Ledden et al., 2004). Benthic communities are richer when the seabed contains fines compared to purely sandy substrates (Van Hoey et al., 2004), because nutrients adhere to fines (van Raaphorst and Malschaert, 1996). On the other hand, fines may have a detrimental effect on the pelagic ecosystem when remobilized from the seabed. Once fines are suspended, they can abrade fish gills, leading to gill damage in several fish species (Au et al., 2004;Sutherland and Meyer, 2007). Furthermore, suspended fines increase the turbidity of the water, thereby attenuating the light climate and thus the growth rate of phytoplankton (e.g., Anthony et al., 2004;Van Duin et al., 2001). Favourable or not, fines influence the ecological functioning of shallow coastal seas. Stephens and Diesing (2015) and Bockelmann et al. (2018) were the first to quantify the spatial distribution of fines for the entire North Sea, based on a large number of seabed samples. They quantified the seabed sediment composition (e.g. mass percentages of fines, sand and gravel) of the entire North Sea by using a geostatistical approach. Because of their large spatial extent, the resolution of these maps is relatively low (Bockelmann et al., 2018;Stephens and Diesing, 2015). Furthermore, these studies did not explicitly include the effect of human activities on seabed sediment composition. This makes them less applicable to coastal areas, where environmental gradients are large and human activities are ubiquitous. To understand how fines are distributed in the seabed of coastal areas, the processes which play a role on smaller spatial scales have to be understood, including the role of human activities.
In this paper, we analyse the fine sediment distribution in the seabed of the Dutch Coastal Zone (DCZ), a coastal area characterized by strong environmental gradients and a variety of human interventions. The aim of this study is to identify the processes and conditions which determine the distribution of fines in the surficial seabed of the DCZ, and assess the response of the system to human interventions. We present a new, extensive sediment dataset, collected by the Port of Rotterdam authority in the DCZ from 2006 to 2014 (Borst and Vellinga, 2012). These data are used to map the distribution of fines in the seabed of the DCZ at unique spatiotemporal scales. To interpret these distributions, we develop and apply a conceptual framework, utilizing hydrodynamic model data and bathymetric data of the DCZ. This framework is used to evaluate the effects of human interventions on the distribution of fines in the DCZ. Additionally, we discuss other areas where the framework could be applied.
The paper is structured as follows: In Section 2 we present our study area and discuss the sediment dataset, bathymetric and hydrodynamic data. Next, we present the results of our analyses on the sediment dataset and introduce our conceptual framework. We then apply this to the study area, evaluating the effect of large-scale human interventions, and also discuss how the framework could be used in other areas.

Study area
The study area is depicted in Fig. 1a and covers part of the Dutch coastal zone (DCZ) (Fettweis and Van Den Eynde, 2003;Visser et al., 1991), which is situated in the southern North Sea (Fig. 1a). The DCZ is a shallow coastal shelf sea with maximum water depths up to 30 m and tidal currents with maximum velocities ranging between 0.7 and 1.1 m/s ( de Kok, 1996;van der Giessen et al., 1990). The progressive tidal wave propagates through the North Sea in a counter-clockwise direction (Kelvin wave). It has an amplitude of 1-2 m along the Dutch coast (van der Hout et al., 2015;Visser et al., 1991), with tidal currents oriented mainly parallel to the shore. Furthermore, the outflow of the River Rhine induces a Region of Freshwater Influence (ROFI), which extends for over 100 km along the coast with an average width of less than 20 km (Pietrzak et al., 2011). This ROFI determines the vertical current structure and resulting suspended matter distribution (de Boer et al., 2009;Pietrzak et al., 2011;Simpson et al., 1993;Souza and Simpson, 1996).
Apart from these physical traits, this area is known for a multitude of human activities taking place. Fig. 1b shows a selection. From multiple offshore platforms, gas and oil is extracted and several major shipping lanes cross the DCZ, where opposing traffic lanes are separated by separation zones. More recently, wind farms have been constructed and are planned. Closer to the shore, sand mining areas and disposal sites are found. At these disposal sites, sediment dredged from harbours is deposited, containing large amounts of fines. Sand from the mining areas serves multiple purposes: it is used for coastal protection, construction and land reclamations.
A major land reclamation realized in the past decade was Maasvlakte 2 (MV2). MV2 is the recent extension of the Port of Rotterdam, constructed mainly in 2009 and 2010. This required a total volume of 220 million m 3 of sand, which was mined from the MV2 sand mining pit, located approximately 10 km west of the River Rhine outflow. During 2009 and 2010, a total volume of 170 million m 3 was mined (Borst and van Tongeren, 2012). After 2010, sand mining for MV2 continued for several years, but at substantially smaller rates (de Jong, 2016).

Sediment samples dataset
We present a new dataset consisting of more than 1700 bed samples. This dataset is established from an extensive monitoring programme carried out between 2006 and 2014 by the Port of Rotterdam authority (Borst et al., 2013;Borst and Vellinga, 2012). The aim of this programme was to monitor the far-field and near-field effects of sand mining for MV2 on the benthic ecosystem. Within this programme, bed samples were collected at pre-defined sampling stations from 2006 until 2014 on a yearly basis, except in 2007. However, the exact sampling locations varied a bit from year to year around the predefined stations. In every sampling year, samples were collected in the period of April to June, which is the post-storm season.
In 2009, 2013 and 2014, only near-field effects of the sand mining were monitored. To establish these near-field effects, 100 to 120 seabed samples were collected within a densely sampled 15 km radius around the MV2 sand mining pit. In 2006In , 2010 to 300 stations were visited to monitor far-field effects of the sand mining. In these years, the sampling domain covered the majority of the study area shown in Fig. 1a. The far-field sampling domain includes the near-field domain, but with a lower sampling density. Still, most dense sampling was done around the MV2 sand mining pit. To the north and south, sampling density decreased. Fig. 1a shows the spatial sampling density in 2012, while Table 1 lists the number of stations visited each year, classifying the sampling years according to sampling domain and MV2 chronology.
A standard protocol was followed to collect sediment samples. A large seabed sample, with a maximum height of 25 cm and 30 cm diameter, was taken with a boxcorer. Three Perspex tubes (length: 15 cm, diameter: 10 mm) were inserted into the boxcorer sample, before the overlying water was siphoned off, not disturbing the sediment interface. These tubes were carefully removed from the mother sample, removing the excess sediment around. Each tube was then split into two parts: an upper part (0-5 cm from the surface) and a lower part (5-10 cm). The lower parts of each tube were combined and stored in one 20 ml vial, i.e., the lower subsample. The same procedure was followed for the upper parts, i.e., the upper subsample. The vials were labelled and stored in a freezer at −20°C. After all stations were visited, the vials were taken to the laboratory ensuring the sediment remained frozen.
The grain size distribution of the samples was determined in the laboratory. First, the subsamples were freeze-dried and passed over a 1 mm sieve. Then, the sieved material was homogenized in local tap water and part of it inserted into a Malvern Mastersizer 2000. The Malvern Mastersizer determines the grain size distribution of a sediment sample by laser diffraction and returns the volume percentage of different size classes. The volume percentage of particles smaller than 63 μm is returned as a separate size class. We refer to this size class as fines, and do not distinguish between the clay and silt fractions.
We assume the volumetric fines percentage measured with the Malvern Mastersizer is close to the gravimetric fines percentage. This is valid if the density of the sediment does not vary considerably, i.e., when the amount of organic matter in the sediment samples is limited (Callesen et al., 2018;Yang et al., 2015). Based on Loss on Ignition (LOI) data of the sediment samples this is an appropriate assumption, as LOI was smaller than 2% for more than 95% of the samples.
As the grain size distribution of the upper and lower subsamples taken from the boxcore is statistically dependant, they are not treated as separate samples. We define the average grain size distribution as the average of the two subsamples. This reflects the grain size distribution of the surficial seabed, i.e. the top 10 cm, for a visited sampling station per sampling year. Furthermore, to account for measuring accuracy, a sample is classified as containing fines if its fines percentage is at least 0.1%. If the fines percentage is smaller, it is classified as a sample with no fines.
To characterize the sediment composition for each station, we aggregated the particle size distribution of the samples collected during the various sampling years. However, sampling at a particular station was not carried out at exactly the same location over the years, while the sampling density and domain also varied. Therefore, we introduce a spatial clustering procedure to assess which data are attributed to which station.
The spatial clustering consists of three subsequent steps. First, a circular buffer is defined around each sampling point in QGIS. As the sampling density varied across the study area, the radius of this buffer depends on the sampling point location. Four main sampling density zones were defined: very sparse, sparse, dense and very dense (Fig.  1a). The corresponding buffer radiuses for each zone are listed in   Table 2. An example is shown in Fig. 2 for an arbitrary sampling station in the sparse sampling zone. Second, sediment composition data are aggregated to form a data cluster if their buffers overlap. When a buffer overlaps any other buffer, its data is added to the cluster. For the sampling station in Fig. 2, the data cluster represents the samples from 2006, 2008, 2010, 2011 and 2012. Third, data clusters were designed such that they do not contain multiple samples from one year, except for the very dense sampling zone (Fig. 1a). In this zone, clusters may contain multiple samples taken during one year.
After aggregating the data, the mean and standard deviation were calculated for the fines percentage percluster. Furthermore, we established the fraction of samples in a cluster which contain fines. This fraction is an estimate for the probability of fines being present in a sample for any cluster.

Bathymetric and hydrodynamic data
We use bathymetric data collected by the Netherlands Hydrographic Office and Rijkswaterstaat, already interpolated to an equidistant grid with 25 × 25 m 2 resolution (Damen et al., 2018;NLHO and Deltares, 2019). Bathymetric data collected during multiple years was merged onto a single grid, as the area of interest was only partially surveyed during subsequent years. For the pre-MV2 bathymetry, we use data collected from 1994 until 2008, and for the post-MV2 bathymetry, data collected from 1994 until 2015. If areas were surveyed multiple times, the bathymetry was based on the latest survey. Missing values on the merged grid were filled by linearly interpolating from surrounding grid points within a 500 m radius, using inverse distance weighting.
Flow velocities and salinity for the study area are extracted from a validated three-dimensional hydrodynamic model, with 10 equidistant vertical layers (Arcadis, 2014;Arcadis and Deltares, 2019). The southern and northern boundaries of this model are located at 51.1°N and 52.8°N, respectively. Its eastern boundary lies along the Dutch shoreline and its western boundary runs parallel to the shoreline, 50 km offshore. The resolution of the curvilinear model grid is lowest at the western boundary with cell sizes of 2500 × 2500 m 2 . It refines in shoreward direction and is highest in the dense and very dense sampling zones (Fig. 1a). In the area of interest, cell sizes range from 250 × 350 m 2 to 500 × 700 m 2 , where the along-shore length of the grid cell is smallest. Note that the computational grid is therefore much coarser than the bathymetric grid, the relevance of which is discussed in Section 4.
With this model, Arcadis (2014) carried out hindcast simulations for the years 2006 to 2014. For each model year, Arcadis (2014) updated the model with the latest bathymetric data. Water levels and salinity at the seaward boundaries of the model were taken from the southern North Sea (ZUNO) model (Gautier and Caires, 2015). River discharges at the landward boundary were based on output from a calibrated 1D model of the fresh water distribution in the Rhine-Meuse delta (SOBEK) and measurements. The model has been validated for water level, temperature and salinity. For these three quantities, model performance was assessed as: -Water level: BIAS = 5 cm, RMSE 0 = 8 cm -Temperature: BIAS = −0,5°C, RMSE 0 = 0,5°C -Salinity: BIAS = 0.5 PSU, RMSE 0 = 1.5 PSU More information on the model setup and validation can be found in Alkyon (2010) and Arcadis (2014). The model output was resampled to a 500 × 500 m 2 grid using inverse distance weighting for all sampling areas shown in Fig. 1a. We define two representative years for pre-MV2 and post-MV2 hydrodynamic conditions: 2008 (pre-MV2) and 2012 (post-MV2).
In our analyses, we relate the distribution of fines to current-as well as wave-induced bed shear stresses. Current-only bed shear stresses (τ b, c ) are calculated based on the flow velocities computed by the Arcadis (2014) model, following Soulsby (1997): Here, ρ denotes the density of seawater (1030 kg/m 3 ), U b is the magnitude of the flow velocity at the lowest model level and C D is the drag coefficient. Thus, computed bed shear stresses are positive scalar values. The value of C D is determined by the bed roughness length z 0 and the height above the bed z, according to: Here, , and d 50 is the median sand grain size, for which we take d 50 = 250 μm.
Wave-induced bed shear stresses were taken from the MoS 2 model (Cronin and Blaas, 2013).

Mapping the distribution of fines in the Dutch coastal zone
An overview of the monitoring programme results is presented here. We focus on the percentages of fines found in the surficial seabed layer (i.e., 0-10 cm from the seabed surface). First, the results of all years (2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014) are discussed, a total of 1790 bed samples. The aggregated data are discussed later. The distribution over the years is indicated in Table 3. The mass percentage of fines in a sample is denoted as φ fines . For the samples containing fines, we compute a conditional mean percentage, 〈φ fines 〉, and a conditional standard deviation φ fines '.
The percentage of samples containing fines in the far-field domain ranged between 17% and 33% ( Table 3). The mean fines percentage increased from 3.7% in 2006 to about 5-6% in the period 2008-2011 and then further to 9.2% in 2012, with a considerably higher standard deviation in 2012.
The percentage of near-field samples containing fines was larger than in the far-field domain (Table 3), with percentages varying between 43% and 80%. In this area, 〈φ fines 〉 was generally higher than in the far-field. Even though local temporal differences exist in fines percentage, the difference between the far-field and near-field sampling years is primarily attributed to the difference in sampling density and domain.
The spatial distribution of fines in the surficial seabed of the DCZ, aggregated from all sampling years, is presented in Fig. 3a. It shows four classes of average fines percentage per sampling station. The location of these stations reflects the sampling coordinates of 2012. Only stations with clusters containing at least three samples are displayed (Section 2.2). The classes for fines percentages are conform van Alphen (1987), to allow for comparison.
The highest percentages of fines are found directly north and south of the River Rhine outflow. In these areas, mean fines percentages may range up to 25%. To the north, most fines are found on the lower shoreface, in a narrow alongshore strip 2 to 3 km wide. This strip extends to the northern boundary of the study domain. To the south, fines are mainly found in the troughs of tidal ridges and in the former tidal channels in front of now-closed estuaries. Both the average fines percentages and fraction of samples containing fines decrease with increasing distance offshore. This fraction is small or equal to zero for the majority of stations beyond 20 km offshore (Fig. 3b).

Comparison with historical data & spatial distribution
The percentage of fines in the surficial seabed of the DCZ has been mapped before by Eisma (1968) and van Alphen (1987). Although the sediment collection method, determination of the grain size distribution, and definition of surficial seabed, i.e. ranging from top 5 cm (Eisma, 1968) to top 5-15 cm (van Alphen, 1987), differ from the present study, a qualitative comparison is possible. Fig. 4 shows that the results of Eisma (1968) and van Alphen (1987) are globally similar to the current results. This is consistent with the fact that no long-term trends in the alongshore flux of fines in the DCZ have been reported (Cronin and Blaas, 2015;Eisma, 1981;Salden, 1998;van Alphen, 1990). Eisma (1968) characterizes the shoreface between Rotterdam and Den Helder as "an area with fine grained deposits about parallel to the coast at roughly 5-15 m depth" (left panel Fig. 4). This area is also reflected by our data in the right panel of Fig. 4. Though Eisma's definition of fines is slightly different (i.e. b50 μm), the agreement between the two datasets is promising. Moreover, this implies that for our analyses, their interpretation and application of our concept elsewhere, the precise definition of fines does not seem crucial.
Later, van Alphen (1987) presented a more detailed analysis of the fines distribution in the DCZ, based on data collected between 1964 and 1987 (middle panel Fig. 4). van Alphen (1987) notes that fines are found in several areas: in the former tidal channels of the Southern Delta, around Loswal Noord (close to Hook of Holland), the site where sediment dredged from the Port of Rotterdam was disposed until 1996, and in a narrow 1-2 km wide strip along the coast. In this strip, fines are mainly found in troughs between the breaker bars and around the 10 m depth contour. The fines distribution map by van Alphen (1987) resembles that of the current study. However, spatial patterns around the River Rhine outflow have changed. Furthermore, in the current study some fines are also found further offshore. In all three studies, the alongshore occurrence of fines coincide with a zone where computed current-only bed shear stresses are low (Fig. 3c). This is caused by a decrease in tidal velocities at smaller water depths. On the other hand, wave-induced shear stresses become larger at smaller water depth, as these stresses scale inversely quadratic with water depth (Fig. 3d). Thus, fines in the DCZ are mainly found where tidal velocities are small and wave stresses are larger than in the majority of the DCZ. At water depths beyond 10 m, our simulations show that, averaged over a year, the role of waves reduces compared to the flow-induced bed shear stresses (detailed results not presented).
To examine the cross-shore variability in fines percentage, we plot all samples containing fines as a function of distance perpendicular from the shore. Fig. 5a shows that the majority of samples with fines is found within 20 km from the shore, with fines percentages generally higher than 1% and up to 60%. Further offshore, the fines percentages decrease, ranging between 0.5 and 2%.
In Fig. 5b, all samples are grouped into four bins, based on their distance from the shore. Beyond 20 km offshore, less than 5% of the samples contain fines. This is consistent with van Alphen (1987), who concluded that fines are virtually absent beyond 20 km offshore. Closer to the shore, fines were found in 30% of all samples in the 10-20 km bin and in 57% of all samples taken within 10 km from the shore. However, within this last area the spatial and temporal variability in the measured fines percentage is significant.

Temporal and spatial variability in Rotterdam area
To investigate the variability in fines percentages within the nearshore zone (b10 km offshore), we zoom in on the area around the River Rhine outflow. Since the River Rhine outflow is also the entrance to the port of Rotterdam, we will refer to this area as the Rotterdam area. As the construction of MV2 was the major human activity in this area during the 2006-2014 period, sampling years were classified relative to the construction of MV2. Hence, the sampling years 2006 and 2008 are pre-MV2 and the sampling years 2011-2014 are post-MV2 (Table 1).
After applying the spatial clustering as described in Section 2.2, we established the pre-MV2 and post-MV2 fines percentage. Since the pre-MV2 and post-MV2 sampling layout differ to some extent, the difference in fines percentage per station was calculated by comparing points within 1000 m distance (Fig. 6). This provided sufficient distance between sampling stations, but also allowed to compare enough stations. The average pre-MV2 fines percentage per station was subtracted from the post-MV2 fines percentage, i.e., a positive value in the right panel of Fig. 6 indicates an increase in fines for the post-MV2 period.
In the left and middle panel of Fig. 6, only stations with clusters containing at least two samples are displayed. We have also indicated where human activities take place.
For the pre-MV2 years (left panel Fig. 6), we observe the highest fines percentages around the disposal sites Verdiepte Loswal (1) and Noordwest (2), and to the northeast of these sites. Fines are also found within several kilometres from the shore, mainly in the vicinity of the River Rhine outflow. However, large spatial gradients in fines percentages are observed everywhere. Sampling stations without fines are adjacent to stations where fines percentages are between 2 and 10% or even exceed 10%.
Overall, fines percentages are higher for the post-MV2 years than for the pre-MV2 years (middle and right panel Fig. 6). The areas northeast of the disposal sites are still characterized by high percentages of fines. Furthermore, fines percentages increase by more than 10% in the area directly south of the MV2 reclamation (indicated with the cross-hatched area). Considerable spatial gradients persist in this area.
High fines percentages, of up to 30%, are observed in and around the sand mining pit approximately 10 km offshore from the MV2 reclamation. Though this area was not sampled in the pre-MV2 period, historical data does not show these high fines percentages in this area, and are not expected on the basis of our analysis below. Therefore, these high fines percentages likely are a recent development.
In the most southern part of the sampling domain, there are several stations with persistent high percentages of fines. These stations are located either in the troughs of tidal ridges, or in former tidal channels where fines have accumulated after closure of the estuaries. In the following section, we will use a conceptual framework to interpret the observed trends.

Conceptual framework
To interpret the presented spatial distribution of fines, we propose a conceptual framework. This can be used to analyse the natural distribution of fines in a coastal zone and the effects of human interventions thereupon. It consists of three components: 1. Source of fines 2. Transport pathways of fines 3. Accumulation potential for fines These components are schematically drawn in Fig. 7. The presence or absence of fines in the seabed depends on all three factors. (1) Multiple local sources of fines exist within a coastal zone. From these sources, fine sediment can follow different transport pathways. (2) The exact pathways are not meaningful, as these are erratic owing to the erratic driving forces. Therefore, we construct envelopes around a large number of potential pathways, representing mean dispersion patterns.
(3) If conditions for accumulation are favourable in an area, fines can deposit and accumulate in and on the bed. Such areas are referred to as potential accumulation areas.
For fines to be present in the seabed, all conditions have to be met. For example, an area can be very calm fulfilling the local conditions of an accumulation area, but if there is no pathway from a source to that area, no fines will accumulate (Fig. 7). Oppositely, there can be a large supply of fines, but local accumulation potential determines whether this yields a temporary (Fig. 7 -II) or a permanent deposit of fines ( Fig. 7 -IV). Furthermore, permanent deposits can interrupt the transport pathway of fines.
The natural spatial distribution of fines only depends on the undisturbed interplay between these three components. Human activities can modify the sources, transport pathways and/or accumulation potential.   Eisma (1968) and van Alphen (1987). The period during which samples were collected in each study is mentioned above the corresponding panel.
The three components indicated in Fig. 7 are elaborated upon below. Sources of fines can have a natural or anthropogenic origin. Natural sources of fines are erosion of geological layers (Adriaens et al., 2018), coastal (cliff) erosion (Eisma, 1981), riverine input (Salomons and Eysink, 1981) and input from other seas or oceans (McManus and Prandle, 1997). Furthermore, fines which were buried within the seabed during calm conditions can be remobilized during storms (van Kessel et al., 2011). Hence, on an annual timescale parts of the seabed may alternately be an accumulation zone and a source.
Anthropogenic sources of fines include disposal of sediment from maintenance dredging (Fettweis et al., 2009) and sand mining overflow (Nichols et al., 1990;Spearman et al., 2011). These sources are represented by a mass flux (ɸ fines ):   Here, T denotes a timescale and m fines denotes the dry mass of fines. Sediment transport pathways in shallow coastal areas are the consequence of a myriad of combinations of barotropic and baroclinic processes governed by tide, wind, waves, and density-driven flows (Otto et al., 1990). To exactly define these pathways, one would either need a high-resolution sediment dataset (McLaren and Powys, 1991) or a complex numerical model (Kim and Lim, 2009). Instead of using the exact pathways, we propose to use the envelope of the pathways. These envelopes are similar to the Depth of Transport concept introduced by Valiente et al. (2019).
The accumulation potential is defined as a parameter reflecting the interaction between the local bathymetry (i.e. the local geomorphological features, such as bedforms) and the prevailing hydrodynamic conditions. It is defined as a binary parameter, which is either high or low. If the accumulation potential is low, fines may be transported into an area, but it is unlikely that they can accumulate on/in the bed. If the accumulation potential is high, sediment deposits on the bed and remains there (Fig. 7).
As hydrodynamic conditions in a coastal zone are determined by tidal currents and waves, they are strongly time-dependent. To make the accumulation potential independent of time, we define a representative parameter for the prevailing hydrodynamic conditions. van Kessel et al. (2011) hypothesize that fines are remobilized from the seabed during storms. Afterwards, during calm conditions, these fines will deposit and are buried in the seabed again. However, these can only accumulate if conditions are calm enough. Hence, in our analyses of post-storm season data, it is not the energetic conditions which determine the accumulation potential, but rather whether the calm conditions are calm enough.
If the wave height over water depth ratio is relatively small, calm conditions are best represented by the magnitude of the tidal current. We assume this is valid for the majority of the DCZ. Calm conditions can then be quantified by selecting an appropriate percentile of the yearly current-only bed shear stress. We use the 90th percentile of the current-only bed shear stress, as it provides a proxy for the maximum tide-induced bed shear stresses during a spring-neap cycle, and denote it with . The current-only bed shear stress has been defined in Section 2.3.
The interaction between bathymetry and hydrodynamics manifests itself on a variety of scales. The larger, regional scale is characterized by geomorphological features such as tidal ridges, sand waves, navigation channels and large sand mining pits, while the smallest scale is determined by the dimensions of ripples. This smallest scale cannot be resolved in any field dataset or hydrodynamic model output. However, the scale-dependency of this interaction is crucial for the local behaviour of fines, and should be explicitly included in the assessment of the accumulation potential.
For the (small-scale) bathymetric contribution to the accumulation potential, we use the DEV parameter proposed by De Reu et al. (2013). DEV expresses the bathymetric level of a central point (z b ) relative to the bathymetry in its direct vicinity. DEV is based on the Bathymetric Position Index (BPI) (Iampietro et al., 2005;Verfaillie et al., 2006;Wilson et al., 2007). BPI measures the difference between the elevation of a point (z b ) and the average elevation () in a circle with radius (R) around it (Wilson and Gallant, 2000): n R indicates the number of observations within the circle. DEV is a modification of BPI, and uses the standard deviation (σ z ) of the bathymetry within radius R to normalize the BPI: The DEV R parameter depends directly on the selected spatial scale, as both σ z and depend on the radius R. Hence, from a proper choice of R, bathymetry-induced sub grid effects, e.g. in bed shear stresses, can be captured. A positive DEV R value means the bed level of that point is relatively high with respect to its surroundings. It therefore experiences larger bed shear stresses than its surroundings. A negative value means a relatively low bed level, with relatively low bed shear stresses.
The bathymetric and hydrodynamic contributions to the accumulation potential are classified through low or high accumulation potential areas, as sketched in Fig. 8. For a bed shear stress lower than a critical value (), accumulation of fines is always expected. With increasing current-only bed shear stress, the relative elevation of an area becomes important. Above a certain bed shear stress (), accumulation is no longer possible as currents are too strong.

Application of framework to Rotterdam Area
In this section, we apply the framework to the area around the River Rhine outflow, to study how human activities in this area have influenced the spatial distribution of fines in the seabed (Fig. 6). We specifically focus on the impact of the construction of MV2. We consider all We start by defining (1) the sources of fines in the DCZ, including the fines that enter the study area from outside the domain. This links to (2) the transport pathway envelopes. Finally, (3) the accumulation potential in the area is determined. Both the pre-MV2 and post-MV2 period are considered.
The main natural sources of fines in and around the DCZ are erosion of fines from geological layers in the Belgian Coastal Zone (e.g. Adriaens et al., 2018), fines entering the North Sea from the Atlantic Ocean (Eisma, 1981;McManus and Prandle, 1997) and a small contribution from the Rhine, Meuse and Scheldt (Laane et al., 1999;Salomons and Eysink, 1981). Hence, the major sources of fines are not located in the DCZ, but south of it. These fines are transported alongshore in a residual north-easterly direction, with the yearly transport flux estimated at 22 ± 10 MT/year (van der Hout et al., 2015).
The transport of fines along the Dutch coast is predominantly determined by tidal currents and their modification by the Rhine ROFI. The freshwater discharge leads to a salinity difference in cross-shore direction, e.g. Souza and Simpson (1996); Pietrzak et al. (2011), inducing a net shoreward near-bed transport of fines. This net transport results from cross-shore density gradients and tidal straining (de Boer et al., 2009;van der Hout et al., 2015). Storms may occasionally transport fines in offshore direction, but this sediment is returned onshore by the previously described processes (Flores et al., 2017).
Once the fines enter the DCZ, about 10% deposits in the estuaries in the southwestern delta. The northward flux at the Rotterdam area is still in the order of 20 ± 10 MT/year (Eisma, 1981). A considerable amount of fines deposits and accumulates in the entrance channels and harbor basins of the Port of Rotterdam. These are dredged regularly and disposed on the disposal sites Verdiepte Loswal and Loswal Noordwest (Fig. 6). From there, the alongshore transport mostly continues in northeasterly direction. From 2000 to 2016, an average of 0.6 and 2.2 MT fines were disposed yearly at Loswal Noordwest and Verdiepte Loswal, respectively (Hendriks and Schuurman, 2017). Though no new fines are introduced into the DCZ, we include these sites in our analyses as they are the main disposal sites in the DCZ and can buffer fines permanently or temporarily.
A major anthropogenic source of fines originates from the overflow during sand mining, thus located at the MV2 sand mining pit. In 2009 and 2010, approximately 2 MT fines were yearly released this way (van Kessel et al., 2011). As this is the only major additional source of fines in the DCZ during the construction of MV2, we investigate whether it has contributed to the post-MV2 distribution of fines. The magnitude of this source strongly decreased after 2010, as extracted sand volumes strongly decreased in subsequent years.
In Fig. 9a and b, we illustrate the location of these sources and their assumed transport pathway envelopes, for the pre-MV2 and the post-MV2 situation, respectively. Though the MV2 sand mining is a temporary source, it is included in Fig. 9b as it may have affected post-MV2 fines percentages. For the yearly natural transport flux, we only draw its envelope, as its major sources lie outside the Rotterdam area. The arrow indicates the residual transport direction.
Several assumptions have been made to establish the envelope of the sediment transport pathways for the anthropogenic sources. Assuming that the majority of fine sediment transport in the DCZ takes place within the Rhine ROFI, the offshore boundary for the transport pathways is determined by the offshore limit of the ROFI. This is assumed to be at the 31 PSU mean surface salinity contour, which was assessed from the hydrodynamic model output. On the nearshore side the transport envelopes areultimatelybounded by the land boundary. Within this area, the transport pathway envelope is expected to develop along the mean salinity contours. The alongshore (north-south) boundaries are determined by the north-easterly residual current along the Dutch coast, with a magnitude of 0.10-0.15 m/s (Simpson, 1997;van der Giessen et al., 1990). As tidal currents are the main alongshore transporting agent, the southern boundary is set at one tidal excursion south of a local source (approximately 10 km). To the north, the envelope extends in north-easterly direction with time.
No significant difference was found between the computed mean salinity contours for the pre-MV2 and post-MV2 periods (Fig. 9). No substantial differences between the pathway envelope of the alongshore flux between the pre-MV2 and post-MV2 period are therefore expected. The magnitude of the natural alongshore flux is not substantially affected by the construction of MV2, nor have flow patterns on the scale of the DCZ changed considerably (Cronin and Blaas, 2015). Therefore, the magnitude of the alongshore flux entering the study domain is not substantially affected by the construction of MV2.
The transport pathway envelopes for the anthropogenic sources (i.e., the disposal sites and the sand mining pit) overlap substantially. Fines from the MV2 sand mining pit which are transported northward, either end up in the pathway envelope of Loswal Noordwest or of the Verdiepte Loswal. It is therefore difficult to discriminate between the effects of different human interventions in the DCZ, as fines from different sources can be transported to the same location. However, in combination with the local accumulation potential, the envelopes may give a good impression on where to expect an increase or decrease of fines percentages in the seabed. The next step is to quantify the conceptual accumulation potential diagram (Fig. 8) for the Rotterdam area. The relationship between DEV and is quantified by applying a logistic regression to a selection of the sediment dataset. This regression predicts a binary response, i.e. the presence or absence of fines in a sediment sample (cf. Section 3.1), to a set of explanatory variables. We only discriminate between presence or absence of fines, as their amount is strongly determined by the non-uniform supply of fines (Fig. 9). All samples taken in the post-MV2 period within 20 km from the shore are included in this regression, as few fines are encountered beyond (Section 3.2).
DEV R and are used as the explanatory variables in this regression. They are calculated for every selected sample, using the highresolution bathymetric and hydrodynamic model data (Section 2.3). The DEV R value of each sampling point is assessed for a radius R of 1000 m, i.e., DEV 1000 . Radiuses of 250, 500 and 1500 m were also tested. However, the 1000 m radius was used in the analysis below as it can represent the larger geomorphological features in the area (i.e. bedforms, navigation channels and sand mining pits), while still contrasting the bathymetric differences adequately. has been calculated for the representative post-MV2 year 2012. Furthermore, the logistic regression is applied between a of 0.15 Pa and a of 0.60 Pa. The results of this regression are shown in Table 4.  Fig. 10, distinguishing between the classes introduced in Fig. 3. Fig. 10 shows that the bed shear stress itself is not discriminative for predicting the presence of fines in the range between 0.15 and 0.6 Pa. DEV 1000 greatly improves the predicted accumulation potential. This likely implies that the presence of fines is strongly influenced by subgrid effects, i.e. morphological features smaller than the model resolution of 500 × 500 m 2 . Thus, the power of the DEV parameter lies in resolving these sub-grid effects, quantifying local bed features on a 25 × 25 m 2 resolution.
The accumulation potential classification of Fig. 10 is applicable for supply-limited systems. In such systems, seabed topography and bed shear stresses mainly determine the distribution of fines. However, the DCZ cannot be regarded entirely as a supply-limited system. There are also areas where there is either no supply of fines, or where this supply is abundant. When there is no supply, the majority of sampling points will be randomly distributed in the 0% category. When supply is abundant, sub-grid hydrodynamic conditions are not discriminative, as eroded fines are replaced continuously (Fig. 7 -II). Then, sampling data will be randomly distributed in the more than 10% fines class.
The accumulation potential in the study domain is visualized in Fig.  11, by combining the hydrodynamic and bathymetric data with the accumulation potential classification. These maps show alternating areas with high and low accumulation potential, both for the pre-and post-MV2 period. The average fines percentages for both periods (Fig. 6) are also plotted in these maps. Fig. 11 shows that the accumulation potential parameter provides a good explanation of the variability in fines percentages measured in the Rotterdam area. Stations where no fines are found, are generally located in low accumulation potential areas. Stations containing fines are generally located in high accumulation potential areas. The accumulation potential classification predicts these sites correctly for 65% of the data. Below, we discuss the areas indicated and numbered in Fig. 11.
First we combine our analysis of accumulation potential with the transport pathway envelopes of the major anthropogenic sources, indicated with a light grey hatch in Fig. 11a and b. Where the Verdiepte Loswal (1) envelope overlaps areas with high accumulation potential, high fines percentages are found for both the pre-MV2 and post-MV2 period. This also holds for Loswal Noordwest (2), although fines percentages are lower. Where the transport pathway envelope of the MV2 sand mining pit (4) overlaps with the envelopes of the two disposal sites, a relatively large increase in fines percentage from the pre-MV2 to post-MV2 period (right panel Figs. 6 and 11b) is observed. The increase must be the result of the MV2 sand mining activities, as the amount of   (3) the fines percentages are mostly zero for both periods. This can be explained by the disposal strategy: fines were disposed here until 1996, but later only sand. As a result, Loswal Noord (3) lies relatively high and the accumulation potential is therefore low.
In the post-MV2 period (Fig. 11b), the MV2 sand mining pit (4) forms a major accumulation. The high accumulation potential in this pit is mainly due to its 20 m larger depth, though the 90th percentile of the bed shear stress is also smaller than during the pre-MV2 period ( Fig. 12a and b). This shows how sand mining influences the distribution of fines in two ways. Sand mining itself acts as a source of fines due to overflow from the hopper. During post-dredging conditions, the resulting pit becomes a sink for fines.
Directly south of MV2, accumulation potential increased, because of a local decrease in bed shear stress. This is caused by a change in the tidal flow pattern, as MV2 protrudes further into the North Sea, deflecting the tidal flow. This leads to a decrease in tidal current magnitude directly to the north and south of MV2, but also to tidal flow contraction directly west of MV2, where the tidal current magnitude increased (Fig. 12b). Indeed, the data show an increase in fines percentage over the pre-to post-MV2 period directly south and north of MV2 (Figs. 6 and 11).
Beyond 15-20 km offshore, the spatial distribution of accumulation potential is mainly determined by large-scale bedforms, leading to alternating accumulation potential patterns. A substantial part of the offshore area is thus classified as a potential accumulation zone. Nevertheless, fines are virtually absent here (Section 3.2). This is explained by a lacking supply of fines, reflecting either the absence of sources or because transport pathways remain closer to shore.
Another zone with high fines percentages, earlier addressed by both Eisma (1968) andvan Alphen (1987), is found within 2-3 km from the shore along the Holland coast (area 5 in Fig. 11). In this area, computed tide-induced bed shear stresses are so low, that it is classified as a potential accumulation zone (Fig. 11), in spite of local high DEV-values. At low , data will mostly fall in the high accumulation potential range (see Figs. 8 and 10).

Accumulation potential and the effect of human interventions
The proposed conceptual framework can be used to assess the effect of different human interventions on the spatial distribution of fines in the seabed. These interventions likely affect the local accumulation potential, which can be illustrated by the MV2 sand mining pit and MV2 reclamation. Due to the reclamation, the area directly south of MV2 became an accumulation zone for fines (Fig. 11), while the sand mining pit also is an accumulation zone.
These interventions affected local accumulation potential in different ways, as visualized in the accumulation potential diagrams of Fig.  13. DEV 1000 values decreased strongly in the sand mining pit because of its 20 m larger depth, accompanied by a small decrease in bed shear stress. Thus, the accumulation potential in the sand mining pit increased (Fig. 13a). South of MV2, DEV 1000 values remain more or less constant, implying small bathymetrical changes. While pre-MV2 data indicate accumulation potential was already high in parts of the area, accumulation conditions became even more favourable because of the decreased bed shear stress in the post-MV2 period (Fig. 13b). Earlier human interventions in the DCZ have also led to a local increase in accumulation potential. Closure of estuaries in the southwestern delta led to a strong decrease in tidal currents, resulting in accumulation of fines in the former tidal channels (van Alphen, 1987). Such an area where bed shear stresses became very low, lies directly west of the Haringvliet mouth, area (6) in Fig. 11. Due to the closure of the Haringvliet estuary and the construction of Maasvlakte 1, a sheltered area was created (Elias et al., 2017). Although our sediment dataset does not provide information, multiple studies have shown that accumulation of fines in this area is significant (Elias et al., 2017;Piekhaar and Kort, 1983;van Alphen, 1987;van Heteren, 2002). The computed accumulation potential is indeed high both in the Haringvliet mouth and in the former tidal channels of the southwestern delta (Fig. 11). This confirms that our conceptual framework is capable of predicting the effect of closing tidal inlets on the distribution of fines in the seabed, and is consistent with previous studies.
The accumulation potential concept can also be used to quantify sub-grid effects in assessing the fines percentage in the seabed from numerical model simulations. The computational grid size in numerical models always exceeds the spatial dimensions of (small) bed forms. We have shown that small-scale elevation differences of the seabed strongly influence the presence of fines in the bed, and therefore explain the large spatial variability in observed fines percentage. Hence, if bathymetrical information is available at scales smaller than the computational grid size, the accumulation potential can be used to obtain a first order estimate of the variability of fines percentage over the computational cells. This may be relevant, for instance, for analyses of the ecological functioning of the system.

Variations in hydrodynamic forcing and response of the seabed
The spatial distribution of fines in the seabed results from the hydrodynamic forcing and the response of the seabed thereupon. The forcing is driven by tidal currents and waves. The seabed response consists of two major parts. Fines which previously accumulated are remobilized by wave action during storms and then transported by tidal currents. Afterwards, during calm conditions, these fines will deposit and are buried in the seabed again. Multiple timescales are associated with the forcing and seabed response, which are crucial to our analysis.
Logically, the fines percentage in the seabed decreases during storms because of remobilization. During subsequent calm conditions, fines percentages may increase again due to deposition and burial. Since storms occur frequently in the DCZ during winter, a seasonal variability in fines percentage is likely. Therefore, sampling was carried out after winter (storm season) to allow for a proper assessment of anthropogenic effects, undisturbed by seasonal variations (Section 2.2).
The accumulation of fines during calm conditions depends on the available time for deposition and burial, which should be sufficiently long. Hence, it depends on the ratio between the timescale for deposition and burial, and the timescale of hydrodynamic forcing variations. This ratio should be reflected by the hydrodynamic contribution to the accumulation potential. The first timescale is subject of ongoing research, but is likely in the order of several days to a week. The second depends on the dominant contribution to the hydrodynamic forcing, as currents and waves vary on different timescales themselves.
Tidal currents vary on diurnal and fortnightly scales, while waves occur more erratic. If the wave height over water depth ratio is small Fig. 11. Accumulation potential map in the Rotterdam area for the pre-MV2 (a) and post-MV2 (b) period. The accumulation potential is either high (green) or low (white), based on the accumulation potential classification. The difference map (c) shows changes in accumulation potential from pre-MV2 to post-MV2 period. We discuss the numbered sites with black contour lines below.  apart from storms, calm conditions are governed by the magnitude of the tidal current. Prolonged calm conditions then relate to the maximum currents that occur during a spring-neap cycle, i.e. the fortnightly timescale. Therefore, we have schematised the hydrodynamic contribution to the accumulation potential through the 90th percentile of the current-only bed shear stress, as it is representative for these maximum currents.
This approach is valid for the majority of the sampling locations, as they lie at water depths larger than 10 m. At such water depths in the DCZ, the seabed is only subject to large wave-induced shear stress during storms. However, with decreasing water depth, waves become more important and will also define calm conditions. Therefore, the accumulation potential is not accurately predicted in the shallow nearshore, i.e. the breaker zone and around (e.g. parts of area 5 in Fig. 11). The differences between observed fines percentage and its prediction by the accumulation potential are thus partly explained by how calm conditions are schematised.

Application of conceptual framework to other coastal shelf seas
The conceptual framework was developed for the southern North Sea, but its concept is generic and it can be used to analyse the dynamics of fines in other coastal shelf seas, provided they have a sandy substrate. The framework allows to study the effect of human activities on these dynamics in a structured way. This enables engineers to effectively assess and thus mitigate the impact of these activities. Here, we provide several examples of how the conceptual framework could be utilized.
In the Seine estuary, France, material dredged from harbours is disposed at sites off the coast of Le Havre. Marmin et al. (2014) demonstrated that relocating this disposal site for, mostly fine, sediment is constrained by both economic and natural restrictions. Apart from the local effects on biota at the dumping site, the far-field effects need to be studied as well. Marmin et al. (2014) mention "these effects depend on a variety of environmental conditions (Essink, 1999), but differ greatly from one site to the other, thus general conclusions are difficult to draw". By estimating the transport pathways from a site and determining accumulation potential in the area, the influence of this newly introduced source can be assessed more specifically. Furthermore, the conceptual framework can be used to assess the differences between a concentrated disposal strategy at one location versus smaller disposal locations across a larger area.
Human interventions in the shelf seas surrounding Australia are subject to strict regulations, as they pose a threat to the present coral reefs (Erftemeijer et al., 2012). In terms of the conceptual framework, accumulation potential is high between and on the coral reefs Jones et al., 2015). Thus, to avoid smothering of the coral reefs, the transport pathway envelope from disposal/dredging sites has to be established. This can form the basis for revised project design or the implementation of mitigating measures.
The North Sea is one of the world's most actively studied shelf seas. However, in many shelf seas the seabed has only been sparsely sampled, such as the Andaman Sea. For example, Kamp-Nielsen et al. (2002) and Feldens et al. (2012) observed patches of fine sediment off the coast of Thailand from local high-resolution surveys. The presence of these patches depended on small-scale topography, with sharp boundaries between a patch and the sandy environment. Combining hydrodynamic model results with bathymetric data provides a valuable first estimate of where to expect fines. This can then be incorporated in the design of seabed sampling campaigns. The framework can be utilized to trace back the fine sediment to its respective sources, once the pathways are established.

Conclusions
In this paper, we studied the processes and conditions which determine the distribution of fines in the surficial seabed of the Dutch coastal zone (DCZ). A new dataset was analysed to determine this spatial distribution and then compared with previous datasets. The large extent and high spatiotemporal resolution of the new dataset enables to study fine sediment dynamics in the North Sea at scales smaller than before.
At mega-scale, the spatial distribution of fines in the DCZ is generally in agreement with previous work by Eisma (1968) and van Alphen (1987). Virtually no fines are found beyond 20 km offshore. Further nearshore, variability in fines percentage is found on both smaller (tens of metres) and medium (kilometres) scales. Locally, fines percentages exceed 10-20%. Highest percentages are found within 2 to 3 km from the shore, north and south of the River Rhine outflow, and in the former tidal channels in front of closed estuaries. Large-scale human interventions invoked local changes of the fines percentage in the seabed in the order of 10% and above.
To analyse the large-scale distribution of fines in the seabed and enable quantitative prediction of the effect of large-scale human interventions, a conceptual framework was developed. This framework consists of three components: (1) sources of fines; (2) transport pathways from these sources; and (3) accumulation potential. It was shown that the large-scale distribution of fines in the DCZ is mainly determined by the first two components, whereas accumulation potential mainly influences the local distribution. Differences in fines distribution in response to the construction of MV2, a seaward extension of the Port of Rotterdam, were caused by an additional source of fines released from overflowing during sand mining, and by local changes in accumulation potentialmost notably in the deep mining pit and the sheltered zone south of MV2.
The new framework enables the assessment of individual human interventions in terms of source, pathway and/or accumulation effects. It further allows for the assessment of cumulative effects due to multiple interventions (and their interactions) in one area. Such analyses can establish a sound basis for Environmental Impact Assessments and may form a starting point for successive analyses on ecological effects. In this way, the framework developed here can help engineers and policy makers to assess how human interventions affect the ecosystems like the North Sea and to limit or mitigate their environmental impact.

Declaration of competing interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.