A Two Decadal (1993–2012) Numerical Assessment of Sediment Dynamics in the Northern Gulf of Mexico

Abstract

We adapted the coupled ocean-sediment transport model to the northern Gulf of Mexico to examine sediment dynamics on seasonal-to-decadal time scales as well as its response to decreased fluvial inputs from the Mississippi-Atchafalaya River. Sediment transport on the shelf exhibited contrasting conditions in a year, with strong westward transport in spring, fall, and winter, and relatively weak eastward transport in summer. Sedimentation rate varied from almost zero on the open shelf to more than 10 cm/year near river mouths. A phase shift in river discharge was detected in 1999 and was associated with the El Niño-Southern Oscillation (ENSO) event, after which, water and sediment fluxes decreased by ~20% and ~40%, respectively. Two sensitivity tests were carried out to examine the response of sediment dynamics to high and low river discharge, respectively. With a decreased fluvial supply, sediment flux and sedimentation rate were largely reduced in areas proximal to the deltas, which might accelerate the land loss in down-coast bays and estuaries. The results of two sensitivity tests indicated the decreased river discharge would largely affect sediment balance in waters around the delta. The impact from decreased fluvial input was minimum on the sandy shoals ~100 km west of the Mississippi Delta, where deposition of fluvial sediments was highly affected by winds.

1. Introduction

The Mississippi-Atchafalaya River system has the third largest drainage basin (3.3 × 10⁶ km²) and seventh largest freshwater discharge (380 km³/year) in the world [1,2,3]. About two thirds of the sediments and water are delivered by the Mississippi River and the rest are diverted to the Atchafalaya River [4,5]. Over the past decades, especially after the 1950s, sediment flux from the Mississippi-Atchafalaya River has decreased dramatically [6,7,8,9]. The deficit of sediment supply, together with the eustatic sea level rise, results in severe coastal erosion and land loss [10,11,12]. The average rate of land loss was 88 km²/year from 1956 to 2000, and an additional loss of 1329 km² is projected by 2050 [13]. Climate change within the Mississippi River watershed has been identified as a significant factor controlling long-term variations of river discharge. Wavelet analysis of the North America annual freshwater discharge indicated a four- to eight-year oscillation, which is correlated with the El Niño-Southern Oscillation event (ENSO hereafter) [14,15,16,17]. A longer temporal scale variation (25-year) is associated with the bi-decadal precipitation oscillation related to the North Atlantic Oscillation (NAO) and Pacific Decadal Oscillation (PDO) [18,19,20]. The most dominant change point of the Mississippi River water discharge was detected around 1970, featured by the eight- to sixteen-year and three- to six-year modes [21]. Anthropogenic activities have been widely discussed in previous studies as another factor influencing river discharge [22]. In the 1950s, dam construction on the Missouri River resulted in substantial sediment flux reduction [7]. The estimated loss of fluvial sediment load at the Mississippi River mouth was ~225 Mt/year in the period of 1950 to 1975 [15].

Hydrodynamics in the northern Gulf of Mexico (nGoM hereafter) exhibit contrasting patterns over a year due to the shift in the direction of prevailing winds [23,24,25,26]. During non-summer months, the westward Louisiana Coastal Current (LCC) prevails because of strong easterly wind. In summer months, intensive westerly winds cease the LCC and the currents reverse to eastward [27]. Waves in the nGoM are introduced by both local winds and remote swell propagation [28]. For tidal schemes, K1, O1 and M2 are the most dominant constituents and the tidal currents maximized to ~9 cm/s near Atchafalaya Bay over the Louisiana–Texas shelf [29]. The maximum tidal range is about 0.6 m [30,31]. Due to the high fluvial discharge and relatively low-energy environment, initial deposition of the fluvial sediments usually happens <30 km off the river mouth [32,33,34]. During episodic events such as hurricanes and winter storms, strong hydrodynamics induced by energetic winds can transport fluvial sediments further offshore [35,36,37,38,39,40,41]. Two depocenters in the nGoM with a deposition >1 cm per year have been identified: one is around the bird-foot delta in the Louisiana Bight and the other is in the Atchafalaya Bay [42]. Radionuclide chronologies of sediment cores around the bird-foot delta indicated a decreasing deposition rate as water became deeper [43,44,45,46]. For the Atchafalaya shelf; however, high deposition rates were found 10–12 km offshore on the clinoform foreset [43]. A possible explanation of this fast deposition is that the fluid mud escaped from the delta topset [47,48,49,50,51,52,53]. Although nourished by the largest river system in North America, the Mississippi Delta and adjacent coast is still suffering from severe erosion [54]. For example, the Barataria Bay, which is adjacent to the Mississippi main channel to the west, has been experiencing substantial land loss (16.9 km²/year) and barrier island retreat over the past decades [55,56]. Many efforts have been made for coastal restoration purposes, and most noticeable examples are the sediment diversion via Davis Pond diversion and sediment emplacement over the barrier islands [57,58,59,60,61]. Nevertheless, the diverted sediments and sand materials dredged from the inner shelf are still insufficient to balance the land loss in the bay.

Submarine shoals over the western Louisiana shelf (e.g., Tiger/Trinity and Ship Shoals; Figure 1) are reworked prograded deltaic headlands formed during low sea level stand [62]. Due to relatively high sandy content and little muddy overburden, the transgressive shoals are treated as potential sand sources for coastal restoration and beach nourishment [59,63,64]. Although the total sand volume of these deposits is massive, recent surveys show that the total dredgeable sands are highly restricted by oil infrastructures, environmental concerns, and cultural resources [65,66]. Given the importance of sandy shoals in coastline protection, understanding long-term sediment dynamics over these shoals and its interaction with hydrodynamics and rivers are essential.

Existing radionuclide studies provide valuable information of deposition rate on seasonal to decadal scales, and the difference between short-term and long-term deposition rates implies the relative importance of episodic events [44,45,46,67,68]. Nevertheless, radionuclide chronology cannot quantitatively evaluate hydrodynamics’ impact on sedimentation, and physical reworking introduced by waves and currents might compromise the temporal resolution and accuracy of such measurements. Moreover, radionuclide chronology can hardly detect the variation of deposition rate in areas with high deposition rate on a decadal scale (e.g., around the Mississippi Delta) due to its low temporal resolution (e.g., ²¹⁰Pb, half-life 22.4 years) or short temporal scale (e.g., ²³⁴Th, half-life 21.4 days; ⁷Be, half-life 53.2 days). There is a substantial knowledge gap of shelf deposition’s response to decreased fluvial sediment flux over the past decades. As an alternative method, numerical model has been widely applied to the nGoM to investigate hydrodynamics and sediment dynamics on different temporal scales [30,36,39,40]. In this study, we used the Coupled Ocean–Atmosphere–Wave–Sediment Transport Modeling system (COAWST) to investigate fluvial sediment dynamics over the continental shelf on seasonal to decadal scales [69,70]. Compared with existing modeling studies in this region, this study is the first effort to investigate shelf sediment dynamics up to a decadal time scale. The objective of this study is to understand: i) the seasonal sedimentation patterns on the continental shelf; ii) the impacts of rivers, winds, and waves on sediment transport; and iii) the bay-shelf sediment exchange and sediment dynamics over submarine shoals.

2. Model Setup

The COAWST model (version 3.2) consists of three state-of-the-art numerical models: The Regional Ocean Modeling System (ROMS; [71,72]) for ocean hydrodynamics; Simulating Waves Nearshore (SWAN; [73]) for waves; and Weather Research and Forecasting Model (WRF-ARW; [74]) for atmospheric simulation. The Community Sediment Transport Modeling System (CSTMS; [70]) is incorporated into the ocean model (ROMS) to simulate sediment transport and deposition. The Model Coupling Toolkit (MCT; [75]) is used for information exchange among different models. In this study, we disabled the atmosphere (WRF) coupling to focus on interactions among wave, ocean, and sediment transport over a 20-year period (1993–2012), which was determined by the availability of the model inputs (Hybrid Coordinate Ocean Model (HYCOM)). Details of our model configuration are described next.

2.1. Ocean and Sediment Models

ROMS (svn 797) is a three-dimensional, free surface, terrain following model that solves Reynolds-Averaged Navier–Stokes (RANS) equations based on the hydrostatic and Boussinesq assumptions [76,77]. We used a “two-step” offline nesting method to reduce the computational cost of the 20-year coupled ocean–wave–sediment simulation. First, we performed a two-way-coupled simulation (wave–ocean) on the Gulf of Mexico (GoM) domain, with 36 weighted vertical layers at a 5 km horizontal resolution. We then utilized the GoM model results as the boundary condition to drive a higher resolution domain covering the nGoM at a 1 km horizontal resolution (Figure 1 and Figure 2; meshes of model domains see Supplementary Materials). Compared to previous sediment transport models in this region [35,37,40], our nesting mesh resolution was high enough to resolve physical and sediment transport processes over the shelf and the structured grid made long-term simulation (20 years) applicable. For the GoM domain, initial conditions of current velocity, sea level, temperature, and salinity were interpolated from the 1/12° data assimilated Hybrid Coordinate Ocean Model (HYCOM/NCODA, GLBu0.08/expt_19.0 and expt_19.1; [78]). The barotropic velocity boundary condition was prescribed following Flather [79]. The baroclinic velocity, temperature, and salinity were specified using the Orlanski-type radiation boundary condition [80]. We extracted the Oregon State University Tidal Inversion Software (OTIS; [81]) regional tidal solution and interpolated it on the model grid as tidal forcing. The 6-hourly, 38 km horizontal resolution Climate Forecast System Reanalysis (CFSR; [82,83]; http://cfs.ncep.noaa.gov) was utilized as meteorological momentum and buoyancy forcing due to its high quality. Monthly average freshwater and suspended sediment inputs from 39 rivers debouching into the GoM were retrieved from United States Geological Survey (USGS) Water Data for the Nation website (http://nwis.waterdata.usgs.gov) and applied as boundary condition. The stations selected for river inputs were the most downstream sites with consecutive available data. Sediment bedload from rivers was not considered in our simulation. The mesh bathymetry was interpolated and smoothed from ETOPO1 dataset [84]. We employed the Mellor-Yamada level-2.5 closure scheme [85] to estimate vertical turbulent mixing. We chose the SSW_BBL module [86,87] for bottom boundary layer parameterization, which calculates both wave- and current-induced bottom shear stress for momentum and sediment resuspension. Model outputs were saved every day for analysis. The time steps for the GoM and nGoM domains were specified as 300 and 120 s, respectively.

The sediment model (CSTMS) integrates several modules to simulate sediment transport, stratigraphy and geomorphology [70]. Suspended sediment transport in the water column is calculated by solving advection-diffusion equation. The other two additional sediment source/sink terms are sediment resuspension from the seabed and sediment vertical settling. As bottom shear stress calculated in SSW_BBL exceeds critical shear stress (details see

Table 1. Sediment characteristics parameterization.

Sediment Type	Grain Diameter (mm)	Settling Velocity (mm/s)	Critical Shear Stress (Pa)	Erosion Rate (10−4 kg/m2/s)
Mud_01(Mississippi River)	0.004	0.1	0.10	5
Mud_02(Mississippi River)	0.03	0.1	0.16	5
Mud_03(Atchafalaya River)	0.004	0.1	0.10	5
Mud_04(Atchafalaya River)	0.03	0.1	0.16	5
Sand_01(Mississippi River)	0.0625	1	0.20	5
Sand_02(Atchafalaya River)	0.0625	1	0.20	5
Sand_03(seabed)	0.14	1	100.0	5

), pre-deposited sediment will be resuspended into the bottom water layer and the resuspension flux is estimated following Ariathurai and Arulanandan [88]. We defined four cohesive and two non-cohesive sediment classes for river inputs. Sediment concentration in the water column was initialized as zero. As this study focused on dynamics of riverine sediments, we prescribed one non-cohesive, resuspension-resistant class as shelf sediments with high critical shear stress (100 Pa) following Harris et al. [89]. To achieve the most reasonable sediment parameterization, we performed a series of sensitivity tests based on the studies by Xu et al. [35,37] and compared our simulation results with ²¹⁰Pb-derived deposition rate ([45,46,68]; core locations see Figure 1). In Table 1 we listed a summary of the sediment model parameterization used in this study, which reproduced the most reasonable deposition rates over the shelf. We prescribed four layers of sediment on the sea floor, each with a thickness of 1.0 m. Seabed erosion–deposition was based on non-cohesive parameterizations [70,90]. Due to the lack of suspended sediment observations at the open boundary, suspended sediment concentration (SSC) at the boundaries of GoM domain was set to zero and we applied the gradient boundary condition to avert unreal artificial sediment plumes along the boundaries. In this study we only simulated suspended sediment transport and bedload transport over the shelf was not considered.

2.2. Wave Model

The Simulating Waves Nearshore model (SWAN, version 41.01) was employed to simulate wind–wave generation and propagation. The SWAN model is based on a Eulerian formulation of the discrete spectral balance of action density that accounts for refractive propagation over bathymetry and current fields [73]. Other incorporated physical processes include wave-wave interaction, white-capping, bottom dissipation, and depth-induced wave breaking. The two SWAN model grids (GoM and nGoM) were the same as those of the ocean models (ROMS). The initial wave spectra were computed from the CFSR wind speed using the deep-water growth curve [91]. The breaker index (certain ratio between wave height and water depth at which wave breaks) and the proportionality coefficient of the dissipation rate were set to 0.73 and 1.0, respectively. The expression of Madsen et al. [92] was applied to estimate the bottom friction. The time step of wave simulations in GoM and nGoM domains were the same as corresponding ocean models (300 and 120 s).

2.3. Model Nesting and Coupling

For both ocean and wave models, we first performed a 20-year, two-way (ROMS-SWAN) coupled simulation on the GoM domain covering the period of 01/01/1993–12/31/2012 (step 1 in Figure 2). The ocean model sent sea surface current velocity, water level, and bathymetry to the wave model, and the wave model sent wave parameters (e.g., significant wave height, wave length, wave direction, etc.) back to the ocean model. The variable exchange interval was specified as 1 h. On completion of the GoM simulation, we interpolated model simulated physics (sea-level, velocity, salinity, temperature, and significant wave height, wave period, wave direction) to the nGoM domains (wave and ocean) as boundary conditions and performed the nGoM 20-year simulation (step 2 in Figure 2). Although COAWST supports a real-time coupling between the ocean and wave models, information interchange between models can greatly slow down the long-term simulations. To speed up the simulations (one benchmark run and two sensitivity tests) in the nGoM domain, we first ran the wave model independently and then utilized model outputs (wave direction, near bottom wave period, and bottom wave orbital velocity) at a 6 h interval to drive the ocean and sediment models.

3. Model Validation

We validated the performance of each nGoM model (wave, ocean, and sediment) using available in-situ measurements. For wave, we gathered monthly-averaged significant wave height at three buoy stations from the National Data Buoy Center (NDBC, https://www.ndbc.noaa.gov; locations see Figure 1). The model–data comparison revealed good agreement between simulated and observed significant wave height (R = 0.94; Figure 3a). To evaluate wave model’s performance on daily scale, we compared time series of observed and model simulated wave height at 42040 station in 2007 and the correlation coefficient was 0.92 (Figure 3b).

To evaluate the model’s skill of resolving long-term coastal hydrodynamics, we retrieved water level records from four NOAA tidal gauges at Calcasieu Pass (station ID: 8768094), Dauphin Island (8735180), Port Fourchon (8762075), and Pilots Station East (8760922). We applied a 36 h low pass filter to both simulated and observed time series. An example is shown in Figure 4, where model simulated sea level anomaly was compared against observations at four stations in 2008. The correlation coefficients were ≥0.80 at all four stations and the two surges brought by hurricanes Gustav and Ike, respectively, were captured. To further evaluate model simulated water level over a longer time period, a statistical assessment is shown in the form of a Taylor diagram (Figure 5), which presents the correlation coefficients, centered root mean square difference (RMSD), and normalized standard deviations of annual sea level anomaly time series [93]. Most correlation coefficients varied from 0.7 to 0.9, and the standard deviation ratios were less than 2. We interpolated simulated salinity to the observation sites at corresponding period and compared it against available measurements from the Southeast Area Monitoring and Assessment Program (SEAMAP; http://seamap.gsmfc.org, data are depth-averaged), which has 2145 data points covering the period from 1993 to 2012. The model-observation comparison in Figure 6 indicates that the ocean model is capable of reproducing the pattern of salinity distribution, with low salinity water embracing coastal Louisiana over the inner shelf and high salinity water further offshore.

For sediment model, we compared simulated sedimentation rate against published estimations based on radionuclide data (²¹⁰Pb; core locations see Figure 1). Our sediment model was capable of capturing the magnitude and variation of the sedimentation rates at these sites. The model–data correlation coefficient was 0.69 and the root mean square error (RMSE) was 0.78 (Figure 7). Moreover, we compared the simulated surface SSC against the map derived from Moderate Resolution Imaging Spectroradiometer (MODIS-aqua; Figure 8). Due to the presence of dense clouds, sun glint, and water vapor in the coastal region, it was a challenge to retrieve a set of consecutive satellite images with satisfying quality. We selected one cloud-free satellite image for each season from December 2009 to November 2010 and applied the nGoM SSC algorithm by Miller and McKee [94]. To highlight the turbid water on the shelf, the region where surface SSC < 1 mg/L was masked out. In spring, the Mississippi River sediment plume with high SSC (>100 mg/L) extended southwest to the 200 m isobath (Figure 8a,b). Turbid water from the Atchafalaya River dominated the entire Atchafalaya Bay and coastal water (water depth < 20 m). Westward sediment transport could be detected over the coastal Chenier Plain, where westward alongshore current was strong. In summer, both SSC and the spatial scale of sediment plume reduced dramatically due to calm weather and low fluvial discharge (Figure 8c,d). The difference between the model result and satellite image in summer were likely due to i) the application of atmospheric correction in the more oligotrophic summer shelf waters, and ii) surface water particle characteristics (e.g., smaller particle size; [95]) during summer. In fall and winter, the shapes of sediment plume were similar to that in spring. The westward transport along the Chenier Plain coast was even stronger due to intensified easterly winds (Figure 8e–h). Although such one-frame comparison might not fully capture the seasonality of sediment plume, which could be easily altered by fluvial discharge and wind condition, our model reproduced the spatial distribution pattern and the magnitude of surface SSC. The above model–data comparisons gave us the confidence that this ocean–wave–sediment model is capable of resolving the major features of the seasonal to decadal scale variability in hydro- and sediment dynamics.

4. Results

In this section we first present the seasonal variations of hydrodynamics and sediment dynamics, followed by an analysis of sedimentation pattern over the 20-year simulation period. We then present the results of two sensitivity tests to assess shelf deposition’s response to the high and low river input scenarios.

4.1. Seasonal Variations of Hydro- and Sediment Dynamics

Hydrodynamics in the nGoM is heavily influenced by the prevailing winds. We plotted wind fields measured at a buoy station east of the bird-foot delta, covering the period of 1995–2012 (station 42040, location see Figure 1, wind roses see Figure 9). The study region was dominated by strong southeast winds in spring (March, April, and May), south to southwest winds in summer (June, July, and August), east and northeast winds in fall (September, October, and November), and north and southwest winds in winter (December, January, and February). Among all seasons, westerly winds only prevailed in summer with relatively low intensity compared with easterly winds in other seasons. A 90th percentile of the westerly winds in summer was at 7.3 m/s, while a 90th percentile of the easterly winds in spring, fall, and winter were at 9.1, 9.9, and 9.9 m/s, respectively.

Figure 10 shows the depth-averaged current fields and bottom shear stress induced by currents averaged over each season. The current direction indicates that westward flow dominated the broad western Louisiana shelf in non-summer seasons (Figure 10a,c,d). Currents between 20 and 50 m isobaths shifted to eastward in summer due to the weak westerly winds (Figure 10b). Over the eastern Louisiana–Mississippi–Alabama shelf, current fields did not show strong seasonality and east- and northeastward currents prevailed throughout the year. Current-induced bottom shear stress (τ_current) maximized in fall and winter, and the highest τ_current was found to the southeast of the bird-foot delta, reaching more than 0.1 Pa. The spatial patterns of wave-induced bottom shear stress (τ_wave) were quite similar among different seasons, with the lowest intensity in summer (Figure 11). High τ_wave was found nearshore (water depth < 20 m), including sandy shoals over the inner shelf, around the bird-foot delta, and to the east of the Chandeleur islands. As a major driving force of resuspension, the maximum τ_wave was estimated above 0.2 Pa, which was 2–3 times higher than the maximum τ_current in most regions except the southeast of bird-foot delta, where τ_wave and τ_current were comparable.

The seasonal mean sedimentation rates of fluvial materials and riverine suspended sediment flux (SSF) based on the 20-year simulation results are shown in Figure 12. High sedimentation rate (>1 cm/season) and strong SSF (>0.1 kg/m/s) were simulated near the Mississippi and Atchafalaya River mouths in spring (Figure 12a). Westward alongshore sediment transport in spring, fall, and winter dominated the Louisiana–Texas shelf (Figure 12a,c,d). In summer, deposition was patchy and the intensity of westward SSF was largely reduced (Figure 12b). Offshore sediment transport in the Atchafalaya Bay was stronger in winter and spring due to high fluvial discharge and strong resuspension. Over the Louisiana–Mississippi–Alabama shelf, eastward sediment transport was dominant in spring, summer, and winter (Figure 12a,b,d). Deposition minimized in fall because of the westward SSF and low fluvial inputs from the bird-foot delta (Figure 12c).

4.2. Interannual Variation of Sedimentation Rate

Spatially-averaged sedimentation rates over the shelf (h < 200 m) and inner shelf (h < 50 m) fluctuated between 1 and 8 mm/year during our simulation period with a similar temporal variation pattern (Figure 13a). Sedimentation rate on the inner shelf was ~0.2–2 mm/year higher than that over the entire shelf due to its proximity to the river mouths. To explore the influence of fluvial discharge to shelf deposition, we plotted sediment flux of the Mississippi and Atchafalaya Rivers (Figure 13b; Data source: http://nwis.waterdata.usgs.gov; station St. Francisville (07373420) for the Mississippi River and station Melville (07381495) for the Atchafalaya River), and found sedimentation rates over the shelf and inner shelf to be highly correlated with fluvial sediment discharge (correlation coefficient: 0.80). The dramatic decrease of both sedimentation rate and fluvial sediment discharge in 1998–2000 can be explained by the shift of ENSO phase from a strong El Niño episode (1997/98) to a strong La Niña episode (1999/2000) [96]. It is noteworthy that we did not include pre-deposited sediment resuspension and other sediment sources (e.g., coastal erosion), which could introduce more uncertainties. To achieve more accurate results, these processes should be taken in account in future studies. Compared with long-term sediment deposition on the shelf, the significance of short-term, event-driven sediment transport and dispersal have been highlighted in the last several decades. During the passage of cold fronts and hurricanes, the magnitude of shelf deposition can reach several cm, which is an order of magnitude higher than annual sediment deposition under calmer hydrodynamic condition [37,40,41]. Therefore, quantitative estimations of intensity and frequency of event-driven sediment dispersal is of importance in studying long-term sediment transport.

4.3. Spatial Pattern of Deposition

We averaged model simulated sedimentation rate from 1993 to 2012 based on the changes in the thickness of sediment layers. As shown in Figure 14, sedimentation rate in the nGoM varied greatly from more than 10 cm/year to almost zero. The highest sedimentation rate was simulated just off the mouths of the Mississippi and Atchafalaya Rivers. High sedimentation rate (>5 cm/year) around the Atchafalaya River estuary and its decreasing trend in an offshore direction suggested that most sediments debouching into the nGoM through the Atchafalaya River were retained in the bay. For the Mississippi River; however, high sedimentation rate distributed on both sides of the bird-foot delta due to bidirectional (eastward and westward) fluvial sediments dispersal with the shift of alongshore currents direction in different seasons. Over the western Louisiana shelf, fluvial sediments were transported westward crossing 93° W and deposited over the shelf. Sedimentation rate over the Louisiana–Mississippi–Alabama shelf was less than 1 cm/year. Over the 20-year simulation period, very little fluvial sediment was deposited in waters deeper than 500 m, indicating limited cross-shelf suspended sediment transport.

4.4. Sensitivity Tests for High and Low Fluvial Discharge

The Mississippi-Atchafalaya River system was the dominant sediment source in the nGoM and both water and suspended sediment fluxes exhibited strong seasonality: most high fluxes started from late winter until the end of spring (Figure 15a,b and Figure 16a,b). The peak flow appeared from February to April, ranging between 60 and 100 km³/month for the Mississippi River. Water flux of the Atchafalaya River was about 55% lower than that of the Mississippi River. In summer and fall, water fluxes of both rivers were only 20%–30% of their maxima. SSF showed similar temporal pattern as that of water. The highest monthly SSF (61.7 Mt/month) of the Mississippi River was in April 1995 due to the concurrence of peak streamflow (74.1 km³/month) and SSC (832 mg/L). The SSF peak of the Atchafalaya River was lower than 20 Mt/month except in 1998 and 1999. Here we applied the non-parametric change-point Pettitt test [97] to monthly water and suspended sediment fluxes to detect the presence of any points of change over the modeling period (1993–2012). As a statistical test used to detect the characteristics of changes, the non-parametric Pettitt test has been widely used in previous hydroclimatic studies [98,99,100]. The non-parametric statistic $U_{t,T}$ is defined as:

$$U_{t,T}={\displaystyle \sum }_{i=1}^t{\displaystyle \sum }_{j=t+1}^{T}sgn\left(X_i-X_j\right)$$

(1)

where $T$ is the length of the time series, $t$ is the time of the shift, $X$ is monthly water flux or SSF, and

$$sgn\left(X_i-X_j\right)=\{\begin{array}{c} 1 \left(X_i-X_j> 0\right) \\ 0 \left(X_i-X_j= 0\right)\\ -1 \left(X_i-X_j< 0\right)\end{array}$$

(2)

If a change point exists, the value of $\left|U_{t,T}\right|$ increases with $t$ to its maximum and then decrease and the most significant change point is established at the time when $\left|U_{t,T}\right|$ is equal to $K_T=\max (|U_{t,T}|)$. The significance probability $p$ of $K_T$ is estimated with

$$p=2\exp \left(\frac{-6K_T^2}{T^3+T^2}\right)$$

(3)

If $p<0.05$, a significant change point is confirmed and time series show different features before and after the change point.

As shown in Figure 15c,d and Figure 16c,d, the most significant change points of water and SSF occurred in 1999. This remarkable change can be ascribed to the phase shift of ENSO from the strong El Niño episode of 1997/98 to the strong La Niña episode of 1999/2000 [96]. The strong decreasing variability of other streamflow-related climate indices (e.g., NAO, PDO) after 1999 were likely related to this climatic/hydrologic regime shift [20]. We, therefore, divided our study period (1993–2012) into two time spans, before (1993–1998) and after (1999–2012) the change point. For the Mississippi River, the annual mean water flux declined from 528.8 km³/year in the first span to 433.5 km³/year in the second. The annual mean SSF almost halved from 122.0 Mt/year (1993–1998) to 69.8 Mt/year (1999–2012). For the Atchafalaya River, the annual mean water flux went down from 248.3 km³/year to 191.9 km³/year and the SSF decreased by 34% (from 69.0 to 45.1 Mt/year). To unravel the difference of river discharge after the change point, we compared the multi-year monthly mean water and sediment fluxes in these two spans (Figure 15e,f and Figure 16e,f). In general, monthly water flux in the first span was higher than that in the second span for both rivers. The major difference for the monthly mean Mississippi River discharge between the two spans was found between January and May, when water and sediment fluxes were high (Figure 15e,f). Unlike the Mississippi River, the decreases of the Atchafalaya River monthly mean water and sediment fluxes from the first to the second span were relatively constant in each month (Figure 16e,f).

To assess the impact of fluvial discharge changes on sediment dispersal over the shelf, we conducted two 20-year sensitivity tests using the 1993–1998 and 1999–2012 monthly mean SSF to represent the high and low river discharge scenarios, respectively. Since few studies quantitatively estimate the Mississippi River channel evolution and its contribution to fluvial sediment flux (proximal sediment supply), fluvial sediment supply variation due to river bed scour and deposition was not considered in this study [32]. Over the 20-year simulation period, the sedimentation rate over the entire shelf turned to be lower in the low fluvial discharge scenario than that under high fluvial discharge scenario (Figure 17). Substantial reduction in sedimentation rate was simulated around the bird-foot delta and in the Atchafalaya Bay. Sedimentation rate difference between the two tests was minimum in waters >200 m deep (Figure 17c), suggesting the impact from reduced river inputs might limit to the shelf water.

5. Discussion

Triggered by both anthropogenic activities and natural forces, Louisiana’s coast has been experiencing severe land loss over the last several decades although efforts have been dedicated to land building through marsh creation, sediment diversion, barrier island restoration and shoreline protection [101]. In this section we use the results from our 20-year simulation and scenario tests to assess the impact of a decreased fluvial inputs on two areas of interests: 1) the exchange of Louisiana Bight with Barataria Bay where serious land loss is undergoing, and the 2) the distal transgressive sandy shoals on the Louisiana Shelf, which provide the sandy dredging materials for coastal restoration purpose.

5.1. Bay-Shelf Exchange

Our sensitivity tests indicate that the sedimentation in nGoM can be greatly affected by the changes in fluvial discharge, especially in areas adjacent to the river mouths and deltas. These results shed light on possible projects from the ongoing Louisiana Coastal Master Plan, which will divert a large amount of water and sediments away from the main channel to coastal bays where new land is expected to be built [101].

Here we focused on the circulation and difference of sedimentation rate between the low and high discharge scenarios in the Louisiana Bight and Barataria Bay region (Figure 18), where 1177 ± 106 km² of land was lost in the period of 1932–2016 [54]. Our model identified two transport pathways of the Mississippi-derived sediments: 1) A direct northwestward alongshore transport from the Southwest Pass; and 2) a gyre-induced clockwise transport, which joins the alongshore transport near Sandy Point (Figure 18). While these two pathways have been previously reported [30,102,103], our sensitivity tests, for the first time, indicated that sharp decrease of sedimentation rate due to the decline of fluvial sediment discharge was expected in waters around the bird-foot delta (>10 cm/year) and within the clockwise gyre (up to 1 cm/year), where sedimentation rate was higher than that over the entire shelf (Figure 17). Since our model used fluvial SSC and discharge measurements from USGS river gages as river input, sediment in-channel storage was treated unchanged although it varies with proximal/distal sediment supply, flow regimes and sediment particle grain size [32]. To better quantify fluvial sediment deposition and sediment flux around the bird-foot delta, the effects of sediment dispersal in the river channel should be taken into account in future study.

It has been estimated that the net sediment transport through tidal inlets between the Barataria Bay and the Louisiana Bight were seaward, and the SSF was ~8800 ton/day with 85% of the flow variability in the pass resulting from tides [30,104]. Fitzgerald et al. [103] related the growth of ebb-tidal deltas outside the Barataria Bay to eroded inlet and alongshore transport. Although few studies quantitatively investigated the contributions from different sediment sources (e.g., coastal erosion; resuspension and fluvial discharge) to the depositions close to the bay mouth, large sediment inputs, associated with low salinity and intensified stratification, was observed through tidal inlets, which suggests that sediments from the Mississippi River can be transported into the Barataria Bay during flood tides [105]. Due to the relatively coarse spatial resolution of our model (1 km) and the lack of information to prescribe sediment inputs from coastal erosion, this study still cannot quantify the importance of the Mississippi fluvial sediments to Barataria Bay’s sediment budget. However, if sediment discharge from the Mississippi River keeps decreasing in the future, we expect less sediment to be transported to the bay via tidal inlets.

5.2. Sediment Dynamics over Submarine Shoals

5.2.1. River Supply

By the end of our 20-year simulation, the fractions of river-derived sediments in the surficial seabed layer for Tiger Shoal, Trinity Shoal, and Ship Shoal were 17.6%, 7.1%, and 10.0%, respectively (

Table 2. Percentage of fluvial sediments over each shoal.

Shoals	Benchmark	High River Scenario	Low River Scenario
Tiger Shoal	17.6%	17.1%	16.5%
Ship Shoal	10.0%	12.1%	10.1%
Trinity Shoal	7.1%	7.6%	6.5%

). The values were comparable with the estimations of previous geotechnical investigations [67,106,107]. Moreover, the variation of shoal-wide fluvial sediment fraction under high and low river discharge scenarios was less than 2.1%, indicating the impact from the changes in riverine inputs was not significant over shoals (Table 2). The low percentages of modern fluvial sediments were due to the limited supply of fluvial sediments and resuspension induced by strong hydrodynamics. Previous studies revealed that modern fluvial sediments transport to the sandy shoals through bedload transport was trivial because most materials over the three shoals were relict coarse sediments [59,107], so we did not incorporate bedload in our simulation. Under calm weather conditions, sediment plume of the Atchafalaya River was mainly confined within the Atchafalaya Bay, and only a small amount of suspended sediments could be transported over the shelf [27,108]. However, during episodic events (i.e., cold fronts and hurricanes), previously deposited riverine sediments were resuspended and transported offshore, suggesting a direct yet intermittent supply of fine riverine sediments to the shoals [59,109]. Besides, wave-supported fluid mud movement is another important mechanism in terms of fluvial sediment across-shelf transport over the muddy Atchafalaya Shelf [48,52,110]. Since most fluid mud observations and modeling studies only focus on short-term period (several days to weeks), the importance of fluid mud transport to sandy shoals on decadal scales is still unclear.

5.2.2. Hydrodynamics

Given strong bottom shear stress induced by shallow water depth, sediment remobilization over sandy shoals is potentially intensive. To investigate the temporal variation of hydrodynamics related to sediment resuspension over these transgressive shoals, we calculated spatially-averaged, monthly-mean bottom shear stress induced by currents and waves (τ_cw) over each shoal. The highest critical shear stress (τ_c) of fluvial sediments (0.2 Pa; see Table 1) was treated as the threshold of “strong resuspension” and the number of days with τ_cw > τ_c was counted to represent the duration of strong resuspension. As shown in Figure 19a–d, both τ_cw and the number of days with excessive bottom shear stress maximized in cold season. About 80% of the days with excessive bottom shear stress was found between October and April (82.2% for Tiger Shoal, 80.2% for Ship Shoal, and 78.8% for Trinity Shoal, respectively). Such unevenly temporal distribution indicates that most resuspension over the shoals happens during cold season when hydrodynamic are stronger. To quantitatively estimate the inter-annual variation of bottom resuspension over the three shoals, we plotted the annual mean bottom shear stress (τ_m), wind speed (data source: CFSR) and the number of days with excessive bottom shear stress in one year over each shoal (right panel of Figure 19). Both τ_m and the number of days with strong resuspension (τ > τ_c) peaked in 1998 and 2008, and the variation of hydrodynamics was highly correlated with wind speed (Figure 19e–h). Although previous investigations found the inter-annual variations of strong meteorological and hydrodynamic conditions can be ascribed to stratospheric ozone depletion, available latent heat, expansion of Hadley cell and large-scale circulation pattern shift [111,112,113,114,115,116,117,118], the balance and interactions between these factors are still less understood and it is still a challenge to directly link their influence with regional hydro- and sediment dynamics. In general, sediment dynamics over the transgressive shoals is mainly impacted by wind-induced hydrodynamics rather than fluvial inputs.

5.3. Limitations and Future Work

Our 20-year simulation reproduced the overall pattern of the transport and dispersal of river-derived sediments in the nGoM. However, it is noteworthy that some important sediment transport processes and mechanisms were not included in our model. First of all, this study only focused on the dynamics of fluvial sediments, thus the resuspension of shelf sediments was not considered. Such simplification will underestimate the SSF over the shallow shelf, where resuspension of shelf sediments can be an important source in addition to river inputs. During intensive events (hurricanes or cold fronts), sediment resuspension could be an order of magnitude higher than fluvial discharge [39,119,120]. Secondly, coastal erosion was not included in our model. The Mississippi River Delta has been experiencing severe land loss over the past decades. Combination of natural processes (e.g., storms, subsidence, and salt water intrusion) and human activities (e.g., artificial channel, oil industry, and urbanization) accelerates the erosion process and a large amount of eroded sediments can be transported to the coastal water [121,122,123]. Thirdly, our model only simulated suspended load. Although model results indicated that deposition of river sediments mainly occurred in waters that are <200 m deep, one should not rule out the possible cross-shelf transport induced by bedload that were not included in our model. For instance, Corbett et al. [44] pointed out that cross-isobath sediment supply to the shelf break could be attributed to subaqueous slides and slumps, where sedimentation rates could be higher than that around the bird-foot delta. Ross et al. [124] found sediment flux to the Mississippi Canyon is more related to the sediment availability rather than current speeds and hurricanes can greatly increase the sediment transport through canyons. Besides, wave-supported gravity flow (i.e., fluid mud), as an important mechanism of the transport of fine sediments, has been reported on the Atchafalaya Shelf [48,50,125]. Although the velocity of fluid mud transport is slower than suspended load, its high concentration (>10 g/L) can substantially increase the sediment flux and change the deposition pattern near river estuaries [126]. In addition, fluid mud can incur bottom turbulence dissipation, which can be a dominant feature over the shelf off the Atchafalaya Bay [52,110]. Bedload transport of non-cohesive sediment and its interaction with hydrodynamics is also important to the formation and geomorphological changes of an erodible bed [127]. Last but not least, baroclinic estuarine circulation, tidal pumping effects and sediment storage in the river channel cannot be resolved in our model since fluvial discharge was treated as point source in ROMS. Such simplification in estuarine dynamics can affect the sediment flux estimation from river estuaries to the shelf [32,128]. Further investigation that accounts for the above-mentioned processes is expected for a more comprehensive analysis of regional sediment dynamics.

6. Conclusions

We adapted the coupled ocean-sediment transport model to the northern Gulf of Mexico to investigate sediment dynamics on seasonal to decadal time scales. Extensive model-data comparisons were carried out to evaluate model performance. Our 20-year model simulation reveals that:

(1)

Strong easterly winds prevailed in non-summer seasons. Relatively weak westerly winds in summer reversed currents between 20 and 50 m isobaths to an eastward direction. Wave- and current-induced bottom shear stresses exhibited similar temporal (strong in winter and weak in summer) and spatial (higher over the inner shelf) patterns. High sedimentation rate (>1 cm/season) and SSF (>0.1 kg/m/s) were found in spring near river mouths. During summer, calm hydrodynamics and reversed coastal currents resulted in weak eastward SSF over the Louisiana–Texas shelf. Deposition on the Louisiana–Mississippi–Alabama shelf became negligible in fall;

(2)

Over the 20-year simulation, sedimentation rate ranged from almost zero to more than 10 cm/year in waters near the river mouth and surrounding the delta. Interannual variation of sedimentation rates over the shelf (h < 200 m) and inner shelf (h < 50 m) were highly correlated with the fluvial sediment flux. Mississippi-derived sediments dispersed on both sides of the bird-foot delta, while the Atchafalaya-derived sediments were mainly confined in the Atchafalaya Bay. Two major pathways for the Mississippi River-derived sediment were identified: A direct westward alongshore transport from the Southwest Pass, and a gyre-induced clockwise transport centered in Louisiana Bight;

(3)

A change point was detected in 1999 in the time series of water and sediment discharge from the Mississippi-Atchafalaya River over the period of 1993–2012. This change point was correlated with the shift of ENSO from a strong warm phase to a strong cold phase. The annual mean water and sediment fluxes decreased sharply from the 1993–1998 period to the 1999–2012 period. Model sensitivity tests indicated that the influence of decreased river inputs on sedimentation rate was limited to waters near the river mouths, which reduced sediment transport into the Barataria Bay during flood tide and potentially worsen the ongoing land loss in the bay;

(4)

Model simulated percentages of fluvial sediments over the Tiger, Trinity, and Ship Shoals were less than 18%, indicating the variation of river sediment flux might have limited impact on local sedimentation. Sediment dynamics over these distal sandy bodies were mostly affected by the strong winds in cold season between October and April.